Logo ROOT   6.14/05
Reference Guide
SpecFuncMathMore.h
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) 2004 ROOT Foundation, CERN/PH-SFT *
9  * *
10  * This library is free software; you can redistribute it and/or *
11  * modify it under the terms of the GNU General Public License *
12  * as published by the Free Software Foundation; either version 2 *
13  * of the License, or (at your option) any later version. *
14  * *
15  * This library is distributed in the hope that it will be useful, *
16  * but WITHOUT ANY WARRANTY; without even the implied warranty of *
17  * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU *
18  * General Public License for more details. *
19  * *
20  * You should have received a copy of the GNU General Public License *
21  * along with this library (see file COPYING); if not, write *
22  * to the Free Software Foundation, Inc., 59 Temple Place, Suite *
23  * 330, Boston, MA 02111-1307 USA, or contact the author. *
24  * *
25  **********************************************************************/
26 
27 /**
28 
29 Special mathematical functions.
30 The naming and numbering of the functions is taken from
31 <A HREF="http://www.open-std.org/jtc1/sc22/wg21/docs/papers/2004/n1687.pdf">
32 Matt Austern,
33 (Draft) Technical Report on Standard Library Extensions,
34 N1687=04-0127, September 10, 2004</A>
35 
36 @author Created by Andras Zsenei on Mon Nov 8 2004
37 
38 @defgroup SpecFunc Special functions
39 
40 */
41 
42 
43 
44 
45 
46 #ifndef ROOT_Math_SpecFuncMathMore
47 #define ROOT_Math_SpecFuncMathMore
48 
49 
50 
51 
52 namespace ROOT {
53 namespace Math {
54 
55  /** @name Special Functions from MathMore */
56 
57 
58  /**
59 
60 
61  Computes the generalized Laguerre polynomials for
62  \f$ n \geq 0, m > -1 \f$.
63  They are defined in terms of the confluent hypergeometric function.
64  For integer values of m they can be defined in terms of the Laguerre polynomials \f$L_n(x)\f$:
65 
66  \f[ L_{n}^{m}(x) = (-1)^{m} \frac{d^m}{dx^m} L_{n+m}(x) \f]
67 
68 
69  For detailed description see
70  <A HREF="http://mathworld.wolfram.com/LaguerrePolynomial.html">Mathworld</A>.
71  The implementation used is that of
72  <A HREF="http://www.gnu.org/software/gsl/manual/html_node/Laguerre-Functions.html">GSL</A>.
73 
74  This function is an extension of C++0x, also consistent in GSL,
75  Abramowitz and Stegun 1972 and MatheMathica that uses non-integer values for m.
76  C++0x calls for 'int m', more restrictive than necessary.
77  The definition for was incorrect in 'n1687.pdf', but fixed in
78  <A HREF="http://www.open-std.org/jtc1/sc22/wg21/docs/papers/2005/n1836.pdf">n1836.pdf</A>,
79  the most recent draft as of 2007-07-01
80 
81 
82  @ingroup SpecFunc
83 
84  */
85  // [5.2.1.1] associated Laguerre polynomials
86 
87  double assoc_laguerre(unsigned n, double m, double x);
88 
89 
90 
91 
92  /**
93 
94  Computes the associated Legendre polynomials.
95 
96  \f[ P_{l}^{m}(x) = (1-x^2)^{m/2} \frac{d^m}{dx^m} P_{l}(x) \f]
97 
98  with \f$m \geq 0\f$, \f$ l \geq m \f$ and \f$ |x|<1 \f$.
99  There are two sign conventions for associated Legendre polynomials.
100  As is the case with the above formula, some authors (e.g., Arfken
101  1985, pp. 668-669) omit the Condon-Shortley phase \f$(-1)^m\f$,
102  while others include it (e.g., Abramowitz and Stegun 1972).
103  One possible way to distinguish the two conventions is due to
104  Abramowitz and Stegun (1972, p. 332), who use the notation
105 
106  \f[ P_{lm} (x) = (-1)^m P_{l}^{m} (x)\f]
107 
108  to distinguish the two. For detailed description see
109  <A HREF="http://mathworld.wolfram.com/LegendrePolynomial.html">
110  Mathworld</A>. The implementation used is that of
111  <A HREF="http://www.gnu.org/software/gsl/manual/html_node/Associated-Legendre-Polynomials-and-Spherical-Harmonics.html">GSL</A>.
112 
113  The definition uses is the one of C++0x, \f$ P_{lm}\f$, while GSL and MatheMatica use the \f$P_{l}^{m}\f$ definition. Note that C++0x and GSL definitions agree instead for the normalized associated Legendre polynomial,
114  sph_legendre(l,m,theta).
115 
116  @ingroup SpecFunc
117 
118  */
119  // [5.2.1.2] associated Legendre functions
120 
121  double assoc_legendre(unsigned l, unsigned m, double x);
122 
123 
124 
125 
126 
127  /**
128 
129  Calculates the complete elliptic integral of the first kind.
130 
131  \f[ K(k) = F(k, \pi / 2) = \int_{0}^{\pi /2} \frac{d \theta}{\sqrt{1 - k^2 \sin^2{\theta}}} \f]
132 
133  with \f$0 \leq k^2 \leq 1\f$. For detailed description see
134  <A HREF="http://mathworld.wolfram.com/CompleteEllipticIntegraloftheFirstKind.html">
135  Mathworld</A>. The implementation used is that of
136  <A HREF="http://www.gnu.org/software/gsl/manual/gsl-ref_7.html#SEC100">GSL</A>.
137 
138  @ingroup SpecFunc
139 
140  */
141  // [5.2.1.4] (complete) elliptic integral of the first kind
142 
143  double comp_ellint_1(double k);
144 
145 
146 
147 
148  /**
149 
150  Calculates the complete elliptic integral of the second kind.
151 
152  \f[ E(k) = E(k , \pi / 2) = \int_{0}^{\pi /2} \sqrt{1 - k^2 \sin^2{\theta}} d \theta \f]
153 
154  with \f$0 \leq k^2 \leq 1\f$. For detailed description see
155  <A HREF="http://mathworld.wolfram.com/CompleteEllipticIntegraloftheSecondKind.html">
156  Mathworld</A>. The implementation used is that of
157  <A HREF="http://www.gnu.org/software/gsl/manual/gsl-ref_7.html#SEC100">GSL</A>.
158 
159  @ingroup SpecFunc
160 
161  */
162  // [5.2.1.5] (complete) elliptic integral of the second kind
163 
164  double comp_ellint_2(double k);
165 
166 
167 
168 
169  /**
170 
171  Calculates the complete elliptic integral of the third kind.
172 
173  \f[ \Pi (n, k, \pi / 2) = \int_{0}^{\pi /2} \frac{d \theta}{(1 - n \sin^2{\theta})\sqrt{1 - k^2 \sin^2{\theta}}} \f]
174 
175  with \f$0 \leq k^2 \leq 1\f$. There are two sign conventions
176  for elliptic integrals of the third kind. Some authors (Abramowitz
177  and Stegun,
178  <A HREF="http://mathworld.wolfram.com/EllipticIntegraloftheThirdKind.html">
179  Mathworld</A>,
180  <A HREF="http://www.open-std.org/jtc1/sc22/wg21/docs/papers/2004/n1687.pdf">
181  C++ standard proposal</A>) use the above formula, while others
182  (<A HREF="http://www.gnu.org/software/gsl/manual/gsl-ref_7.html#SEC95">
183  GSL</A>, <A HREF="http://planetmath.org/encyclopedia/EllipticIntegralsAndJacobiEllipticFunctions.html">
184  Planetmath</A> and
185  <A HREF="https://cern-tex.web.cern.ch/cern-tex/shortwrupsdir/c346/top.html">
186  CERNLIB</A>) use the + sign in front of n in the denominator. In
187  order to be C++ compliant, the present library uses the former
188  convention. The implementation used is that of
189  <A HREF="http://www.gnu.org/software/gsl/manual/gsl-ref_7.html#SEC101">GSL</A>.
190 
191  @ingroup SpecFunc
192 
193  */
194  // [5.2.1.6] (complete) elliptic integral of the third kind
195  double comp_ellint_3(double n, double k);
196 
197 
198 
199 
200  /**
201 
202  Calculates the confluent hypergeometric functions of the first kind.
203 
204  \f[ _{1}F_{1}(a;b;z) = \frac{\Gamma(b)}{\Gamma(a)} \sum_{n=0}^{\infty} \frac{\Gamma(a+n)}{\Gamma(b+n)} \frac{z^n}{n!} \f]
205 
206  For detailed description see
207  <A HREF="http://mathworld.wolfram.com/ConfluentHypergeometricFunctionoftheFirstKind.html">
208  Mathworld</A>. The implementation used is that of
209  <A HREF="http://www.gnu.org/software/gsl/manual/gsl-ref_7.html#SEC125">GSL</A>.
210 
211  @ingroup SpecFunc
212 
213  */
214  // [5.2.1.7] confluent hypergeometric functions
215 
216  double conf_hyperg(double a, double b, double z);
217 
218 
219  /**
220 
221  Calculates the confluent hypergeometric functions of the second kind, known also as Kummer function of the second kind,
222  it is related to the confluent hypergeometric functions of the first kind.
223 
224  \f[ U(a,b,z) = \frac{ \pi}{ \sin{\pi b } } \left[ \frac{ _{1}F_{1}(a,b,z) } {\Gamma(a-b+1) }
225  - \frac{ z^{1-b} { _{1}F_{1}}(a-b+1,2-b,z)}{\Gamma(a)} \right] \f]
226 
227  For detailed description see
228  <A HREF="http://mathworld.wolfram.com/ConfluentHypergeometricFunctionoftheSecondKind.html">
229  Mathworld</A>. The implementation used is that of
230  <A HREF="http://www.gnu.org/software/gsl/manual/gsl-ref_7.html#SEC125">GSL</A>.
231  This function is not part of the C++ standard proposal
232 
233  @ingroup SpecFunc
234 
235  */
236  // confluent hypergeometric functions of second type
237 
238  double conf_hypergU(double a, double b, double z);
239 
240 
241 
242  /**
243 
244  Calculates the modified Bessel function of the first kind
245  (also called regular modified (cylindrical) Bessel function).
246 
247  \f[ I_{\nu} (x) = i^{-\nu} J_{\nu}(ix) = \sum_{k=0}^{\infty} \frac{(\frac{1}{2}x)^{\nu + 2k}}{k! \Gamma(\nu + k + 1)} \f]
248 
249  for \f$x>0, \nu > 0\f$. For detailed description see
250  <A HREF="http://mathworld.wolfram.com/ModifiedBesselFunctionoftheFirstKind.html">
251  Mathworld</A>. The implementation used is that of
252  <A HREF="http://www.gnu.org/software/gsl/manual/gsl-ref_7.html#SEC71">GSL</A>.
253 
254  @ingroup SpecFunc
255 
256  */
257  // [5.2.1.8] regular modified cylindrical Bessel functions
258 
259  double cyl_bessel_i(double nu, double x);
260 
261 
262 
263 
264  /**
265 
266  Calculates the (cylindrical) Bessel functions of the first kind (also
267  called regular (cylindrical) Bessel functions).
268 
269  \f[ J_{\nu} (x) = \sum_{k=0}^{\infty} \frac{(-1)^k(\frac{1}{2}x)^{\nu + 2k}}{k! \Gamma(\nu + k + 1)} \f]
270 
271  For detailed description see
272  <A HREF="http://mathworld.wolfram.com/BesselFunctionoftheFirstKind.html">
273  Mathworld</A>. The implementation used is that of
274  <A HREF="http://www.gnu.org/software/gsl/manual/gsl-ref_7.html#SEC69">GSL</A>.
275 
276  @ingroup SpecFunc
277 
278  */
279  // [5.2.1.9] cylindrical Bessel functions (of the first kind)
280 
281  double cyl_bessel_j(double nu, double x);
282 
283 
284 
285 
286 
287  /**
288 
289  Calculates the modified Bessel functions of the second kind
290  (also called irregular modified (cylindrical) Bessel functions).
291 
292  \f[ K_{\nu} (x) = \frac{\pi}{2} i^{\nu + 1} (J_{\nu} (ix) + iN(ix)) = \left\{ \begin{array}{cl} \frac{\pi}{2} \frac{I_{-\nu}(x) - I_{\nu}(x)}{\sin{\nu \pi}} & \mbox{for non-integral $\nu$} \\ \frac{\pi}{2} \lim{\mu \to \nu} \frac{I_{-\mu}(x) - I_{\mu}(x)}{\sin{\mu \pi}}
293 & \mbox{for integral $\nu$} \end{array} \right. \f]
294 
295  for \f$x>0, \nu > 0\f$. For detailed description see
296  <A HREF="http://mathworld.wolfram.com/ModifiedBesselFunctionoftheSecondKind.html">
297  Mathworld</A>. The implementation used is that of
298  <A HREF="http://www.gnu.org/software/gsl/manual/gsl-ref_7.html#SEC72">GSL</A>.
299 
300  @ingroup SpecFunc
301 
302  */
303  // [5.2.1.10] irregular modified cylindrical Bessel functions
304 
305  double cyl_bessel_k(double nu, double x);
306 
307 
308 
309 
310  /**
311 
312  Calculates the (cylindrical) Bessel functions of the second kind
313  (also called irregular (cylindrical) Bessel functions or
314  (cylindrical) Neumann functions).
315 
316  \f[ N_{\nu} (x) = Y_{\nu} (x) = \left\{ \begin{array}{cl} \frac{J_{\nu} \cos{\nu \pi}-J_{-\nu}(x)}{\sin{\nu \pi}} & \mbox{for non-integral $\nu$} \\ \lim{\mu \to \nu} \frac{J_{\mu} \cos{\mu \pi}-J_{-\mu}(x)}{\sin{\mu \pi}} & \mbox{for integral $\nu$} \end{array} \right. \f]
317 
318  For detailed description see
319  <A HREF="http://mathworld.wolfram.com/BesselFunctionoftheSecondKind.html">
320  Mathworld</A>. The implementation used is that of
321  <A HREF="http://www.gnu.org/software/gsl/manual/gsl-ref_7.html#SEC70">GSL</A>.
322 
323  @ingroup SpecFunc
324 
325  */
326  // [5.2.1.11] cylindrical Neumann functions;
327  // cylindrical Bessel functions (of the second kind)
328 
329  double cyl_neumann(double nu, double x);
330 
331 
332 
333 
334  /**
335 
336  Calculates the incomplete elliptic integral of the first kind.
337 
338  \f[ F(k, \phi) = \int_{0}^{\phi} \frac{d \theta}{\sqrt{1 - k^2 \sin^2{\theta}}} \f]
339 
340  with \f$0 \leq k^2 \leq 1\f$. For detailed description see
341  <A HREF="http://mathworld.wolfram.com/EllipticIntegraloftheFirstKind.html">
342  Mathworld</A>. The implementation used is that of
343  <A HREF="http://www.gnu.org/software/gsl/manual/gsl-ref_7.html#SEC101">GSL</A>.
344 
345  @param k
346  @param phi angle in radians
347 
348  @ingroup SpecFunc
349 
350  */
351  // [5.2.1.12] (incomplete) elliptic integral of the first kind
352  // phi in radians
353 
354  double ellint_1(double k, double phi);
355 
356 
357 
358 
359  /**
360 
361  Calculates the complete elliptic integral of the second kind.
362 
363  \f[ E(k , \phi) = \int_{0}^{\phi} \sqrt{1 - k^2 \sin^2{\theta}} d \theta \f]
364 
365  with \f$0 \leq k^2 \leq 1\f$. For detailed description see
366  <A HREF="http://mathworld.wolfram.com/EllipticIntegraloftheSecondKind.html">
367  Mathworld</A>. The implementation used is that of
368  <A HREF="http://www.gnu.org/software/gsl/manual/gsl-ref_7.html#SEC101">GSL</A>.
369 
370  @param k
371  @param phi angle in radians
372 
373  @ingroup SpecFunc
374 
375  */
376  // [5.2.1.13] (incomplete) elliptic integral of the second kind
377  // phi in radians
378 
379  double ellint_2(double k, double phi);
380 
381 
382 
383 
384  /**
385 
386  Calculates the complete elliptic integral of the third kind.
387 
388  \f[ \Pi (n, k, \phi) = \int_{0}^{\phi} \frac{d \theta}{(1 - n \sin^2{\theta})\sqrt{1 - k^2 \sin^2{\theta}}} \f]
389 
390  with \f$0 \leq k^2 \leq 1\f$. There are two sign conventions
391  for elliptic integrals of the third kind. Some authors (Abramowitz
392  and Stegun,
393  <A HREF="http://mathworld.wolfram.com/EllipticIntegraloftheThirdKind.html">
394  Mathworld</A>,
395  <A HREF="http://www.open-std.org/jtc1/sc22/wg21/docs/papers/2004/n1687.pdf">
396  C++ standard proposal</A>) use the above formula, while others
397  (<A HREF="http://www.gnu.org/software/gsl/manual/gsl-ref_7.html#SEC95">
398  GSL</A>, <A HREF="http://planetmath.org/encyclopedia/EllipticIntegralsAndJacobiEllipticFunctions.html">
399  Planetmath</A> and
400  <A HREF="https://cern-tex.web.cern.ch/cern-tex/shortwrupsdir/c346/top.html">
401  CERNLIB</A>) use the + sign in front of n in the denominator. In
402  order to be C++ compliant, the present library uses the former
403  convention. The implementation used is that of
404  <A HREF="http://www.gnu.org/software/gsl/manual/gsl-ref_7.html#SEC101">GSL</A>.
405 
406  @param n
407  @param k
408  @param phi angle in radians
409 
410  @ingroup SpecFunc
411 
412  */
413  // [5.2.1.14] (incomplete) elliptic integral of the third kind
414  // phi in radians
415 
416  double ellint_3(double n, double k, double phi);
417 
418 
419 
420 
421  /**
422 
423  Calculates the exponential integral.
424 
425  \f[ Ei(x) = - \int_{-x}^{\infty} \frac{e^{-t}}{t} dt \f]
426 
427  For detailed description see
428  <A HREF="http://mathworld.wolfram.com/ExponentialIntegral.html">
429  Mathworld</A>. The implementation used is that of
430  <A HREF="http://www.gnu.org/software/gsl/manual/gsl-ref_7.html#SEC115">GSL</A>.
431 
432  @ingroup SpecFunc
433 
434  */
435  // [5.2.1.15] exponential integral
436 
437  double expint(double x);
438  double expint_n(int n, double x);
439 
440 
441 
442  // [5.2.1.16] Hermite polynomials
443 
444  //double hermite(unsigned n, double x);
445 
446 
447 
448 
449 
450  /**
451 
452  Calculates Gauss' hypergeometric function.
453 
454  \f[ _{2}F_{1}(a,b;c;x) = \frac{\Gamma(c)}{\Gamma(a) \Gamma(b)} \sum_{n=0}^{\infty} \frac{\Gamma(a+n)\Gamma(b+n)}{\Gamma(c+n)} \frac{x^n}{n!} \f]
455 
456  For detailed description see
457  <A HREF="http://mathworld.wolfram.com/HypergeometricFunction.html">
458  Mathworld</A>. The implementation used is that of
459  <A HREF="http://www.gnu.org/software/gsl/manual/gsl-ref_7.html#SEC125">GSL</A>.
460 
461  @ingroup SpecFunc
462 
463  */
464  // [5.2.1.17] hypergeometric functions
465 
466  double hyperg(double a, double b, double c, double x);
467 
468 
469 
470  /**
471 
472  Calculates the Laguerre polynomials
473 
474  \f[ P_{l}(x) = \frac{ e^x}{n!} \frac{d^n}{dx^n} (x^n - e^{-x}) \f]
475 
476  for \f$x \geq 0 \f$ in the Rodrigues representation.
477  They corresponds to the associated Laguerre polynomial of order m=0.
478  See Abramowitz and Stegun, (22.5.16)
479  For detailed description see
480  <A HREF="http://mathworld.wolfram.com/LaguerrePolynomial.html">
481  Mathworld</A>.
482  The are implemented using the associated Laguerre polynomial of order m=0.
483 
484  @ingroup SpecFunc
485 
486  */
487  // [5.2.1.18] Laguerre polynomials
488 
489  double laguerre(unsigned n, double x);
490 
491 
492  /**
493 
494  Calculates the Legendre polynomials.
495 
496  \f[ P_{l}(x) = \frac{1}{2^l l!} \frac{d^l}{dx^l} (x^2 - 1)^l \f]
497 
498  for \f$l \geq 0, |x|\leq1\f$ in the Rodrigues representation.
499  For detailed description see
500  <A HREF="http://mathworld.wolfram.com/LegendrePolynomial.html">
501  Mathworld</A>. The implementation used is that of
502  <A HREF="http://www.gnu.org/software/gsl/manual/gsl-ref_7.html#SEC129">GSL</A>.
503 
504  @ingroup SpecFunc
505 
506  */
507  // [5.2.1.19] Legendre polynomials
508 
509  double legendre(unsigned l, double x);
510 
511 
512 
513 
514  /**
515 
516  Calculates the Riemann zeta function.
517 
518  \f[ \zeta (x) = \left\{ \begin{array}{cl} \sum_{k=1}^{\infty}k^{-x} & \mbox{for $x > 1$} \\ 2^x \pi^{x-1} \sin{(\frac{1}{2}\pi x)} \Gamma(1-x) \zeta (1-x) & \mbox{for $x < 1$} \end{array} \right. \f]
519 
520  For detailed description see
521  <A HREF="http://mathworld.wolfram.com/RiemannZetaFunction.html">
522  Mathworld</A>. The implementation used is that of
523  <A HREF="http://www.gnu.org/software/gsl/manual/gsl-ref_7.html#SEC149">GSL</A>.
524 
525  CHECK WHETHER THE IMPLEMENTATION CALCULATES X<1
526 
527  @ingroup SpecFunc
528 
529  */
530  // [5.2.1.20] Riemann zeta function
531 
532  double riemann_zeta(double x);
533 
534 
535  /**
536 
537  Calculates the spherical Bessel functions of the first kind
538  (also called regular spherical Bessel functions).
539 
540  \f[ j_{n}(x) = \sqrt{\frac{\pi}{2x}} J_{n+1/2}(x) \f]
541 
542  For detailed description see
543  <A HREF="http://mathworld.wolfram.com/SphericalBesselFunctionoftheFirstKind.html">
544  Mathworld</A>. The implementation used is that of
545  <A HREF="http://www.gnu.org/software/gsl/manual/gsl-ref_7.html#SEC73">GSL</A>.
546 
547  @ingroup SpecFunc
548 
549  */
550  // [5.2.1.21] spherical Bessel functions of the first kind
551 
552  double sph_bessel(unsigned n, double x);
553 
554 
555  /**
556 
557  Computes the spherical (normalized) associated Legendre polynomials,
558  or spherical harmonic without azimuthal dependence (\f$e^(im\phi)\f$).
559 
560  \f[ Y_l^m(theta,0) = \sqrt{(2l+1)/(4\pi)} \sqrt{(l-m)!/(l+m)!} P_l^m(cos \theta) \f]
561 
562  for \f$m \geq 0, l \geq m\f$,
563  where the Condon-Shortley phase \f$(-1)^m\f$ is included in P_l^m(x)
564  This function is consistent with both C++0x and GSL,
565  even though there is a discrepancy in where to include the phase.
566  There is no reference in Abramowitz and Stegun.
567 
568 
569  @ingroup SpecFunc
570 
571  */
572 
573  // [5.2.1.22] spherical associated Legendre functions
574 
575  double sph_legendre(unsigned l, unsigned m, double theta);
576 
577 
578  /**
579 
580  Calculates the spherical Bessel functions of the second kind
581  (also called irregular spherical Bessel functions or
582  spherical Neumann functions).
583 
584  \f[ n_n(x) = y_n(x) = \sqrt{\frac{\pi}{2x}} N_{n+1/2}(x) \f]
585 
586  For detailed description see
587  <A HREF="http://mathworld.wolfram.com/SphericalBesselFunctionoftheSecondKind.html">
588  Mathworld</A>. The implementation used is that of
589  <A HREF="http://www.gnu.org/software/gsl/manual/gsl-ref_7.html#SEC74">GSL</A>.
590 
591  @ingroup SpecFunc
592 
593  */
594  // [5.2.1.23] spherical Neumann functions
595 
596  double sph_neumann(unsigned n, double x);
597 
598  /**
599 
600  Calculates the Airy function Ai
601 
602  \f[ Ai(x) = \frac{1}{\pi} \int\limits_{0}^{\infty} \cos(xt + t^3/3) dt \f]
603 
604  For detailed description see
605  <A HREF="http://mathworld.wolfram.com/AiryFunctions.html">
606  Mathworld</A>
607  and <A HREF="http://www.nrbook.com/abramowitz_and_stegun/page_446.htm">Abramowitz&Stegun, Sect. 10.4</A>.
608  The implementation used is that of
609  <A HREF="http://www.gnu.org/software/gsl/manual/html_node/Airy-Functions.html">GSL</A>.
610 
611  @ingroup SpecFunc
612 
613  */
614  // Airy function Ai
615 
616  double airy_Ai(double x);
617 
618  /**
619 
620  Calculates the Airy function Bi
621 
622  \f[ Bi(x) = \frac{1}{\pi} \int\limits_{0}^{\infty} [\exp(xt-t^3/3) + \cos(xt + t^3/3)] dt \f]
623 
624  For detailed description see
625  <A HREF="http://mathworld.wolfram.com/AiryFunctions.html">
626  Mathworld</A>
627  and <A HREF="http://www.nrbook.com/abramowitz_and_stegun/page_446.htm">Abramowitz&Stegun, Sect. 10.4</A>.
628  The implementation used is that of
629  <A HREF="http://www.gnu.org/software/gsl/manual/html_node/Airy-Functions.html">GSL</A>.
630 
631  @ingroup SpecFunc
632 
633  */
634  // Airy function Bi
635 
636  double airy_Bi(double x);
637 
638  /**
639 
640  Calculates the derivative of the Airy function Ai
641 
642  \f[ Ai(x) = \frac{1}{\pi} \int\limits_{0}^{\infty} \cos(xt + t^3/3) dt \f]
643 
644  For detailed description see
645  <A HREF="http://mathworld.wolfram.com/AiryFunctions.html">
646  Mathworld</A>
647  and <A HREF="http://www.nrbook.com/abramowitz_and_stegun/page_446.htm">Abramowitz&Stegun, Sect. 10.4</A>.
648  The implementation used is that of
649  <A HREF="http://www.gnu.org/software/gsl/manual/html_node/Derivatives-of-Airy-Functions.html">GSL</A>.
650 
651  @ingroup SpecFunc
652 
653  */
654  // Derivative of the Airy function Ai
655 
656  double airy_Ai_deriv(double x);
657 
658  /**
659 
660  Calculates the derivative of the Airy function Bi
661 
662  \f[ Bi(x) = \frac{1}{\pi} \int\limits_{0}^{\infty} [\exp(xt-t^3/3) + \cos(xt + t^3/3)] dt \f]
663 
664  For detailed description see
665  <A HREF="http://mathworld.wolfram.com/AiryFunctions.html">
666  Mathworld</A>
667  and <A HREF="http://www.nrbook.com/abramowitz_and_stegun/page_446.htm">Abramowitz&Stegun, Sect. 10.4</A>.
668  The implementation used is that of
669  <A HREF="http://www.gnu.org/software/gsl/manual/html_node/Derivatives-of-Airy-Functions.html">GSL</A>.
670 
671  @ingroup SpecFunc
672 
673  */
674  // Derivative of the Airy function Bi
675 
676  double airy_Bi_deriv(double x);
677 
678  /**
679 
680  Calculates the zeroes of the Airy function Ai
681 
682  \f[ Ai(x) = \frac{1}{\pi} \int\limits_{0}^{\infty} \cos(xt + t^3/3) dt \f]
683 
684  For detailed description see
685  <A HREF="http://mathworld.wolfram.com/AiryFunctionZeros.html">
686  Mathworld</A>
687  and <A HREF="http://www.nrbook.com/abramowitz_and_stegun/page_446.htm">Abramowitz&Stegun, Sect. 10.4</A>.
688  The implementation used is that of
689  <A HREF="http://www.gnu.org/software/gsl/manual/html_node/Zeros-of-Airy-Functions.html">GSL</A>.
690 
691  @ingroup SpecFunc
692 
693  */
694  // s-th zero of the Airy function Ai
695 
696  double airy_zero_Ai(unsigned int s);
697 
698  /**
699 
700  Calculates the zeroes of the Airy function Bi
701 
702  \f[ Bi(x) = \frac{1}{\pi} \int\limits_{0}^{\infty} [\exp(xt-t^3/3) + \cos(xt + t^3/3)] dt \f]
703 
704  For detailed description see
705  <A HREF="http://mathworld.wolfram.com/AiryFunctionZeros.html">
706  Mathworld</A>
707  and <A HREF="http://www.nrbook.com/abramowitz_and_stegun/page_446.htm">Abramowitz&Stegun, Sect. 10.4</A>.
708  The implementation used is that of
709  <A HREF="http://www.gnu.org/software/gsl/manual/html_node/Zeros-of-Airy-Functions.html">GSL</A>.
710 
711  @ingroup SpecFunc
712 
713  */
714  // s-th zero of the Airy function Bi
715 
716  double airy_zero_Bi(unsigned int s);
717 
718  /**
719 
720  Calculates the zeroes of the derivative of the Airy function Ai
721 
722  \f[ Ai(x) = \frac{1}{\pi} \int\limits_{0}^{\infty} \cos(xt + t^3/3) dt \f]
723 
724  For detailed description see
725  <A HREF="http://mathworld.wolfram.com/AiryFunctionZeros.html">
726  Mathworld</A>
727  and <A HREF="http://www.nrbook.com/abramowitz_and_stegun/page_446.htm">Abramowitz&Stegun, Sect. 10.4</A>.
728  The implementation used is that of
729  <A HREF="http://www.gnu.org/software/gsl/manual/html_node/Zeros-of-Derivatives-of-Airy-Functions.html">GSL</A>.
730 
731  @ingroup SpecFunc
732 
733  */
734  // s-th zero of the derivative of the Airy function Ai
735 
736  double airy_zero_Ai_deriv(unsigned int s);
737 
738  /**
739 
740  Calculates the zeroes of the derivative of the Airy function Bi
741 
742  \f[ Bi(x) = \frac{1}{\pi} \int\limits_{0}^{\infty} [\exp(xt-t^3/3) + \cos(xt + t^3/3)] dt \f]
743 
744  For detailed description see
745  <A HREF="http://mathworld.wolfram.com/AiryFunctionZeros.html">
746  Mathworld</A>
747  and <A HREF="http://www.nrbook.com/abramowitz_and_stegun/page_446.htm">Abramowitz&Stegun, Sect. 10.4</A>.
748  The implementation used is that of
749  <A HREF="http://www.gnu.org/software/gsl/manual/html_node/Zeros-of-Derivatives-of-Airy-Functions.html">GSL</A>.
750 
751  @ingroup SpecFunc
752 
753  */
754  // s-th zero of the derivative of the Airy function Bi
755 
756  double airy_zero_Bi_deriv(unsigned int s);
757 
758  /**
759 
760  Calculates the Wigner 3j coupling coefficients
761 
762  (ja jb jc
763  ma mb mc)
764 
765  where ja,ma,...etc are integers or half integers.
766  The function takes as input arguments only integers which corresponds
767  to half integer units, e.g two_ja = 2 * ja
768 
769  For detailed description see
770  <A HREF="http://mathworld.wolfram.com/Wigner3j-Symbol.html.html">
771  Mathworld</A>.
772  The implementation used is that of
773  <A HREF="http://www.gnu.org/software/gsl/manual/html_node/3_002dj-Symbols.html#g_t3_002dj-Symbols">GSL</A>.
774 
775  @ingroup SpecFunc
776 
777  */
778 
779  double wigner_3j(int two_ja, int two_jb, int two_jc, int two_ma, int two_mb, int two_mc);
780 
781  /**
782 
783  Calculates the Wigner 6j coupling coefficients
784 
785  (ja jb jc
786  jd je jf)
787 
788  where ja,jb,...etc are integers or half integers.
789  The function takes as input arguments only integers which corresponds
790  to half integer units, e.g two_ja = 2 * ja
791 
792  For detailed description see
793  <A HREF="http://mathworld.wolfram.com/Wigner6j-Symbol.html">
794  Mathworld</A>.
795  The implementation used is that of
796  <A HREF="http://www.gnu.org/software/gsl/manual/html_node/6_002dj-Symbols.html#g_t6_002dj-Symbols">GSL</A>.
797 
798  @ingroup SpecFunc
799 
800  */
801 
802  double wigner_6j(int two_ja, int two_jb, int two_jc, int two_jd, int two_je, int two_jf);
803 
804  /**
805 
806  Calculates the Wigner 9j coupling coefficients
807 
808  (ja jb jc
809  jd je jf
810  jg jh ji)
811 
812  where ja,jb...etc are integers or half integers.
813  The function takes as input arguments only integers which corresponds
814  to half integer units, e.g two_ja = 2 * ja
815 
816 
817  For detailed description see
818  <A HREF="http://mathworld.wolfram.com/Wigner9j-Symbol.html">
819  Mathworld</A>.
820  The implementation used is that of
821  <A HREF="http://www.gnu.org/software/gsl/manual/html_node/9_002dj-Symbols.html#g_t9_002dj-Symbols">GSL</A>.
822 
823  @ingroup SpecFunc
824 
825  */
826 
827  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);
828 
829 
830 
831 } // namespace Math
832 } // namespace ROOT
833 
834 
835 #endif //ROOT_Math_SpecFuncMathMore
double airy_zero_Bi(unsigned int s)
Calculates the zeroes of the Airy function Bi.
auto * m
Definition: textangle.C:8
Namespace for new ROOT classes and functions.
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...
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 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.
auto * a
Definition: textangle.C:12
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.
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.
static constexpr double s
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.
auto * l
Definition: textangle.C:4
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 ...
#define c(i)
Definition: RSha256.hxx:101
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
double cyl_bessel_k(double nu, double x)
Calculates the modified Bessel functions of the second kind (also called irregular modified (cylindri...