Make chart example more better
[Math-GSL.git] / SF.i
blobcb02021a83f6203f1ffa630ea629d973fa45313a
1 %module "Math::GSL::SF"
2 %include "typemaps.i"
4 %apply double *OUTPUT { double * sn, double * cn, double * dn, double * sgn };
6 %{
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"
35 #endif
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"
72 #endif
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"
82 %perlcode %{
84 @EXPORT_airy = qw/
85 gsl_sf_airy_Ai_e
86 gsl_sf_airy_Ai
87 gsl_sf_airy_Bi_e
88 gsl_sf_airy_Bi
89 gsl_sf_airy_Ai_scaled_e
90 gsl_sf_airy_Ai_scaled
91 gsl_sf_airy_Bi_scaled_e
92 gsl_sf_airy_Bi_scaled
93 gsl_sf_airy_Ai_deriv_e
94 gsl_sf_airy_Ai_deriv
95 gsl_sf_airy_Bi_deriv_e
96 gsl_sf_airy_Bi_deriv
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
102 gsl_sf_airy_zero_Ai
103 gsl_sf_airy_zero_Bi_e
104 gsl_sf_airy_zero_Bi
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
110 @EXPORT_bessel =qw/
111 gsl_sf_bessel_J0_e
112 gsl_sf_bessel_J0
113 gsl_sf_bessel_J1_e
114 gsl_sf_bessel_J1
115 gsl_sf_bessel_Jn_e
116 gsl_sf_bessel_Jn
117 gsl_sf_bessel_Jn_array
118 gsl_sf_bessel_Y0_e
119 gsl_sf_bessel_Y0
120 gsl_sf_bessel_Y1_e
121 gsl_sf_bessel_Y1
122 gsl_sf_bessel_Yn_e
123 gsl_sf_bessel_Yn
124 gsl_sf_bessel_Yn_array
125 gsl_sf_bessel_I0_e
126 gsl_sf_bessel_I0
127 gsl_sf_bessel_I1_e
128 gsl_sf_bessel_I1
129 gsl_sf_bessel_In_e
130 gsl_sf_bessel_In
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
139 gsl_sf_bessel_K0_e
140 gsl_sf_bessel_K0
141 gsl_sf_bessel_K1_e
142 gsl_sf_bessel_K1
143 gsl_sf_bessel_Kn_e
144 gsl_sf_bessel_Kn
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
153 gsl_sf_bessel_j0_e
154 gsl_sf_bessel_j0
155 gsl_sf_bessel_j1_e
156 gsl_sf_bessel_j1
157 gsl_sf_bessel_j2_e
158 gsl_sf_bessel_j2
159 gsl_sf_bessel_jl_e
160 gsl_sf_bessel_jl
161 gsl_sf_bessel_jl_array
162 gsl_sf_bessel_jl_steed_array
163 gsl_sf_bessel_y0_e
164 gsl_sf_bessel_y0
165 gsl_sf_bessel_y1_e
166 gsl_sf_bessel_y1
167 gsl_sf_bessel_y2_e
168 gsl_sf_bessel_y2
169 gsl_sf_bessel_yl_e
170 gsl_sf_bessel_yl
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
190 gsl_sf_bessel_Jnu_e
191 gsl_sf_bessel_Jnu
192 gsl_sf_bessel_Ynu_e
193 gsl_sf_bessel_Ynu
194 gsl_sf_bessel_sequence_Jnu_e
195 gsl_sf_bessel_Inu_scaled_e
196 gsl_sf_bessel_Inu_scaled
197 gsl_sf_bessel_Inu_e
198 gsl_sf_bessel_Inu
199 gsl_sf_bessel_Knu_scaled_e
200 gsl_sf_bessel_Knu_scaled
201 gsl_sf_bessel_Knu_e
202 gsl_sf_bessel_Knu
203 gsl_sf_bessel_lnKnu_e
204 gsl_sf_bessel_lnKnu
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/
213 gsl_sf_clausen_e
214 gsl_sf_clausen
216 @EXPORT_hydrogenic = qw/
217 gsl_sf_hydrogenicR_1_e
218 gsl_sf_hydrogenicR_1
219 gsl_sf_hydrogenicR_e
220 gsl_sf_hydrogenicR
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
228 gsl_sf_coulomb_CL_e
229 gsl_sf_coulomb_CL_array
231 @EXPORT_coupling = qw/
232 gsl_sf_coupling_3j_e
233 gsl_sf_coupling_3j
234 gsl_sf_coupling_6j_e
235 gsl_sf_coupling_6j
236 gsl_sf_coupling_RacahW_e
237 gsl_sf_coupling_RacahW
238 gsl_sf_coupling_9j_e
239 gsl_sf_coupling_9j
240 gsl_sf_coupling_6j_INCORRECT_e
241 gsl_sf_coupling_6j_INCORRECT
243 @EXPORT_dawson = qw/
244 gsl_sf_dawson_e
245 gsl_sf_dawson
247 @EXPORT_debye = qw/
248 gsl_sf_debye_1_e
249 gsl_sf_debye_1
250 gsl_sf_debye_2_e
251 gsl_sf_debye_2
252 gsl_sf_debye_3_e
253 gsl_sf_debye_3
254 gsl_sf_debye_4_e
255 gsl_sf_debye_4
256 gsl_sf_debye_5_e
257 gsl_sf_debye_5
258 gsl_sf_debye_6_e
259 gsl_sf_debye_6
261 @EXPORT_dilog = qw/
262 gsl_sf_dilog_e
263 gsl_sf_dilog
264 gsl_sf_complex_dilog_xy_e
265 gsl_sf_complex_dilog_e
268 @EXPORT_misc = qw/
269 gsl_sf_complex_spence_xy_e
270 gsl_sf_multiply_e
271 gsl_sf_multiply
272 gsl_sf_multiply_err_e
274 @EXPORT_elliptic = qw/
275 gsl_sf_ellint_Kcomp_e
276 gsl_sf_ellint_Kcomp
277 gsl_sf_ellint_Ecomp_e
278 gsl_sf_ellint_Ecomp
279 gsl_sf_ellint_Pcomp_e
280 gsl_sf_ellint_Pcomp
281 gsl_sf_ellint_Dcomp_e
282 gsl_sf_ellint_Dcomp
283 gsl_sf_ellint_F_e
284 gsl_sf_ellint_F
285 gsl_sf_ellint_E_e
286 gsl_sf_ellint_E
287 gsl_sf_ellint_P_e
288 gsl_sf_ellint_P
289 gsl_sf_ellint_D_e
290 gsl_sf_ellint_D
291 gsl_sf_ellint_RC_e
292 gsl_sf_ellint_RC
293 gsl_sf_ellint_RD_e
294 gsl_sf_ellint_RD
295 gsl_sf_ellint_RF_e
296 gsl_sf_ellint_RF
297 gsl_sf_ellint_RJ_e
298 gsl_sf_ellint_RJ
299 gsl_sf_elljac_e
301 @EXPORT_error = qw/
302 gsl_sf_erfc_e
303 gsl_sf_erfc
304 gsl_sf_log_erfc_e
305 gsl_sf_log_erfc
306 gsl_sf_erf_e
307 gsl_sf_erf
308 gsl_sf_erf_Z_e
309 gsl_sf_erf_Q_e
310 gsl_sf_erf_Z
311 gsl_sf_erf_Q
312 gsl_sf_hazard_e
313 gsl_sf_hazard
315 push @EXPORT_misc, qw/
316 gsl_sf_exp_e
317 gsl_sf_exp
318 gsl_sf_exp_e10_e
319 gsl_sf_exp_mult_e
320 gsl_sf_exp_mult
321 gsl_sf_exp_mult_e10_e
322 gsl_sf_expm1_e
323 gsl_sf_expm1
324 gsl_sf_exprel_e
325 gsl_sf_exprel
326 gsl_sf_exprel_2_e
327 gsl_sf_exprel_2
328 gsl_sf_exprel_n_e
329 gsl_sf_exprel_n
330 gsl_sf_exp_err_e
331 gsl_sf_exp_err_e10_e
332 gsl_sf_exp_mult_err_e
333 gsl_sf_exp_mult_err_e10_e
334 gsl_sf_expint_E1_e
335 gsl_sf_expint_E1
336 gsl_sf_expint_E2_e
337 gsl_sf_expint_E2
338 gsl_sf_expint_En_e
339 gsl_sf_expint_En
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
346 gsl_sf_expint_Ei_e
347 gsl_sf_expint_Ei
348 gsl_sf_expint_Ei_scaled_e
349 gsl_sf_expint_Ei_scaled
350 gsl_sf_Shi_e
351 gsl_sf_Shi
352 gsl_sf_Chi_e
353 gsl_sf_Chi
354 gsl_sf_expint_3_e
355 gsl_sf_expint_3
356 gsl_sf_Si_e
357 gsl_sf_Si
358 gsl_sf_Ci_e
359 gsl_sf_Ci
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
365 gsl_sf_fermi_dirac_0
366 gsl_sf_fermi_dirac_1_e
367 gsl_sf_fermi_dirac_1
368 gsl_sf_fermi_dirac_2_e
369 gsl_sf_fermi_dirac_2
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/
382 gsl_sf_legendre_Pl_e
383 gsl_sf_legendre_Pl
384 gsl_sf_legendre_Pl_array
385 gsl_sf_legendre_Pl_deriv_array
386 gsl_sf_legendre_P1_e
387 gsl_sf_legendre_P2_e
388 gsl_sf_legendre_P3_e
389 gsl_sf_legendre_P1
390 gsl_sf_legendre_P2
391 gsl_sf_legendre_P3
392 gsl_sf_legendre_Q0_e
393 gsl_sf_legendre_Q0
394 gsl_sf_legendre_Q1_e
395 gsl_sf_legendre_Q1
396 gsl_sf_legendre_Ql_e
397 gsl_sf_legendre_Ql
398 gsl_sf_legendre_Plm_e
399 gsl_sf_legendre_Plm
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
412 gsl_sf_legendre_H3d
413 gsl_sf_legendre_H3d_array
415 @EXPORT_gamma = qw/
416 gsl_sf_lngamma_e
417 gsl_sf_lngamma
418 gsl_sf_lngamma_sgn_e
419 gsl_sf_gamma_e
420 gsl_sf_gamma
421 gsl_sf_gammastar_e
422 gsl_sf_gammastar
423 gsl_sf_gammainv_e
424 gsl_sf_gammainv
425 gsl_sf_lngamma_complex_e
426 gsl_sf_gamma_inc_Q_e
427 gsl_sf_gamma_inc_Q
428 gsl_sf_gamma_inc_P_e
429 gsl_sf_gamma_inc_P
430 gsl_sf_gamma_inc_e
431 gsl_sf_gamma_inc
433 @EXPORT_factorial = qw/
434 gsl_sf_fact_e
435 gsl_sf_fact
436 gsl_sf_doublefact_e
437 gsl_sf_doublefact
438 gsl_sf_lnfact_e
439 gsl_sf_lnfact
440 gsl_sf_lndoublefact_e
441 gsl_sf_lndoublefact
443 @EXPORT_hypergeometric = qw/
444 gsl_sf_hyperg_0F1_e
445 gsl_sf_hyperg_0F1
446 gsl_sf_hyperg_1F1_int_e
447 gsl_sf_hyperg_1F1_int
448 gsl_sf_hyperg_1F1_e
449 gsl_sf_hyperg_1F1
450 gsl_sf_hyperg_U_int_e
451 gsl_sf_hyperg_U_int
452 gsl_sf_hyperg_U_int_e10_e
453 gsl_sf_hyperg_U_e
454 gsl_sf_hyperg_U
455 gsl_sf_hyperg_U_e10_e
456 gsl_sf_hyperg_2F1_e
457 gsl_sf_hyperg_2F1
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
464 gsl_sf_hyperg_2F0_e
465 gsl_sf_hyperg_2F0
467 @EXPORT_laguerre = qw/
468 gsl_sf_laguerre_1_e
469 gsl_sf_laguerre_2_e
470 gsl_sf_laguerre_3_e
471 gsl_sf_laguerre_1
472 gsl_sf_laguerre_2
473 gsl_sf_laguerre_3
474 gsl_sf_laguerre_n_e
475 gsl_sf_laguerre_n
477 push @EXPORT_misc, qw/
478 gsl_sf_taylorcoeff_e
479 gsl_sf_taylorcoeff
480 gsl_sf_lnchoose_e
481 gsl_sf_lnchoose
482 gsl_sf_choose_e
483 gsl_sf_choose
484 gsl_sf_lnpoch_e
485 gsl_sf_lnpoch
486 gsl_sf_lnpoch_sgn_e
487 gsl_sf_poch_e
488 gsl_sf_poch
489 gsl_sf_pochrel_e
490 gsl_sf_pochrel
491 gsl_sf_lnbeta_e
492 gsl_sf_lnbeta
493 gsl_sf_lnbeta_sgn_e
494 gsl_sf_beta_e
495 gsl_sf_beta
496 gsl_sf_beta_inc_e
497 gsl_sf_beta_inc
498 gsl_sf_gegenpoly_1_e
499 gsl_sf_gegenpoly_2_e
500 gsl_sf_gegenpoly_3_e
501 gsl_sf_gegenpoly_1
502 gsl_sf_gegenpoly_2
503 gsl_sf_gegenpoly_3
504 gsl_sf_gegenpoly_n_e
505 gsl_sf_gegenpoly_n
506 gsl_sf_gegenpoly_array
507 gsl_sf_lambert_W0_e
508 gsl_sf_lambert_W0
509 gsl_sf_lambert_Wm1_e
510 gsl_sf_lambert_Wm1
511 gsl_sf_conicalP_half_e
512 gsl_sf_conicalP_half
513 gsl_sf_conicalP_mhalf_e
514 gsl_sf_conicalP_mhalf
515 gsl_sf_conicalP_0_e
516 gsl_sf_conicalP_0
517 gsl_sf_conicalP_1_e
518 gsl_sf_conicalP_1
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
523 gsl_sf_log_e
524 gsl_sf_log
525 gsl_sf_log_abs_e
526 gsl_sf_log_abs
527 gsl_sf_complex_log_e
528 gsl_sf_log_1plusx_e
529 gsl_sf_log_1plusx
530 gsl_sf_log_1plusx_mx_e
531 gsl_sf_log_1plusx_mx
532 gsl_sf_pow_int_e
533 gsl_sf_pow_int
534 gsl_sf_psi_int_e
535 gsl_sf_psi_int
536 gsl_sf_psi_e
537 gsl_sf_psi
538 gsl_sf_psi_1piy_e
539 gsl_sf_psi_1piy
540 gsl_sf_complex_psi_e
541 gsl_sf_psi_1_int_e
542 gsl_sf_psi_1_int
543 gsl_sf_psi_1_e
544 gsl_sf_psi_1
545 gsl_sf_psi_n_e
546 gsl_sf_psi_n
547 gsl_sf_result_smash_e
548 gsl_sf_synchrotron_1_e
549 gsl_sf_synchrotron_1
550 gsl_sf_synchrotron_2_e
551 gsl_sf_synchrotron_2
553 @EXPORT_mathieu = qw/
554 gsl_sf_mathieu_a_array
555 gsl_sf_mathieu_b_array
556 gsl_sf_mathieu_a
557 gsl_sf_mathieu_b
558 gsl_sf_mathieu_a_coeff
559 gsl_sf_mathieu_b_coeff
560 gsl_sf_mathieu_alloc
561 gsl_sf_mathieu_free
562 gsl_sf_mathieu_ce
563 gsl_sf_mathieu_se
564 gsl_sf_mathieu_ce_array
565 gsl_sf_mathieu_se_array
566 gsl_sf_mathieu_Mc
567 gsl_sf_mathieu_Ms
568 gsl_sf_mathieu_Mc_array
569 gsl_sf_mathieu_Ms_array
571 @EXPORT_transport = qw/
572 gsl_sf_transport_2_e
573 gsl_sf_transport_2
574 gsl_sf_transport_3_e
575 gsl_sf_transport_3
576 gsl_sf_transport_4_e
577 gsl_sf_transport_4
578 gsl_sf_transport_5_e
579 gsl_sf_transport_5
581 @EXPORT_trig = qw/
582 gsl_sf_sin_e
583 gsl_sf_sin
584 gsl_sf_sin_pi_x_e
585 gsl_sf_cos_e
586 gsl_sf_cos_pi_x_e
587 gsl_sf_cos
588 gsl_sf_hypot_e
589 gsl_sf_hypot
590 gsl_sf_complex_sin_e
591 gsl_sf_complex_cos_e
592 gsl_sf_complex_logsin_e
593 gsl_sf_sinc_e
594 gsl_sf_sinc
595 gsl_sf_lnsinh_e
596 gsl_sf_lnsinh
597 gsl_sf_lncosh_e
598 gsl_sf_lncosh
599 gsl_sf_polar_to_rect
600 gsl_sf_rect_to_polar
601 gsl_sf_sin_err_e
602 gsl_sf_cos_err_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
609 gsl_sf_atanint_e
610 gsl_sf_atanint
612 @EXPORT_zeta = qw/
613 gsl_sf_zeta_int_e
614 gsl_sf_zeta_int
615 gsl_sf_zeta_e
616 gsl_sf_zeta
617 gsl_sf_zetam1_e
618 gsl_sf_zetam1
619 gsl_sf_zetam1_int_e
620 gsl_sf_zetam1_int
621 gsl_sf_hzeta_e
622 gsl_sf_hzeta
624 @EXPORT_eta = qw/
625 gsl_sf_eta_int_e
626 gsl_sf_eta_int
627 gsl_sf_eta_e
628 gsl_sf_eta
630 @EXPORT_vars = qw/
631 GSL_SF_GAMMA_XMAX
632 GSL_SF_FACT_NMAX
633 GSL_SF_DOUBLEFACT_NMAX
634 GSL_SF_MATHIEU_COEFF
637 @EXPORT_OK = (
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,
642 @EXPORT_vars
645 %EXPORT_TAGS = (
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 ],
672 __END__
674 =head1 NAME
676 Math::GSL::SF - Special Functions
678 =head1 SYNOPSIS
680 use Math::GSL::SF qw /:all/;
682 =head1 DESCRIPTION
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 :
690 my $val = $r->{val};
692 my $error = $r->{err};
694 Here is a list of all included functions:
696 =over
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.
704 =back
706 =over
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.
714 =back
716 =over
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.
724 =back
726 =over
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.
734 =back
736 =over
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.
744 =back
746 =over
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.
754 =back
756 =over
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.
764 =back
766 =over
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.
774 =back
776 =over
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.
784 =back
786 =over
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.
794 =back
796 =over
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.
804 =back
806 =over
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.
814 =back
816 =over
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.
824 =back
826 =over
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.
834 =back
836 =over
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).
844 =back
846 =over
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.
850 =back
852 =over
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.
860 =back
862 =over
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.
870 =back
872 =over
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.
880 =back
882 =over
884 =item C<gsl_sf_bessel_Yn_array>
888 =back
890 =over
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).
898 =back
900 =over
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).
908 =back
910 =over
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).
918 =back
920 =over
922 =item C<gsl_sf_bessel_In_array>
926 =back
928 =over
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).
936 =back
938 =over
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).
946 =back
948 =over
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)
956 =back
958 =over
960 =item C<gsl_sf_bessel_In_scaled_array>
964 =back
966 =over
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.
974 =back
976 =over
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.
984 =back
986 =over
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.
994 =back
996 =over
998 =item C<gsl_sf_bessel_Kn_array>
1002 =back
1004 =over
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.
1012 =back
1014 =over
1016 =item C<gsl_sf_bessel_K1_scaled_e($x, $result)>
1018 =item C<gsl_sf_bessel_K1_scaled($x)>
1022 =back
1024 =over
1026 =item C<gsl_sf_bessel_Kn_scaled_e($n, $x, $result)>
1028 =item C<gsl_sf_bessel_Kn_scaled($n, $x)>
1032 =back
1034 =over
1036 =item C<gsl_sf_bessel_Kn_scaled_array >
1040 =back
1042 =over
1044 =item C<gsl_sf_bessel_j0_e($x, $result)>
1046 =item C<gsl_sf_bessel_j0($x)>
1050 =back
1052 =over
1054 =item C<gsl_sf_bessel_j1_e($x, $result)>
1056 =item C<gsl_sf_bessel_j1($x)>
1060 =back
1062 =over
1064 =item C<gsl_sf_bessel_j2_e($x, $result)>
1066 =item C<gsl_sf_bessel_j2($x)>
1070 =back
1072 =over
1074 =item C<gsl_sf_bessel_jl_e($l, $x, $result)>
1076 =item C<gsl_sf_bessel_jl($l, $x)>
1080 =back
1082 =over
1084 =item C<gsl_sf_bessel_jl_array>
1088 =back
1090 =over
1092 =item C<gsl_sf_bessel_jl_steed_array>
1096 =back
1098 =over
1100 =item C<gsl_sf_bessel_y0_e($x, $result)>
1102 =item C<gsl_sf_bessel_y0($x)>
1106 =back
1108 =over
1110 =item C<gsl_sf_bessel_y1_e($x, $result)>
1112 =item C<gsl_sf_bessel_y1($x)>
1116 =back
1118 =over
1120 =item C<gsl_sf_bessel_y2_e($x, $result)>
1122 =item C<gsl_sf_bessel_y2($x)>
1126 =back
1128 =over
1130 =item C<gsl_sf_bessel_yl_e($l, $x, $result)>
1132 =item C<gsl_sf_bessel_yl($l, $x)>
1136 =back
1138 =over
1140 =item C<gsl_sf_bessel_yl_array>
1144 =back
1146 =over
1148 =item C<gsl_sf_bessel_i0_scaled_e($x, $result)>
1150 =item C<gsl_sf_bessel_i0_scaled($x)>
1154 =back
1156 =over
1158 =item C<gsl_sf_bessel_i1_scaled_e($x, $result)>
1160 =item C<gsl_sf_bessel_i1_scaled($x)>
1164 =back
1166 =over
1168 =item C<gsl_sf_bessel_i2_scaled_e($x, $result)>
1170 =item C<gsl_sf_bessel_i2_scaled($x)>
1174 =back
1176 =over
1178 =item C<gsl_sf_bessel_il_scaled_e($l, $x, $result)>
1180 =item C<gsl_sf_bessel_il_scaled($x)>
1184 =back
1186 =over
1188 =item C<gsl_sf_bessel_il_scaled_array>
1192 =back
1194 =over
1196 =item C<gsl_sf_bessel_k0_scaled_e($x, $result)>
1198 =item C<gsl_sf_bessel_k0_scale($x)>
1202 =back
1204 =over
1206 =item C<gsl_sf_bessel_k1_scaled_e($x, $result)>
1208 =item C<gsl_sf_bessel_k1_scaled($x)>
1212 =back
1214 =over
1216 =item C<gsl_sf_bessel_k2_scaled_e($x, $result) >
1218 =item C<gsl_sf_bessel_k2_scaled($x)>
1222 =back
1224 =over
1226 =item C<gsl_sf_bessel_kl_scaled_e($l, $x, $result)>
1228 =item C<gsl_sf_bessel_kl_scaled($l, $x)>
1232 =back
1234 =over
1236 =item C<gsl_sf_bessel_kl_scaled_array>
1240 =back
1242 =over
1244 =item C<gsl_sf_bessel_Jnu_e($nu, $x, $result)>
1246 =item C<gsl_sf_bessel_Jnu($nu, $x)>
1250 =back
1252 =over
1254 =item C<gsl_sf_bessel_sequence_Jnu_e >
1258 =back
1260 =over
1262 =item C<gsl_sf_bessel_Ynu_e($nu, $x, $result)>
1264 =item C<gsl_sf_bessel_Ynu($nu, $x)>
1268 =back
1270 =over
1272 =item C<gsl_sf_bessel_Inu_scaled_e($nu, $x, $result)>
1274 =item C<gsl_sf_bessel_Inu_scaled($nu, $x)>
1278 =back
1280 =over
1282 =item C<gsl_sf_bessel_Inu_e($nu, $x, $result)>
1284 =item C<gsl_sf_bessel_Inu($nu, $x)>
1288 =back
1290 =over
1292 =item C<gsl_sf_bessel_Knu_scaled_e($nu, $x, $result)>
1294 =item C<gsl_sf_bessel_Knu_scaled($nu, $x)>
1298 =back
1300 =over
1302 =item C<gsl_sf_bessel_Knu_e($nu, $x, $result)>
1304 =item C<gsl_sf_bessel_Knu($nu, $x)>
1308 =back
1310 =over
1312 =item C<gsl_sf_bessel_lnKnu_e($nu, $x, $result)>
1314 =item C<gsl_sf_bessel_lnKnu($nu, $x)>
1318 =back
1320 =over
1322 =item C<gsl_sf_bessel_zero_J0_e($s, $result)>
1324 =item C<gsl_sf_bessel_zero_J0($s)>
1328 =back
1330 =over
1332 =item C<gsl_sf_bessel_zero_J1_e($s, $result)>
1334 =item C<gsl_sf_bessel_zero_J1($s)>
1338 =back
1340 =over
1342 =item C<gsl_sf_bessel_zero_Jnu_e($nu, $s, $result)>
1344 =item C<gsl_sf_bessel_zero_Jnu($nu, $s)>
1348 =back
1350 =over
1352 =item C<gsl_sf_clausen_e($x, $result)>
1354 =item C<gsl_sf_clausen($x)>
1358 =back
1360 =over
1362 =item C<gsl_sf_hydrogenicR_1_e($Z, $r, $result)>
1364 =item C<gsl_sf_hydrogenicR_1($Z, $r)>
1368 =back
1370 =over
1372 =item C<gsl_sf_hydrogenicR_e($n, $l, $Z, $r, $result)>
1374 =item C<gsl_sf_hydrogenicR($n, $l, $Z, $r)>
1378 =back
1380 =over
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> -
1396 =back
1398 =over
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,
1405 (ja jb jc
1406 ma mb mc)
1407 where the arguments are given in half-integer units, ja = $two_ja/2, ma = $two_ma/2, etc.
1409 =back
1411 =over
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,
1418 {ja jb jc
1419 jd je jf}
1420 where the arguments are given in half-integer units, ja = $two_ja/2, ma = $two_ma/2, etc.
1422 =back
1424 =over
1426 =item C<gsl_sf_coupling_RacahW_e>
1428 =item C<gsl_sf_coupling_RacahW>
1432 =back
1434 =over
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,
1442 {ja jb jc
1443 jd je jf
1444 jg jh ji}
1445 where the arguments are given in half-integer units, ja = two_ja/2, ma = two_ma/2, etc.
1447 =back
1449 =over
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.
1457 =back
1459 =over
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)).
1467 =back
1469 =over
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)).
1477 =back
1479 =over
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)).
1487 =back
1489 =over
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)).
1497 =back
1499 =over
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)).
1507 =back
1509 =over
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)).
1517 =back
1519 =over
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).
1527 =back
1529 =over
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> -
1537 =back
1539 =over
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.
1549 =back
1551 =over
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.
1560 =back
1562 =over
1564 =item C<gsl_sf_ellint_Ecomp_e($k, $mode, $result)>
1566 =item C<gsl_sf_ellint_Ecomp($k, $mode)>
1570 =back
1572 =over
1574 =item C<gsl_sf_ellint_Pcomp_e($k, $n, $mode, $result)>
1576 =item C<gsl_sf_ellint_Pcomp($k, $n, $mode)>
1580 =back
1582 =over
1584 =item C<gsl_sf_ellint_Dcomp_e>
1586 =item C<gsl_sf_ellint_Dcomp >
1590 =back
1592 =over
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.
1600 =back
1602 =over
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.
1610 =back
1612 =over
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.
1620 =back
1622 =over
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.
1630 =back
1632 =over
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.
1640 =back
1642 =over
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.
1650 =back
1652 =over
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.
1660 =back
1662 =over
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.
1670 =back
1672 =over
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).
1682 =back
1684 =over
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)).
1692 =back
1694 =over
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).
1702 =back
1704 =over
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).
1712 =back
1714 =over
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.
1722 =back
1724 =over
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.
1732 =back
1734 =over
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.
1742 =back
1744 =over
1746 =item C<gsl_sf_exp_e10_e> -
1748 =back
1750 =over
1752 =item C<gsl_sf_exp_mult_e >
1754 =item C<gsl_sf_exp_mult>
1758 =back
1760 =over
1762 =item C<gsl_sf_exp_mult_e10_e> -
1764 =back
1766 =over
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.
1774 =back
1776 =over
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.
1784 =back
1786 =over
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.
1794 =back
1796 =over
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)) + ...
1805 = 1F1 (1,1+N,x)
1807 =back
1809 =over
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> -
1819 =back
1821 =over
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.
1829 =back
1831 =over
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.
1841 =back
1843 =over
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.
1853 =back
1855 =over
1857 =item C<gsl_sf_expint_E1_scaled_e >
1859 =item C<gsl_sf_expint_E1_scaled>
1863 =back
1865 =over
1867 =item C<gsl_sf_expint_E2_scaled_e>
1869 =item C<gsl_sf_expint_E2_scaled >
1873 =back
1875 =over
1877 =item C<gsl_sf_expint_En_scaled_e>
1879 =item C<gsl_sf_expint_En_scaled>
1883 =back
1885 =over
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.
1893 =back
1895 =over
1897 =item C<gsl_sf_expint_Ei_scaled_e>
1899 =item C<gsl_sf_expint_Ei_scaled >
1903 =back
1905 =over
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.
1913 =back
1915 =over
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).
1923 =back
1925 =over
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.
1933 =back
1935 =over
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.
1943 =back
1945 =over
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.
1953 =back
1955 =over
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).
1963 =back
1965 =over
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).
1973 =back
1975 =over
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)).
1983 =back
1985 =over
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)).
1993 =back
1995 =over
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)).
2003 =back
2005 =over
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).
2013 =back
2015 =over
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).
2023 =back
2025 =over
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).
2033 =back
2035 =over
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).
2043 =back
2045 =over
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
2053 =back
2055 =over
2057 =item C<gsl_sf_legendre_Pl_array>
2059 =item C<gsl_sf_legendre_Pl_deriv_array>
2063 =back
2065 =over
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.
2081 =back
2083 =over
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.
2091 =back
2093 =over
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.
2101 =back
2103 =over
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.
2111 =back
2113 =over
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.
2121 =back
2123 =over
2125 =item C<gsl_sf_legendre_Plm_array>
2127 =item C<gsl_sf_legendre_Plm_deriv_array >
2131 =back
2133 =over
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).
2141 =back
2143 =over
2145 =item C<gsl_sf_legendre_sphPlm_array >
2147 =item C<gsl_sf_legendre_sphPlm_deriv_array>
2151 =back
2153 =over
2155 =item C<gsl_sf_legendre_array_size> -
2157 =back
2159 =over
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.
2167 =back
2169 =over
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).
2173 =back
2175 =over
2177 =item C<gsl_sf_gamma_e >
2179 =item C<gsl_sf_gamma>
2183 =back
2185 =over
2187 =item C<gsl_sf_gammastar_e>
2189 =item C<gsl_sf_gammastar >
2193 =back
2195 =over
2197 =item C<gsl_sf_gammainv_e>
2199 =item C<gsl_sf_gammainv>
2203 =back
2205 =over
2207 =item C<gsl_sf_lngamma_complex_e >
2211 =back
2213 =over
2215 =item C<gsl_sf_gamma_inc_Q_e>
2217 =item C<gsl_sf_gamma_inc_Q>
2221 =back
2223 =over
2225 =item C<gsl_sf_gamma_inc_P_e >
2227 =item C<gsl_sf_gamma_inc_P>
2231 =back
2233 =over
2235 =item C<gsl_sf_gamma_inc_e>
2237 =item C<gsl_sf_gamma_inc >
2241 =back
2243 =over
2245 =item C<gsl_sf_taylorcoeff_e>
2247 =item C<gsl_sf_taylorcoeff>
2251 =back
2253 =over
2255 =item C<gsl_sf_fact_e >
2257 =item C<gsl_sf_fact>
2261 =back
2263 =over
2265 =item C<gsl_sf_doublefact_e>
2267 =item C<gsl_sf_doublefact >
2271 =back
2273 =over
2275 =item C<gsl_sf_lnfact_e>
2277 =item C<gsl_sf_lnfact>
2281 =back
2283 =over
2285 =item C<gsl_sf_lndoublefact_e >
2287 =item C<gsl_sf_lndoublefact>
2291 =back
2293 =over
2295 =item C<gsl_sf_lnchoose_e>
2297 =item C<gsl_sf_lnchoose >
2301 =back
2303 =over
2305 =item C<gsl_sf_choose_e>
2307 =item C<gsl_sf_choose>
2311 =back
2313 =over
2315 =item C<gsl_sf_lnpoch_e >
2317 =item C<gsl_sf_lnpoch>
2321 =back
2323 =over
2325 =item C<gsl_sf_lnpoch_sgn_e>
2329 =back
2331 =over
2333 =item C<gsl_sf_poch_e >
2335 =item C<gsl_sf_poch>
2339 =back
2341 =over
2343 =item C<gsl_sf_pochrel_e>
2345 =item C<gsl_sf_pochrel >
2349 =back
2351 =over
2353 =item C<gsl_sf_lnbeta_e>
2355 =item C<gsl_sf_lnbeta>
2359 =back
2361 =over
2363 =item C<gsl_sf_lnbeta_sgn_e >
2367 =back
2369 =over
2371 =item C<gsl_sf_beta_e>
2373 =item C<gsl_sf_beta>
2377 =back
2379 =over
2381 =item C<gsl_sf_beta_inc_e >
2383 =item C<gsl_sf_beta_inc>
2387 =back
2389 =over
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>
2405 =back
2407 =over
2409 =item C<gsl_sf_gegenpoly_n_e>
2411 =item C<gsl_sf_gegenpoly_n >
2415 =back
2417 =over
2419 =item C<gsl_sf_gegenpoly_array>
2421 =item C<gsl_sf_hyperg_0F1_e>
2423 =item C<gsl_sf_hyperg_0F1 >
2427 =back
2429 =over
2431 =item C<gsl_sf_hyperg_1F1_int_e>
2433 =item C<gsl_sf_hyperg_1F1_int>
2437 =back
2439 =over
2441 =item C<gsl_sf_hyperg_1F1_e >
2443 =item C<gsl_sf_hyperg_1F1>
2447 =back
2449 =over
2451 =item C<gsl_sf_hyperg_U_int_e>
2453 =item C<gsl_sf_hyperg_U_int >
2457 =back
2459 =over
2461 =item C<gsl_sf_hyperg_U_int_e10_e>
2465 =back
2467 =over
2469 =item C<gsl_sf_hyperg_U_e>
2471 =item C<gsl_sf_hyperg_U >
2475 =back
2477 =over
2479 =item C<gsl_sf_hyperg_U_e10_e>
2483 =back
2485 =over
2487 =item C<gsl_sf_hyperg_2F1_e>
2489 =item C<gsl_sf_hyperg_2F1 >
2493 =back
2495 =over
2497 =item C<gsl_sf_hyperg_2F1_conj_e>
2499 =item C<gsl_sf_hyperg_2F1_conj>
2503 =back
2505 =over
2507 =item C<gsl_sf_hyperg_2F1_renorm_e >
2509 =item C<gsl_sf_hyperg_2F1_renorm>
2513 =back
2515 =over
2517 =item C<gsl_sf_hyperg_2F1_conj_renorm_e>
2519 =item C<gsl_sf_hyperg_2F1_conj_renorm >
2523 =back
2525 =over
2527 =item C<gsl_sf_hyperg_2F0_e>
2529 =item C<gsl_sf_hyperg_2F0>
2533 =back
2535 =over
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>
2551 =back
2553 =over
2555 =item C<gsl_sf_laguerre_n_e >
2557 =item C<gsl_sf_laguerre_n>
2561 =back
2563 =over
2565 =item C<gsl_sf_lambert_W0_e>
2567 =item C<gsl_sf_lambert_W0 >
2571 =back
2573 =over
2575 =item C<gsl_sf_lambert_Wm1_e>
2577 =item C<gsl_sf_lambert_Wm1>
2581 =back
2583 =over
2585 =item C<gsl_sf_conicalP_half_e >
2587 =item C<gsl_sf_conicalP_half>
2591 =back
2593 =over
2595 =item C<gsl_sf_conicalP_mhalf_e>
2597 =item C<gsl_sf_conicalP_mhalf >
2601 =back
2603 =over
2605 =item C<gsl_sf_conicalP_0_e>
2607 =item C<gsl_sf_conicalP_0>
2611 =back
2613 =over
2615 =item C<gsl_sf_conicalP_1_e >
2617 =item C<gsl_sf_conicalP_1>
2621 =back
2623 =over
2625 =item C<gsl_sf_conicalP_sph_reg_e>
2627 =item C<gsl_sf_conicalP_sph_reg >
2631 =back
2633 =over
2635 =item C<gsl_sf_conicalP_cyl_reg_e>
2637 =item C<gsl_sf_conicalP_cyl_reg>
2641 =back
2643 =over
2645 =item C<gsl_sf_legendre_H3d_0_e >
2647 =item C<gsl_sf_legendre_H3d_0>
2651 =back
2653 =over
2655 =item C<gsl_sf_legendre_H3d_1_e>
2657 =item C<gsl_sf_legendre_H3d_1 >
2661 =back
2663 =over
2665 =item C<gsl_sf_legendre_H3d_e>
2667 =item C<gsl_sf_legendre_H3d>
2671 =back
2673 =over
2675 =item C<gsl_sf_legendre_H3d_array >
2679 =back
2681 =over
2683 =item C<gsl_sf_log_e>
2685 =item C<gsl_sf_log>
2689 =back
2691 =over
2693 =item C<gsl_sf_log_abs_e >
2695 =item C<gsl_sf_log_abs>
2699 =back
2701 =over
2703 =item C<gsl_sf_complex_log_e>
2707 =back
2709 =over
2711 =item C<gsl_sf_log_1plusx_e >
2713 =item C<gsl_sf_log_1plusx>
2717 =back
2719 =over
2721 =item C<gsl_sf_log_1plusx_mx_e>
2723 =item C<gsl_sf_log_1plusx_mx >
2727 =back
2729 =over
2731 =item C<gsl_sf_mathieu_a_array>
2733 =item C<gsl_sf_mathieu_b_array>
2737 =back
2739 =over
2741 =item C<gsl_sf_mathieu_a >
2743 =item C<gsl_sf_mathieu_b>
2747 =back
2749 =over
2751 =item C<gsl_sf_mathieu_a_coeff>
2753 =item C<gsl_sf_mathieu_b_coeff >
2757 =back
2759 =over
2761 =item C<gsl_sf_mathieu_alloc>
2765 =back
2767 =over
2769 =item C<gsl_sf_mathieu_free>
2773 =back
2775 =over
2777 =item C<gsl_sf_mathieu_ce >
2779 =item C<gsl_sf_mathieu_se>
2783 =back
2785 =over
2787 =item C<gsl_sf_mathieu_ce_array>
2789 =item C<gsl_sf_mathieu_se_array >
2793 =back
2795 =over
2797 =item C<gsl_sf_mathieu_Mc>
2799 =item C<gsl_sf_mathieu_Ms>
2803 =back
2805 =over
2807 =item C<gsl_sf_mathieu_Mc_array >
2809 =item C<gsl_sf_mathieu_Ms_array>
2813 =back
2815 =over
2817 =item C<gsl_sf_pow_int_e>
2819 =item C<gsl_sf_pow_int >
2823 =back
2825 =over
2827 =item C<gsl_sf_psi_int_e>
2829 =item C<gsl_sf_psi_int>
2833 =back
2835 =over
2837 =item C<gsl_sf_psi_e >
2839 =item C<gsl_sf_psi>
2843 =back
2845 =over
2847 =item C<gsl_sf_psi_1piy_e>
2849 =item C<gsl_sf_psi_1piy >
2853 =back
2855 =over
2857 =item C<gsl_sf_complex_psi_e>
2861 =back
2863 =over
2865 =item C<gsl_sf_psi_1_int_e>
2867 =item C<gsl_sf_psi_1_int >
2871 =back
2873 =over
2875 =item C<gsl_sf_psi_1_e >
2877 =item C<gsl_sf_psi_1>
2881 =back
2883 =over
2885 =item C<gsl_sf_psi_n_e >
2887 =item C<gsl_sf_psi_n>
2891 =back
2893 =over
2895 =item C<gsl_sf_result_smash_e>
2899 =back
2901 =over
2903 =item C<gsl_sf_synchrotron_1_e >
2905 =item C<gsl_sf_synchrotron_1>
2909 =back
2911 =over
2913 =item C<gsl_sf_synchrotron_2_e>
2915 =item C<gsl_sf_synchrotron_2 >
2919 =back
2921 =over
2923 =item C<gsl_sf_transport_2_e>
2925 =item C<gsl_sf_transport_2>
2929 =back
2931 =over
2933 =item C<gsl_sf_transport_3_e >
2935 =item C<gsl_sf_transport_3>
2939 =back
2941 =over
2943 =item C<gsl_sf_transport_4_e>
2945 =item C<gsl_sf_transport_4 >
2949 =back
2951 =over
2953 =item C<gsl_sf_transport_5_e>
2955 =item C<gsl_sf_transport_5>
2959 =back
2961 =over
2963 =item C<gsl_sf_sin_e >
2965 =item C<gsl_sf_sin>
2969 =back
2971 =over
2973 =item C<gsl_sf_cos_e>
2975 =item C<gsl_sf_cos >
2979 =back
2981 =over
2983 =item C<gsl_sf_hypot_e>
2985 =item C<gsl_sf_hypot>
2989 =back
2991 =over
2993 =item C<gsl_sf_complex_sin_e >
2997 =back
2999 =over
3001 =item C<gsl_sf_complex_cos_e>
3005 =back
3007 =over
3009 =item C<gsl_sf_complex_logsin_e>
3013 =back
3015 =over
3017 =item C<gsl_sf_sinc_e >
3019 =item C<gsl_sf_sinc>
3023 =back
3025 =over
3027 =item C<gsl_sf_lnsinh_e>
3029 =item C<gsl_sf_lnsinh >
3033 =back
3035 =over
3037 =item C<gsl_sf_lncosh_e>
3039 =item C<gsl_sf_lncosh>
3043 =back
3045 =over
3047 =item C<gsl_sf_polar_to_rect >
3051 =back
3053 =over
3055 =item C<gsl_sf_rect_to_polar>
3059 =back
3061 =over
3063 =item C<gsl_sf_sin_err_e>
3065 =item C<gsl_sf_cos_err_e >
3069 =back
3071 =over
3073 =item C<gsl_sf_angle_restrict_symm_e>
3075 =item C<gsl_sf_angle_restrict_symm>
3079 =back
3081 =over
3083 =item C<gsl_sf_angle_restrict_pos_e >
3085 =item C<gsl_sf_angle_restrict_pos>
3089 =back
3091 =over
3093 =item C<gsl_sf_angle_restrict_symm_err_e>
3095 =item C<gsl_sf_angle_restrict_pos_err_e >
3097 =over
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.
3105 =back
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 >
3133 =back
3135 This module also contains the following constants used as mode in various of those functions :
3137 =over
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.
3145 =back
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 :
3154 =over
3156 =item C<airy>
3158 =item C<bessel>
3160 =item C<clausen>
3162 =item C<hydrogenic>
3164 =item C<coulumb>
3166 =item C<coupling>
3168 =item C<dawson>
3170 =item C<debye>
3172 =item C<dilog>
3174 =item C<factorial>
3176 =item C<misc>
3178 =item C<elliptic>
3180 =item C<error>
3182 =item C<hypergeometric>
3184 =item C<laguerre>
3186 =item C<legendre>
3188 =item C<gamma>
3190 =item C<transport>
3192 =item C<trig>
3194 =item C<zeta>
3196 =item C<eta>
3198 =item C<vars>
3200 =back
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
3207 =head1 EXAMPLES
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.
3217 =head1 AUTHORS
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.
3228 =cut