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