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

Fix offset incrementation on the diagonal. 

Solves the translation matrix projection discrepancy.


Former-commit-id: 7bb7dbba47aa5ac3df32177098c096a86d2879fd
This commit is contained in:
Marek Nečada 2019-03-12 15:43:28 +02:00
parent b97f3e6e2e
commit 8efd582daf
1 changed files with 15 additions and 13 deletions
Show Stats Download Patch File Download Diff File Expand all files Collapse all files

28
qpms/scatsystem.c
View File

@ -1120,13 +1120,14 @@ complex double *qpms_scatsys_build_translation_matrix_full(
const cart3_t posR = ss->p[piR].pos; const cart3_t posR = ss->p[piR].pos;
size_t fullvec_offsetC = 0; size_t fullvec_offsetC = 0;
for(qpms_ss_pi_t piC = 0; piC < ss->p_count; ++piC) { for(qpms_ss_pi_t piC = 0; piC < ss->p_count; ++piC) {
if(piC == piR) continue; // The diagonal will be dealt with later.
const qpms_vswf_set_spec_t *bspecC = ss->tm[ss->p[piC].tmatrix_id]->spec; const qpms_vswf_set_spec_t *bspecC = ss->tm[ss->p[piC].tmatrix_id]->spec;
const cart3_t posC = ss->p[piC].pos; if(piC != piR) { // The diagonal will be dealt with later.
QPMS_ENSURE_SUCCESS(qpms_trans_calculator_get_trans_array_lc3p(ss->c, const cart3_t posC = ss->p[piC].pos;
target + fullvec_offsetR*full_len + fullvec_offsetC, QPMS_ENSURE_SUCCESS(qpms_trans_calculator_get_trans_array_lc3p(ss->c,
bspecR, full_len, bspecC, 1, target + fullvec_offsetR*full_len + fullvec_offsetC,
k, posR, posC)); bspecR, full_len, bspecC, 1,
k, posR, posC));
}
fullvec_offsetC += bspecC->n; fullvec_offsetC += bspecC->n;
} }
assert(fullvec_offsetC = full_len); assert(fullvec_offsetC = full_len);
@ -1162,19 +1163,20 @@ complex double *qpms_scatsys_build_modeproblem_matrix_full(
// dest particle T-matrix // dest particle T-matrix
const complex double *tmmR = ss->tm[ss->p[piR].tmatrix_id]->m; const complex double *tmmR = ss->tm[ss->p[piR].tmatrix_id]->m;
for(qpms_ss_pi_t piC = 0; piC < ss->p_count; ++piC) { for(qpms_ss_pi_t piC = 0; piC < ss->p_count; ++piC) {
if(piC == piR) continue; // The diagonal will be dealt with later.
const qpms_vswf_set_spec_t *bspecC = ss->tm[ss->p[piC].tmatrix_id]->spec; const qpms_vswf_set_spec_t *bspecC = ss->tm[ss->p[piC].tmatrix_id]->spec;
const cart3_t posC = ss->p[piC].pos; if(piC != piR) { // The diagonal will be dealt with later.
QPMS_ENSURE_SUCCESS(qpms_trans_calculator_get_trans_array_lc3p(ss->c, const cart3_t posC = ss->p[piC].pos;
tmp, // tmp is S(piR<-piC) QPMS_ENSURE_SUCCESS(qpms_trans_calculator_get_trans_array_lc3p(ss->c,
bspecR, bspecC->n, bspecC, 1, tmp, // tmp is S(piR<-piC)
k, posR, posC)); bspecR, bspecC->n, bspecC, 1,
cblas_zgemm(CblasRowMajor, CblasNoTrans, CblasNoTrans, k, posR, posC));
cblas_zgemm(CblasRowMajor, CblasNoTrans, CblasNoTrans,
bspecR->n /*m*/, bspecC->n /*n*/, bspecR->n /*k*/, bspecR->n /*m*/, bspecC->n /*n*/, bspecR->n /*k*/,
&one/*alpha*/, tmmR/*a*/, bspecR->n/*lda*/, &one/*alpha*/, tmmR/*a*/, bspecR->n/*lda*/,
tmp/*b*/, bspecC->n/*ldb*/, &zero/*beta*/, tmp/*b*/, bspecC->n/*ldb*/, &zero/*beta*/,
target + fullvec_offsetR*full_len + fullvec_offsetC /*c*/, target + fullvec_offsetR*full_len + fullvec_offsetC /*c*/,
full_len /*ldc*/); full_len /*ldc*/);
}
fullvec_offsetC += bspecC->n; fullvec_offsetC += bspecC->n;
} }
fullvec_offsetR += bspecR->n; fullvec_offsetR += bspecR->n;