|
template<class T > |
T | special::integral::trapz (double h, NumberArray< T > const &y) |
| Uniform trapezoidal integration. More...
|
|
template<class T > |
T | special::integral::trapz (NumberArray< T > const &x, NumberArray< T > const &y) |
| Non-uniform trapezoidal integration. More...
|
|
template<class T > |
T | special::integral::simpson (double h, NumberArray< T > const &y) |
| Uniform Simpson integration. More...
|
|
template<class T > |
T | special::integral::pow_exp_hyperg1F1 (T a, T b, T c, T u, T v, T x, double epsrel=1e-10, unsigned maxiter=1000) |
| Compute integral of the confluent hypergeometric function. More...
|
|
template<class T > |
NumberArray< T > | special::integral::romberg (const ArrayView< T > y) |
| Romberg integration. More...
|
|
int | special::coulomb_zeros (double eta, int L, int nzeros, double *zeros, double epsrel=1e-8) |
| Get zeros of the Coulomb wave function \( F_L(-1/k,kr) \). More...
|
|
std::vector< std::vector< int > > | special::FdB_partition (int n) |
| Faa di Bruno partitioning. More...
|
|
template<class T , class OuterFunctionDerivative , class InnerFunctionDerivative > |
T | special::chain_rule (OuterFunctionDerivative Df, InnerFunctionDerivative Dg, int n, double x) |
| Chain rule for n-th derivative. More...
|
|
Complex | hydro_P (unsigned n, unsigned l, Complex z) |
|
Complex | dhydro_P (unsigned n, unsigned l, Complex z) |
|
Complex | ric_j (int n, Complex z) |
|
cArray | ric_jv (int lmax, Complex z) |
|
Complex | dric_j (int n, Complex z) |
|
cArray | dric_jv (int lmax, Complex z) |
|
double | ric_j (int n, double x) |
| Ricatti-Bessel function of the first kind, \( \hat{j}_n(x) \). More...
|
|
double | ric_n (int n, double x) |
| Ricatti-Bessel function of the second kind, \( \hat{y}_n(x) \). More...
|
|
double | sph_i_scaled (int n, double x) |
| Scaled modified spherical Bessel function of the first kind, \( \mathrm{e}^{-x} i_n(x) \). More...
|
|
double | ric_i_scaled (int n, double x) |
| Scaled modified Ricatti-Bessel function of the first kind, \( \mathrm{e}^{-x} \hat{i}_n(x) \). More...
|
|
double | sph_i (int n, double x) |
| Modified spherical Bessel function of the first kind, \( i_n(x) \). More...
|
|
double | ric_i (int n, double x) |
| Modified Ricatti-Bessel function of the first kind, \( \hat{i}_n(x) \). More...
|
|
double | sph_k_scaled (int n, double x) |
| Scaled modified spherical Bessel function of the second kind, \( \mathrm{e}^x k_n(x) \). More...
|
|
double | ric_k_scaled (int n, double x) |
| Scaled modified Ricatti-Bessel function of the second kind, \( \mathrm{e}^x \hat{k}_n(x) \). More...
|
|
double | sph_k (int n, double x) |
| Modified spherical Bessel function of the second kind, \( k_n(x) \). More...
|
|
double | ric_k (int n, double x) |
| Modified Ricatti-Bessel function of the second kind, \( \hat{k}_n(x) \). More...
|
|
Complex | ric_h_plus (int n, double x) |
| Ricatti-Hankel function of the first kind, \( \hat{h}_n^{(+)}(x) \). More...
|
|
Complex | sphY (int l, int m, double theta, double phi) |
| Spherical harmonic function. More...
|
|
Complex | sphBiY (int l1, int l2, int L, int M, double theta1, double phi1, double theta2, double phi2) |
| Bi-polar spherical harmonic function. More...
|
|
double | dric_j (int n, double x) |
| Derivative of Ricatti-Bessel function of the first kind, \( \hat{j}_n'(x) \). More...
|
|
double | dric_n (int n, double x) |
| Derivative of Ricatti-Bessel function of the second kind, \( \hat{y}_n'(x) \). More...
|
|
double | dric_i (int n, double x) |
| Derivative of the modified Ricatti-Bessel function of the second kind, \( \hat{i}_n'(x) \). More...
|
|
double | dric_k (int n, double x) |
| Derivative of the modified Ricatti-Bessel function of the second kind, \( \hat{k}_n'(x) \). More...
|
|
double | dric_i_scaled (int n, double x) |
| Derivative of scaled modified Ricatti-Bessel function of the second kind, \( (\mathrm{e}^{-x} \hat{i}_n(x))' \). More...
|
|
double | dric_k_scaled (int n, double x) |
| Scaled derivative of modified Ricatti-Bessel function of the second kind, \( (\mathrm{e}^x \hat{k}_n(x)) \). More...
|
|
Complex | dric_h_plus (int n, double x) |
| Derivative of Ricatti-Hankel function of the first kind, \( {\hat{h}_n^{(+)}}'(x) \). More...
|
|
int | coul_F_michel (int l, double k, double r, double &F, double &Fp) |
| Uniform approximation to the Coulomb wave function. More...
|
|
int | coul_F (int l, double k, double r, double &F, double &Fp) |
| Evaluate Coulomb wave function (and its derivative). More...
|
|
double | coul_F_asy (int l, double k, double r, double sigma=special::constant::Nan) |
| Asymptotic form of the regular Coulomb wave. More...
|
|
double | coul_F_sigma (int l, double k) |
| Coulomb phase shift. More...
|
|
bool | makes_triangle (int two_j1, int two_j2, int two_j3) |
| Check triangle inequality. More...
|
|
double | logdelta (int two_j1, int two_j2, int two_j3) |
| Auxiliary function used in coupling coefficient computations. More...
|
|
double | computef (int lambda, int l1, int l2, int l1p, int l2p, int L) |
|
double | ClebschGordan (int l1, int m1, int l2, int m2, int L, int M) |
| Clebsch-Gordan coefficient. More...
|
|
double | Gaunt (int l1, int m1, int l2, int m2, int l, int m) |
| Gaunt's integral. More...
|
|
int | triangle_count (int L, int maxell) |
|
double | Hyper2F1 (double a, double b, double c, double x) |
| The hypergeometric function \( {}_2F_1 \). More...
|
|
|
const double | special::constant::e = M_E |
|
const double | special::constant::pi = M_PI |
|
const double | special::constant::pi_half = M_PI_2 |
|
const double | special::constant::pi_quart = M_PI_4 |
|
const double | special::constant::pi_sqrt = M_SQRTPI |
|
const double | special::constant::pi_inv = M_1_PI |
|
const double | special::constant::pi_two_inv = M_2_PI |
|
const double | special::constant::pi_two_invsqrt = M_2_SQRTPI |
|
const double | special::constant::sqrt_two = M_SQRT2 |
|
const double | special::constant::sqrt_three = M_SQRT3 |
|
const double | special::constant::sqrt_half = M_SQRT1_2 |
|
const double | special::constant::log2e = M_LOG2E |
|
const double | special::constant::log10e = M_LOG10E |
|
const double | special::constant::ln2 = M_LN2 |
|
const double | special::constant::ln10 = M_LN10 |
|
const double | special::constant::lnpi = M_LNPI |
|
const double | special::constant::euler = M_EULER |
|
const double | special::constant::Inf = std::numeric_limits<double>::infinity() |
|
const double | special::constant::Nan = std::numeric_limits<double>::quiet_NaN() |
|
class RadialFunction | exception |
|
#define | Wigner3j(a, b, c, d, e, f) Wigner3j_2(2*(a),2*(b),2*(c),2*(d),2*(e),2*(f)) |
|
double | Wigner3j_2 (int two_j1, int two_j2, int two_j3, int two_m1, int two_m2, int two_m3) |
| Wigner 3j coefficient. More...
|
|
#define | Wigner6j(a, b, c, d, e, f) Wigner6j_2(2*(a),2*(b),2*(c),2*(d),2*(e),2*(f)) |
|
double | Wigner6j_2 (int two_j1, int two_j2, int two_j3, int two_j4, int two_j5, int two_j6) |
| Wigner 6j coefficient. More...
|
|
double Wigner3j_2 |
( |
int |
two_j1, |
|
|
int |
two_j2, |
|
|
int |
two_j3, |
|
|
int |
two_m1, |
|
|
int |
two_m2, |
|
|
int |
two_m3 |
|
) |
| |
Compute the Wigner 3j coefficient. Even though there is a routine "gsl_sf_coupling_3j" in GSL library, it has to be implemented anew, because of factorial overflows. This routine computes all factorials only in logarithms using the standard function "lgamma". The formula is taken from Edmonds, A. R.: Angular momentum in quantum mechanics, Princeton 1968.
\[ \left( \matrix{j_1 7 j_2 & j_3 \cr m_1 & m_2 & m_3} \right) = \epsilon(j_1,j_2,j_3)\Delta(j_1,j_2,j_3) \delta_{m_1+m_2+m_3}^0 (-1)^{j_1-j_2-m_3} \sqrt{(j_1+m_1)! (j_1-m_1)! (j_2+m_2)! (j_2-m_2)! (j_3+m_3)! (j_3-m_3)!} \sum_k \frac{(-1)^2}{k! (j_1+j_2-j_3-k)! (j_1-m_1-k)! (j_2+m_2-k)! (j_3-j_2+m_1+k)! (j_3-j_1-m_2+k)!} \ , \]
\[ \epsilon(j_1,j_2,j_3) = \cases{1 & triangle inequality satisfied \cr 0 & otherwise} \ . \]
See logdelta for definition of the triangle function \( \Delta(a,b,c) \). Note that the arguments \( j_1, j_2, j_3, m_1, m_2, m_3 \) need to be supplied doubled, as \( 2j_1, 2j_2, 2j_3, 2m_1, 2m_2, 2m_3 \) so that the parameters can be considered integral even though the half-integral angular momentum is allowed.
double Wigner6j_2 |
( |
int |
two_j1, |
|
|
int |
two_j2, |
|
|
int |
two_j3, |
|
|
int |
two_j4, |
|
|
int |
two_j5, |
|
|
int |
two_j6 |
|
) |
| |
Compute the Wigner 6j coefficient. Even though there is a routine "gsl_sf_coupling_6j" in GSL library, it has to be implemented anew, because of factorial overflows. This routine computes all factorials only in logarithms using the standard function "lgamma". The formula is taken from Edmonds, A. R.: Angular momentum in quantum mechanics, Princeton 1968.
\[ \left\{ \matrix{j_1 & j_2 & j_3 \cr j_4 & j_5 & j_6} \right\} = \Delta(j_1,j_2,j_3) \Delta(j_3,j_4,j_5) \Delta(j_1,j_5,j_6) \Delta(j_2,j_4,j_6) \sum_k \frac{(-1)^k (k+1)!}{(k-j_1-j_2-j_3)! (k-j_1-j_5-j_6)! (k-j_2-j_4-j_6)! (k-j_3-j_4-j_5)! (j_1+j_2+j_4+j_5-k)! (j_1+j_3+j_4+j_6-k)! (j_2+j_3+j_5+j_6-k)!} \ . \]
See logdelta for definition of the triangle function \( \Delta(a,b,c) \). Note that the arguments \( j_1, j_2, j_3, j_4, j_5, j_6 \) need to be supplied doubled, as \( 2j_1, 2j_2, 2j_3, 2j_4, 2j_5, 2j_6 \) so that the parameters can be considered integral even though the half-integral angular momentum is allowed.