Logo ROOT   6.08/07
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 
289 // [5.2.1.16] Hermite polynomials
290 // (26.x.10)
291 
292 //double hermite(unsigned n, double x) {
293 //}
294 
295 
296 
297 
298 // [5.2.1.17] hypergeometric functions
299 // (26.x.13)
300 
301 double hyperg(double a, double b, double c, double x) {
302 
303  return gsl_sf_hyperg_2F1(a, b, c, x);
304 
305 }
306 
307 
308 
309 // [5.2.1.18] Laguerre polynomials
310 // (26.x.11)
311 
312 double laguerre(unsigned n, double x) {
313  return gsl_sf_laguerre_n(n, 0, x);
314 }
315 
316 
317 
318 
319 // [5.2.1.19] Legendre polynomials
320 // (26.x.7)
321 
322 double legendre(unsigned l, double x) {
323 
324  return gsl_sf_legendre_Pl(l, x);
325 
326 }
327 
328 
329 
330 // [5.2.1.20] Riemann zeta function
331 // (26.x.22)
332 
333 double riemann_zeta(double x) {
334 
335  return gsl_sf_zeta(x);
336 
337 }
338 
339 
340 
341 // [5.2.1.21] spherical Bessel functions of the first kind
342 // (26.x.5)
343 
344 double sph_bessel(unsigned n, double x) {
345 
346  return gsl_sf_bessel_jl(n, x);
347 
348 }
349 
350 
351 
352 // [5.2.1.22] spherical associated Legendre functions
353 // (26.x.9) ??????????
354 
355 double sph_legendre(unsigned l, unsigned m, double theta) {
356 
357  return gsl_sf_legendre_sphPlm(l, m, std::cos(theta));
358 
359 }
360 
361 
362 
363 
364 // [5.2.1.23] spherical Neumann functions
365 // (26.x.6)
366 
367 double sph_neumann(unsigned n, double x) {
368 
369  return gsl_sf_bessel_yl(n, x);
370 
371 }
372 
373 // Airy function Ai
374 
375 double airy_Ai(double x) {
376 
377  return gsl_sf_airy_Ai(x, GSL_PREC_DOUBLE);
378 
379 }
380 
381 // Airy function Bi
382 
383 double airy_Bi(double x) {
384 
385  return gsl_sf_airy_Bi(x, GSL_PREC_DOUBLE);
386 
387 }
388 
389 // Derivative of the Airy function Ai
390 
391 double airy_Ai_deriv(double x) {
392 
393  return gsl_sf_airy_Ai_deriv(x, GSL_PREC_DOUBLE);
394 
395 }
396 
397 // Derivative of the Airy function Bi
398 
399 double airy_Bi_deriv(double x) {
400 
401  return gsl_sf_airy_Bi_deriv(x, GSL_PREC_DOUBLE);
402 
403 }
404 
405 // s-th zero of the Airy function Ai
406 
407 double airy_zero_Ai(unsigned int s) {
408 
409  return gsl_sf_airy_zero_Ai(s);
410 
411 }
412 
413 // s-th zero of the Airy function Bi
414 
415 double airy_zero_Bi(unsigned int s) {
416 
417  return gsl_sf_airy_zero_Bi(s);
418 
419 }
420 
421 // s-th zero of the derivative of the Airy function Ai
422 
423 double airy_zero_Ai_deriv(unsigned int s) {
424 
425  return gsl_sf_airy_zero_Ai_deriv(s);
426 
427 }
428 
429 // s-th zero of the derivative of the Airy function Bi
430 
431 double airy_zero_Bi_deriv(unsigned int s) {
432 
433  return gsl_sf_airy_zero_Bi_deriv(s);
434 
435 }
436 
437 // wigner coefficient
438 
439 double wigner_3j(int ja, int jb, int jc, int ma, int mb, int mc) {
440  return gsl_sf_coupling_3j(ja,jb,jc,ma,mb,mc);
441 }
442 
443 double wigner_6j(int ja, int jb, int jc, int jd, int je, int jf) {
444  return gsl_sf_coupling_6j(ja,jb,jc,jd,je,jf);
445 }
446 
447 double wigner_9j(int ja, int jb, int jc, int jd, int je, int jf, int jg, int jh, int ji) {
448  return gsl_sf_coupling_9j(ja,jb,jc,jd,je,jf,jg,jh,ji);
449 }
450 
451 } // namespace Math
452 } // namespace ROOT
double airy_zero_Bi(unsigned int s)
Calculates the zeroes of the Airy function Bi.
This namespace contains pre-defined functions to be used in conjuction with TExecutor::Map and TExecu...
Definition: StringConv.hxx:21
double sph_legendre(unsigned l, unsigned m, double theta)
Computes the spherical (normalized) associated Legendre polynomials, or spherical harmonic without az...
return c
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.
TArc * a
Definition: textangle.C:12
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_t x[n]
Definition: legend1.C:17
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.
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.
TMarker * m
Definition: textangle.C:8
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.
TLine * l
Definition: textangle.C:4
double airy_zero_Ai(unsigned int s)
Calculates the zeroes of the Airy function Ai.
double airy_Bi(double x)
Calculates the Airy function Bi.
double laguerre(unsigned n, double x)
Calculates the Laguerre polynomials.
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.
const Int_t n
Definition: legend1.C:16
#define PI
double cyl_bessel_k(double nu, double x)
Calculates the modified Bessel functions of the second kind (also called irregular modified (cylindri...