qpms/vswf.h

   1 /*! \file vswf.h
   2  * \brief Vector spherical wavefunctions.
   3  *
   4  * N.B. for the Legendre polynomial norm definitions, see
   5  * <a href="https://www.gnu.org/software/gsl/doc/html/specfunc.html#associated-legendre-polynomials-and-spherical-harmonics">the corresponding section of GSL docs</a>
   6  * or <a href="http://git.savannah.gnu.org/cgit/gsl.git/tree/specfunc/legendre_source.c">gsl/specfunc/legendre_source.c</a>.
   7  */
   8 #ifndef QPMS_VSWF_H
   9 #define QPMS_VSWF_H
  10 #include <unistd.h> // ssize_t
  11 #include "qpms_types.h"
  12 #include <gsl/gsl_sf_legendre.h>
  13
  14 // -------------- Typedefs (function prototypes) for qpms_vswf_spec_t ----------------------
  15
  16 /// Calculates the (regular VSWF) expansion coefficients of an external incident field.
  17 typedef qpms_errno_t (*qpms_incfield_t)(
  18         /// Target non-NULL array of the regular VSWF expansion coefficients of length bspec->n.
  19         complex double *target,
  20         const qpms_vswf_set_spec_t *bspec,
  21         const cart3_t evalpoint, ///< Point at which the VSWF expansion is made.
  22         const void *args, ///< Pointer to additional function-specific arguments.
  23         bool add ///< If true, add to target; rewrite target if false.
  24 );
  25
  26 // ---------------Methods for qpms_vswf_set_spec_t-----------------------
  27 //
  28 /// Creates a qpms_vswf_set_spec_t structure with an empty list of wave indices.
  29 qpms_vswf_set_spec_t *qpms_vswf_set_spec_init(void);
  30 /// Appends a VSWF index to a \ref qpms_vswf_set_spec_t, also updating metadata.
  31 qpms_errno_t qpms_vswf_set_spec_append(qpms_vswf_set_spec_t *self, qpms_uvswfi_t u);
  32 /// Destroys a \ref qpms_vswf_set_spec_t.
  33 void qpms_vswf_set_spec_free(qpms_vswf_set_spec_t *);
  34 /// Compares two vswf basis specs.
  35 /**
  36  * Checks whether ilist is the same and of the same length.
  37  * If yes, returns true, else returns false.
  38  */
  39 bool qpms_vswf_set_spec_isidentical(const qpms_vswf_set_spec_t *a,
  40                 const qpms_vswf_set_spec_t *b);
  41 /// Copies an instance of qpms_vswf_set_spec_t
  42 qpms_vswf_set_spec_t *qpms_vswf_set_spec_copy(const qpms_vswf_set_spec_t *orig);
  43 /// Creates an instance of qpms_vswf_set_spec_t in the 'traditional' layout.
  44 qpms_vswf_set_spec_t *qpms_vswf_set_spec_from_lMax(qpms_l_t lMax,
  45                 qpms_normalisation_t norm);
  46
  47 /// Finds the position of a given index in the bspec's ilist.
  48 /** If not found, returns -1. */
  49 // TODO more consistency in types (here size_t vs. ptrdiff_t).
  50 static inline ssize_t qpms_vswf_set_spec_find_uvswfi(const qpms_vswf_set_spec_t *bspec,
  51                 const qpms_uvswfi_t index) {
  52         for(size_t i = 0; i < bspec->n; ++i)
  53                 if (bspec->ilist[i] == index)
  54                         return i;
  55         return -1;
  56 }
  57
  58 /// Creates an index mapping between two bspecs.
  59 /**
  60  * Creates an array r such that small->ilist[i] == big->ilist[r[i]].
  61  * It's not lossless if the two bspecs contain different combinations of waves.
  62  *
  63  * Preferably, big->ilist contains everything small->ilist does.
  64  * If small->ilist[i] is not found in big->ilist, r[i] will be set to ~(size_t)0.
  65  *
  66  * Discard with free() after use.
  67  */
  68 size_t *qpms_vswf_set_reindex(const qpms_vswf_set_spec_t *small, const qpms_vswf_set_spec_t *big);
  69
  70
  71 /// Evaluates a set of VSWF basis functions at a given point.
  72 /** The list of basis wave indices is specified in \a setspec;
  73  *  \a setspec->norm must be set as well.
  74  */
  75 qpms_errno_t qpms_uvswf_fill(
  76                 csphvec_t *const target, ///< Target array of size at least setspec->n.
  77                 const qpms_vswf_set_spec_t *setspec,
  78                 csph_t kr, ///< Evaluation point.
  79                 qpms_bessel_t btyp);
  80
  81 /// Evaluates field specified by SVWF coefficients at a given point.
  82 /** SVWF coefficients in \a coeffs must be ordered according to \a setspec->ilist.
  83  */
  84 csphvec_t qpms_eval_uvswf(const qpms_vswf_set_spec_t *setspec,
  85                 const complex double *coeffs, ///< SVWF coefficient vector of size setspec->n.
  86                 csph_t kr, ///< Evaluation point.
  87                 qpms_bessel_t btyp);
  88
  89
  90 // --- qpms_incfield_t instances and their arguments
  91
  92 /// Parameter structure for qpms_incfield_planewave()
  93 typedef struct qpms_incfield_planewave_params_t {
  94         bool use_cartesian; ///< If true, wave direction k and amplitude E are specified in cartesian coordinates (via k.cart, E.cart). If false, k is specified in spherical coordinates and E are specified in the corresponding geographical coordinates (via k.sph, E.sph).
  95         union {
  96                 ccart3_t cart;
  97                 csph_t sph;
  98         } k; ///< Wave vector.
  99         union {
 100                 ccart3_t cart;
 101                 csphvec_t sph;
 102         } E; ///< Electric field amplitude at origin.
 103 } qpms_incfield_planewave_params_t;
 104
 105
 106 /// Calculates the (regular VSWF) expansion coefficients of a plane wave.
 107 /** The wave amplitude and wave vector is defined by struct qpms_incfield_planewave_params_t.
 108  *
 109  *  If the wave vector and amplitude are not orthogonal (i.e. the plane wave is not
 110  *  fully transversal) and bspec->lMaxL is non-negative,
 111  *  the corresponding longitudinal components are calculated as well.
 112  *
 113  *  For complex k vectors, the implementation is not completely correct right now.
 114  *  Locally, it corresponds to decomposition of a plane wave with a real \a k
 115  *  (using the real part of the \a k supplied), just the whole decomposition
 116  *  is modulated by the origin-dependent factor \f$ \vect E e^{i \vect k \cdot \vect r} \f$
 117  *  with \f$ \vect k \f$ complex.
 118  */
 119 qpms_errno_t qpms_incfield_planewave(
 120         /// Target non-NULL array of the regular VSWF expansion coefficients of length bspec->n.
 121         complex double *target,
 122         const qpms_vswf_set_spec_t *bspec,
 123         const cart3_t evalpoint, ///< Point at which the VSWF expansion is made.
 124         const void *args, ///< Pointer to additional function-specific arguments (converted to (const qpms_incfield_planewave_params_t *)).
 125         bool add ///< If true, add to target; rewrite target if false.
 126 );
 127
 128 // -----------------------------------------------------------------------
 129
 130 /// Electric wave N.
 131 csphvec_t qpms_vswf_single_el(qpms_m_t m, qpms_l_t n, sph_t kdlj,
 132                 qpms_bessel_t btyp, qpms_normalisation_t norm);
 133 /// Magnetic wave M.
 134 csphvec_t qpms_vswf_single_mg(qpms_m_t m, qpms_l_t n, sph_t kdlj,
 135                 qpms_bessel_t btyp, qpms_normalisation_t norm);
 136 /// Electric wave N, complex wave number version.
 137 csphvec_t qpms_vswf_single_el_csph(qpms_m_t m, qpms_l_t n, csph_t kdlj,
 138                 qpms_bessel_t btyp, qpms_normalisation_t norm);
 139 /// Magnetic wave M, complex wave number version..
 140 csphvec_t qpms_vswf_single_mg_csph(qpms_m_t m, qpms_l_t n, csph_t kdlj,
 141                 qpms_bessel_t btyp, qpms_normalisation_t norm);
 142
 143 /// Set of electric and magnetic VSWF values in spherical coordinate basis.
 144 /** This is supposed to contain all the waves up to $l = lMax$.
 145  *
 146  *  For a completely custom set of waves, use \ref qpms_uvswfset_sph_t instead.
 147  */
 148 typedef struct qpms_vswfset_sph_t {
 149         //qpms_normalisation_t norm;
 150         qpms_l_t lMax;
 151         //qpms_y_t nelem;
 152         //sph_t kdlj
 153         csphvec_t *el, *mg;
 154 } qpms_vswfset_sph_t;
 155
 156 qpms_errno_t qpms_legendre_deriv_y_get(double **result, double **result_deriv, double x, qpms_l_t lMax,
 157                 gsl_sf_legendre_t lnorm, double csphase); // free() result and result_deriv yourself!
 158 qpms_errno_t qpms_legendre_deriv_y_fill(double *where, double *where_deriv, double x,
 159                 qpms_l_t lMax, gsl_sf_legendre_t lnorm, double csphase);
 160
 161
 162 /// Evaluate the zeroth-degree longitudinal VSWF \f$ \mathbf{L}_0^0 \f$.
 163 /**
 164  * Any `norm` is being ignored right now.
 165  */
 166 csphvec_t qpms_vswf_L00(
 167                 csph_t kdrj, //< VSWF evaluation point.
 168                 qpms_bessel_t btyp,
 169                 qpms_normalisation_t norm //< Ignored!
 170                 );
 171
 172 /// Evaluate VSWFs at a given point from \a l = 1 up to a given degree \a lMax.
 173 /**
 174  * The target arrays \a resultL, \a resultM, \a resultN have to be large enough to contain
 175  * \a lMax * (\a lMax + 2) elements. If NULL is passed instead, the corresponding SVWF type
 176  * is not evaluated.
 177  *
 178  * Does not evaluate the zeroth-order wave \f$ \mathbf{L}_0^0 \f$.
 179  * If you need that, use qpms_vswf_L00().
 180  *
 181  * If \a kdrj.r == 0 and \a btyp != QPMS_BESSEL_REGULAR, returns QPMS_ESING and fills
 182  * targets with NaNs.
 183  *
 184  * \return 0 (QPMS_SUCCESS) on success, otherwise.
 185  */
 186 qpms_errno_t qpms_vswf_fill(
 187         csphvec_t *resultL, //< Target array for longitudinal VSWFs.
 188         csphvec_t *resultM, //< Target array for magnetic VSWFs.
 189         csphvec_t *resultN, //< Target array for electric VSWFs.
 190         qpms_l_t lMax, //< Maximum multipole degree to be calculated.
 191         sph_t kdrj, //< VSWF evaluation point.
 192         qpms_bessel_t btyp, qpms_normalisation_t norm);
 193
 194 // Should give the same results: for consistency checks
 195 qpms_errno_t qpms_vswf_fill_alternative(csphvec_t *resultL, csphvec_t *resultM, csphvec_t *resultN, qpms_l_t lMax, sph_t kdrj,
 196                 qpms_bessel_t btyp, qpms_normalisation_t norm);
 197
 198 /// Evaluate VSWFs at a given point from \a l = 1 up to a given degree \a lMax (complex \a kr version).
 199 /**
 200  * The target arrays \a resultL, \a resultM, \a resultN have to be large enough to contain
 201  * \a lMax * (\a lMax + 2) elements. If NULL is passed instead, the corresponding SVWF type
 202  * is not evaluated.
 203  *
 204  * Does not evaluate the zeroth-order wave \f$ \mathbf{L}_0^0 \f$.
 205  * If you need that, use qpms_vswf_L00().
 206  */
 207 qpms_errno_t qpms_vswf_fill_csph(
 208         csphvec_t *resultL, //< Target array for longitudinal VSWFs.
 209         csphvec_t *resultM, //< Target array for magnetic VSWFs.
 210         csphvec_t *resultN, //< Target array for electric VSWFs.
 211         qpms_l_t lMax, //< Maximum multipole degree to be calculated.
 212         csph_t kdrj, //< VSWF evaluation point.
 213         qpms_bessel_t btyp, qpms_normalisation_t norm);
 214
 215 /// Evaluate vector spherical harmonics at \a dir.
 216 /**
 217  * The length of each of the target arrays shall be `lMax * (lMax + 2)`.
 218  * If any of the target arrays pointers is NULL, the corresponding VSH will not be evaluated.
 219  * The "zeroth" radial VSH \f$ \vshrad_0^0 \f$ is not evaluated.
 220  */
 221 qpms_errno_t qpms_vecspharm_fill(
 222                 csphvec_t *const a1target, ///< Target array for radial VSH \f$ \vshrad \f$.
 223                 csphvec_t *const a2target, ///< Target array for "rotational" VSH \f$ \vshrot \f$.
 224                 csphvec_t *const a3target, ///< Target array for "gradiental" VSH \f$ \vshgrad \f$.
 225                 qpms_l_t lMax, sph_t dir, qpms_normalisation_t norm);
 226 qpms_errno_t qpms_vecspharm_dual_fill(csphvec_t *const a1target, csphvec_t *const a2target, csphvec_t *const a3target,
 227                 qpms_l_t lMax, sph_t dir, qpms_normalisation_t norm);
 228
 229 qpms_errno_t qpms_planewave2vswf_fill_cart(cart3_t wavedir, ccart3_t amplitude,
 230                 complex double *targt_longcoeff, complex double *target_mgcoeff, complex double *target_elcoeff,
 231                 qpms_l_t lMax, qpms_normalisation_t norm);
 232 qpms_errno_t qpms_planewave2vswf_fill_sph(sph_t wavedir, csphvec_t amplitude,
 233                 complex double *targt_longcoeff, complex double *target_mgcoeff, complex double *target_elcoeff,
 234                 qpms_l_t lMax, qpms_normalisation_t norm);
 235
 236
 237 csphvec_t qpms_eval_vswf(sph_t where,
 238                 complex double *longcoeffs, complex double *mgcoeffs, complex double *elcoeffs,
 239                 qpms_l_t lMax, qpms_bessel_t btyp, qpms_normalisation_t norm);
 240
 241 csphvec_t qpms_eval_vswf_csph(csph_t where,
 242                 complex double *longcoeffs, complex double *mgcoeffs, complex double *elcoeffs,
 243                 qpms_l_t lMax, qpms_bessel_t btyp, qpms_normalisation_t norm);
 244
 245 qpms_vswfset_sph_t *qpms_vswfset_make(qpms_l_t lMax, sph_t kdlj,
 246                 qpms_bessel_t btyp, qpms_normalisation_t norm);//NI
 247 void qpms_vswfset_sph_pfree(qpms_vswfset_sph_t *);//NI
 248
 249 #endif // QPMS_VSWF_H