QPMS
Electromagnetic multiple scattering library and toolkit.
Data Structures | Macros | Typedefs | Functions
scatsystem.h File Reference

Modern interface for finite lattice calculations, including symmetries. More...

#include "qpms_types.h"
#include "vswf.h"
#include "tmatrices.h"
#include <stdbool.h>
Include dependency graph for scatsystem.h:

Go to the source code of this file.

Data Structures

struct  qpms_particle_t
 A particle, defined by its T-matrix and position. More...
 
struct  qpms_particle_tid_t
 A particle, defined by its T-matrix INDEX and position, to be used in qpms_scatsys_t. More...
 
struct  qpms_ss_orbit_type_t
 Structure describing a particle's "orbit type" under symmetry group actions in a system. More...
 
struct  qpms_ss_particle_orbitinfo
 Auxillary type used in qpms_scatsys_t that identifies the particle's orbit and its id inside that orbit. More...
 
struct  qpms_ss_derived_tmatrix_t
 Auxillary type used in qpms_scatsys_t: A recepy to create another T-matrices by symmetry operations. More...
 
struct  qpms_scatsys_periodic_info_t
 
struct  qpms_scatsys_t
 Common "class" for system of scatterers, both periodic and non-periodic. More...
 
struct  qpms_scatsys_at_omega_t
 
struct  qpms_ss_LU
 LU factorisation (LAPACKE_zgetrf) result holder. More...
 
struct  qpms_scatsys_at_omega_k_t
 Scattering system at a given frequency and k-vector. Used only with periodic systems. More...
 

Macros

#define QPMS_SS_P_ORBITINFO_UNDEF   (-1)
 This labels that the particle has not yet been assigned to an orbit.
 

Typedefs

typedef struct qpms_particle_t qpms_particle_t
 A particle, defined by its T-matrix and position. More...
 
typedef struct qpms_finite_group_t qpms_finite_group_t
 
typedef struct qpms_particle_tid_t qpms_particle_tid_t
 A particle, defined by its T-matrix INDEX and position, to be used in qpms_scatsys_t.
 
typedef qpms_gmi_t qpms_ss_orbit_pi_t
 Auxilliary type used in qpms_ss_orbit_type_t for labeling particles inside orbits.
 
typedef qpms_ss_tmi_t qpms_ss_oti_t
 Auxilliary type used for labeling orbit types.
 
typedef struct qpms_ss_orbit_type_t qpms_ss_orbit_type_t
 Structure describing a particle's "orbit type" under symmetry group actions in a system. More...
 
typedef ptrdiff_t qpms_ss_osn_t
 "serial number" of av orbit in a given type.
 
typedef struct qpms_ss_particle_orbitinfo qpms_ss_particle_orbitinfo_t
 Auxillary type used in qpms_scatsys_t that identifies the particle's orbit and its id inside that orbit.
 
typedef struct qpms_ss_derived_tmatrix_t qpms_ss_derived_tmatrix_t
 Auxillary type used in qpms_scatsys_t: A recepy to create another T-matrices by symmetry operations.
 
typedef struct qpms_scatsys_periodic_info_t qpms_scatsys_periodic_info_t
 
typedef struct qpms_scatsys_t qpms_scatsys_t
 Common "class" for system of scatterers, both periodic and non-periodic. More...
 
typedef struct qpms_scatsys_at_omega_t qpms_scatsys_at_omega_t
 
typedef struct qpms_ss_LU qpms_ss_LU
 LU factorisation (LAPACKE_zgetrf) result holder.
 
typedef struct qpms_scatsys_at_omega_k_t qpms_scatsys_at_omega_k_t
 Scattering system at a given frequency and k-vector. Used only with periodic systems. More...
 

Functions

void qpms_scatsystem_set_nthreads (long n)
 Overrides the number of threads spawned by the paralellized functions. More...
 
static const qpms_vswf_set_spec_tqpms_ss_bspec_tmi (const qpms_scatsys_t *ss, qpms_ss_tmi_t tmi)
 Retrieve the bspec of tmi'th element of ss->tm.
 
static const qpms_vswf_set_spec_tqpms_ss_bspec_pi (const qpms_scatsys_t *ss, qpms_ss_pi_t pi)
 Retrieve the bspec of pi'th particle in ss->p.
 
qpms_scatsys_at_omega_tqpms_scatsys_apply_symmetry (const qpms_scatsys_t *orig, const struct qpms_finite_group_t *sym, complex double omega, const struct qpms_tolerance_spec_t *tol)
 Creates a new scatsys by applying a symmetry group onto a "proto-scatsys", copying particles if needed. More...
 
void qpms_scatsys_free (qpms_scatsys_t *s)
 Destroys the result of qpms_scatsys_apply_symmetry or qpms_scatsys_load.
 
void qpms_scatsys_at_omega_free (qpms_scatsys_at_omega_t *ssw)
 Destroys a qpms_scatsys_at_omega_t. More...
 
qpms_scatsys_at_omega_tqpms_scatsys_at_omega (const qpms_scatsys_t *ss, complex double omega)
 Evaluates scattering system T-matrices at a given frequency. More...
 
complex double * qpms_scatsys_irrep_transform_matrix (complex double *target_U, const qpms_scatsys_t *ss, qpms_iri_t iri)
 Creates a "full" transformation matrix U that takes a full vector and projects it onto an symmetry adapted basis. More...
 
complex double * qpms_scatsys_irrep_pack_matrix_stupid (complex double *target_packed, const complex double *orig_full, const qpms_scatsys_t *ss, qpms_iri_t iri)
 Projects a "big" matrix onto an irrep (slow reference implementation). More...
 
complex double * qpms_scatsys_irrep_unpack_matrix_stupid (complex double *target_full, const complex double *orig_packed, const qpms_scatsys_t *ss, qpms_iri_t iri, bool add)
 Transforms a big "packed" matrix into the full basis (slow reference implementation). More...
 
complex double * qpms_scatsys_irrep_pack_matrix (complex double *target_packed, const complex double *orig_full, const qpms_scatsys_t *ss, qpms_iri_t iri)
 Projects a "big" matrix onto an irrep. More...
 
complex double * qpms_scatsys_irrep_unpack_matrix (complex double *target_full, const complex double *orig_packed, const qpms_scatsys_t *ss, qpms_iri_t iri, bool add)
 Transforms a big "packed" matrix into the full basis. More...
 
complex double * qpms_scatsys_irrep_pack_vector (complex double *target_packed, const complex double *orig_full, const qpms_scatsys_t *ss, qpms_iri_t iri)
 Projects a "big" vector onto an irrep. More...
 
complex double * qpms_scatsys_irrep_unpack_vector (complex double *target_full, const complex double *orig_packed, const qpms_scatsys_t *ss, qpms_iri_t iri, bool add)
 Transforms a big "packed" vector into the full basis. More...
 
complex double * qpms_scatsys_build_translation_matrix_full (complex double *target, const qpms_scatsys_t *ss, complex double k)
 Global translation matrix. More...
 
complex double * qpms_scatsysw_build_modeproblem_matrix_full (complex double *target, const qpms_scatsys_at_omega_t *ssw)
 Creates the full \( (I - WS) \) matrix of the periodic scattering system. More...
 
complex double * qpms_scatsys_build_translation_matrix_e_full (complex double *target, const qpms_scatsys_t *ss, complex double k, qpms_bessel_t J)
 As qpms_scatsys_build_translation_full() but with choice of Bessel function type. More...
 
complex double * qpms_scatsys_build_translation_matrix_e_irrep_packed (complex double *target, const qpms_scatsys_t *ss, qpms_iri_t iri, complex double k, qpms_bessel_t J)
 Global translation matrix with selectable Bessel function, projected on an irrep. More...
 
complex double * qpms_scatsysw_build_modeproblem_matrix_irrep_packed (complex double *target, const qpms_scatsys_at_omega_t *ssw, qpms_iri_t iri)
 Creates the mode problem matrix \( (I - TS) \) directly in the irrep-packed form. More...
 
complex double * qpms_scatsysw_build_modeproblem_matrix_irrep_packed_orbitorderR (complex double *target, const qpms_scatsys_at_omega_t *ssw, qpms_iri_t iri)
 Alternative implementation of qpms_scatsysw_build_modeproblem_matrix_irrep_packed(). More...
 
complex double * qpms_scatsysw_build_modeproblem_matrix_irrep_packed_serial (complex double *target, const qpms_scatsys_at_omega_t *ssw, qpms_iri_t iri)
 Alternative (serial reference) implementation of qpms_scatsysw_build_modeproblem_matrix_irrep_packed(). More...
 
void qpms_ss_LU_free (qpms_ss_LU)
 
qpms_ss_LU qpms_scatsysw_build_modeproblem_matrix_full_LU (complex double *target, int *target_piv, const qpms_scatsys_at_omega_t *ssw)
 Builds an LU-factorised mode/scattering problem \( (I - TS) \) matrix from scratch. Nonperiodic systems only. More...
 
qpms_ss_LU qpms_scatsysw_build_modeproblem_matrix_irrep_packed_LU (complex double *target, int *target_piv, const qpms_scatsys_at_omega_t *ssw, qpms_iri_t iri)
 Builds an irrep-packed LU-factorised mode/scattering problem matrix from scratch. More...
 
qpms_ss_LU qpms_scatsysw_modeproblem_matrix_full_factorise (complex double *modeproblem_matrix_full, int *target_piv, const qpms_scatsys_at_omega_t *ssw, const struct qpms_scatsys_at_omega_k_t *sswk)
 Computes LU factorisation of a pre-calculated mode/scattering problem matrix, replacing its contents. More...
 
qpms_ss_LU qpms_scatsysw_modeproblem_matrix_irrep_packed_factorise (complex double *modeproblem_matrix_irrep_packed, int *target_piv, const qpms_scatsys_at_omega_t *ssw, qpms_iri_t iri)
 Computes LU factorisation of a pre-calculated irrep-packed mode/scattering problem matrix, replacing its contents. More...
 
complex double * qpms_scatsys_scatter_solve (complex double *target_f, const complex double *a_inc, qpms_ss_LU ludata)
 Solves a (possibly partial, irrep-packed) scattering problem \( (I-TS)f = Ta_\mathrm{inc} \) using a pre-factorised \( (I-TS) \). More...
 
complex double * qpms_scatsyswk_build_modeproblem_matrix_full (complex double *target, const qpms_scatsys_at_omega_k_t *sswk)
 Creates the full \( (I - WS) \) matrix of the periodic scattering system. More...
 
complex double * qpms_scatsys_periodic_build_translation_matrix_full (complex double *target, const qpms_scatsys_t *ss, complex double wavenumber, const cart3_t *wavevector, double eta)
 Global translation matrix. More...
 
complex double * qpms_scatsyswk_build_translation_matrix_full (complex double *target, const qpms_scatsys_at_omega_k_t *sswk)
 Global translation matrix. More...
 
qpms_ss_LU qpms_scatsyswk_build_modeproblem_matrix_full_LU (complex double *target, int *target_piv, const qpms_scatsys_at_omega_k_t *sswk)
 Builds an LU-factorised mode/scattering problem \( (I - TS) \) matrix from scratch. Periodic systems only. More...
 
struct beyn_result_tqpms_scatsys_periodic_find_eigenmodes (const qpms_scatsys_t *ss, const double k[3], complex double omega_centre, double omega_rr, double omega_ri, size_t contour_npoints, double rank_tol, size_t rank_sel_min, double res_tol)
 Searches for periodic scattering system's eigenmodes using Beyn's algorithm. More...
 
qpms_errno_t qpms_scatsys_dump (qpms_scatsys_t *ss, char *path)
 NOT IMPLEMENTED Dumps a qpms_scatsys_t structure to a file.
 
qpms_scatsys_tqpms_scatsys_load (char *path)
 NOT IMPLEMENTED Reads a qpms_scatsys_t structure from a file.
 
complex double * qpms_orbit_action_matrix (complex double *target, const qpms_ss_orbit_type_t *orbit, const qpms_vswf_set_spec_t *bspec, const struct qpms_finite_group_t *sym, const qpms_gmi_t g)
 Constructs a "full matrix action" of a point group element for an orbit type. More...
 
complex double * qpms_orbit_irrep_projector_matrix (complex double *target, const qpms_ss_orbit_type_t *orbit, const qpms_vswf_set_spec_t *bspec, const struct qpms_finite_group_t *sym, const qpms_iri_t iri)
 Constructs a dense matrix representation of a irrep projector for an orbit type. More...
 
complex double * qpms_orbit_irrep_basis (size_t *basis_size, complex double *target, const qpms_ss_orbit_type_t *orbit, const qpms_vswf_set_spec_t *bspec, const struct qpms_finite_group_t *sym, const qpms_iri_t iri)
 TODO DOC!!!!! More...
 
complex double * qpms_scatsys_incident_field_vector_full (complex double *target_full, const qpms_scatsys_t *ss, qpms_incfield_t field_at_point, const void *args, bool add)
 Creates an incident field vector in the full basis, given a function that evaluates the field expansions at points. More...
 
complex double * qpms_scatsysw_apply_Tmatrices_full (complex double *target_full, const complex double *inc_full, const qpms_scatsys_at_omega_t *ssw)
 Applies T-matrices onto an incident field vector in the full basis. More...
 
struct beyn_result_tqpms_scatsys_finite_find_eigenmodes (const qpms_scatsys_t *ss, qpms_iri_t iri, complex double omega_centre, double omega_rr, double omega_ri, size_t contour_npoints, double rank_tol, size_t rank_sel_min, double res_tol)
 Searches for finite scattering system's eigenmodes using Beyn's algorithm. More...
 
ccart3_t qpms_scatsys_scattered_E (const qpms_scatsys_t *ss, qpms_bessel_t typ, complex double wavenumber, const complex double *scatcoeff_full, cart3_t evalpoint)
 Evaluates scattered electric field at a point, given a full vector of scattered field coefficients. More...
 
ccart3_t qpms_scatsysw_scattered_E (const qpms_scatsys_at_omega_t *ssw, qpms_bessel_t typ, const complex double *scatcoeff_full, cart3_t evalpoint)
 Evaluates scattered electric field at a point, given a full vector of scattered field coefficients. More...
 
ccart3_t qpms_scatsys_scattered_E__alt (const qpms_scatsys_t *ss, qpms_bessel_t typ, complex double wavenumber, const complex double *scatcoeff_full, cart3_t evalpoint)
 Evaluates scattered electric field at a point, given a full vector of scattered field coefficients. More...
 
ccart3_t qpms_scatsysw_scattered_E__alt (const qpms_scatsys_at_omega_t *ssw, qpms_bessel_t typ, const complex double *scatcoeff_full, cart3_t evalpoint)
 Evaluates scattered electric field at a point, given a full vector of scattered field coefficients. More...
 
ccart3_t qpms_scatsyswk_scattered_E (const qpms_scatsys_at_omega_k_t *sswk, qpms_bessel_t typ, const complex double *scatcoeff_full, cart3_t evalpoint)
 Evaluates scattered electric field at a point, given a full vector of scattered field coefficients. (Periodic system.) More...
 
qpms_errno_t qpms_scatsyswk_scattered_field_basis (ccart3_t *target, const qpms_scatsys_at_omega_k_t *sswk, qpms_bessel_t typ, cart3_t evalpoint)
 Evaluates "scattered" field basis functions in a periodic system. More...
 
double qpms_ss_adjusted_eta (const qpms_scatsys_t *ss, complex double wavenumber, const double wavevector[3])
 Adjusted Ewadl parameter to avoid high-frequency breakdown.
 

Detailed Description

Modern interface for finite lattice calculations, including symmetries.

N.B. Only "reasonably normalised" waves are supported now in most of the functions defined here, i.e. those that can be rotated by the usual Wigner matrices, i.e. the "power" or "spharm" -normalised ones.

TODO FIXME check whether Condon-Shortley phase can have some nasty influence here; I fear that yes.

Typedef Documentation

◆ qpms_particle_t

A particle, defined by its T-matrix and position.

This is rather only an auxillary intermediate structure to ultimately build an qpms_scatsys_t instance

◆ qpms_scatsys_at_omega_k_t

Scattering system at a given frequency and k-vector. Used only with periodic systems.

N.B. use as a stack variable now, but this might become heap-allocated in the future (with own con- and destructor)

◆ qpms_scatsys_t

Common "class" for system of scatterers, both periodic and non-periodic.

Infinite periodic structures (those with lattice_dimension > 0) have the per filled. These are ignored for finite systems (lattice_dimension == 0).

◆ qpms_ss_orbit_type_t

Structure describing a particle's "orbit type" under symmetry group actions in a system.

Each particle has its orbit with respect to a symmetry group of the system in which the particle lies, i.e. a set of particles onto which the group operations map the original particle.

(TODO DOC improve the following explanation:) Typically, there will be only a small number of different (T-matrix, particle stabiliser) pairs in the system. We can group the particles accordingly, into the same "orbit types" that will allow to do certain operations only once for each "orbit type", saving memory and (perhaps) time.

Each particle will then have assigned:

  1. an orbit type,
  2. an ID inside that orbit.

TODO DOC how the process of assigning the particle IDs actually work, orbit type (non-)uniqueness.

Memory is managed by qpms_scatspec_t; qpms_ss_orbit_type_t does not own anything.

Function Documentation

◆ qpms_orbit_action_matrix()

complex double* qpms_orbit_action_matrix ( complex double *  target,
const qpms_ss_orbit_type_t orbit,
const qpms_vswf_set_spec_t bspec,
const struct qpms_finite_group_t sym,
const qpms_gmi_t  g 
)

Constructs a "full matrix action" of a point group element for an orbit type.

TODO detailed doc

Parameters
targetTarget array. If NULL, a new one is allocated. The size of the array is (orbit->size * bspec->n)**2 (it makes sense to assume all the T-matrices share their spec).
orbitThe orbit (type).
bspecBase spec of the t-matrices (we don't know it from orbit, as it has only T-matrix indices.
symThe symmetry group used to generate the orbit (must have rep3d filled).
gThe index of the operation in sym to represent.

◆ qpms_orbit_irrep_basis()

complex double* qpms_orbit_irrep_basis ( size_t *  basis_size,
complex double *  target,
const qpms_ss_orbit_type_t orbit,
const qpms_vswf_set_spec_t bspec,
const struct qpms_finite_group_t sym,
const qpms_iri_t  iri 
)

TODO DOC!!!!!

Parameters
basis_sizeHere theh size of theh basis shall be saved,
targetTarget array. If NULL, a new one is allocated. The size of the array is basis_size * (orbit->size * bspec->n) (it makes sense to assume all the T-matrices share their spec).
orbitThe orbit (type).
bspecBase spec of the t-matrices (we don't know it from orbit, as it has only T-matrix indices.
symThe symmetry group used to generate the orbit (must have rep3d filled).
iriThe index of the irreducible representation of sym.

◆ qpms_orbit_irrep_projector_matrix()

complex double* qpms_orbit_irrep_projector_matrix ( complex double *  target,
const qpms_ss_orbit_type_t orbit,
const qpms_vswf_set_spec_t bspec,
const struct qpms_finite_group_t sym,
const qpms_iri_t  iri 
)

Constructs a dense matrix representation of a irrep projector for an orbit type.

TODO detailed doc

Parameters
targetTarget array. If NULL, a new one is allocated. The size of the array is (orbit->size * bspec->n)**2 (it makes sense to assume all the T-matrices share their spec).
orbitThe orbit (type).
bspecBase spec of the t-matrices (we don't know it from orbit, as it has only T-matrix indices.
symThe symmetry group used to generate the orbit (must have rep3d filled).
iriThe index of the irreducible representation of sym.

◆ qpms_scatsys_apply_symmetry()

qpms_scatsys_at_omega_t* qpms_scatsys_apply_symmetry ( const qpms_scatsys_t orig,
const struct qpms_finite_group_t sym,
complex double  omega,
const struct qpms_tolerance_spec_t tol 
)

Creates a new scatsys by applying a symmetry group onto a "proto-scatsys", copying particles if needed.

In fact, it copies everything except the vswf set specs and qpms_abstract_tmatrix_t instances, so keep them alive until scatsys is destroyed.

The following fields must be filled in the "proto- scattering system" orig:

  • orig->lattice_dimension
  • orig->medium – The pointers are copied to the new qpms_scatsys_t instance; the target qpms_abstract_tmatrix_t objects must be kept alive before all the resulting qpms_scatsys_t instances are properly destroyed.
  • orig->tmg – The pointers are copied to the new qpms_scatsys_t instance; the target qpms_abstract_tmatrix_t objects must be kept alive before all the resulting qpms_scatsys_t instances are properly destroyed. The pointers from orig->tmg, however, are copied.
  • orig->tmg_count
  • orig->tm – Must be filled, although the operations will typically be identities (QPMS_TMATRIX_OPERATION_NOOP). N.B. these NOOPs might be replaced with some symmetrisation operation in the resulting "full" qpms_scatsys_t instance.
  • orig->tm_count
  • orig->p
  • orig->p_count

For periodic systems, the corresponding number of orig->per->lattice_basis[] elements must be filled as well.

For periodic systems, only trivial group is currently supported. Non-trivial groups will cause undefined behaviour.

The resulting qpms_scatsys_t is obtained by actually evaluating the T-matrices at the given frequency omega and where applicable, these are compared by their values with given tolerances. The T-matrix generators are expected to preserve the point group symmetries for all frequencies.

This particular function will likely change in the future.

Returns
An instance sso of qpms_scatsys_omega_t. Note that sso->ss must be saved by the caller before destroying sso (with qpms_scatsys_at_omega_free(), and destroyed only afterwards with qpms_scatsys_free() when not needed anymore. sso->ss can be reused for different frequency by a qpms_scatsys_at_omega() call.

◆ qpms_scatsys_at_omega()

qpms_scatsys_at_omega_t* qpms_scatsys_at_omega ( const qpms_scatsys_t ss,
complex double  omega 
)

Evaluates scattering system T-matrices at a given frequency.

Free the result using qpms_scatsys_at_omega_free() when done.

◆ qpms_scatsys_at_omega_free()

void qpms_scatsys_at_omega_free ( qpms_scatsys_at_omega_t ssw)

◆ qpms_scatsys_build_translation_matrix_e_full()

complex double* qpms_scatsys_build_translation_matrix_e_full ( complex double *  target,
const qpms_scatsys_t ss,
complex double  k,
qpms_bessel_t  J 
)

As qpms_scatsys_build_translation_full() but with choice of Bessel function type.

Might be useful for evaluation of cross sections and testing.

Parameters
targetTarget memory with capacity for ss->fecv_size**2 elements. If NULL, new will be allocated.
kWave number to use in the translation matrix.
JBessel function type.

◆ qpms_scatsys_build_translation_matrix_e_irrep_packed()

complex double* qpms_scatsys_build_translation_matrix_e_irrep_packed ( complex double *  target,
const qpms_scatsys_t ss,
qpms_iri_t  iri,
complex double  k,
qpms_bessel_t  J 
)

Global translation matrix with selectable Bessel function, projected on an irrep.

The diagonal (particle self-) blocks are currently filled with zeros. This may change in the future.

Parameters
targetTarget memory with capacity for ss->fecv_size**2 elements. If NULL, new will be allocated.
kWave number to use in the translation matrix.

◆ qpms_scatsys_build_translation_matrix_full()

complex double* qpms_scatsys_build_translation_matrix_full ( complex double *  target,
const qpms_scatsys_t ss,
complex double  k 
)

Global translation matrix.

The diagonal (particle self-) block are filled with zeros (even for regular Bessel waves). This may change in the future.

Parameters
targetTarget memory with capacity for ss->fecv_size**2 elements. If NULL, new will be allocated.
kWave number to use in the translation matrix.

◆ qpms_scatsys_finite_find_eigenmodes()

struct beyn_result_t* qpms_scatsys_finite_find_eigenmodes ( const qpms_scatsys_t ss,
qpms_iri_t  iri,
complex double  omega_centre,
double  omega_rr,
double  omega_ri,
size_t  contour_npoints,
double  rank_tol,
size_t  rank_sel_min,
double  res_tol 
)

Searches for finite scattering system's eigenmodes using Beyn's algorithm.

Currently, elliptical contour is used.

TODO In the future, this will probably support irrep decomposition as well, but it does not make much sense for periodic / small systems, as in their case the bottleneck is the T-matrix and translation matrix evaluation rather than the linear algebra.

Parameters
iriA valid irrep index to search only in one irrep, or QPMS_NO_IRREP for solving the full system.
omega_centreCenter of the ellipse inside which the eigenfreqs are searched for.
omega_rrReal half-axis of the ellipse inside which the eigenfreqs are searched for.
omega_riImaginary half-axis of the ellipse inside which the eigenfreqs are searched for.
contour_npointsNumber of elliptic contour discretisation points (preferably even number)
rank_tol(default: `1e-4`) TODO DOC.
rank_sel_minMinimum number of eigenvalue candidates, even if they don't pass rank_tol.
res_tol(default: `0.0`) TODO DOC.

◆ qpms_scatsys_incident_field_vector_full()

complex double* qpms_scatsys_incident_field_vector_full ( complex double *  target_full,
const qpms_scatsys_t ss,
qpms_incfield_t  field_at_point,
const void *  args,
bool  add 
)

Creates an incident field vector in the full basis, given a function that evaluates the field expansions at points.

TODO detailed doc!

Returns
target_full if target_full was not NULL, otherwise the newly allocated array.
Parameters
target_fullTarget array. If NULL, a new one is allocated. The length of the array is ss->fecv_size.
argsPointer passed as the last argument to (*field_at_point)()
addIf true, add to target_full; rewrite target_full if false.

◆ qpms_scatsys_irrep_pack_matrix()

complex double* qpms_scatsys_irrep_pack_matrix ( complex double *  target_packed,
const complex double *  orig_full,
const qpms_scatsys_t ss,
qpms_iri_t  iri 
)

Projects a "big" matrix onto an irrep.

TODO doc

◆ qpms_scatsys_irrep_pack_matrix_stupid()

complex double* qpms_scatsys_irrep_pack_matrix_stupid ( complex double *  target_packed,
const complex double *  orig_full,
const qpms_scatsys_t ss,
qpms_iri_t  iri 
)

Projects a "big" matrix onto an irrep (slow reference implementation).

TODO doc

◆ qpms_scatsys_irrep_pack_vector()

complex double* qpms_scatsys_irrep_pack_vector ( complex double *  target_packed,
const complex double *  orig_full,
const qpms_scatsys_t ss,
const qpms_iri_t  iri 
)

Projects a "big" vector onto an irrep.

TODO doc

◆ qpms_scatsys_irrep_transform_matrix()

complex double* qpms_scatsys_irrep_transform_matrix ( complex double *  target_U,
const qpms_scatsys_t ss,
qpms_iri_t  iri 
)

Creates a "full" transformation matrix U that takes a full vector and projects it onto an symmetry adapted basis.

Mostly as a reference and a debugging tool, as multiplicating these big matrices would be inefficient.

TODO doc about shape etc.

◆ qpms_scatsys_irrep_unpack_matrix()

complex double* qpms_scatsys_irrep_unpack_matrix ( complex double *  target_full,
const complex double *  orig_packed,
const qpms_scatsys_t ss,
qpms_iri_t  iri,
bool  add 
)

Transforms a big "packed" matrix into the full basis.

TODO doc

◆ qpms_scatsys_irrep_unpack_matrix_stupid()

complex double* qpms_scatsys_irrep_unpack_matrix_stupid ( complex double *  target_full,
const complex double *  orig_packed,
const qpms_scatsys_t ss,
qpms_iri_t  iri,
bool  add 
)

Transforms a big "packed" matrix into the full basis (slow reference implementation).

TODO doc

Transforms a big "packed" matrix into the full basis (slow reference implementation).

◆ qpms_scatsys_irrep_unpack_vector()

complex double* qpms_scatsys_irrep_unpack_vector ( complex double *  target_full,
const complex double *  orig_packed,
const qpms_scatsys_t ss,
const qpms_iri_t  iri,
bool  add 
)

Transforms a big "packed" vector into the full basis.

TODO doc

◆ qpms_scatsys_periodic_build_translation_matrix_full()

complex double* qpms_scatsys_periodic_build_translation_matrix_full ( complex double *  target,
const qpms_scatsys_t ss,
complex double  wavenumber,
const cart3_t wavevector,
double  eta 
)

Global translation matrix.

Parameters
targetTarget memory with capacity for ss->fecv_size**2 elements. If NULL, new will be allocated.
wavenumberWave number to use in the translation matrix.
wavevectorWavevector / pseudomomentum in cartesian coordinates.
etaEwald parameter eta. Pass 0 or NaN to use the default value in ss.

◆ qpms_scatsys_periodic_find_eigenmodes()

struct beyn_result_t* qpms_scatsys_periodic_find_eigenmodes ( const qpms_scatsys_t ss,
const double  k[3],
complex double  omega_centre,
double  omega_rr,
double  omega_ri,
size_t  contour_npoints,
double  rank_tol,
size_t  rank_sel_min,
double  res_tol 
)

Searches for periodic scattering system's eigenmodes using Beyn's algorithm.

Currently, elliptical contour is used.

TODO In the future, this will probably support irrep decomposition as well, but for the case of periodic / small systems, the bottleneck is the T-matrix and translation matrix evaluation rather than the linear algebra.

Parameters
kWavevector in cartesian coordinates (must lie in the lattice plane).
omega_centreCenter of the ellipse inside which the eigenfreqs are searched for.
omega_rrReal half-axis of the ellipse inside which the eigenfreqs are searched for.
omega_riImaginary half-axis of the ellipse inside which the eigenfreqs are searched for.
contour_npointsNumber of elliptic contour discretisation points (preferably even number)
rank_tol(default: `1e-4`) TODO DOC.
rank_sel_minMinimum number of eigenvalue candidates, even if they don't pass rank_tol.
res_tol(default: `0.0`) TODO DOC.

◆ qpms_scatsys_scatter_solve()

complex double* qpms_scatsys_scatter_solve ( complex double *  target_f,
const complex double *  a_inc,
qpms_ss_LU  ludata 
)

Solves a (possibly partial, irrep-packed) scattering problem \( (I-TS)f = Ta_\mathrm{inc} \) using a pre-factorised \( (I-TS) \).

Parameters
target_fTarget (full or irrep-packed, depending on `ludata.full`) array for f. If NULL, a new one is allocated.
a_incIncident field expansion coefficient vector a (full or irrep-packed, depending on `ludata.full`).
ludataPre-factorised \( I - TS \) matrix data.

◆ qpms_scatsys_scattered_E()

ccart3_t qpms_scatsys_scattered_E ( const qpms_scatsys_t ss,
qpms_bessel_t  typ,
complex double  wavenumber,
const complex double *  scatcoeff_full,
cart3_t  evalpoint 
)

Evaluates scattered electric field at a point, given a full vector of scattered field coefficients.

This function evaluates the field \( \vect E (\vect r ) \), with any given wavenumber of the background medium and any given vector of the excitation coefficients \( \wckcout \).

Returns
Complex electric field at the point defined by evalpoint.
See also
qpms_scatsysw_scattered_E()
qpms_scatsyswk_scattered_E() for periodic systems.
Parameters
typBessel function kind to use (for scattered fields, use QPMS_HANKEL_PLUS).
wavenumberWavenumber of the background medium.
scatcoeff_fullFull vector of the scattered field coefficients \( \wckcout \).
evalpointA point \( \vect r \), at which the field is evaluated.

◆ qpms_scatsys_scattered_E__alt()

ccart3_t qpms_scatsys_scattered_E__alt ( const qpms_scatsys_t ss,
qpms_bessel_t  typ,
complex double  wavenumber,
const complex double *  scatcoeff_full,
cart3_t  evalpoint 
)

Evaluates scattered electric field at a point, given a full vector of scattered field coefficients.

This is an alternative implementation of qpms_scatsys_scattered_E(), and should give the same results up to rounding errors.

Returns
Complex electric field at the point defined by evalpoint.
See also
qpms_scatsys_scattered_E()
Parameters
typBessel function kind to use (for scattered fields, use QPMS_HANKEL_PLUS).
wavenumberWavenumber of the background medium.
scatcoeff_fullFull vector of the scattered field coefficients \( \wckcout \).
evalpointA point \( \vect r \), at which the field is evaluated.

◆ qpms_scatsystem_set_nthreads()

void qpms_scatsystem_set_nthreads ( long  n)

Overrides the number of threads spawned by the paralellized functions.

TODO MORE DOC which are those?

◆ qpms_scatsysw_apply_Tmatrices_full()

complex double* qpms_scatsysw_apply_Tmatrices_full ( complex double *  target_full,
const complex double *  inc_full,
const qpms_scatsys_at_omega_t ssw 
)

Applies T-matrices onto an incident field vector in the full basis.

Parameters
inc_fullTarget vector array. If NULL, a new one is allocated.
sswIncident field coefficient vector. Must not be NULL.

◆ qpms_scatsysw_build_modeproblem_matrix_full()

complex double* qpms_scatsysw_build_modeproblem_matrix_full ( complex double *  target,
const qpms_scatsys_at_omega_t ssw 
)

Creates the full \( (I - WS) \) matrix of the periodic scattering system.

Returns
target on success, NULL on error.
Parameters
targetTarget memory with capacity for ss->fecv_size**2 elements. If NULL, new will be allocated.

◆ qpms_scatsysw_build_modeproblem_matrix_full_LU()

qpms_ss_LU qpms_scatsysw_build_modeproblem_matrix_full_LU ( complex double *  target,
int *  target_piv,
const qpms_scatsys_at_omega_t ssw 
)

Builds an LU-factorised mode/scattering problem \( (I - TS) \) matrix from scratch. Nonperiodic systems only.

Parameters
targetPre-allocated target array. Optional (if NULL, new one is allocated).
target_pivPre-allocated pivot array. Optional (if NULL, new one is allocated).

◆ qpms_scatsysw_build_modeproblem_matrix_irrep_packed()

complex double* qpms_scatsysw_build_modeproblem_matrix_irrep_packed ( complex double *  target,
const qpms_scatsys_at_omega_t ssw,
qpms_iri_t  iri 
)

Creates the mode problem matrix \( (I - TS) \) directly in the irrep-packed form.

Returns
target on success, NULL on error.
Parameters
targetTarget memory with capacity for ss->fecv_size**2 elements. If NULL, new will be allocated.
iriIndex of the irreducible representation in ssw->ss->sym

◆ qpms_scatsysw_build_modeproblem_matrix_irrep_packed_LU()

qpms_ss_LU qpms_scatsysw_build_modeproblem_matrix_irrep_packed_LU ( complex double *  target,
int *  target_piv,
const qpms_scatsys_at_omega_t ssw,
qpms_iri_t  iri 
)

Builds an irrep-packed LU-factorised mode/scattering problem matrix from scratch.

Parameters
targetPre-allocated target array. Optional (if NULL, new one is allocated).
target_pivPre-allocated pivot array. Optional (if NULL, new one is allocated).

◆ qpms_scatsysw_build_modeproblem_matrix_irrep_packed_orbitorderR()

complex double* qpms_scatsysw_build_modeproblem_matrix_irrep_packed_orbitorderR ( complex double *  target,
const qpms_scatsys_at_omega_t ssw,
qpms_iri_t  iri 
)

Alternative implementation of qpms_scatsysw_build_modeproblem_matrix_irrep_packed().

Returns
target on success, NULL on error.
Parameters
targetTarget memory with capacity for ss->fecv_size**2 elements. If NULL, new will be allocated.
iriIndex of the irreducible representation in ssw->ss->sym

◆ qpms_scatsysw_build_modeproblem_matrix_irrep_packed_serial()

complex double* qpms_scatsysw_build_modeproblem_matrix_irrep_packed_serial ( complex double *  target,
const qpms_scatsys_at_omega_t ssw,
qpms_iri_t  iri 
)

Alternative (serial reference) implementation of qpms_scatsysw_build_modeproblem_matrix_irrep_packed().

Returns
target on success, NULL on error.
Parameters
targetTarget memory with capacity for ss->fecv_size**2 elements. If NULL, new will be allocated.
iriIndex of the irreducible representation in ssw->ss->sym

◆ qpms_scatsysw_modeproblem_matrix_full_factorise()

qpms_ss_LU qpms_scatsysw_modeproblem_matrix_full_factorise ( complex double *  modeproblem_matrix_full,
int *  target_piv,
const qpms_scatsys_at_omega_t ssw,
const struct qpms_scatsys_at_omega_k_t sswk 
)

Computes LU factorisation of a pre-calculated mode/scattering problem matrix, replacing its contents.

Parameters
modeproblem_matrix_fullPre-calculated mode problem matrix (I-TS). Mandatory.
target_pivPre-allocated pivot array. Optional (if NULL, new one is allocated).
sswMust be filled for non-periodic systems.
sswkMust be filled for periodic systems, otherwise must be NULL.

◆ qpms_scatsysw_modeproblem_matrix_irrep_packed_factorise()

qpms_ss_LU qpms_scatsysw_modeproblem_matrix_irrep_packed_factorise ( complex double *  modeproblem_matrix_irrep_packed,
int *  target_piv,
const qpms_scatsys_at_omega_t ssw,
qpms_iri_t  iri 
)

Computes LU factorisation of a pre-calculated irrep-packed mode/scattering problem matrix, replacing its contents.

Parameters
modeproblem_matrix_irrep_packedPre-calculated mode problem matrix (I-TS). Mandatory.
target_pivPre-allocated pivot array. Optional (if NULL, new one is allocated).

◆ qpms_scatsysw_scattered_E()

ccart3_t qpms_scatsysw_scattered_E ( const qpms_scatsys_at_omega_t ssw,
qpms_bessel_t  typ,
const complex double *  scatcoeff_full,
cart3_t  evalpoint 
)

Evaluates scattered electric field at a point, given a full vector of scattered field coefficients.

This function evaluates the field \( \vect E (\vect r ) \), with any given vector of the excitation coefficients \( \wckcout \).

Returns
Complex electric field at the point defined by evalpoint.
See also
qpms_scatsys_scattered_E()
qpms_scatsyswk_scattered_E() for periodic systems.
Parameters
typBessel function kind to use (for scattered fields, use QPMS_HANKEL_PLUS).
scatcoeff_fullFull vector of the scattered field coefficients \( \wckcout \).
evalpointA point \( \vect r \), at which the field is evaluated.

◆ qpms_scatsysw_scattered_E__alt()

ccart3_t qpms_scatsysw_scattered_E__alt ( const qpms_scatsys_at_omega_t ssw,
qpms_bessel_t  typ,
const complex double *  scatcoeff_full,
cart3_t  evalpoint 
)

Evaluates scattered electric field at a point, given a full vector of scattered field coefficients.

This is an alternative implementation of qpms_scatsys_scattered_E(), and should give the same results up to rounding errors.

Returns
Complex electric field at the point defined by evalpoint.
See also
qpms_scatsysw_scattered_E()
Parameters
typBessel function kind to use (for scattered fields, use QPMS_HANKEL_PLUS).
scatcoeff_fullFull vector of the scattered field coefficients \( \wckcout \).
evalpointA point \( \vect r \), at which the field is evaluated.

◆ qpms_scatsyswk_build_modeproblem_matrix_full()

complex double* qpms_scatsyswk_build_modeproblem_matrix_full ( complex double *  target,
const qpms_scatsys_at_omega_k_t sswk 
)

Creates the full \( (I - WS) \) matrix of the periodic scattering system.

Returns
target on success, NULL on error.
Parameters
targetTarget memory with capacity for ss->fecv_size**2 elements. If NULL, new will be allocated.

◆ qpms_scatsyswk_build_modeproblem_matrix_full_LU()

qpms_ss_LU qpms_scatsyswk_build_modeproblem_matrix_full_LU ( complex double *  target,
int *  target_piv,
const qpms_scatsys_at_omega_k_t sswk 
)

Builds an LU-factorised mode/scattering problem \( (I - TS) \) matrix from scratch. Periodic systems only.

Parameters
targetPre-allocated target array. Optional (if NULL, new one is allocated).
target_pivPre-allocated pivot array. Optional (if NULL, new one is allocated).

◆ qpms_scatsyswk_build_translation_matrix_full()

complex double* qpms_scatsyswk_build_translation_matrix_full ( complex double *  target,
const qpms_scatsys_at_omega_k_t sswk 
)

Global translation matrix.

Parameters
targetTarget memory with capacity for ss->fecv_size**2 elements. If NULL, new will be allocated.

◆ qpms_scatsyswk_scattered_E()

ccart3_t qpms_scatsyswk_scattered_E ( const qpms_scatsys_at_omega_k_t sswk,
qpms_bessel_t  typ,
const complex double *  scatcoeff_full,
cart3_t  evalpoint 
)

Evaluates scattered electric field at a point, given a full vector of scattered field coefficients. (Periodic system.)

This function evaluates the field \( \vect E (\vect r ) \), with any given vector of the excitation coefficients \( \wckcout \).

Returns
Complex electric field at the point defined by evalpoint.
Bug:
Currently implemented only for btyp == QPMS_HANKEL_PLUS.
See also
qpms_scatsys_scattered_E(), qpms_scatsysw_scattered_E() for finite systems.
Parameters
typBessel function kind to use (only QPMS_HANKEL_PLUS is currently supported).
scatcoeff_fullFull vector of the scattered field coefficients \( \wckcout \).
evalpointA point \( \vect r \), at which the field is evaluated.

◆ qpms_scatsyswk_scattered_field_basis()

qpms_errno_t qpms_scatsyswk_scattered_field_basis ( ccart3_t target,
const qpms_scatsys_at_omega_k_t sswk,
qpms_bessel_t  typ,
cart3_t  evalpoint 
)

Evaluates "scattered" field basis functions in a periodic system.

See also
qpms_uvswf_fill() evaluates a set of VSWF basis functions (for finite systems).
Parameters
targetTarget array of size sswk->ssw->ss->fecv_size
typBessel function kind to use (only QPMS_HANKEL_PLUS is currently supponted).
evalpointA point \( \vect r \), at which the basis is evaluated.