From c445b83593bab4b6317ade014e260c9f32513f51 Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?Marek=20Ne=C4=8Dada?= Date: Thu, 16 Jan 2020 07:52:50 +0200 Subject: [PATCH] Implement previously forgotten functions. Former-commit-id: df5215ad0349eb75bd2c7bee60f1fac50d23eb99 --- qpms/scatsystem.c | 53 ++++++++++++++++++++++++++++++++++++++++++++++- qpms/tmatrices.c | 18 ++++++++++++++++ qpms/tmatrices.h | 11 ++++++++++ 3 files changed, 81 insertions(+), 1 deletion(-) diff --git a/qpms/scatsystem.c b/qpms/scatsystem.c index ef29643..1af486f 100644 --- a/qpms/scatsystem.c +++ b/qpms/scatsystem.c @@ -476,7 +476,58 @@ void qpms_scatsys_free(qpms_scatsys_t *ss) { free(ss); } +void qpms_scatsys_at_omega_refill(qpms_scatsys_at_omega_t *ssw, + complex double omega) { + const qpms_scatsys_t * const ss = ssw->ss; + ssw->omega = omega; + ssw->medium = qpms_epsmu_generator_eval(ss->medium, omega); + ssw->wavenumber = qpms_wavenumber(omega, ssw->medium); + qpms_tmatrix_t **tmatrices_preop; + QPMS_CRASHING_CALLOC(tmatrices_preop, ss->tmg_count, sizeof(*tmatrices_preop)); + for (qpms_ss_tmgi_t tmgi = 0; tmgi < ss->tmg_count; ++tmgi) + tmatrices_preop[tmgi] = qpms_tmatrix_init_from_function(ss->tmg[tmgi], omega); + for (qpms_ss_tmi_t tmi = 0; tmi < ss->tm_count; ++tmi) + qpms_tmatrix_apply_operation_replace(ssw->tm[tmi], &ss->tm[tmi].op, + tmatrices_preop[ss->tm[tmi].tmgi]); + for (qpms_ss_tmgi_t tmgi = 0; tmgi < ss->tmg_count; ++tmgi) + qpms_tmatrix_free(tmatrices_preop[tmgi]); + free(tmatrices_preop); +} +qpms_scatsys_at_omega_t *qpms_scatsys_at_omega(const qpms_scatsys_t *ss, + complex double omega) { + // TODO + qpms_scatsys_at_omega_t *ssw; + QPMS_CRASHING_MALLOC(ssw, sizeof(*ssw)); + ssw->omega = omega; + ssw->ss = ss; + ssw->medium = qpms_epsmu_generator_eval(ss->medium, omega); + ssw->wavenumber = qpms_wavenumber(omega, ssw->medium); + qpms_tmatrix_t **tmatrices_preop; + QPMS_CRASHING_CALLOC(tmatrices_preop, ss->tmg_count, sizeof(*tmatrices_preop)); + for (qpms_ss_tmgi_t tmgi = 0; tmgi < ss->tmg_count; ++tmgi) + tmatrices_preop[tmgi] = qpms_tmatrix_init_from_function(ss->tmg[tmgi], omega); + for (qpms_ss_tmi_t tmi = 0; tmi < ss->tm_count; ++tmi) { + ssw->tm[tmi] = qpms_tmatrix_apply_operation(&ss->tm[tmi].op, + tmatrices_preop[ss->tm[tmi].tmgi]); //<- main difference to .._refill() + QPMS_ENSURE(ssw->tm[tmi], + "Got NULL pointer from qpms_tmatrix_apply_operation"); + } + for (qpms_ss_tmgi_t tmgi = 0; tmgi < ss->tmg_count; ++tmgi) + qpms_tmatrix_free(tmatrices_preop[tmgi]); + free(tmatrices_preop); + return ssw; +} + +void qpms_scatsys_at_omega_free(qpms_scatsys_at_omega_t *ssw) { + if (ssw) { + if(ssw->tm) + for(qpms_ss_tmi_t i = 0; i < ssw->ss->tm_count; ++i) + qpms_tmatrix_free(ssw->tm[i]); + free(ssw->tm); + } + free(ssw); +} // (copypasta from symmetries.c) // TODO at some point, maybe support also other norms. @@ -1050,7 +1101,7 @@ complex double *qpms_scatsys_build_translation_matrix_e_full( -complex double *qpms_scatsys_at_omega_build_modeproblem_matrix_full( +complex double *qpms_scatsysw_build_modeproblem_matrix_full( /// Target memory with capacity for ss->fecv_size**2 elements. If NULL, new will be allocated. complex double *target, const qpms_scatsys_at_omega_t *ssw diff --git a/qpms/tmatrices.c b/qpms/tmatrices.c index 96516ba..c60f146 100644 --- a/qpms/tmatrices.c +++ b/qpms/tmatrices.c @@ -1140,7 +1140,25 @@ void qpms_tmatrix_operation_compose_chain_init(qpms_tmatrix_operation_t *dest, o->opmem_size = opmem_size; } +qpms_tmatrix_t *qpms_tmatrix_mv(qpms_tmatrix_t *dest, + const qpms_tmatrix_t *orig) { + QPMS_ENSURE(qpms_vswf_set_spec_isidentical(dest->spec, orig->spec), + "Basis specifications must be identical!"); + memcpy(dest->m, orig->m, SQ(orig->spec->n)*sizeof(complex double)); + return dest; +} +qpms_tmatrix_t *qpms_tmatrix_apply_operation_replace(qpms_tmatrix_t *dest, + const qpms_tmatrix_operation_t *op, const qpms_tmatrix_t *orig) { + //QPMS_ENSURE(qpms_vswf_set_spec_isidentical(dest->spec, orig->spec), + // "Basis specifications must be identical!"); + // TODO should I check also for dest->owns_m? + // FIXME this is highly inoptimal; the hierarchy should be such + // that _operation() and operation_inplace() call this, and not the + // other way around + qpms_tmatrix_mv(dest, orig); + return qpms_tmatrix_apply_operation_inplace(op, dest); +} qpms_tmatrix_t *qpms_tmatrix_apply_operation( const qpms_tmatrix_operation_t *f, const qpms_tmatrix_t *orig) { diff --git a/qpms/tmatrices.h b/qpms/tmatrices.h index 8b4d2e6..117e3f9 100644 --- a/qpms/tmatrices.h +++ b/qpms/tmatrices.h @@ -27,6 +27,13 @@ qpms_tmatrix_t *qpms_tmatrix_init(const qpms_vswf_set_spec_t *bspec); /// Copies a T-matrix, allocating new array for the T-matrix data. qpms_tmatrix_t *qpms_tmatrix_copy(const qpms_tmatrix_t *T); +/// Copies a T-matrix contents between two pre-allocated T-matrices. +/** orig->spec and dest->spec must be identical. + * + * \returns \a dest + */ +qpms_tmatrix_t *qpms_tmatrix_mv(qpms_tmatrix_t *dest, const qpms_tmatrix_t *orig); + /// Destroys a T-matrix. void qpms_tmatrix_free(qpms_tmatrix_t *t); @@ -623,6 +630,10 @@ qpms_tmatrix_t *qpms_tmatrix_apply_operation(const qpms_tmatrix_operation_t *op, /// Apply an operation on a T-matrix and replace it with the result. qpms_tmatrix_t *qpms_tmatrix_apply_operation_inplace(const qpms_tmatrix_operation_t *op, qpms_tmatrix_t *orig); +/// Apply an operation on a T-matrix and replace another one it with the result. +qpms_tmatrix_t *qpms_tmatrix_apply_operation_replace(qpms_tmatrix_t *dest, + const qpms_tmatrix_operation_t *op, const qpms_tmatrix_t *orig); + /// (Recursively) deallocates internal data of qpms_tmatrix_operation_t. /** It does NOT clear the memory pointed to by op it self, but only * heap-allocated data of certain operations, if labeled as owned by it.