QPMS
/
qpms
1
0
Fork 0
Code Issues Pull Requests Projects Releases Wiki Activity

Correct plane wave decomposition at origin (wrong elsewhere for longitudinal waves 

Former-commit-id: 4dc70fde23e31828fb959da429d716bd39e84e47
This commit is contained in:
Marek Nečada 2018-02-04 15:05:41 +00:00
parent dbc118e712
commit 4dc151d5b9
1 changed files with 13 additions and 3 deletions
Show Stats Download Patch File Download Diff File Expand all files Collapse all files

16
qpms/vswf.c
View File

@ -230,13 +230,23 @@ qpms_errno_t qpms_vswf_fill(csphvec_t *const longtarget, csphvec_t * const mgtar
double const *ptau = pt.tau; double const *ptau = pt.tau;
csphvec_t *plong = longtarget, *pmg = mgtarget, *pel = eltarget; csphvec_t *plong = longtarget, *pmg = mgtarget, *pel = eltarget;
for(qpms_l_t l = 1; l <= lMax; ++l) { for(qpms_l_t l = 1; l <= lMax; ++l) {
complex double besfac = *pbes / kr.r; complex double besfac;
complex double besderfac = *(pbes-1) - l * besfac; complex double besderfac;
if (kr.r) {
besfac = *pbes / kr.r;
} else {
besfac = (1 == l) ? 1/3. : 0;
}
besderfac = *(pbes-1) - l * besfac;
for(qpms_m_t m = -l; m <= l; ++m) { for(qpms_m_t m = -l; m <= l; ++m) {
complex double eimf = cexp(m * kr.phi * I); complex double eimf = cexp(m * kr.phi * I);
if (longtarget) { if (longtarget) {
complex double longfac = sqrt(l*(l+1)) * eimf; complex double longfac = sqrt(l*(l+1)) * eimf;
plong->rc = (besderfac-besfac) * (*pleg) * longfac; plong->rc = // FATAL FIXME: I get wrong result here for plane wave re-expansion
// whenever kr.r > 0 (for waves with longitudinal component, ofcoz)
/*(*(pbes-1) - (l+1)/kr.r* *pbes)*/
(besderfac-besfac)
* (*pleg) * longfac;
plong->thetac = *ptau * besfac * longfac; plong->thetac = *ptau * besfac * longfac;
plong->phic = *ppi * I * besfac * longfac; plong->phic = *ppi * I * besfac * longfac;
++plong; ++plong;