1 %module
"Math::GSL::SF"
4 %apply double
*OUTPUT { double
* sn
, double
* cn
, double
* dn
, double
* sgn
};
7 #include
"gsl/gsl_types.h"
8 #include
"gsl/gsl_version.h"
9 #include
"gsl/gsl_mode.h"
10 #include
"gsl/gsl_sf.h"
11 #include
"gsl/gsl_sf_airy.h"
12 #include
"gsl/gsl_sf_bessel.h"
13 #include
"gsl/gsl_sf_clausen.h"
14 #include
"gsl/gsl_sf_coulomb.h"
15 #include
"gsl/gsl_sf_coupling.h"
16 #include
"gsl/gsl_sf_dawson.h"
17 #include
"gsl/gsl_sf_debye.h"
18 #include
"gsl/gsl_sf_dilog.h"
19 #include
"gsl/gsl_sf_elementary.h"
20 #include
"gsl/gsl_sf_ellint.h"
21 #include
"gsl/gsl_sf_elljac.h"
22 #include
"gsl/gsl_sf_erf.h"
23 #include
"gsl/gsl_sf_exp.h"
24 #include
"gsl/gsl_sf_expint.h"
25 #include
"gsl/gsl_sf_fermi_dirac.h"
26 #include
"gsl/gsl_sf_gamma.h"
27 #include
"gsl/gsl_sf_gegenbauer.h"
28 #include
"gsl/gsl_sf_hyperg.h"
29 #include
"gsl/gsl_sf_laguerre.h"
30 #include
"gsl/gsl_sf_lambert.h"
31 #include
"gsl/gsl_sf_legendre.h"
32 #include
"gsl/gsl_sf_log.h"
33 #ifdef GSL_VERSION
&& GSL_VERSION == "1.11"
34 #include
"gsl/gsl_sf_mathieu.h"
36 #include
"gsl/gsl_sf_pow_int.h"
37 #include
"gsl/gsl_sf_psi.h"
38 #include
"gsl/gsl_sf_result.h"
39 #include
"gsl/gsl_sf_synchrotron.h"
40 #include
"gsl/gsl_sf_transport.h"
41 #include
"gsl/gsl_sf_trig.h"
42 #include
"gsl/gsl_sf_zeta.h"
44 %include
"gsl/gsl_types.h"
45 %include
"gsl/gsl_version.h"
46 %include
"gsl/gsl_mode.h"
47 %include
"gsl/gsl_sf.h"
48 %include
"gsl/gsl_sf_airy.h"
49 %include
"gsl/gsl_sf_bessel.h"
50 %include
"gsl/gsl_sf_clausen.h"
51 %include
"gsl/gsl_sf_coulomb.h"
52 %include
"gsl/gsl_sf_coupling.h"
53 %include
"gsl/gsl_sf_dawson.h"
54 %include
"gsl/gsl_sf_debye.h"
55 %include
"gsl/gsl_sf_dilog.h"
56 %include
"gsl/gsl_sf_elementary.h"
57 %include
"gsl/gsl_sf_ellint.h"
58 %include
"gsl/gsl_sf_elljac.h"
59 %include
"gsl/gsl_sf_erf.h"
60 %include
"gsl/gsl_sf_exp.h"
61 %include
"gsl/gsl_sf_expint.h"
62 %include
"gsl/gsl_sf_fermi_dirac.h"
63 %include
"gsl/gsl_sf_gamma.h"
64 %include
"gsl/gsl_sf_gegenbauer.h"
65 %include
"gsl/gsl_sf_hyperg.h"
66 %include
"gsl/gsl_sf_laguerre.h"
67 %include
"gsl/gsl_sf_lambert.h"
68 %include
"gsl/gsl_sf_legendre.h"
69 %include
"gsl/gsl_sf_log.h"
70 #ifdef GSL_VERSION
&& GSL_VERSION == '1.11'
71 %include
"gsl/gsl_sf_mathieu.h"
73 %include
"gsl/gsl_sf_pow_int.h"
74 %include
"gsl/gsl_sf_psi.h"
75 %include
"gsl/gsl_sf_result.h"
76 %include
"gsl/gsl_sf_synchrotron.h"
77 %include
"gsl/gsl_sf_transport.h"
78 %include
"gsl/gsl_sf_trig.h"
79 %include
"gsl/gsl_sf_zeta.h"
89 gsl_sf_airy_Ai_scaled_e
91 gsl_sf_airy_Bi_scaled_e
93 gsl_sf_airy_Ai_deriv_e
95 gsl_sf_airy_Bi_deriv_e
97 gsl_sf_airy_Ai_deriv_scaled_e
98 gsl_sf_airy_Ai_deriv_scaled
99 gsl_sf_airy_Bi_deriv_scaled_e
100 gsl_sf_airy_Bi_deriv_scaled
101 gsl_sf_airy_zero_Ai_e
103 gsl_sf_airy_zero_Bi_e
105 gsl_sf_airy_zero_Ai_deriv_e
106 gsl_sf_airy_zero_Ai_deriv
107 gsl_sf_airy_zero_Bi_deriv_e
108 gsl_sf_airy_zero_Bi_deriv
117 gsl_sf_bessel_Jn_array
124 gsl_sf_bessel_Yn_array
131 gsl_sf_bessel_In_array
132 gsl_sf_bessel_I0_scaled_e
133 gsl_sf_bessel_I0_scaled
134 gsl_sf_bessel_I1_scaled_e
135 gsl_sf_bessel_I1_scaled
136 gsl_sf_bessel_In_scaled_e
137 gsl_sf_bessel_In_scaled
138 gsl_sf_bessel_In_scaled_array
145 gsl_sf_bessel_Kn_array
146 gsl_sf_bessel_K0_scaled_e
147 gsl_sf_bessel_K0_scaled
148 gsl_sf_bessel_K1_scaled_e
149 gsl_sf_bessel_K1_scaled
150 gsl_sf_bessel_Kn_scaled_e
151 gsl_sf_bessel_Kn_scaled
152 gsl_sf_bessel_Kn_scaled_array
161 gsl_sf_bessel_jl_array
162 gsl_sf_bessel_jl_steed_array
171 gsl_sf_bessel_yl_array
172 gsl_sf_bessel_i0_scaled_e
173 gsl_sf_bessel_i0_scaled
174 gsl_sf_bessel_i1_scaled_e
175 gsl_sf_bessel_i1_scaled
176 gsl_sf_bessel_i2_scaled_e
177 gsl_sf_bessel_i2_scaled
178 gsl_sf_bessel_il_scaled_e
179 gsl_sf_bessel_il_scaled
180 gsl_sf_bessel_il_scaled_array
181 gsl_sf_bessel_k0_scaled_e
182 gsl_sf_bessel_k0_scaled
183 gsl_sf_bessel_k1_scaled_e
184 gsl_sf_bessel_k1_scaled
185 gsl_sf_bessel_k2_scaled_e
186 gsl_sf_bessel_k2_scaled
187 gsl_sf_bessel_kl_scaled_e
188 gsl_sf_bessel_kl_scaled
189 gsl_sf_bessel_kl_scaled_array
194 gsl_sf_bessel_sequence_Jnu_e
195 gsl_sf_bessel_Inu_scaled_e
196 gsl_sf_bessel_Inu_scaled
199 gsl_sf_bessel_Knu_scaled_e
200 gsl_sf_bessel_Knu_scaled
203 gsl_sf_bessel_lnKnu_e
205 gsl_sf_bessel_zero_J0_e
206 gsl_sf_bessel_zero_J0
207 gsl_sf_bessel_zero_J1_e
208 gsl_sf_bessel_zero_J1
209 gsl_sf_bessel_zero_Jnu_e
210 gsl_sf_bessel_zero_Jnu
212 @EXPORT_clausen
= qw
/
216 @EXPORT_hydrogenic
= qw
/
217 gsl_sf_hydrogenicR_1_e
222 @EXPORT_coulumb
= qw
/
223 gsl_sf_coulomb_wave_FG_e
224 gsl_sf_coulomb_wave_F_array
225 gsl_sf_coulomb_wave_FG_array
226 gsl_sf_coulomb_wave_FGp_array
227 gsl_sf_coulomb_wave_sphF_array
229 gsl_sf_coulomb_CL_array
231 @EXPORT_coupling
= qw
/
236 gsl_sf_coupling_RacahW_e
237 gsl_sf_coupling_RacahW
240 gsl_sf_coupling_6j_INCORRECT_e
241 gsl_sf_coupling_6j_INCORRECT
264 gsl_sf_complex_dilog_xy_e
265 gsl_sf_complex_dilog_e
269 gsl_sf_complex_spence_xy_e
272 gsl_sf_multiply_err_e
274 @EXPORT_elliptic
= qw
/
275 gsl_sf_ellint_Kcomp_e
277 gsl_sf_ellint_Ecomp_e
279 gsl_sf_ellint_Pcomp_e
281 gsl_sf_ellint_Dcomp_e
315 push @EXPORT_misc
, qw
/
321 gsl_sf_exp_mult_e10_e
332 gsl_sf_exp_mult_err_e
333 gsl_sf_exp_mult_err_e10_e
340 gsl_sf_expint_E1_scaled_e
341 gsl_sf_expint_E1_scaled
342 gsl_sf_expint_E2_scaled_e
343 gsl_sf_expint_E2_scaled
344 gsl_sf_expint_En_scaled_e
345 gsl_sf_expint_En_scaled
348 gsl_sf_expint_Ei_scaled_e
349 gsl_sf_expint_Ei_scaled
361 @EXPORT_fermi_dirac
= qw
/
362 gsl_sf_fermi_dirac_m1_e
363 gsl_sf_fermi_dirac_m1
364 gsl_sf_fermi_dirac_0_e
366 gsl_sf_fermi_dirac_1_e
368 gsl_sf_fermi_dirac_2_e
370 gsl_sf_fermi_dirac_int_e
371 gsl_sf_fermi_dirac_int
372 gsl_sf_fermi_dirac_mhalf_e
373 gsl_sf_fermi_dirac_mhalf
374 gsl_sf_fermi_dirac_half_e
375 gsl_sf_fermi_dirac_half
376 gsl_sf_fermi_dirac_3half_e
377 gsl_sf_fermi_dirac_3half
378 gsl_sf_fermi_dirac_inc_0_e
379 gsl_sf_fermi_dirac_inc_0
381 @EXPORT_legendre
= qw
/
384 gsl_sf_legendre_Pl_array
385 gsl_sf_legendre_Pl_deriv_array
398 gsl_sf_legendre_Plm_e
400 gsl_sf_legendre_Plm_array
401 gsl_sf_legendre_Plm_deriv_array
402 gsl_sf_legendre_sphPlm_e
403 gsl_sf_legendre_sphPlm
404 gsl_sf_legendre_sphPlm_array
405 gsl_sf_legendre_sphPlm_deriv_array
406 gsl_sf_legendre_array_size
407 gsl_sf_legendre_H3d_0_e
408 gsl_sf_legendre_H3d_0
409 gsl_sf_legendre_H3d_1_e
410 gsl_sf_legendre_H3d_1
411 gsl_sf_legendre_H3d_e
413 gsl_sf_legendre_H3d_array
425 gsl_sf_lngamma_complex_e
433 @EXPORT_factorial
= qw
/
440 gsl_sf_lndoublefact_e
443 @EXPORT_hypergeometric
= qw
/
446 gsl_sf_hyperg_1F1_int_e
447 gsl_sf_hyperg_1F1_int
450 gsl_sf_hyperg_U_int_e
452 gsl_sf_hyperg_U_int_e10_e
455 gsl_sf_hyperg_U_e10_e
458 gsl_sf_hyperg_2F1_conj_e
459 gsl_sf_hyperg_2F1_conj
460 gsl_sf_hyperg_2F1_renorm_e
461 gsl_sf_hyperg_2F1_renorm
462 gsl_sf_hyperg_2F1_conj_renorm_e
463 gsl_sf_hyperg_2F1_conj_renorm
467 @EXPORT_laguerre
= qw
/
477 push @EXPORT_misc
, qw
/
506 gsl_sf_gegenpoly_array
511 gsl_sf_conicalP_half_e
513 gsl_sf_conicalP_mhalf_e
514 gsl_sf_conicalP_mhalf
519 gsl_sf_conicalP_sph_reg_e
520 gsl_sf_conicalP_sph_reg
521 gsl_sf_conicalP_cyl_reg_e
522 gsl_sf_conicalP_cyl_reg
530 gsl_sf_log_1plusx_mx_e
547 gsl_sf_result_smash_e
548 gsl_sf_synchrotron_1_e
550 gsl_sf_synchrotron_2_e
553 @EXPORT_mathieu
= qw
/
554 gsl_sf_mathieu_a_array
555 gsl_sf_mathieu_b_array
558 gsl_sf_mathieu_a_coeff
559 gsl_sf_mathieu_b_coeff
564 gsl_sf_mathieu_ce_array
565 gsl_sf_mathieu_se_array
568 gsl_sf_mathieu_Mc_array
569 gsl_sf_mathieu_Ms_array
571 @EXPORT_transport
= qw
/
592 gsl_sf_complex_logsin_e
603 gsl_sf_angle_restrict_symm_e
604 gsl_sf_angle_restrict_symm
605 gsl_sf_angle_restrict_pos_e
606 gsl_sf_angle_restrict_pos
607 gsl_sf_angle_restrict_symm_err_e
608 gsl_sf_angle_restrict_pos_err_e
633 GSL_SF_DOUBLEFACT_NMAX
638 @EXPORT_airy
, @EXPORT_bessel
, @EXPORT_clausen
, @EXPORT_hydrogenic
,
639 @EXPORT_coulumb
, @EXPORT_coupling
, @EXPORT_dawson
, @EXPORT_debye
,
640 @EXPORT_dilog
, @EXPORT_misc
, @EXPORT_elliptic
, @EXPORT_error
, @EXPORT_legendre
,
641 @EXPORT_gamma
, @EXPORT_transport
, @EXPORT_trig
, @EXPORT_zeta
, @EXPORT_eta
,
646 all
=> [ @EXPORT_OK
],
647 airy
=> [ @EXPORT_airy
],
648 bessel
=> [ @EXPORT_bessel
],
649 clausen
=> [ @EXPORT_clausen
],
650 coulumb
=> [ @EXPORT_coulumb
],
651 coupling
=> [ @EXPORT_coupling
],
652 dawson
=> [ @EXPORT_dawson
],
653 debye
=> [ @EXPORT_debye
],
654 dilog
=> [ @EXPORT_dilog
],
655 eta
=> [ @EXPORT_eta
],
656 elliptic
=> [ @EXPORT_elliptic
],
657 error
=> [ @EXPORT_error
],
658 factorial
=> [ @EXPORT_factorial
],
659 gamma
=> [ @EXPORT_gamma
],
660 hydrogenic
=> [ @EXPORT_hydrogenic
],
661 hypergeometric
=> [ @EXPORT_hypergeometric
],
662 laguerre
=> [ @EXPORT_laguerre
],
663 legendre
=> [ @EXPORT_legendre
],
664 mathieu
=> [ @EXPORT_mathieu
],
665 misc
=> [ @EXPORT_misc
],
666 transport
=> [ @EXPORT_transport
],
667 trig
=> [ @EXPORT_trig
],
668 vars
=> [ @EXPORT_vars
],
669 zeta
=> [ @EXPORT_zeta
],
676 Math
::GSL
::SF
- Special Functions
680 use Math
::GSL
::SF qw
/:all
/;
684 This module contains a data structure named gsl_sf_result. To create a new one use
686 $r
= Math
::GSL
::SF
::gsl_sf_result_struct-
>new
;
688 You can then access the elements of the structure in this way
:
692 my $error
= $r-
>{err
};
694 Here is a list of all included functions
:
698 =item C
<gsl_sf_airy_Ai_e
($x
, $mode
)>
700 =item C
<gsl_sf_airy_Ai
($x
, $mode
, $result
)>
702 - These routines compute the Airy function Ai
($x
) with an accuracy specified by $mode. $mode should be $GSL_PREC_DOUBLE
, $GSL_PREC_SINGLE or $GSL_PREC_APPROX. $result is a gsl_sf_result structure.
708 =item C
<gsl_sf_airy_Bi_e
($x
, $mode
, $result
)>
710 =item C
<gsl_sf_airy_Bi
($x
, $mode
)>
712 - These routines compute the Airy function Bi
($x
) with an accuracy specified by $mode. $mode should be $GSL_PREC_DOUBLE
, $GSL_PREC_SINGLE or $GSL_PREC_APPROX. $result is a gsl_sf_result structure.
718 =item C
<gsl_sf_airy_Ai_scaled_e
($x
, $mode
, $result
)>
720 =item C
<gsl_sf_airy_Ai_scaled
($x
, $mode
)>
722 - These routines compute a scaled version of the Airy function S_A
($x
) Ai
($x
). For $x
>0 the scaling factor S_A
($x
) is \exp
(+(2/3) $x
**(3/2)), and is
1 for $x
<0. $result is a gsl_sf_result structure.
728 =item C
<gsl_sf_airy_Bi_scaled_e
($x
, $mode
, $result
)>
730 =item C
<gsl_sf_airy_Bi_scaled
($x
, $mode
)>
732 - These routines compute a scaled version of the Airy function S_B
($x
) Bi
($x
). For $x
>0 the scaling factor S_B
($x
) is exp
(-(2/3) $x
**(3/2)), and is
1 for $x
<0. $result is a gsl_sf_result structure.
738 =item C
<gsl_sf_airy_Ai_deriv_e
($x
, $mode
, $result
)>
740 =item C
<gsl_sf_airy_Ai_deriv
($x
, $mode
)>
742 - These routines compute the Airy function derivative Ai'
($x
) with an accuracy specified by $mode. $result is a gsl_sf_result structure.
748 =item C
<gsl_sf_airy_Bi_deriv_e
($x
, $mode
, $result
)>
750 =item C
<gsl_sf_airy_Bi_deriv
($x
, $mode
)>
752 -These routines compute the Airy function derivative Bi'
($x
) with an accuracy specified by $mode. $result is a gsl_sf_result structure.
758 =item C
<gsl_sf_airy_Ai_deriv_scaled_e
($x
, $mode
, $result
)>
760 =item C
<gsl_sf_airy_Ai_deriv_scaled
($x
, $mode
)>
762 -These routines compute the scaled Airy function derivative S_A
(x
) Ai'
(x
). For x
>0 the scaling factor S_A
(x
) is \exp
(+(2/3) x^
(3/2)), and is
1 for x
<0. $result is a gsl_sf_result structure.
768 =item C
<gsl_sf_airy_Bi_deriv_scaled_e
($x
, $mode
, $result
)>
770 =item C
<gsl_sf_airy_Bi_deriv_scaled
($x
, $mode
)>
772 -These routines compute the scaled Airy function derivative S_B
(x
) Bi'
(x
). For x
>0 the scaling factor S_B
(x
) is exp
(-(2/3) x^
(3/2)), and is
1 for x
<0. $result is a gsl_sf_result structure.
778 =item C
<gsl_sf_airy_zero_Ai_e
($s
, $result
)>
780 =item C
<gsl_sf_airy_zero_Ai
($s
)>
782 -These routines compute the location of the s-th zero of the Airy function Ai
($x
). $result is a gsl_sf_result structure.
788 =item C
<gsl_sf_airy_zero_Bi_e
($s
, $result
)>
790 =item C
<gsl_sf_airy_zero_Bi
($s
)>
792 -These routines compute the location of the s-th zero of the Airy function Bi
($x
). $result is a gsl_sf_result structure.
798 =item C
<gsl_sf_airy_zero_Ai_deriv_e
($s
, $result
)>
800 =item C
<gsl_sf_airy_zero_Ai_deriv
($s
)>
802 -These routines compute the location of the s-th zero of the Airy function derivative Ai'
(x
). $result is a gsl_sf_result structure.
808 =item C
<gsl_sf_airy_zero_Bi_deriv_e
($s
, $result
)>
810 =item C
<gsl_sf_airy_zero_Bi_deriv
($s
)>
812 - These routines compute the location of the s-th zero of the Airy function derivative Bi'
(x
). $result is a gsl_sf_result structure.
818 =item C
<gsl_sf_bessel_J0_e
($x
, $result
)>
820 =item C
<gsl_sf_bessel_J0
($x
)>
822 -These routines compute the regular cylindrical Bessel function of zeroth order
, J_0
($x
). $result is a gsl_sf_result structure.
828 =item C
<gsl_sf_bessel_J1_e
($x
, $result
)>
830 =item C
<gsl_sf_bessel_J1
($x
)>
832 - These routines compute the regular cylindrical Bessel function of first order
, J_1
($x
). $result is a gsl_sf_result structure.
838 =item C
<gsl_sf_bessel_Jn_e
($n
, $x
, $result
)>
840 =item C
<gsl_sf_bessel_Jn
($n
, $x
)>
842 -These routines compute the regular cylindrical Bessel function of order n
, J_n
($x
).
848 =item C
<gsl_sf_bessel_Jn_array
($nmin
, $nmax
, $x
, $result_array
)> - This routine computes the values of the regular cylindrical Bessel functions J_n
($x
) for n from $nmin to $nmax inclusive
, storing the results in the array $result_array. The values are computed using recurrence relations for efficiency
, and therefore may differ slightly from the exact values.
854 =item C
<gsl_sf_bessel_Y0_e
($x
, $result
)>
856 =item C
<gsl_sf_bessel_Y0
($x
)>
858 - These routines compute the irregular spherical Bessel function of zeroth order
, y_0
(x
) = -\cos
(x
)/x.
864 =item C
<gsl_sf_bessel_Y1_e
($x
, $result
)>
866 =item C
<gsl_sf_bessel_Y1
($x
)>
868 -These routines compute the irregular spherical Bessel function of first order
, y_1
(x
) = -(\cos
(x
)/x
+ \sin
(x
))/x.
874 =item C
<gsl_sf_bessel_Yn_e
>($n
, $x
, $result
)
876 =item C
<gsl_sf_bessel_Yn
($n
, $x
)>
878 -These routines compute the irregular cylindrical Bessel function of order $n
, Y_n
(x
), for x
>0.
884 =item C
<gsl_sf_bessel_Yn_array
>
892 =item C
<gsl_sf_bessel_I0_e
($x
, $result
)>
894 =item C
<gsl_sf_bessel_I0
($x
)>
896 -These routines compute the regular modified cylindrical Bessel function of zeroth order
, I_0
(x
).
902 =item C
<gsl_sf_bessel_I1_e
($x
, $result
)>
904 =item C
<gsl_sf_bessel_I1
($x
)>
906 -These routines compute the regular modified cylindrical Bessel function of first order
, I_1
(x
).
912 =item C
<gsl_sf_bessel_In_e
($n
, $x
, $result
)>
914 =item C
<gsl_sf_bessel_In
($n
, $x
)>
916 -These routines compute the regular modified cylindrical Bessel function of order $n
, I_n
(x
).
922 =item C
<gsl_sf_bessel_In_array
>
930 =item C
<gsl_sf_bessel_I0_scaled_e
($x
, $result
)>
932 =item C
<gsl_sf_bessel_I0_scaled
($x
)>
934 -These routines compute the scaled regular modified cylindrical Bessel function of zeroth order \exp
(-|x|
) I_0
(x
).
940 =item C
<gsl_sf_bessel_I1_scaled_e
($x
, $result
)>
942 =item C
<gsl_sf_bessel_I1_scaled
($x
)>
944 -These routines compute the scaled regular modified cylindrical Bessel function of first order \exp
(-|x|
) I_1
(x
).
950 =item C
<gsl_sf_bessel_In_scaled_e
($n
, $x
, $result
)>
952 =item C
<gsl_sf_bessel_In_scaled
($n
, $x
)>
954 -These routines compute the scaled regular modified cylindrical Bessel function of order $n
, \exp
(-|x|
) I_n
(x
)
960 =item C
<gsl_sf_bessel_In_scaled_array
>
968 =item C
<gsl_sf_bessel_K0_e
($x
, $result
)>
970 =item C
<gsl_sf_bessel_K0
($x
)>
972 -These routines compute the irregular modified cylindrical Bessel function of zeroth order
, K_0
(x
), for x
> 0.
978 =item C
<gsl_sf_bessel_K1_e
($x
, $result
)>
980 =item C
<gsl_sf_bessel_K1
($x
)>
982 -These routines compute the irregular modified cylindrical Bessel function of first order
, K_1
(x
), for x
> 0.
988 =item C
<gsl_sf_bessel_Kn_e
($n
, $x
, $result
)>
990 =item C
<gsl_sf_bessel_Kn
($n
, $x
)>
992 -These routines compute the irregular modified cylindrical Bessel function of order $n
, K_n
(x
), for x
> 0.
998 =item C
<gsl_sf_bessel_Kn_array
>
1006 =item C
<gsl_sf_bessel_K0_scaled_e
($x
, $result
)>
1008 =item C
<gsl_sf_bessel_K0_scaled
($x
)>
1010 -These routines compute the scaled irregular modified cylindrical Bessel function of zeroth order \exp
(x
) K_0
(x
) for x
>0.
1016 =item C
<gsl_sf_bessel_K1_scaled_e
($x
, $result
)>
1018 =item C
<gsl_sf_bessel_K1_scaled
($x
)>
1026 =item C
<gsl_sf_bessel_Kn_scaled_e
($n
, $x
, $result
)>
1028 =item C
<gsl_sf_bessel_Kn_scaled
($n
, $x
)>
1036 =item C
<gsl_sf_bessel_Kn_scaled_array
>
1044 =item C
<gsl_sf_bessel_j0_e
($x
, $result
)>
1046 =item C
<gsl_sf_bessel_j0
($x
)>
1054 =item C
<gsl_sf_bessel_j1_e
($x
, $result
)>
1056 =item C
<gsl_sf_bessel_j1
($x
)>
1064 =item C
<gsl_sf_bessel_j2_e
($x
, $result
)>
1066 =item C
<gsl_sf_bessel_j2
($x
)>
1074 =item C
<gsl_sf_bessel_jl_e
($l
, $x
, $result
)>
1076 =item C
<gsl_sf_bessel_jl
($l
, $x
)>
1084 =item C
<gsl_sf_bessel_jl_array
>
1092 =item C
<gsl_sf_bessel_jl_steed_array
>
1100 =item C
<gsl_sf_bessel_y0_e
($x
, $result
)>
1102 =item C
<gsl_sf_bessel_y0
($x
)>
1110 =item C
<gsl_sf_bessel_y1_e
($x
, $result
)>
1112 =item C
<gsl_sf_bessel_y1
($x
)>
1120 =item C
<gsl_sf_bessel_y2_e
($x
, $result
)>
1122 =item C
<gsl_sf_bessel_y2
($x
)>
1130 =item C
<gsl_sf_bessel_yl_e
($l
, $x
, $result
)>
1132 =item C
<gsl_sf_bessel_yl
($l
, $x
)>
1140 =item C
<gsl_sf_bessel_yl_array
>
1148 =item C
<gsl_sf_bessel_i0_scaled_e
($x
, $result
)>
1150 =item C
<gsl_sf_bessel_i0_scaled
($x
)>
1158 =item C
<gsl_sf_bessel_i1_scaled_e
($x
, $result
)>
1160 =item C
<gsl_sf_bessel_i1_scaled
($x
)>
1168 =item C
<gsl_sf_bessel_i2_scaled_e
($x
, $result
)>
1170 =item C
<gsl_sf_bessel_i2_scaled
($x
)>
1178 =item C
<gsl_sf_bessel_il_scaled_e
($l
, $x
, $result
)>
1180 =item C
<gsl_sf_bessel_il_scaled
($x
)>
1188 =item C
<gsl_sf_bessel_il_scaled_array
>
1196 =item C
<gsl_sf_bessel_k0_scaled_e
($x
, $result
)>
1198 =item C
<gsl_sf_bessel_k0_scale
($x
)>
1206 =item C
<gsl_sf_bessel_k1_scaled_e
($x
, $result
)>
1208 =item C
<gsl_sf_bessel_k1_scaled
($x
)>
1216 =item C
<gsl_sf_bessel_k2_scaled_e
($x
, $result
) >
1218 =item C
<gsl_sf_bessel_k2_scaled
($x
)>
1226 =item C
<gsl_sf_bessel_kl_scaled_e
($l
, $x
, $result
)>
1228 =item C
<gsl_sf_bessel_kl_scaled
($l
, $x
)>
1236 =item C
<gsl_sf_bessel_kl_scaled_array
>
1244 =item C
<gsl_sf_bessel_Jnu_e
($nu
, $x
, $result
)>
1246 =item C
<gsl_sf_bessel_Jnu
($nu
, $x
)>
1254 =item C
<gsl_sf_bessel_sequence_Jnu_e
>
1262 =item C
<gsl_sf_bessel_Ynu_e
($nu
, $x
, $result
)>
1264 =item C
<gsl_sf_bessel_Ynu
($nu
, $x
)>
1272 =item C
<gsl_sf_bessel_Inu_scaled_e
($nu
, $x
, $result
)>
1274 =item C
<gsl_sf_bessel_Inu_scaled
($nu
, $x
)>
1282 =item C
<gsl_sf_bessel_Inu_e
($nu
, $x
, $result
)>
1284 =item C
<gsl_sf_bessel_Inu
($nu
, $x
)>
1292 =item C
<gsl_sf_bessel_Knu_scaled_e
($nu
, $x
, $result
)>
1294 =item C
<gsl_sf_bessel_Knu_scaled
($nu
, $x
)>
1302 =item C
<gsl_sf_bessel_Knu_e
($nu
, $x
, $result
)>
1304 =item C
<gsl_sf_bessel_Knu
($nu
, $x
)>
1312 =item C
<gsl_sf_bessel_lnKnu_e
($nu
, $x
, $result
)>
1314 =item C
<gsl_sf_bessel_lnKnu
($nu
, $x
)>
1322 =item C
<gsl_sf_bessel_zero_J0_e
($s
, $result
)>
1324 =item C
<gsl_sf_bessel_zero_J0
($s
)>
1332 =item C
<gsl_sf_bessel_zero_J1_e
($s
, $result
)>
1334 =item C
<gsl_sf_bessel_zero_J1
($s
)>
1342 =item C
<gsl_sf_bessel_zero_Jnu_e
($nu
, $s
, $result
)>
1344 =item C
<gsl_sf_bessel_zero_Jnu
($nu
, $s
)>
1352 =item C
<gsl_sf_clausen_e
($x
, $result
)>
1354 =item C
<gsl_sf_clausen
($x
)>
1362 =item C
<gsl_sf_hydrogenicR_1_e
($Z
, $r
, $result
)>
1364 =item C
<gsl_sf_hydrogenicR_1
($Z
, $r
)>
1372 =item C
<gsl_sf_hydrogenicR_e
($n
, $l
, $Z
, $r
, $result
)>
1374 =item C
<gsl_sf_hydrogenicR
($n
, $l
, $Z
, $r
)>
1382 =item C
<gsl_sf_coulomb_wave_FG_e
($eta
, $x
, $L_F
, $k
, $F
, gsl_sf_result
* Fp
, gsl_sf_result
* G
, $Gp
)> - This function computes the Coulomb wave functions F_L
(\eta
,x
), G_
{L-k
}(\eta
,x
) and their derivatives F'_L
(\eta
,x
), G'_
{L-k
}(\eta
,x
) with respect to $x. The parameters are restricted to L
, L-k
> -1/2, x
> 0 and integer $k. Note that L itself is not restricted to being an integer. The results are stored in the parameters $F
, $G for the function values and $Fp
, $Gp for the derivative values. $F
, $G
, $Fp
, $Gp are all gsl_result structs. If an overflow occurs
, $GSL_EOVRFLW is returned and scaling exponents are returned as second and third values.
1384 =item C
<gsl_sf_coulomb_wave_F_array
> -
1386 =item C
<gsl_sf_coulomb_wave_FG_array
> -
1388 =item C
<gsl_sf_coulomb_wave_FGp_array
> -
1390 =item C
<gsl_sf_coulomb_wave_sphF_array
> -
1392 =item C
<gsl_sf_coulomb_CL_e
($L
, $eta
, $result
)> - This function computes the Coulomb wave function normalization constant C_L
($eta
) for $L
> -1.
1394 =item C
<gsl_sf_coulomb_CL_arrayi
> -
1400 =item C
<gsl_sf_coupling_3j_e
($two_ja
, $two_jb
, $two_jc
, $two_ma
, $two_mb
, $two_mc
, $result
)>
1402 =item C
<gsl_sf_coupling_3j
($two_ja
, $two_jb
, $two_jc
, $two_ma
, $two_mb
, $two_mc
)>
1404 - These routines compute the Wigner
3-j coefficient
,
1407 where the arguments are given in half-integer units
, ja
= $two_ja
/2, ma
= $two_ma
/2, etc.
1413 =item C
<gsl_sf_coupling_6j_e
($two_ja
, $two_jb
, $two_jc
, $two_jd
, $two_je
, $two_jf
, $result
)>
1415 =item C
<gsl_sf_coupling_6j
($two_ja
, $two_jb
, $two_jc
, $two_jd
, $two_je
, $two_jf
)>
1417 - These routines compute the Wigner
6-j coefficient
,
1420 where the arguments are given in half-integer units
, ja
= $two_ja
/2, ma
= $two_ma
/2, etc.
1426 =item C
<gsl_sf_coupling_RacahW_e
>
1428 =item C
<gsl_sf_coupling_RacahW
>
1436 =item C
<gsl_sf_coupling_9j_e
($two_ja
, $two_jb
, $two_jc
, $two_jd
, $two_je
, $two_jf
, $two_jg
, $two_jh
, $two_ji
, $result
)>
1438 =item C
<gsl_sf_coupling_9j
($two_ja
, $two_jb
, $two_jc
, $two_jd
, $two_je
, $two_jf
, $two_jg
, $two_jh
, $two_ji
)>
1440 -These routines compute the Wigner
9-j coefficient
,
1445 where the arguments are given in half-integer units
, ja
= two_ja
/2, ma
= two_ma
/2, etc.
1451 =item C
<gsl_sf_dawson_e
($x
, $result
)>
1453 =item C
<gsl_sf_dawson
($x
)>
1455 -These routines compute the value of Dawson's integral for $x.
1461 =item C
<gsl_sf_debye_1_e
($x
, $result
)>
1463 =item C
<gsl_sf_debye_1
($x
)>
1465 -These routines compute the first-order Debye function D_1
(x
) = (1/x
) \int_0^x dt
(t
/(e^t
- 1)).
1471 =item C
<gsl_sf_debye_2_e
($x
, $result
)>
1473 =item C
<gsl_sf_debye_2
($x
)>
1475 -These routines compute the second-order Debye function D_2
(x
) = (2/x^
2) \int_0^x dt
(t^
2/(e^t
- 1)).
1481 =item C
<gsl_sf_debye_3_e
($x
, $result
)>
1483 =item C
<gsl_sf_debye_3
($x
)>
1485 -These routines compute the third-order Debye function D_3
(x
) = (3/x^
3) \int_0^x dt
(t^
3/(e^t
- 1)).
1491 =item C
<gsl_sf_debye_4_e
($x
, $result
)>
1493 =item C
<gsl_sf_debye_4
($x
)>
1495 -These routines compute the fourth-order Debye function D_4
(x
) = (4/x^
4) \int_0^x dt
(t^
4/(e^t
- 1)).
1501 =item C
<gsl_sf_debye_5_e
($x
, $result
)>
1503 =item C
<gsl_sf_debye_5
($x
)>
1505 -These routines compute the fifth-order Debye function D_5
(x
) = (5/x^
5) \int_0^x dt
(t^
5/(e^t
- 1)).
1511 =item C
<gsl_sf_debye_6_e
($x
, $result
)>
1513 =item C
<gsl_sf_debye_6
($x
)>
1515 -These routines compute the sixth-order Debye function D_6
(x
) = (6/x^
6) \int_0^x dt
(t^
6/(e^t
- 1)).
1521 =item C
<gsl_sf_dilog_e
($x
, $result
)>
1523 =item C
<gsl_sf_dilog
($x
)>
1525 - These routines compute the dilogarithm for a real argument. In Lewin's notation this is Li_2
(x
), the real part of the dilogarithm of a real x. It is defined by the integral representation Li_2
(x
) = - \Re \int_0^x ds \log
(1-s
) / s. Note that \Im
(Li_2
(x
)) = 0 for x
<= 1, and
-\pi\log
(x
) for x
> 1. Note that Abramowitz
& Stegun refer to the Spence integral S(x)=Li_2(1-x) as the dilogarithm rather than Li_2(x).
1531 =item C
<gsl_sf_complex_dilog_xy_e
> -
1533 =item C
<gsl_sf_complex_dilog_e
($r
, $theta
, $result_re
, $result_im
)> - This function computes the full complex-valued dilogarithm for the complex argument z
= r \exp
(i \theta
). The real and imaginary parts of the result are returned in the $result_re and $result_im gsl_result structs.
1535 =item C
<gsl_sf_complex_spence_xy_e
> -
1541 =item C
<gsl_sf_multiply
>
1543 =item C
<gsl_sf_multiply_e
($x
, $y
, $result
)> - This function multiplies $x and $y storing the product and its associated error in $result.
1545 =item C
<gsl_sf_multiply_err_e
($x
, $dx
, $y
, $dy
, $result
)> - This function multiplies $x and $y with associated absolute errors $dx and $dy. The product xy
+/- xy \sqrt
((dx
/x
)^
2 +(dy
/y
)^
2) is stored in $result.
1554 =item C
<gsl_sf_ellint_Kcomp_e
($k
, $mode
, $result
)>
1556 =item C
<gsl_sf_ellint_Kcomp
($k
, $mode
)>
1558 -These routines compute the complete elliptic integral K
($k
) to the accuracy specified by the mode variable mode. Note that Abramowitz
& Stegun define this function in terms of the parameter m = k^2.
1564 =item C
<gsl_sf_ellint_Ecomp_e
($k
, $mode
, $result
)>
1566 =item C
<gsl_sf_ellint_Ecomp
($k
, $mode
)>
1574 =item C
<gsl_sf_ellint_Pcomp_e
($k
, $n
, $mode
, $result
)>
1576 =item C
<gsl_sf_ellint_Pcomp
($k
, $n
, $mode
)>
1584 =item C
<gsl_sf_ellint_Dcomp_e
>
1586 =item C
<gsl_sf_ellint_Dcomp
>
1594 =item C
<gsl_sf_ellint_F_e
($phi
, $k
, $mode
, $result
)>
1596 =item C
<gsl_sf_ellint_F
($phi
, $k
, $mode
)>
1598 -These routines compute the incomplete elliptic integral F
($phi
,$k
) to the accuracy specified by the mode variable mode. Note that Abramowitz
& Stegun define this function in terms of the parameter m = k^2.
1604 =item C
<gsl_sf_ellint_E_e
($phi
, $k
, $mode
, $result
)>
1606 =item C
<gsl_sf_ellint_E
($phi
, $k
, $mode
)>
1608 -These routines compute the incomplete elliptic integral E
($phi
,$k
) to the accuracy specified by the mode variable mode. Note that Abramowitz
& Stegun define this function in terms of the parameter m = k^2.
1614 =item C
<gsl_sf_ellint_P_e
($phi
, $k
, $n
, $mode
, $result
)>
1616 =item C
<gsl_sf_ellint_P
($phi
, $k
, $n
, $mode
)>
1618 -These routines compute the incomplete elliptic integral \Pi
(\phi
,k
,n
) to the accuracy specified by the mode variable mode. Note that Abramowitz
& Stegun define this function in terms of the parameters m = k^2 and \sin^2(\alpha) = k^2, with the change of sign n \to -n.
1624 =item C
<gsl_sf_ellint_D_e
($phi
, $k
, $n
, $mode
, $result
)>
1626 =item C
<gsl_sf_ellint_D
($phi
, $k
, $n
, $mode
)>
1628 -These functions compute the incomplete elliptic integral D
(\phi
,k
) which is defined through the Carlson form RD
(x
,y
,z
) by the following relation
, D
(\phi
,k
,n
) = (1/3)(\sin
(\phi
))^
3 RD
(1-\sin^
2(\phi
), 1-k^
2 \sin^
2(\phi
), 1). The argument $n is not used and will be removed in a future release.
1634 =item C
<gsl_sf_ellint_RC_e
($x
, $y
, $mode
, $result
)>
1636 =item C
<gsl_sf_ellint_RC
($x
, $y
, $mode
)>
1638 - These routines compute the incomplete elliptic integral RC
($x
,$y
) to the accuracy specified by the mode variable $mode.
1644 =item C
<gsl_sf_ellint_RD_e
($x
, $y
, $z
, $mode
, $result
)>
1646 =item C
<gsl_sf_ellint_RD
($x
, $y
, $z
, $mode
)>
1648 - These routines compute the incomplete elliptic integral RD
($x
,$y
,$z
) to the accuracy specified by the mode variable $mode.
1654 =item C
<gsl_sf_ellint_RF_e
($x
, $y
, $z
, $mode
, $result
)>
1656 =item C
<gsl_sf_ellint_RF
($x
, $y
, $z
, $mode
)>
1658 - These routines compute the incomplete elliptic integral RF
($x
,$y
,$z
) to the accuracy specified by the mode variable $mode.
1664 =item C
<gsl_sf_ellint_RJ_e
($x
, $y
, $z
, $p
, $mode
, $result
)>
1666 =item C
<gsl_sf_ellint_RJ
($x
, $y
, $z
, $p
, $mode
)>
1668 - These routines compute the incomplete elliptic integral RJ
($x
,$y
,$z
,$p
) to the accuracy specified by the mode variable $mode.
1674 =item C
<gsl_sf_elljac_e
($u
, $m
)> - This function computes the Jacobian elliptic functions sn
(u|m
), cn
(u|m
), dn
(u|m
) by descending Landen transformations. The function returns
0 if the operation succeded
, 1 otherwise and then returns the result of sn
, cn and dn in this order.
1676 =item C
<gsl_sf_erfc_e
($x
, $result
)>
1678 =item C
<gsl_sf_erfc
($x
)>
1680 -These routines compute the complementary error function erfc
(x
) = 1 - erf
(x
) = (2/\sqrt
(\pi
)) \int_x^\infty \exp
(-t^
2).
1686 =item C
<gsl_sf_log_erfc_e
($x
, $result
)>
1688 =item C
<gsl_sf_log_erfc
($x
)>
1690 -These routines compute the logarithm of the complementary error function \log
(\erfc
(x
)).
1696 =item C
<gsl_sf_erf_e
($x
, $result
)>
1698 =item C
<gsl_sf_erf
($x
)>
1700 -These routines compute the error function erf
(x
), where erf
(x
) = (2/\sqrt
(\pi
)) \int_0^x dt \exp
(-t^
2).
1706 =item C
<gsl_sf_erf_Z_e
($x
, $result
)>
1708 =item C
<gsl_sf_erf_Z
($x
)>
1710 -These routines compute the Gaussian probability density function Z
(x
) = (1/\sqrt
{2\pi
}) \exp
(-x^
2/2).
1716 =item C
<gsl_sf_erf_Q_e
($x
, $result
)>
1718 =item C
<gsl_sf_erf_Q
($x
)>
1720 - These routines compute the upper tail of the Gaussian probability function Q
(x
) = (1/\sqrt
{2\pi
}) \int_x^\infty dt \exp
(-t^
2/2). The hazard function for the normal distribution
, also known as the inverse Mill's ratio
, is defined as
, h
(x
) = Z
(x
)/Q
(x
) = \sqrt
{2/\pi
} \exp
(-x^
2 / 2) / \erfc
(x
/\sqrt
2) It decreases rapidly as x approaches
-\infty and asymptotes to h
(x
) \sim x as x approaches
+\infty.
1726 =item C
<gsl_sf_hazard_e
($x
, $result
)>
1728 =item C
<gsl_sf_hazard
($x
)>
1730 - These routines compute the hazard function for the normal distribution.
1736 =item C
<gsl_sf_exp_e
($x
, $result
)>
1738 =item C
<gsl_sf_exp
($x
)>
1740 - These routines provide an exponential function \exp
(x
) using GSL semantics and error checking.
1746 =item C
<gsl_sf_exp_e10_e
> -
1752 =item C
<gsl_sf_exp_mult_e
>
1754 =item C
<gsl_sf_exp_mult
>
1762 =item C
<gsl_sf_exp_mult_e10_e
> -
1768 =item C
<gsl_sf_expm1_e
($x
, $result
)>
1770 =item C
<gsl_sf_expm1
($x
)>
1772 -These routines compute the quantity \exp
(x
)-1 using an algorithm that is accurate for small x.
1778 =item C
<gsl_sf_exprel_e
($x
, $result
)>
1780 =item C
<gsl_sf_exprel
($x
)>
1782 -These routines compute the quantity
(\exp
(x
)-1)/x using an algorithm that is accurate for small x. For small x the algorithm is based on the expansion
(\exp
(x
)-1)/x
= 1 + x
/2 + x^
2/(2*3) + x^
3/(2*3*4) + \dots.
1788 =item C
<gsl_sf_exprel_2_e
($x
, $result
)>
1790 =item C
<gsl_sf_exprel_2
($x
)>
1792 -These routines compute the quantity
2(\exp
(x
)-1-x
)/x^
2 using an algorithm that is accurate for small x. For small x the algorithm is based on the expansion
2(\exp
(x
)-1-x
)/x^
2 = 1 + x
/3 + x^
2/(3*4) + x^
3/(3*4*5) + \dots.
1798 =item C
<gsl_sf_exprel_n_e
($x
, $result
)>
1800 =item C
<gsl_sf_exprel_n
($x
)>
1802 -These routines compute the N-relative exponential
, which is the n-th generalization of the functions gsl_sf_exprel and gsl_sf_exprel2. The N-relative exponential is given by
,
1803 exprel_N
(x
) = N
!/x^N
(\exp
(x
) - \sum_
{k
=0}^
{N-1
} x^k
/k
!)
1804 = 1 + x
/(N
+1) + x^
2/((N
+1)(N
+2)) + ...
1811 =item C
<gsl_sf_exp_err_e
($x
, $dx
, $result
)> - This function exponentiates $x with an associated absolute error $dx.
1813 =item C
<gsl_sf_exp_err_e10_e
> -
1815 =item C
<gsl_sf_exp_mult_err_e
($x
, $dx
, $y
, $dy
, $result
)> -
1817 =item C
<gsl_sf_exp_mult_err_e10_e
> -
1823 =item C
<gsl_sf_expint_E1_e
($x
, $result
)>
1825 =item C
<gsl_sf_expint_E1
($x
)>
1827 -These routines compute the exponential integral E_1
(x
), E_1
(x
) := \Re \int_1^\infty dt \exp
(-xt
)/t.
1833 =item C
<gsl_sf_expint_E2_e
($x
, $result
)>
1835 =item C
<gsl_sf_expint_E2
($x
)>
1837 -These routines compute the second-order exponential integral E_2
(x
),
1838 E_2
(x
) := \Re \int_1^\infty dt \exp
(-xt
)/t^
2.
1845 =item C
<gsl_sf_expint_En_e
($n
, $x
, $result
)>
1847 =item C
<gsl_sf_expint_En
($n
, $x
)>
1849 -These routines compute the exponential integral E_n
(x
) of order n
,
1850 E_n
(x
) := \Re \int_1^\infty dt \exp
(-xt
)/t^n.
1857 =item C
<gsl_sf_expint_E1_scaled_e
>
1859 =item C
<gsl_sf_expint_E1_scaled
>
1867 =item C
<gsl_sf_expint_E2_scaled_e
>
1869 =item C
<gsl_sf_expint_E2_scaled
>
1877 =item C
<gsl_sf_expint_En_scaled_e
>
1879 =item C
<gsl_sf_expint_En_scaled
>
1887 =item C
<gsl_sf_expint_Ei_e
($x
, $result
)>
1889 =item C
<gsl_sf_expint_Ei
($x
)>
1891 -These routines compute the exponential integral Ei
(x
), Ei
(x
) := - PV
(\int_
{-x
}^\infty dt \exp
(-t
)/t
) where PV denotes the principal value of the integral.
1897 =item C
<gsl_sf_expint_Ei_scaled_e
>
1899 =item C
<gsl_sf_expint_Ei_scaled
>
1907 =item C
<gsl_sf_Shi_e
($x
, $result
)>
1909 =item C
<gsl_sf_Shi
($x
)>
1911 -These routines compute the integral Shi
(x
) = \int_0^x dt \sinh
(t
)/t.
1917 =item C
<gsl_sf_Chi_e
($x
, $result
)>
1919 =item C
<gsl_sf_Chi
($x
)>
1921 -These routines compute the integral Chi
(x
) := \Re
[ \gamma_E
+ \log
(x
) + \int_0^x dt
(\cosh
[t
]-1)/t
] , where \gamma_E is the Euler constant
(available as $M_EULER from the Math
::GSL
::Const module
).
1927 =item C
<gsl_sf_expint_3_e
($x
, $result
)>
1929 =item C
<gsl_sf_expint_3
($x
)>
1931 -These routines compute the third-order exponential integral Ei_3
(x
) = \int_0^xdt \exp
(-t^
3) for x
>= 0.
1937 =item C
<gsl_sf_Si_e
($x
, $result
)>
1939 =item C
<gsl_sf_Si
($x
)>
1941 -These routines compute the Sine integral Si
(x
) = \int_0^x dt \sin
(t
)/t.
1947 =item C
<gsl_sf_Ci_e
($x
, $result
)>
1949 =item C
<gsl_sf_Ci
($x
)>
1951 -These routines compute the Cosine integral Ci
(x
) = -\int_x^\infty dt \cos
(t
)/t for x
> 0.
1957 =item C
<gsl_sf_fermi_dirac_m1_e
($x
, $result
)>
1959 =item C
<gsl_sf_fermi_dirac_m1
($x
)>
1961 -These routines compute the complete Fermi-Dirac integral with an index of
-1. This integral is given by F_
{-1}(x
) = e^x
/ (1 + e^x
).
1967 =item C
<gsl_sf_fermi_dirac_0_e
($x
, $result
)>
1969 =item C
<gsl_sf_fermi_dirac_0
($x
)>
1971 -These routines compute the complete Fermi-Dirac integral with an index of
0. This integral is given by F_0
(x
) = \ln
(1 + e^x
).
1977 =item C
<gsl_sf_fermi_dirac_1_e
($x
, $result
)>
1979 =item C
<gsl_sf_fermi_dirac_1
($x
)>
1981 -These routines compute the complete Fermi-Dirac integral with an index of
1, F_1
(x
) = \int_0^\infty dt
(t
/(\exp
(t-x
)+1)).
1987 =item C
<gsl_sf_fermi_dirac_2_e
($x
, $result
)>
1989 =item C
<gsl_sf_fermi_dirac_2
($x
)>
1991 -These routines compute the complete Fermi-Dirac integral with an index of
2, F_2
(x
) = (1/2) \int_0^\infty dt
(t^
2 /(\exp
(t-x
)+1)).
1997 =item C
<gsl_sf_fermi_dirac_int_e
($j
, $x
, $result
)>
1999 =item C
<gsl_sf_fermi_dirac_int
($j
, $x
)>
2001 -These routines compute the complete Fermi-Dirac integral with an integer index of j
, F_j
(x
) = (1/\Gamma
(j
+1)) \int_0^\infty dt
(t^j
/(\exp
(t-x
)+1)).
2007 =item C
<gsl_sf_fermi_dirac_mhalf_e
($x
, $result
)>
2009 =item C
<gsl_sf_fermi_dirac_mhalf
($x
)>
2011 -These routines compute the complete Fermi-Dirac integral F_
{-1/2}(x
).
2017 =item C
<gsl_sf_fermi_dirac_half_e
($x
, $result
)>
2019 =item C
<gsl_sf_fermi_dirac_half
($x
)>
2021 -These routines compute the complete Fermi-Dirac integral F_
{1/2}(x
).
2027 =item C
<gsl_sf_fermi_dirac_3half_e
($x
, $result
)>
2029 =item C
<gsl_sf_fermi_dirac_3half
($x
)>
2031 -These routines compute the complete Fermi-Dirac integral F_
{3/2}(x
).
2037 =item C
<gsl_sf_fermi_dirac_inc_0_e
($x
, $b
, $result
)>
2039 =item C
<gsl_sf_fermi_dirac_inc_0
($x
, $b
, $result
)>
2041 -These routines compute the incomplete Fermi-Dirac integral with an index of zero
, F_0
(x
,b
) = \ln
(1 + e^
{b-x
}) - (b-x
).
2047 =item C
<gsl_sf_legendre_Pl_e
($l
, $x
, $result
)>
2049 =item C
<gsl_sf_legendre_Pl
($l
, $x
)>
2051 -These functions evaluate the Legendre polynomial P_l
(x
) for a specific value of l
, x subject to l
>= 0, |x|
<= 1
2057 =item C
<gsl_sf_legendre_Pl_array
>
2059 =item C
<gsl_sf_legendre_Pl_deriv_array
>
2067 =item C
<gsl_sf_legendre_P1_e
($x
, $result
)>
2069 =item C
<gsl_sf_legendre_P2_e
($x
, $result
)>
2071 =item C
<gsl_sf_legendre_P3_e
($x
, $result
)>
2073 =item C
<gsl_sf_legendre_P1
($x
)>
2075 =item C
<gsl_sf_legendre_P2
($x
)>
2077 =item C
<gsl_sf_legendre_P3
($x
)>
2079 -These functions evaluate the Legendre polynomials P_l
(x
) using explicit representations for l
=1, 2, 3.
2085 =item C
<gsl_sf_legendre_Q0_e
($x
, $result
)>
2087 =item C
<gsl_sf_legendre_Q0
($x
)>
2089 -These routines compute the Legendre function Q_0
(x
) for x
> -1, x
!= 1.
2095 =item C
<gsl_sf_legendre_Q1_e
($x
, $result
)>
2097 =item C
<gsl_sf_legendre_Q1
($x
)>
2099 -These routines compute the Legendre function Q_1
(x
) for x
> -1, x
!= 1.
2105 =item C
<gsl_sf_legendre_Ql_e
($l
, $x
, $result
)>
2107 =item C
<gsl_sf_legendre_Ql
($l
, $x
)>
2109 -These routines compute the Legendre function Q_l
(x
) for x
> -1, x
!= 1 and l
>= 0.
2115 =item C
<gsl_sf_legendre_Plm_e
($l
, $m
, $x
, $result
)>
2117 =item C
<gsl_sf_legendre_Plm
($l
, $m
, $x
)>
2119 -These routines compute the associated Legendre polynomial P_l^m
(x
) for m
>= 0, l
>= m
, |x|
<= 1.
2125 =item C
<gsl_sf_legendre_Plm_array
>
2127 =item C
<gsl_sf_legendre_Plm_deriv_array
>
2135 =item C
<gsl_sf_legendre_sphPlm_e
($l
, $m
, $x
, $result
)>
2137 =item C
<gsl_sf_legendre_sphPlm
($l
, $m
, $x
)>
2139 -These routines compute the normalized associated Legendre polynomial $\sqrt
{(2l+1)/(4\pi
)} \sqrt
{(l-m
)!/(l
+m
)!} P_l^m
(x
)$ suitable for use in spherical harmonics. The parameters must satisfy m
>= 0, l
>= m
, |x|
<= 1. Theses routines avoid the overflows that occur for the standard normalization of P_l^m
(x
).
2145 =item C
<gsl_sf_legendre_sphPlm_array
>
2147 =item C
<gsl_sf_legendre_sphPlm_deriv_array
>
2155 =item C
<gsl_sf_legendre_array_size
> -
2161 =item C
<gsl_sf_lngamma_e
($x
, $result
)>
2163 =item C
<gsl_sf_lngamma
($x
)>
2165 -These routines compute the logarithm of the Gamma function
, \log
(\Gamma
(x
)), subject to x not being a negative integer or zero. For x
<0 the real part of \log
(\Gamma
(x
)) is returned
, which is equivalent to \log
(|\Gamma
(x
)|
). The function is computed using the real Lanczos method.
2171 =item C
<gsl_sf_lngamma_sgn_e
($x
, $result_lg
)> - This routine returns the sign of the gamma function and the logarithm of its magnitude into this order
, subject to $x not being a negative integer or zero. The function is computed using the real Lanczos method. The value of the gamma function can be reconstructed using the relation \Gamma
(x
) = sgn
* \exp
(resultlg
).
2177 =item C
<gsl_sf_gamma_e
>
2179 =item C
<gsl_sf_gamma
>
2187 =item C
<gsl_sf_gammastar_e
>
2189 =item C
<gsl_sf_gammastar
>
2197 =item C
<gsl_sf_gammainv_e
>
2199 =item C
<gsl_sf_gammainv
>
2207 =item C
<gsl_sf_lngamma_complex_e
>
2215 =item C
<gsl_sf_gamma_inc_Q_e
>
2217 =item C
<gsl_sf_gamma_inc_Q
>
2225 =item C
<gsl_sf_gamma_inc_P_e
>
2227 =item C
<gsl_sf_gamma_inc_P
>
2235 =item C
<gsl_sf_gamma_inc_e
>
2237 =item C
<gsl_sf_gamma_inc
>
2245 =item C
<gsl_sf_taylorcoeff_e
>
2247 =item C
<gsl_sf_taylorcoeff
>
2255 =item C
<gsl_sf_fact_e
>
2257 =item C
<gsl_sf_fact
>
2265 =item C
<gsl_sf_doublefact_e
>
2267 =item C
<gsl_sf_doublefact
>
2275 =item C
<gsl_sf_lnfact_e
>
2277 =item C
<gsl_sf_lnfact
>
2285 =item C
<gsl_sf_lndoublefact_e
>
2287 =item C
<gsl_sf_lndoublefact
>
2295 =item C
<gsl_sf_lnchoose_e
>
2297 =item C
<gsl_sf_lnchoose
>
2305 =item C
<gsl_sf_choose_e
>
2307 =item C
<gsl_sf_choose
>
2315 =item C
<gsl_sf_lnpoch_e
>
2317 =item C
<gsl_sf_lnpoch
>
2325 =item C
<gsl_sf_lnpoch_sgn_e
>
2333 =item C
<gsl_sf_poch_e
>
2335 =item C
<gsl_sf_poch
>
2343 =item C
<gsl_sf_pochrel_e
>
2345 =item C
<gsl_sf_pochrel
>
2353 =item C
<gsl_sf_lnbeta_e
>
2355 =item C
<gsl_sf_lnbeta
>
2363 =item C
<gsl_sf_lnbeta_sgn_e
>
2371 =item C
<gsl_sf_beta_e
>
2373 =item C
<gsl_sf_beta
>
2381 =item C
<gsl_sf_beta_inc_e
>
2383 =item C
<gsl_sf_beta_inc
>
2391 =item C
<gsl_sf_gegenpoly_1_e
>
2393 =item C
<gsl_sf_gegenpoly_2_e
>
2395 =item C
<gsl_sf_gegenpoly_3_e
>
2397 =item C
<gsl_sf_gegenpoly_1
>
2399 =item C
<gsl_sf_gegenpoly_2
>
2401 =item C
<gsl_sf_gegenpoly_3
>
2409 =item C
<gsl_sf_gegenpoly_n_e
>
2411 =item C
<gsl_sf_gegenpoly_n
>
2419 =item C
<gsl_sf_gegenpoly_array
>
2421 =item C
<gsl_sf_hyperg_0F1_e
>
2423 =item C
<gsl_sf_hyperg_0F1
>
2431 =item C
<gsl_sf_hyperg_1F1_int_e
>
2433 =item C
<gsl_sf_hyperg_1F1_int
>
2441 =item C
<gsl_sf_hyperg_1F1_e
>
2443 =item C
<gsl_sf_hyperg_1F1
>
2451 =item C
<gsl_sf_hyperg_U_int_e
>
2453 =item C
<gsl_sf_hyperg_U_int
>
2461 =item C
<gsl_sf_hyperg_U_int_e10_e
>
2469 =item C
<gsl_sf_hyperg_U_e
>
2471 =item C
<gsl_sf_hyperg_U
>
2479 =item C
<gsl_sf_hyperg_U_e10_e
>
2487 =item C
<gsl_sf_hyperg_2F1_e
>
2489 =item C
<gsl_sf_hyperg_2F1
>
2497 =item C
<gsl_sf_hyperg_2F1_conj_e
>
2499 =item C
<gsl_sf_hyperg_2F1_conj
>
2507 =item C
<gsl_sf_hyperg_2F1_renorm_e
>
2509 =item C
<gsl_sf_hyperg_2F1_renorm
>
2517 =item C
<gsl_sf_hyperg_2F1_conj_renorm_e
>
2519 =item C
<gsl_sf_hyperg_2F1_conj_renorm
>
2527 =item C
<gsl_sf_hyperg_2F0_e
>
2529 =item C
<gsl_sf_hyperg_2F0
>
2537 =item C
<gsl_sf_laguerre_1_e
>
2539 =item C
<gsl_sf_laguerre_2_e
>
2541 =item C
<gsl_sf_laguerre_3_e
>
2543 =item C
<gsl_sf_laguerre_1
>
2545 =item C
<gsl_sf_laguerre_2
>
2547 =item C
<gsl_sf_laguerre_3
>
2555 =item C
<gsl_sf_laguerre_n_e
>
2557 =item C
<gsl_sf_laguerre_n
>
2565 =item C
<gsl_sf_lambert_W0_e
>
2567 =item C
<gsl_sf_lambert_W0
>
2575 =item C
<gsl_sf_lambert_Wm1_e
>
2577 =item C
<gsl_sf_lambert_Wm1
>
2585 =item C
<gsl_sf_conicalP_half_e
>
2587 =item C
<gsl_sf_conicalP_half
>
2595 =item C
<gsl_sf_conicalP_mhalf_e
>
2597 =item C
<gsl_sf_conicalP_mhalf
>
2605 =item C
<gsl_sf_conicalP_0_e
>
2607 =item C
<gsl_sf_conicalP_0
>
2615 =item C
<gsl_sf_conicalP_1_e
>
2617 =item C
<gsl_sf_conicalP_1
>
2625 =item C
<gsl_sf_conicalP_sph_reg_e
>
2627 =item C
<gsl_sf_conicalP_sph_reg
>
2635 =item C
<gsl_sf_conicalP_cyl_reg_e
>
2637 =item C
<gsl_sf_conicalP_cyl_reg
>
2645 =item C
<gsl_sf_legendre_H3d_0_e
>
2647 =item C
<gsl_sf_legendre_H3d_0
>
2655 =item C
<gsl_sf_legendre_H3d_1_e
>
2657 =item C
<gsl_sf_legendre_H3d_1
>
2665 =item C
<gsl_sf_legendre_H3d_e
>
2667 =item C
<gsl_sf_legendre_H3d
>
2675 =item C
<gsl_sf_legendre_H3d_array
>
2683 =item C
<gsl_sf_log_e
>
2693 =item C
<gsl_sf_log_abs_e
>
2695 =item C
<gsl_sf_log_abs
>
2703 =item C
<gsl_sf_complex_log_e
>
2711 =item C
<gsl_sf_log_1plusx_e
>
2713 =item C
<gsl_sf_log_1plusx
>
2721 =item C
<gsl_sf_log_1plusx_mx_e
>
2723 =item C
<gsl_sf_log_1plusx_mx
>
2731 =item C
<gsl_sf_mathieu_a_array
>
2733 =item C
<gsl_sf_mathieu_b_array
>
2741 =item C
<gsl_sf_mathieu_a
>
2743 =item C
<gsl_sf_mathieu_b
>
2751 =item C
<gsl_sf_mathieu_a_coeff
>
2753 =item C
<gsl_sf_mathieu_b_coeff
>
2761 =item C
<gsl_sf_mathieu_alloc
>
2769 =item C
<gsl_sf_mathieu_free
>
2777 =item C
<gsl_sf_mathieu_ce
>
2779 =item C
<gsl_sf_mathieu_se
>
2787 =item C
<gsl_sf_mathieu_ce_array
>
2789 =item C
<gsl_sf_mathieu_se_array
>
2797 =item C
<gsl_sf_mathieu_Mc
>
2799 =item C
<gsl_sf_mathieu_Ms
>
2807 =item C
<gsl_sf_mathieu_Mc_array
>
2809 =item C
<gsl_sf_mathieu_Ms_array
>
2817 =item C
<gsl_sf_pow_int_e
>
2819 =item C
<gsl_sf_pow_int
>
2827 =item C
<gsl_sf_psi_int_e
>
2829 =item C
<gsl_sf_psi_int
>
2837 =item C
<gsl_sf_psi_e
>
2847 =item C
<gsl_sf_psi_1piy_e
>
2849 =item C
<gsl_sf_psi_1piy
>
2857 =item C
<gsl_sf_complex_psi_e
>
2865 =item C
<gsl_sf_psi_1_int_e
>
2867 =item C
<gsl_sf_psi_1_int
>
2875 =item C
<gsl_sf_psi_1_e
>
2877 =item C
<gsl_sf_psi_1
>
2885 =item C
<gsl_sf_psi_n_e
>
2887 =item C
<gsl_sf_psi_n
>
2895 =item C
<gsl_sf_result_smash_e
>
2903 =item C
<gsl_sf_synchrotron_1_e
>
2905 =item C
<gsl_sf_synchrotron_1
>
2913 =item C
<gsl_sf_synchrotron_2_e
>
2915 =item C
<gsl_sf_synchrotron_2
>
2923 =item C
<gsl_sf_transport_2_e
>
2925 =item C
<gsl_sf_transport_2
>
2933 =item C
<gsl_sf_transport_3_e
>
2935 =item C
<gsl_sf_transport_3
>
2943 =item C
<gsl_sf_transport_4_e
>
2945 =item C
<gsl_sf_transport_4
>
2953 =item C
<gsl_sf_transport_5_e
>
2955 =item C
<gsl_sf_transport_5
>
2963 =item C
<gsl_sf_sin_e
>
2973 =item C
<gsl_sf_cos_e
>
2975 =item C
<gsl_sf_cos
>
2983 =item C
<gsl_sf_hypot_e
>
2985 =item C
<gsl_sf_hypot
>
2993 =item C
<gsl_sf_complex_sin_e
>
3001 =item C
<gsl_sf_complex_cos_e
>
3009 =item C
<gsl_sf_complex_logsin_e
>
3017 =item C
<gsl_sf_sinc_e
>
3019 =item C
<gsl_sf_sinc
>
3027 =item C
<gsl_sf_lnsinh_e
>
3029 =item C
<gsl_sf_lnsinh
>
3037 =item C
<gsl_sf_lncosh_e
>
3039 =item C
<gsl_sf_lncosh
>
3047 =item C
<gsl_sf_polar_to_rect
>
3055 =item C
<gsl_sf_rect_to_polar
>
3063 =item C
<gsl_sf_sin_err_e
>
3065 =item C
<gsl_sf_cos_err_e
>
3073 =item C
<gsl_sf_angle_restrict_symm_e
>
3075 =item C
<gsl_sf_angle_restrict_symm
>
3083 =item C
<gsl_sf_angle_restrict_pos_e
>
3085 =item C
<gsl_sf_angle_restrict_pos
>
3093 =item C
<gsl_sf_angle_restrict_symm_err_e
>
3095 =item C
<gsl_sf_angle_restrict_pos_err_e
>
3099 =item C
<gsl_sf_atanint_e
>
3101 =item C
<gsl_sf_atanint
>
3103 -These routines compute the Arctangent integral
, which is defined as AtanInt
(x
) = \int_0^x dt \arctan
(t
)/t.
3107 =item C
<gsl_sf_zeta_int_e
>
3109 =item C
<gsl_sf_zeta_int
>
3111 =item C
<gsl_sf_zeta_e gsl_sf_zeta
>
3113 =item C
<gsl_sf_zetam1_e
>
3115 =item C
<gsl_sf_zetam1
>
3117 =item C
<gsl_sf_zetam1_int_e
>
3119 =item C
<gsl_sf_zetam1_int
>
3121 =item C
<gsl_sf_hzeta_e
>
3123 =item C
<gsl_sf_hzeta
>
3125 =item C
<gsl_sf_eta_int_e
>
3127 =item C
<gsl_sf_eta_int
>
3129 =item C
<gsl_sf_eta_e
>
3131 =item C
<gsl_sf_eta
>
3135 This module also contains the following constants used as mode in various of those functions
:
3139 =item
* GSL_PREC_DOUBLE
- Double-precision
, a relative accuracy of approximately
2 * 10^
-16.
3141 =item
* GSL_PREC_SINGLE
- Single-precision
, a relative accuracy of approximately
10^
-7.
3143 =item
* GSL_PREC_APPROX
- Approximate values
, a relative accuracy of approximately
5 * 10^
-4.
3147 You can import the functions that you want to use by giving a space separated
3148 list to Math
::GSL
::SF when you use the package. You can also write
3149 use Math
::GSL
::SF qw
/:all
/
3150 to use all avaible functions of the module. Note that
3151 the tag names begin with a colon. Other tags are also available
, here is a
3152 complete list of all tags for this module
:
3182 =item C
<hypergeometric
>
3202 For more informations on the functions
, we refer you to the GSL offcial
3203 documentation
: L
<http
://www.gnu.org
/software
/gsl
/manual
/html_node
/>
3205 Tip
: search on google
: site
:http
://www.gnu.org
/software
/gsl
/manual
/html_node
/name_of_the_function_you_want
3209 This example computes the dilogarithm of
1/10 :
3211 use Math
::GSL
::SF qw
/dilog
/;
3212 my $x
= gsl_sf_dilog
(0.1);
3213 print
"gsl_sf_dilog(0.1) = $x\n";
3215 An example using Math
::GSL
::SF and gnuplot is in the B
<examples
/sf
> folder of the source code.
3219 Jonathan Leto
<jonathan@leto.net
> and Thierry Moisan
<thierry.moisan@gmail.com
>
3221 =head1 COPYRIGHT
AND LICENSE
3223 Copyright
(C
) 2008 Jonathan Leto and Thierry Moisan
3225 This program is free software
; you can redistribute it and
/or modify it
3226 under the same terms as Perl itself.