QPMS
Electromagnetic multiple scattering library and toolkit.
qpms_specfunc.h
Go to the documentation of this file.
1 
4 #ifndef QPMS_SPECFUNC_H
5 #define QPMS_SPECFUNC_H
6 #include "qpms_types.h"
7 #include <gsl/gsl_sf_legendre.h>
8 
9 /******************************************************************************
10  * Spherical Bessel functions *
11  ******************************************************************************/
12 
13 // TODO unify types
14 qpms_errno_t qpms_sph_bessel_fill(qpms_bessel_t typ, qpms_l_t lmax, complex double x, complex double *result_array);
15 
16 
17 
18 typedef struct {
19  qpms_l_t lMax;
20  double *akn; // coefficients as in DLMF 10.49.1
21  //complex double *bkn; // coefficients of the derivatives
23 
24 qpms_sbessel_calculator_t *qpms_sbessel_calculator_init(void);
25 void qpms_sbessel_calculator_pfree(qpms_sbessel_calculator_t *c);
26 
27 qpms_errno_t qpms_sbessel_calc_fill(qpms_sbessel_calculator_t *c, qpms_bessel_t typ, qpms_l_t lmax,
28  double x, complex double *result_array);
29 
30 complex double qpms_sbessel_calc_h1(qpms_sbessel_calculator_t *c, qpms_l_t n, complex double x);
31 qpms_errno_t qpms_sbessel_calc_h1_fill(qpms_sbessel_calculator_t *c, qpms_l_t lmax,
32  complex double x, complex double *result_array);
33 
34 
35 /******************************************************************************
36  * Legendre functions and their "angular derivatives" *
37  ******************************************************************************/
38 
39 /*
40  * N.B. for the norm definitions, see
41  * https://www.gnu.org/software/gsl/manual/html_node/Associated-Legendre-Polynomials-and-Spherical-Harmonics.html
42  * ( gsl/specfunc/legendre_source.c and 7.24.2 of gsl docs
43  */
44 
45 qpms_errno_t qpms_legendre_deriv_y_get(double **result, double **result_deriv, double x, qpms_l_t lMax,
46  gsl_sf_legendre_t lnorm, double csphase); // free() result and result_deriv yourself!
47 qpms_errno_t qpms_legendre_deriv_y_fill(double *where, double *where_deriv, double x,
48  qpms_l_t lMax, gsl_sf_legendre_t lnorm, double csphase);
49 
50 
51 double *qpms_legendre_y_get(double x, qpms_l_t lMax, qpms_normalisation_t norm);//NI
52 double *qpms_legendre0d_y_get(qpms_l_t lMax, qpms_normalisation_t norm); //NI
53 double *qpms_legendre_plus1d_y_get(qpms_l_t lMax, qpms_normalisation_t norm); //NI
54 double *qpms_legendre_minus1d_y_get(qpms_l_t lMax, qpms_normalisation_t norm); //NI
55 
56 
57 
59 
64 typedef struct {
65  //qpms_normalisation_t norm;
66  qpms_l_t lMax;
67  //qpms_y_t nelem;
68  double *leg, *pi, *tau;
69 } qpms_pitau_t;
70 
72 
90 qpms_pitau_t qpms_pitau_get(double theta, qpms_l_t lMax, double csphase);
91 
93 
97 qpms_errno_t qpms_pitau_fill(double *target_leg, double *target_pi, double *target_tau,
98  double theta, qpms_l_t lMax, double csphase);
99 
102 
103 //void qpms_pitau_pfree(qpms_pitau_t*);//NI
104 
105 // Associated Legendre polynomial at zero argument (DLMF 14.5.1) DEPRECATED?
106 double qpms_legendre0(int m, int n);
107 // Associated Legendre polynomial derivative at zero argument (DLMF 14.5.2)
108 double qpms_legendred0(int m, int n);
109 
110 #endif // QPMS_SPECFUNC_H
void qpms_pitau_free(qpms_pitau_t)
Frees the dynamically allocated arrays from qpms_pitau_t.
Definition: legendre.c:119
qpms_errno_t qpms_pitau_fill(double *target_leg, double *target_pi, double *target_tau, double theta, qpms_l_t lMax, double csphase)
Directly fills (pre-allocated) arrays of normalised Legendre and auxillary functions.
Definition: legendre.c:65
qpms_pitau_t qpms_pitau_get(double theta, qpms_l_t lMax, double csphase)
Returns an array of normalised Legendre and auxillary functions.
Definition: legendre.c:54
Common qpms types.
qpms_normalisation_t
Vector spherical wavefuction normalisation and phase convention codes.
Definition: qpms_types.h:104
int qpms_l_t
Type for spherical harmonic degree l.
Definition: qpms_types.h:27
qpms_errno_t
Error codes / return values for certain numerical functions.
Definition: qpms_types.h:85
qpms_bessel_t
Bessel function kinds.
Definition: qpms_types.h:176
Array of Legendre and and auxillary functions.
Definition: qpms_specfunc.h:64
Definition: qpms_specfunc.h:18