Logo ROOT   6.13/01
Reference Guide
SpecFuncMathMore.cxx
Go to the documentation of this file.
1 // @(#)root/mathmore:$Id$
2 // Authors: L. Moneta, A. Zsenei 08/2005
3 
4 // Authors: Andras Zsenei & Lorenzo Moneta 06/2005
5 
6 /**********************************************************************
7  * *
8  * Copyright (c) 2005 , LCG ROOT MathLib Team *
9  * *
10  * *
11  **********************************************************************/
12 
13 #include <cmath>
14 
15 #ifndef PI
16 #define PI 3.14159265358979323846264338328 /* pi */
17 #endif
18 
19 
20 #include "gsl/gsl_sf_bessel.h"
21 #include "gsl/gsl_sf_legendre.h"
22 #include "gsl/gsl_sf_laguerre.h"
23 #include "gsl/gsl_sf_hyperg.h"
24 #include "gsl/gsl_sf_ellint.h"
25 //#include "gsl/gsl_sf_gamma.h"
26 #include "gsl/gsl_sf_expint.h"
27 #include "gsl/gsl_sf_zeta.h"
28 #include "gsl/gsl_sf_airy.h"
29 #include "gsl/gsl_sf_coupling.h"
30 
31 
32 namespace ROOT {
33 namespace Math {
34 
35 
36 
37 
38 // [5.2.1.1] associated Laguerre polynomials
39 // (26.x.12)
40 
41 double assoc_laguerre(unsigned n, double m, double x) {
42 
43  return gsl_sf_laguerre_n(n, m, x);
44 
45 }
46 
47 
48 
49 // [5.2.1.2] associated Legendre functions
50 // (26.x.8)
51 
52 double assoc_legendre(unsigned l, unsigned m, double x) {
53  // add (-1)^m
54  return (m%2 == 0) ? gsl_sf_legendre_Plm(l, m, x) : -gsl_sf_legendre_Plm(l, m, x);
55 
56 }
57 
58 
59 
60 
61 
62 // [5.2.1.4] (complete) elliptic integral of the first kind
63 // (26.x.15.2)
64 
65 double comp_ellint_1(double k) {
66 
67  return gsl_sf_ellint_Kcomp(k, GSL_PREC_DOUBLE);
68 
69 }
70 
71 
72 
73 // [5.2.1.5] (complete) elliptic integral of the second kind
74 // (26.x.16.2)
75 
76 double comp_ellint_2(double k) {
77 
78  return gsl_sf_ellint_Ecomp(k, GSL_PREC_DOUBLE);
79 
80 }
81 
82 
83 
84 // [5.2.1.6] (complete) elliptic integral of the third kind
85 // (26.x.17.2)
86 /**
87 Complete elliptic integral of the third kind.
88 
89 There are two different definitions used for the elliptic
90 integral of the third kind:
91 
92 \f[
93 P(\phi,k,n) = \int_0^\phi \frac{dt}{(1 + n \sin^2{t})\sqrt{1 - k^2 \sin^2{t}}}
94 \f]
95 
96 and
97 
98 \f[
99 P(\phi,k,n) = \int_0^\phi \frac{dt}{(1 - n \sin^2{t})\sqrt{1 - k^2 \sin^2{t}}}
100 \f]
101 
102 the former is adopted by
103 
104 - GSL
105  http://www.gnu.org/software/gsl/manual/gsl-ref_7.html#SEC95
106 
107 - Planetmath
108  http://planetmath.org/encyclopedia/EllipticIntegralsAndJacobiEllipticFunctions.html
109 
110 - CERNLIB
111  https://cern-tex.web.cern.ch/cern-tex/shortwrupsdir/c346/top.html
112 
113  while the latter is used by
114 
115 - Abramowitz and Stegun
116 
117 - Mathematica
118  http://mathworld.wolfram.com/EllipticIntegraloftheThirdKind.html
119 
120 - C++ standard
121  http://www.open-std.org/jtc1/sc22/wg21/docs/papers/2004/n1687.pdf
122 
123  in order to be C++ compliant, we decided to use the latter, hence the
124  change of the sign in the function call to GSL.
125 
126  */
127 
128 double comp_ellint_3(double n, double k) {
129 
130  return gsl_sf_ellint_P(PI/2.0, k, -n, GSL_PREC_DOUBLE);
131 
132 }
133 
134 
135 
136 // [5.2.1.7] confluent hypergeometric functions
137 // (26.x.14)
138 
139 double conf_hyperg(double a, double b, double z) {
140 
141  return gsl_sf_hyperg_1F1(a, b, z);
142 
143 }
144 
145 // confluent hypergeometric functions of second type
146 
147 double conf_hypergU(double a, double b, double z) {
148 
149  return gsl_sf_hyperg_U(a, b, z);
150 
151 }
152 
153 
154 
155 // [5.2.1.8] regular modified cylindrical Bessel functions
156 // (26.x.4.1)
157 
158 double cyl_bessel_i(double nu, double x) {
159 
160  return gsl_sf_bessel_Inu(nu, x);
161 
162 }
163 
164 
165 
166 // [5.2.1.9] cylindrical Bessel functions (of the first kind)
167 // (26.x.2)
168 
169 double cyl_bessel_j(double nu, double x) {
170 
171  return gsl_sf_bessel_Jnu(nu, x);
172 
173 }
174 
175 
176 
177 // [5.2.1.10] irregular modified cylindrical Bessel functions
178 // (26.x.4.2)
179 
180 double cyl_bessel_k(double nu, double x) {
181 
182  return gsl_sf_bessel_Knu(nu, x);
183 
184 }
185 
186 
187 
188 // [5.2.1.11] cylindrical Neumann functions;
189 // cylindrical Bessel functions (of the second kind)
190 // (26.x.3)
191 
192 double cyl_neumann(double nu, double x) {
193 
194  return gsl_sf_bessel_Ynu(nu, x);
195 
196 }
197 
198 
199 
200 // [5.2.1.12] (incomplete) elliptic integral of the first kind
201 // phi in radians
202 // (26.x.15.1)
203 
204 double ellint_1(double k, double phi) {
205 
206  return gsl_sf_ellint_F(phi, k, GSL_PREC_DOUBLE);
207 
208 }
209 
210 
211 
212 // [5.2.1.13] (incomplete) elliptic integral of the second kind
213 // phi in radians
214 // (26.x.16.1)
215 
216 double ellint_2(double k, double phi) {
217 
218  return gsl_sf_ellint_E(phi, k, GSL_PREC_DOUBLE);
219 
220 }
221 
222 
223 
224 // [5.2.1.14] (incomplete) elliptic integral of the third kind
225 // phi in radians
226 // (26.x.17.1)
227 /**
228 
229 Incomplete elliptic integral of the third kind.
230 
231 There are two different definitions used for the elliptic
232 integral of the third kind:
233 
234 \f[
235 P(\phi,k,n) = \int_0^\phi \frac{dt}{(1 + n \sin^2{t})\sqrt{1 - k^2 \sin^2{t}}}
236 \f]
237 
238 and
239 
240 \f[
241 P(\phi,k,n) = \int_0^\phi \frac{dt}{(1 - n \sin^2{t})\sqrt{1 - k^2 \sin^2{t}}}
242 \f]
243 
244 the former is adopted by
245 
246 - GSL
247  http://www.gnu.org/software/gsl/manual/gsl-ref_7.html#SEC95
248 
249 - Planetmath
250  http://planetmath.org/encyclopedia/EllipticIntegralsAndJacobiEllipticFunctions.html
251 
252 - CERNLIB
253  https://cern-tex.web.cern.ch/cern-tex/shortwrupsdir/c346/top.html
254 
255  while the latter is used by
256 
257 - Abramowitz and Stegun
258 
259 - Mathematica
260  http://mathworld.wolfram.com/EllipticIntegraloftheThirdKind.html
261 
262 - C++ standard
263  http://www.open-std.org/jtc1/sc22/wg21/docs/papers/2004/n1687.pdf
264 
265  in order to be C++ compliant, we decided to use the latter, hence the
266  change of the sign in the function call to GSL.
267 
268  */
269 
270 double ellint_3(double n, double k, double phi) {
271 
272  return gsl_sf_ellint_P(phi, k, -n, GSL_PREC_DOUBLE);
273 
274 }
275 
276 
277 
278 // [5.2.1.15] exponential integral
279 // (26.x.20)
280 
281 double expint(double x) {
282 
283  return gsl_sf_expint_Ei(x);
284 
285 }
286 
287 
288 // Generalization of expint(x)
289 //
290 double expint_n(int n, double x) {
291 
292  return gsl_sf_expint_En(n, x);
293 
294 }
295 
296 
297 
298 // [5.2.1.16] Hermite polynomials
299 // (26.x.10)
300 
301 //double hermite(unsigned n, double x) {
302 //}
303 
304 
305 
306 
307 // [5.2.1.17] hypergeometric functions
308 // (26.x.13)
309 
310 double hyperg(double a, double b, double c, double x) {
311 
312  return gsl_sf_hyperg_2F1(a, b, c, x);
313 
314 }
315 
316 
317 
318 // [5.2.1.18] Laguerre polynomials
319 // (26.x.11)
320 
321 double laguerre(unsigned n, double x) {
322  return gsl_sf_laguerre_n(n, 0, x);
323 }
324 
325 
326 
327 
328 // [5.2.1.19] Legendre polynomials
329 // (26.x.7)
330 
331 double legendre(unsigned l, double x) {
332 
333  return gsl_sf_legendre_Pl(l, x);
334 
335 }
336 
337 
338 
339 // [5.2.1.20] Riemann zeta function
340 // (26.x.22)
341 
342 double riemann_zeta(double x) {
343 
344  return gsl_sf_zeta(x);
345 
346 }
347 
348 
349 
350 // [5.2.1.21] spherical Bessel functions of the first kind
351 // (26.x.5)
352 
353 double sph_bessel(unsigned n, double x) {
354 
355  return gsl_sf_bessel_jl(n, x);
356 
357 }
358 
359 
360 
361 // [5.2.1.22] spherical associated Legendre functions
362 // (26.x.9) ??????????
363 
364 double sph_legendre(unsigned l, unsigned m, double theta) {
365 
366  return gsl_sf_legendre_sphPlm(l, m, std::cos(theta));
367 
368 }
369 
370 
371 
372 
373 // [5.2.1.23] spherical Neumann functions
374 // (26.x.6)
375 
376 double sph_neumann(unsigned n, double x) {
377 
378  return gsl_sf_bessel_yl(n, x);
379 
380 }
381 
382 // Airy function Ai
383 
384 double airy_Ai(double x) {
385 
386  return gsl_sf_airy_Ai(x, GSL_PREC_DOUBLE);
387 
388 }
389 
390 // Airy function Bi
391 
392 double airy_Bi(double x) {
393 
394  return gsl_sf_airy_Bi(x, GSL_PREC_DOUBLE);
395 
396 }
397 
398 // Derivative of the Airy function Ai
399 
400 double airy_Ai_deriv(double x) {
401 
402  return gsl_sf_airy_Ai_deriv(x, GSL_PREC_DOUBLE);
403 
404 }
405 
406 // Derivative of the Airy function Bi
407 
408 double airy_Bi_deriv(double x) {
409 
410  return gsl_sf_airy_Bi_deriv(x, GSL_PREC_DOUBLE);
411 
412 }
413 
414 // s-th zero of the Airy function Ai
415 
416 double airy_zero_Ai(unsigned int s) {
417 
418  return gsl_sf_airy_zero_Ai(s);
419 
420 }
421 
422 // s-th zero of the Airy function Bi
423 
424 double airy_zero_Bi(unsigned int s) {
425 
426  return gsl_sf_airy_zero_Bi(s);
427 
428 }
429 
430 // s-th zero of the derivative of the Airy function Ai
431 
432 double airy_zero_Ai_deriv(unsigned int s) {
433 
434  return gsl_sf_airy_zero_Ai_deriv(s);
435 
436 }
437 
438 // s-th zero of the derivative of the Airy function Bi
439 
440 double airy_zero_Bi_deriv(unsigned int s) {
441 
442  return gsl_sf_airy_zero_Bi_deriv(s);
443 
444 }
445 
446 // wigner coefficient
447 
448 double wigner_3j(int ja, int jb, int jc, int ma, int mb, int mc) {
449  return gsl_sf_coupling_3j(ja,jb,jc,ma,mb,mc);
450 }
451 
452 double wigner_6j(int ja, int jb, int jc, int jd, int je, int jf) {
453  return gsl_sf_coupling_6j(ja,jb,jc,jd,je,jf);
454 }
455 
456 double wigner_9j(int ja, int jb, int jc, int jd, int je, int jf, int jg, int jh, int ji) {
457  return gsl_sf_coupling_9j(ja,jb,jc,jd,je,jf,jg,jh,ji);
458 }
459 
460 } // namespace Math
461 } // namespace ROOT
double airy_zero_Bi(unsigned int s)
Calculates the zeroes of the Airy function Bi.
Namespace for new ROOT classes and functions.
Definition: TFoamSampler.h:19
double sph_legendre(unsigned l, unsigned m, double theta)
Computes the spherical (normalized) associated Legendre polynomials, or spherical harmonic without az...
double cyl_bessel_j(double nu, double x)
Calculates the (cylindrical) Bessel functions of the first kind (also called regular (cylindrical) Be...
double airy_zero_Ai_deriv(unsigned int s)
Calculates the zeroes of the derivative of the Airy function Ai.
double ellint_2(double k, double phi)
Calculates the complete elliptic integral of the second kind.
double wigner_6j(int two_ja, int two_jb, int two_jc, int two_jd, int two_je, int two_jf)
Calculates the Wigner 6j coupling coefficients.
double cos(double)
double legendre(unsigned l, double x)
Calculates the Legendre polynomials.
double comp_ellint_1(double k)
Calculates the complete elliptic integral of the first kind.
double sph_bessel(unsigned n, double x)
Calculates the spherical Bessel functions of the first kind (also called regular spherical Bessel fun...
double conf_hyperg(double a, double b, double z)
Calculates the confluent hypergeometric functions of the first kind.
double hyperg(double a, double b, double c, double x)
Calculates Gauss&#39; hypergeometric function.
double wigner_9j(int two_ja, int two_jb, int two_jc, int two_jd, int two_je, int two_jf, int two_jg, int two_jh, int two_ji)
Calculates the Wigner 9j coupling coefficients.
double airy_zero_Bi_deriv(unsigned int s)
Calculates the zeroes of the derivative of the Airy function Bi.
you should not use this method at all Int_t Int_t Double_t Double_t Double_t Int_t Double_t Double_t Double_t Double_t Int_t m
Definition: TRolke.cxx:637
double expint(double x)
Calculates the exponential integral.
double ellint_3(double n, double k, double phi)
Calculates the complete elliptic integral of the third kind.
double riemann_zeta(double x)
Calculates the Riemann zeta function.
double comp_ellint_2(double k)
Calculates the complete elliptic integral of the second kind.
double assoc_laguerre(unsigned n, double m, double x)
Computes the generalized Laguerre polynomials for .
double assoc_legendre(unsigned l, unsigned m, double x)
Computes the associated Legendre polynomials.
double airy_zero_Ai(unsigned int s)
Calculates the zeroes of the Airy function Ai.
double expint_n(int n, double x)
double airy_Bi(double x)
Calculates the Airy function Bi.
double laguerre(unsigned n, double x)
Calculates the Laguerre polynomials.
* x
Deprecated and error prone model selection interface.
Definition: TRolke.cxx:630
double airy_Bi_deriv(double x)
Calculates the derivative of the Airy function Bi.
double comp_ellint_3(double n, double k)
Calculates the complete elliptic integral of the third kind.
double airy_Ai_deriv(double x)
Calculates the derivative of the Airy function Ai.
double ellint_1(double k, double phi)
Calculates the incomplete elliptic integral of the first kind.
double cyl_neumann(double nu, double x)
Calculates the (cylindrical) Bessel functions of the second kind (also called irregular (cylindrical)...
Namespace for new Math classes and functions.
you should not use this method at all Int_t Int_t z
Definition: TRolke.cxx:630
double wigner_3j(int two_ja, int two_jb, int two_jc, int two_ma, int two_mb, int two_mc)
Calculates the Wigner 3j coupling coefficients.
you should not use this method at all Int_t Int_t Double_t Double_t Double_t Int_t Double_t Double_t Double_t Double_t b
Definition: TRolke.cxx:630
double sph_neumann(unsigned n, double x)
Calculates the spherical Bessel functions of the second kind (also called irregular spherical Bessel ...
double conf_hypergU(double a, double b, double z)
Calculates the confluent hypergeometric functions of the second kind, known also as Kummer function o...
double cyl_bessel_i(double nu, double x)
Calculates the modified Bessel function of the first kind (also called regular modified (cylindrical)...
double airy_Ai(double x)
Calculates the Airy function Ai.
#define PI
double cyl_bessel_k(double nu, double x)
Calculates the modified Bessel functions of the second kind (also called irregular modified (cylindri...