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