NAME
PDL::GSL::SF - PDL interface to GSL Special Functions
DESCRIPTION
This is an interface to the Special Function package present in the GNU Scientific Library.
AUTHOR
This file copyright (C) 1999 Christian Pellegrin <chri@infis.univ.trieste.it> All rights reserved. There is no warranty. You are allowed to redistribute this software / documentation under certain conditions. For details, see the file COPYING in the PDL distribution. If this file is separated from the PDL distribution, the copyright notice should be included in the file.
The GSL SF modules were written by G. Jungman.
Calculate all normalized associated Legendre polynomials.
$Plm = gsl_sf_legendre_array($x,'P',4,-1);
The calculation is done for degree 0 <= l <= lmax and order 0 <= m <= l on the range abs(x)<=1.
The parameter norm should be:
- 'S' for Schmidt semi-normalized associated Legendre polynomials S_l^m(x),
- 'Y' for spherical harmonic associated Legendre polynomials Y_l^m(x), or
- 'N' for fully normalized associated Legendre polynomials N_l^m(x).
- 'P' (or any other) for unnormalized associated Legendre polynomials P_l^m(x),
lmax is the maximum degree l. csphase should be (-1) to INCLUDE the Condon-Shortley phase factor (-1)^m, or (+1) to EXCLUDE it.
See "gsl_sf_legendre_array_index" to get the value of l and m in the returned vector. EOD );
pp_def('gsl_sf_legendre_array_index', OtherPars => 'int lmax', Pars => 'int [o]l(n=CALC($COMP(lmax)*($COMP(lmax)+1)/2+$COMP(lmax)+1)); int [o]m(n)', Code => q/ int ell, em, index; for (ell=0; ell<=$COMP(lmax); ell++){ for (em=0; em<=ell; em++){ index = gsl_sf_legendre_array_index(ell,em); $l(n=>index)=ell; $m(n=>index)=em; } }/, Doc =>'=for ref
Calculate the relation between gsl_sf_legendre_arrays index and l and m values.
Note that this function is called differently than the corresponding GSL function, to make it more useful for PDL: here you just input the maximum l (lmax) that was used in gsl_sf_legendre_array and it calculates all l and m values.' );
pp_def('gsl_sf_legendre_sphPlm', GenericTypes => ['D'], OtherPars =>'int l; int m', Pars=>'double x(); double [o]y(); double [o]e()', Code =>' gsl_sf_result r; GSLERR(gsl_sf_legendre_sphPlm_e,($COMP(l),$COMP(m),$x(),&r)) $y() = r.val; $e() = r.err; ', Doc =>'P_lm(x), normalized properly for use in spherical harmonics' );
pp_def('gsl_sf_conicalP_half', GenericTypes => ['D'], OtherPars =>'double lambda', Pars=>'double x(); double [o]y(); double [o]e()', Code =>' gsl_sf_result r; GSLERR(gsl_sf_conicalP_half_e,($COMP(lambda),$x(),&r)) $y() = r.val; $e() = r.err; ', Doc =>'Irregular Spherical Conical Function P^{1/2}_{-1/2 + I lambda}(x)' );
pp_def('gsl_sf_conicalP_mhalf', GenericTypes => ['D'], OtherPars =>'double lambda', Pars=>'double x(); double [o]y(); double [o]e()', Code =>' gsl_sf_result r; GSLERR(gsl_sf_conicalP_mhalf_e,($COMP(lambda),$x(),&r)) $y() = r.val; $e() = r.err; ', Doc =>'Regular Spherical Conical Function P^{-1/2}_{-1/2 + I lambda}(x)' );
pp_def('gsl_sf_conicalP_0', GenericTypes => ['D'], OtherPars =>'double lambda', Pars=>'double x(); double [o]y(); double [o]e()', Code =>' gsl_sf_result r; GSLERR(gsl_sf_conicalP_0_e,($COMP(lambda),$x(),&r)) $y() = r.val; $e() = r.err; ', Doc =>'Conical Function P^{0}_{-1/2 + I lambda}(x)' );
pp_def('gsl_sf_conicalP_1', GenericTypes => ['D'], OtherPars =>'double lambda', Pars=>'double x(); double [o]y(); double [o]e()', Code =>' gsl_sf_result r; GSLERR(gsl_sf_conicalP_1_e,($COMP(lambda),$x(),&r)) $y() = r.val; $e() = r.err; ', Doc =>'Conical Function P^{1}_{-1/2 + I lambda}(x)' );
pp_def('gsl_sf_conicalP_sph_reg', GenericTypes => ['D'], OtherPars =>'int l; double lambda', Pars=>'double x(); double [o]y(); double [o]e()', Code =>' gsl_sf_result r; GSLERR(gsl_sf_conicalP_sph_reg_e,($COMP(l),$COMP(lambda),$x(),&r)) $y() = r.val; $e() = r.err; ', Doc =>'Regular Spherical Conical Function P^{-1/2-l}_{-1/2 + I lambda}(x)' );
pp_def('gsl_sf_conicalP_cyl_reg_e', GenericTypes => ['D'], OtherPars =>'int m; double lambda', Pars=>'double x(); double [o]y(); double [o]e()', Code =>' gsl_sf_result r; GSLERR(gsl_sf_conicalP_cyl_reg_e,($COMP(m),$COMP(lambda),$x(),&r)) $y() = r.val; $e() = r.err; ', Doc =>'Regular Cylindrical Conical Function P^{-m}_{-1/2 + I lambda}(x)' );
pp_def('gsl_sf_legendre_H3d', GenericTypes => ['D'], OtherPars =>'int l; double lambda; double eta', Pars=>'double [o]y(); double [o]e()', Code =>' gsl_sf_result r; GSLERR(gsl_sf_legendre_H3d_e,($COMP(l),$COMP(lambda),$COMP(eta),&r)) $y() = r.val; $e() = r.err; ', Doc =>'lth radial eigenfunction of the Laplacian on the 3-dimensional hyperbolic space.' );
pp_def('gsl_sf_legendre_H3d_array', GenericTypes => ['D'], OtherPars =>'int l=>num; double lambda; double eta', Pars=>'double [o]y(num)', Code =>' GSLERR(gsl_sf_legendre_H3d_array,($COMP(l)-1,$COMP(lambda),$COMP(eta),$P(y))) ', Doc =>'Array of H3d(ell), for l from 0 to n-1.' );
pp_def('gsl_sf_log', GenericTypes => ['D'], Pars=>'double x(); double [o]y(); double [o]e()', Code =>' gsl_sf_result r; GSLERR(gsl_sf_log_e,($x(),&r)) $y() = r.val; $e() = r.err; ', Doc =>'Provide a logarithm function with GSL semantics.' );
pp_def('gsl_sf_complex_log', GenericTypes => ['D'], Pars=>'double zr(); double zi(); double [o]x(); double [o]y(); double [o]xe(); double [o]ye()', Code =>' gsl_sf_result r; gsl_sf_result ri; GSLERR(gsl_sf_complex_log_e,($zr(),$zi(),&r,&ri)) $x() = r.val; $xe() = r.err; $y() = ri.val; $ye() = ri.err; ', Doc =>'Complex Logarithm exp(lnr + I theta) = zr + I zi Returns argument in [-pi,pi].' );
pp_def('gsl_poly_eval', GenericTypes => ['D'], Pars=>'double x(); double c(m); double [o]y()', Code =>' $y() = gsl_poly_eval($P(c),$SIZE(m),$x()); ', Doc =>'c[0] + c[1] x + c[2] x^2 + ... + c[m-1] x^(m-1)' );
pp_def('gsl_sf_pow_int', GenericTypes => ['D'], OtherPars => 'int n', Pars=>'double x(); double [o]y(); double [o]e()', Code =>' gsl_sf_result r; GSLERR(gsl_sf_pow_int_e,($x(),$COMP(n),&r)) $y() = r.val; $e() = r.err; ', Doc =>'Calculate x^n.' );
pp_def('gsl_sf_psi', GenericTypes => ['D'], Pars=>'double x(); double [o]y(); double [o]e()', Code =>' gsl_sf_result r; GSLERR(gsl_sf_psi_e,($x(),&r)) $y() = r.val; $e() = r.err; ', Doc =>'Di-Gamma Function psi(x).' );
pp_def('gsl_sf_psi_1piy', GenericTypes => ['D'], Pars=>'double x(); double [o]y(); double [o]e()', Code =>' gsl_sf_result r; GSLERR(gsl_sf_psi_1piy_e,($x(),&r)) $y() = r.val; $e() = r.err; ', Doc =>'Di-Gamma Function Re[psi(1 + I y)]' );
pp_def('gsl_sf_psi_n', GenericTypes => ['D'], OtherPars => 'int n', Pars=>'double x(); double [o]y(); double [o]e()', Code =>' gsl_sf_result r; GSLERR(gsl_sf_psi_n_e,($COMP(n),$x(),&r)) $y() = r.val; $e() = r.err; ', Doc =>'Poly-Gamma Function psi^(n)(x)' );
pp_def('gsl_sf_synchrotron_1', GenericTypes => ['D'], Pars=>'double x(); double [o]y(); double [o]e()', Code =>' gsl_sf_result r; GSLERR(gsl_sf_synchrotron_1_e,($x(),&r)) $y() = r.val; $e() = r.err; ', Doc =>'First synchrotron function: synchrotron_1(x) = x Integral[ K_{5/3}(t), {t, x, Infinity}]' );
pp_def('gsl_sf_synchrotron_2', GenericTypes => ['D'], Pars=>'double x(); double [o]y(); double [o]e()', Code =>' gsl_sf_result r; GSLERR(gsl_sf_synchrotron_2_e,($x(),&r)) $y() = r.val; $e() = r.err; ', Doc =>'Second synchroton function: synchrotron_2(x) = x * K_{2/3}(x)' );
pp_def('gsl_sf_transport_2', GenericTypes => ['D'], Pars=>'double x(); double [o]y(); double [o]e()', Code =>' gsl_sf_result r; GSLERR(gsl_sf_transport_2_e,($x(),&r)) $y() = r.val; $e() = r.err; ', Doc =>'J(2,x)' );
pp_def('gsl_sf_transport_3', GenericTypes => ['D'], Pars=>'double x(); double [o]y(); double [o]e()', Code =>' gsl_sf_result r; GSLERR(gsl_sf_transport_3_e,($x(),&r)) $y() = r.val; $e() = r.err; ', Doc =>'J(3,x)' );
pp_def('gsl_sf_transport_4', GenericTypes => ['D'], Pars=>'double x(); double [o]y(); double [o]e()', Code =>' gsl_sf_result r; GSLERR(gsl_sf_transport_4_e,($x(),&r)) $y() = r.val; $e() = r.err; ', Doc =>'J(4,x)' );
pp_def('gsl_sf_transport_5', GenericTypes => ['D'], Pars=>'double x(); double [o]y(); double [o]e()', Code =>' gsl_sf_result r; GSLERR(gsl_sf_transport_5_e,($x(),&r)) $y() = r.val; $e() = r.err; ', Doc =>'J(5,x)' );
pp_def('gsl_sf_sin', GenericTypes => ['D'], Pars=>'double x(); double [o]y(); double [o]e()', Code =>' gsl_sf_result r; GSLERR(gsl_sf_sin_e,($x(),&r)) $y() = r.val; $e() = r.err; ', Doc =>'Sin(x) with GSL semantics.' );
pp_def('gsl_sf_cos', GenericTypes => ['D'], Pars=>'double x(); double [o]y(); double [o]e()', Code =>' gsl_sf_result r; GSLERR(gsl_sf_cos_e,($x(),&r)) $y() = r.val; $e() = r.err; ', Doc =>'Cos(x) with GSL semantics.' );
pp_def('gsl_sf_hypot', GenericTypes => ['D'], Pars=>'double x(); double xx(); double [o]y(); double [o]e()', Code =>' gsl_sf_result r; GSLERR(gsl_sf_hypot_e,($x(),$xx(),&r)) $y() = r.val; $e() = r.err; ', Doc =>'Hypot(x,xx) with GSL semantics.' );
pp_def('gsl_sf_complex_sin', GenericTypes => ['D'], Pars=>'double zr(); double zi(); double [o]x(); double [o]y(); double [o]xe(); double [o]ye()', Code =>' gsl_sf_result r; gsl_sf_result ri; GSLERR(gsl_sf_complex_sin_e,($zr(),$zi(),&r,&ri)) $x() = r.val; $xe() = r.err; $y() = ri.val; $ye() = ri.err; ', Doc =>'Sin(z) for complex z' );
pp_def('gsl_sf_complex_cos', GenericTypes => ['D'], Pars=>'double zr(); double zi(); double [o]x(); double [o]y(); double [o]xe(); double [o]ye()', Code =>' gsl_sf_result r; gsl_sf_result ri; GSLERR(gsl_sf_complex_cos_e,($zr(),$zi(),&r,&ri)) $x() = r.val; $xe() = r.err; $y() = ri.val; $ye() = ri.err; ', Doc =>'Cos(z) for complex z' );
pp_def('gsl_sf_complex_logsin', GenericTypes => ['D'], Pars=>'double zr(); double zi(); double [o]x(); double [o]y(); double [o]xe(); double [o]ye()', Code =>' gsl_sf_result r; gsl_sf_result ri; GSLERR(gsl_sf_complex_logsin_e,($zr(),$zi(),&r,&ri)) $x() = r.val; $xe() = r.err; $y() = ri.val; $ye() = ri.err; ', Doc =>'Log(Sin(z)) for complex z' );
pp_def('gsl_sf_lnsinh', GenericTypes => ['D'], Pars=>'double x(); double [o]y(); double [o]e()', Code =>' gsl_sf_result r; GSLERR(gsl_sf_lnsinh_e,($x(),&r)) $y() = r.val; $e() = r.err; ', Doc =>'Log(Sinh(x)) with GSL semantics.' );
pp_def('gsl_sf_lncosh', GenericTypes => ['D'], Pars=>'double x(); double [o]y(); double [o]e()', Code =>' gsl_sf_result r; GSLERR(gsl_sf_lncosh_e,($x(),&r)) $y() = r.val; $e() = r.err; ', Doc =>'Log(Cos(x)) with GSL semantics.' );
pp_def('gsl_sf_polar_to_rect', GenericTypes => ['D'], Pars=>'double r(); double t(); double [o]x(); double [o]y(); double [o]xe(); double [o]ye()', Code =>' gsl_sf_result r; gsl_sf_result ri; GSLERR(gsl_sf_polar_to_rect,($r(),$t(),&r,&ri)) $x() = r.val; $xe() = r.err; $y() = ri.val; $ye() = ri.err; ', Doc =>'Convert polar to rectlinear coordinates.' );
pp_def('gsl_sf_rect_to_polar', GenericTypes => ['D'], Pars=>'double x(); double y(); double [o]r(); double [o]t(); double [o]re(); double [o]te()', Code =>' gsl_sf_result r; gsl_sf_result ri; GSLERR(gsl_sf_rect_to_polar,($x(),$y(),&r,&ri)) $r() = r.val; $re() = r.err; $t() = ri.val; $te() = ri.err; ', Doc =>'Convert rectlinear to polar coordinates. return argument in range [-pi, pi].' );
pp_def('gsl_sf_angle_restrict_symm', GenericTypes => ['D'], Pars=>'double [o]y();', Code =>' GSLERR(gsl_sf_angle_restrict_symm_e,($P(y))) ', Doc =>'Force an angle to lie in the range (-pi,pi].' );
pp_def('gsl_sf_angle_restrict_pos', GenericTypes => ['D'], Pars=>'double [o]y();', Code =>' GSLERR(gsl_sf_angle_restrict_pos_e,($P(y))) ', Doc =>'Force an angle to lie in the range [0,2 pi).' );
pp_def('gsl_sf_sin_err', GenericTypes => ['D'], Pars=>'double x(); double dx(); double [o]y(); double [o]e()', Code =>' gsl_sf_result r; GSLERR(gsl_sf_sin_err_e,($x(),$dx(),&r)) $y() = r.val; $e() = r.err; ', Doc =>'Sin(x) for quantity with an associated error.' );
pp_def('gsl_sf_cos_err', GenericTypes => ['D'], Pars=>'double x(); double dx(); double [o]y(); double [o]e()', Code =>' gsl_sf_result r; GSLERR(gsl_sf_cos_err_e,($x(),$dx(),&r)) $y() = r.val; $e() = r.err; ', Doc =>'Cos(x) for quantity with an associated error.' );
pp_def('gsl_sf_zeta', GenericTypes => ['D'], Pars=>'double x(); double [o]y(); double [o]e()', Code =>' gsl_sf_result r; GSLERR(gsl_sf_zeta_e,($x(),&r)) $y() = r.val; $e() = r.err; ', Doc =>'Riemann Zeta Function zeta(x) = Sum[ k^(-s), {k,1,Infinity} ], s != 1.0' );
pp_def('gsl_sf_hzeta', GenericTypes => ['D'], OtherPars =>'double q', Pars=>'double s(); double [o]y(); double [o]e()', Code =>' gsl_sf_result r; GSLERR(gsl_sf_hzeta_e,($s(), $COMP(q), &r)) $y() = r.val; $e() = r.err; ', Doc =>'Hurwicz Zeta Function zeta(s,q) = Sum[ (k+q)^(-s), {k,0,Infinity} ]' );
pp_def('gsl_sf_eta', GenericTypes => ['D'], Pars=>'double x(); double [o]y(); double [o]e()', Code =>' gsl_sf_result r; GSLERR(gsl_sf_eta_e,($x(),&r)) $y() = r.val; $e() = r.err; ', Doc =>'Eta Function eta(s) = (1-2^(1-s)) zeta(s)' );
pp_done();