1 | // Mathematical Special Functions for -*- C++ -*- |
2 | |
3 | // Copyright (C) 2006-2024 Free Software Foundation, Inc. |
4 | // |
5 | // This file is part of the GNU ISO C++ Library. This library is free |
6 | // software; you can redistribute it and/or modify it under the |
7 | // terms of the GNU General Public License as published by the |
8 | // Free Software Foundation; either version 3, or (at your option) |
9 | // any later version. |
10 | |
11 | // This library is distributed in the hope that it will be useful, |
12 | // but WITHOUT ANY WARRANTY; without even the implied warranty of |
13 | // MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the |
14 | // GNU General Public License for more details. |
15 | |
16 | // Under Section 7 of GPL version 3, you are granted additional |
17 | // permissions described in the GCC Runtime Library Exception, version |
18 | // 3.1, as published by the Free Software Foundation. |
19 | |
20 | // You should have received a copy of the GNU General Public License and |
21 | // a copy of the GCC Runtime Library Exception along with this program; |
22 | // see the files COPYING3 and COPYING.RUNTIME respectively. If not, see |
23 | // <http://www.gnu.org/licenses/>. |
24 | |
25 | /** @file bits/specfun.h |
26 | * This is an internal header file, included by other library headers. |
27 | * Do not attempt to use it directly. @headername{cmath} |
28 | */ |
29 | |
30 | #ifndef _GLIBCXX_BITS_SPECFUN_H |
31 | #define _GLIBCXX_BITS_SPECFUN_H 1 |
32 | |
33 | #include <bits/c++config.h> |
34 | |
35 | #define __glibcxx_want_math_spec_funcs |
36 | #define __glibcxx_want_math_special_functions |
37 | #include <bits/version.h> |
38 | |
39 | #if __cplusplus <= 201403L && __STDCPP_WANT_MATH_SPEC_FUNCS__ == 0 |
40 | # error include <cmath> and define __STDCPP_WANT_MATH_SPEC_FUNCS__ |
41 | #endif |
42 | |
43 | #include <bits/stl_algobase.h> |
44 | #include <limits> |
45 | #include <type_traits> |
46 | |
47 | #include <tr1/gamma.tcc> |
48 | #include <tr1/bessel_function.tcc> |
49 | #include <tr1/beta_function.tcc> |
50 | #include <tr1/ell_integral.tcc> |
51 | #include <tr1/exp_integral.tcc> |
52 | #include <tr1/hypergeometric.tcc> |
53 | #include <tr1/legendre_function.tcc> |
54 | #include <tr1/modified_bessel_func.tcc> |
55 | #include <tr1/poly_hermite.tcc> |
56 | #include <tr1/poly_laguerre.tcc> |
57 | #include <tr1/riemann_zeta.tcc> |
58 | |
59 | namespace std _GLIBCXX_VISIBILITY(default) |
60 | { |
61 | _GLIBCXX_BEGIN_NAMESPACE_VERSION |
62 | |
63 | /** |
64 | * @defgroup mathsf Mathematical Special Functions |
65 | * @ingroup numerics |
66 | * |
67 | * @section mathsf_desc Mathematical Special Functions |
68 | * |
69 | * A collection of advanced mathematical special functions, |
70 | * defined by ISO/IEC IS 29124 and then added to ISO C++ 2017. |
71 | * |
72 | * |
73 | * @subsection mathsf_intro Introduction and History |
74 | * The first significant library upgrade on the road to C++2011, |
75 | * <a href="http://www.open-std.org/JTC1/SC22/WG21/docs/papers/2005/n1836.pdf"> |
76 | * TR1</a>, included a set of 23 mathematical functions that significantly |
77 | * extended the standard transcendental functions inherited from C and declared |
78 | * in @<cmath@>. |
79 | * |
80 | * Although most components from TR1 were eventually adopted for C++11 these |
81 | * math functions were left behind out of concern for implementability. |
82 | * The math functions were published as a separate international standard |
83 | * <a href="http://www.open-std.org/JTC1/SC22/WG21/docs/papers/2010/n3060.pdf"> |
84 | * IS 29124 - Extensions to the C++ Library to Support Mathematical Special |
85 | * Functions</a>. |
86 | * |
87 | * For C++17 these functions were incorporated into the main standard. |
88 | * |
89 | * @subsection mathsf_contents Contents |
90 | * The following functions are implemented in namespace @c std: |
91 | * - @ref assoc_laguerre "assoc_laguerre - Associated Laguerre functions" |
92 | * - @ref assoc_legendre "assoc_legendre - Associated Legendre functions" |
93 | * - @ref beta "beta - Beta functions" |
94 | * - @ref comp_ellint_1 "comp_ellint_1 - Complete elliptic functions of the first kind" |
95 | * - @ref comp_ellint_2 "comp_ellint_2 - Complete elliptic functions of the second kind" |
96 | * - @ref comp_ellint_3 "comp_ellint_3 - Complete elliptic functions of the third kind" |
97 | * - @ref cyl_bessel_i "cyl_bessel_i - Regular modified cylindrical Bessel functions" |
98 | * - @ref cyl_bessel_j "cyl_bessel_j - Cylindrical Bessel functions of the first kind" |
99 | * - @ref cyl_bessel_k "cyl_bessel_k - Irregular modified cylindrical Bessel functions" |
100 | * - @ref cyl_neumann "cyl_neumann - Cylindrical Neumann functions or Cylindrical Bessel functions of the second kind" |
101 | * - @ref ellint_1 "ellint_1 - Incomplete elliptic functions of the first kind" |
102 | * - @ref ellint_2 "ellint_2 - Incomplete elliptic functions of the second kind" |
103 | * - @ref ellint_3 "ellint_3 - Incomplete elliptic functions of the third kind" |
104 | * - @ref expint "expint - The exponential integral" |
105 | * - @ref hermite "hermite - Hermite polynomials" |
106 | * - @ref laguerre "laguerre - Laguerre functions" |
107 | * - @ref legendre "legendre - Legendre polynomials" |
108 | * - @ref riemann_zeta "riemann_zeta - The Riemann zeta function" |
109 | * - @ref sph_bessel "sph_bessel - Spherical Bessel functions" |
110 | * - @ref sph_legendre "sph_legendre - Spherical Legendre functions" |
111 | * - @ref sph_neumann "sph_neumann - Spherical Neumann functions" |
112 | * |
113 | * The hypergeometric functions were stricken from the TR29124 and C++17 |
114 | * versions of this math library because of implementation concerns. |
115 | * However, since they were in the TR1 version and since they are popular |
116 | * we kept them as an extension in namespace @c __gnu_cxx: |
117 | * - @ref __gnu_cxx::conf_hyperg "conf_hyperg - Confluent hypergeometric functions" |
118 | * - @ref __gnu_cxx::hyperg "hyperg - Hypergeometric functions" |
119 | * |
120 | * <!-- @subsection mathsf_general General Features --> |
121 | * |
122 | * @subsection mathsf_promotion Argument Promotion |
123 | * The arguments suppled to the non-suffixed functions will be promoted |
124 | * according to the following rules: |
125 | * 1. If any argument intended to be floating point is given an integral value |
126 | * That integral value is promoted to double. |
127 | * 2. All floating point arguments are promoted up to the largest floating |
128 | * point precision among them. |
129 | * |
130 | * @subsection mathsf_NaN NaN Arguments |
131 | * If any of the floating point arguments supplied to these functions is |
132 | * invalid or NaN (std::numeric_limits<Tp>::quiet_NaN), |
133 | * the value NaN is returned. |
134 | * |
135 | * @subsection mathsf_impl Implementation |
136 | * |
137 | * We strive to implement the underlying math with type generic algorithms |
138 | * to the greatest extent possible. In practice, the functions are thin |
139 | * wrappers that dispatch to function templates. Type dependence is |
140 | * controlled with std::numeric_limits and functions thereof. |
141 | * |
142 | * We don't promote @c float to @c double or @c double to <tt>long double</tt> |
143 | * reflexively. The goal is for @c float functions to operate more quickly, |
144 | * at the cost of @c float accuracy and possibly a smaller domain of validity. |
145 | * Similaryly, <tt>long double</tt> should give you more dynamic range |
146 | * and slightly more pecision than @c double on many systems. |
147 | * |
148 | * @subsection mathsf_testing Testing |
149 | * |
150 | * These functions have been tested against equivalent implementations |
151 | * from the <a href="http://www.gnu.org/software/gsl"> |
152 | * Gnu Scientific Library, GSL</a> and |
153 | * <a href="http://www.boost.org/doc/libs/1_60_0/libs/math/doc/html/index.html">Boost</a> |
154 | * and the ratio |
155 | * @f[ |
156 | * \frac{|f - f_{test}|}{|f_{test}|} |
157 | * @f] |
158 | * is generally found to be within 10<sup>-15</sup> for 64-bit double on |
159 | * linux-x86_64 systems over most of the ranges of validity. |
160 | * |
161 | * @todo Provide accuracy comparisons on a per-function basis for a small |
162 | * number of targets. |
163 | * |
164 | * @subsection mathsf_bibliography General Bibliography |
165 | * |
166 | * @see Abramowitz and Stegun: Handbook of Mathematical Functions, |
167 | * with Formulas, Graphs, and Mathematical Tables |
168 | * Edited by Milton Abramowitz and Irene A. Stegun, |
169 | * National Bureau of Standards Applied Mathematics Series - 55 |
170 | * Issued June 1964, Tenth Printing, December 1972, with corrections |
171 | * Electronic versions of A&S abound including both pdf and navigable html. |
172 | * @see for example http://people.math.sfu.ca/~cbm/aands/ |
173 | * |
174 | * @see The old A&S has been redone as the |
175 | * NIST Digital Library of Mathematical Functions: http://dlmf.nist.gov/ |
176 | * This version is far more navigable and includes more recent work. |
177 | * |
178 | * @see An Atlas of Functions: with Equator, the Atlas Function Calculator |
179 | * 2nd Edition, by Oldham, Keith B., Myland, Jan, Spanier, Jerome |
180 | * |
181 | * @see Asymptotics and Special Functions by Frank W. J. Olver, |
182 | * Academic Press, 1974 |
183 | * |
184 | * @see Numerical Recipes in C, The Art of Scientific Computing, |
185 | * by William H. Press, Second Ed., Saul A. Teukolsky, |
186 | * William T. Vetterling, and Brian P. Flannery, |
187 | * Cambridge University Press, 1992 |
188 | * |
189 | * @see The Special Functions and Their Approximations: Volumes 1 and 2, |
190 | * by Yudell L. Luke, Academic Press, 1969 |
191 | * |
192 | * @{ |
193 | */ |
194 | |
195 | // Associated Laguerre polynomials |
196 | |
197 | /** |
198 | * Return the associated Laguerre polynomial of order @c n, |
199 | * degree @c m: @f$ L_n^m(x) @f$ for @c float argument. |
200 | * |
201 | * @see assoc_laguerre for more details. |
202 | */ |
203 | inline float |
204 | assoc_laguerref(unsigned int __n, unsigned int __m, float __x) |
205 | { return __detail::__assoc_laguerre<float>(__n, __m, __x); } |
206 | |
207 | /** |
208 | * Return the associated Laguerre polynomial of order @c n, |
209 | * degree @c m: @f$ L_n^m(x) @f$. |
210 | * |
211 | * @see assoc_laguerre for more details. |
212 | */ |
213 | inline long double |
214 | assoc_laguerrel(unsigned int __n, unsigned int __m, long double __x) |
215 | { return __detail::__assoc_laguerre<long double>(__n, __m, __x); } |
216 | |
217 | /** |
218 | * Return the associated Laguerre polynomial of nonnegative order @c n, |
219 | * nonnegative degree @c m and real argument @c x: @f$ L_n^m(x) @f$. |
220 | * |
221 | * The associated Laguerre function of real degree @f$ \alpha @f$, |
222 | * @f$ L_n^\alpha(x) @f$, is defined by |
223 | * @f[ |
224 | * L_n^\alpha(x) = \frac{(\alpha + 1)_n}{n!} |
225 | * {}_1F_1(-n; \alpha + 1; x) |
226 | * @f] |
227 | * where @f$ (\alpha)_n @f$ is the Pochhammer symbol and |
228 | * @f$ {}_1F_1(a; c; x) @f$ is the confluent hypergeometric function. |
229 | * |
230 | * The associated Laguerre polynomial is defined for integral |
231 | * degree @f$ \alpha = m @f$ by: |
232 | * @f[ |
233 | * L_n^m(x) = (-1)^m \frac{d^m}{dx^m} L_{n + m}(x) |
234 | * @f] |
235 | * where the Laguerre polynomial is defined by: |
236 | * @f[ |
237 | * L_n(x) = \frac{e^x}{n!} \frac{d^n}{dx^n} (x^ne^{-x}) |
238 | * @f] |
239 | * and @f$ x >= 0 @f$. |
240 | * @see laguerre for details of the Laguerre function of degree @c n |
241 | * |
242 | * @tparam _Tp The floating-point type of the argument @c __x. |
243 | * @param __n The order of the Laguerre function, <tt>__n >= 0</tt>. |
244 | * @param __m The degree of the Laguerre function, <tt>__m >= 0</tt>. |
245 | * @param __x The argument of the Laguerre function, <tt>__x >= 0</tt>. |
246 | * @throw std::domain_error if <tt>__x < 0</tt>. |
247 | */ |
248 | template<typename _Tp> |
249 | inline typename __gnu_cxx::__promote<_Tp>::__type |
250 | assoc_laguerre(unsigned int __n, unsigned int __m, _Tp __x) |
251 | { |
252 | typedef typename __gnu_cxx::__promote<_Tp>::__type __type; |
253 | return __detail::__assoc_laguerre<__type>(__n, __m, __x); |
254 | } |
255 | |
256 | // Associated Legendre functions |
257 | |
258 | /** |
259 | * Return the associated Legendre function of degree @c l and order @c m |
260 | * for @c float argument. |
261 | * |
262 | * @see assoc_legendre for more details. |
263 | */ |
264 | inline float |
265 | assoc_legendref(unsigned int __l, unsigned int __m, float __x) |
266 | { return __detail::__assoc_legendre_p<float>(__l, __m, __x); } |
267 | |
268 | /** |
269 | * Return the associated Legendre function of degree @c l and order @c m. |
270 | * |
271 | * @see assoc_legendre for more details. |
272 | */ |
273 | inline long double |
274 | assoc_legendrel(unsigned int __l, unsigned int __m, long double __x) |
275 | { return __detail::__assoc_legendre_p<long double>(__l, __m, __x); } |
276 | |
277 | |
278 | /** |
279 | * Return the associated Legendre function of degree @c l and order @c m. |
280 | * |
281 | * The associated Legendre function is derived from the Legendre function |
282 | * @f$ P_l(x) @f$ by the Rodrigues formula: |
283 | * @f[ |
284 | * P_l^m(x) = (1 - x^2)^{m/2}\frac{d^m}{dx^m}P_l(x) |
285 | * @f] |
286 | * @see legendre for details of the Legendre function of degree @c l |
287 | * |
288 | * @tparam _Tp The floating-point type of the argument @c __x. |
289 | * @param __l The degree <tt>__l >= 0</tt>. |
290 | * @param __m The order <tt>__m <= l</tt>. |
291 | * @param __x The argument, <tt>abs(__x) <= 1</tt>. |
292 | * @throw std::domain_error if <tt>abs(__x) > 1</tt>. |
293 | */ |
294 | template<typename _Tp> |
295 | inline typename __gnu_cxx::__promote<_Tp>::__type |
296 | assoc_legendre(unsigned int __l, unsigned int __m, _Tp __x) |
297 | { |
298 | typedef typename __gnu_cxx::__promote<_Tp>::__type __type; |
299 | return __detail::__assoc_legendre_p<__type>(__l, __m, __x); |
300 | } |
301 | |
302 | // Beta functions |
303 | |
304 | /** |
305 | * Return the beta function, @f$ B(a,b) @f$, for @c float parameters @c a, @c b. |
306 | * |
307 | * @see beta for more details. |
308 | */ |
309 | inline float |
310 | betaf(float __a, float __b) |
311 | { return __detail::__beta<float>(x: __a, y: __b); } |
312 | |
313 | /** |
314 | * Return the beta function, @f$B(a,b)@f$, for long double |
315 | * parameters @c a, @c b. |
316 | * |
317 | * @see beta for more details. |
318 | */ |
319 | inline long double |
320 | betal(long double __a, long double __b) |
321 | { return __detail::__beta<long double>(x: __a, y: __b); } |
322 | |
323 | /** |
324 | * Return the beta function, @f$B(a,b)@f$, for real parameters @c a, @c b. |
325 | * |
326 | * The beta function is defined by |
327 | * @f[ |
328 | * B(a,b) = \int_0^1 t^{a - 1} (1 - t)^{b - 1} dt |
329 | * = \frac{\Gamma(a)\Gamma(b)}{\Gamma(a+b)} |
330 | * @f] |
331 | * where @f$ a > 0 @f$ and @f$ b > 0 @f$ |
332 | * |
333 | * @tparam _Tpa The floating-point type of the parameter @c __a. |
334 | * @tparam _Tpb The floating-point type of the parameter @c __b. |
335 | * @param __a The first argument of the beta function, <tt> __a > 0 </tt>. |
336 | * @param __b The second argument of the beta function, <tt> __b > 0 </tt>. |
337 | * @throw std::domain_error if <tt> __a < 0 </tt> or <tt> __b < 0 </tt>. |
338 | */ |
339 | template<typename _Tpa, typename _Tpb> |
340 | inline typename __gnu_cxx::__promote_2<_Tpa, _Tpb>::__type |
341 | beta(_Tpa __a, _Tpb __b) |
342 | { |
343 | typedef typename __gnu_cxx::__promote_2<_Tpa, _Tpb>::__type __type; |
344 | return __detail::__beta<__type>(__a, __b); |
345 | } |
346 | |
347 | // Complete elliptic integrals of the first kind |
348 | |
349 | /** |
350 | * Return the complete elliptic integral of the first kind @f$ E(k) @f$ |
351 | * for @c float modulus @c k. |
352 | * |
353 | * @see comp_ellint_1 for details. |
354 | */ |
355 | inline float |
356 | comp_ellint_1f(float __k) |
357 | { return __detail::__comp_ellint_1<float>(__k); } |
358 | |
359 | /** |
360 | * Return the complete elliptic integral of the first kind @f$ E(k) @f$ |
361 | * for long double modulus @c k. |
362 | * |
363 | * @see comp_ellint_1 for details. |
364 | */ |
365 | inline long double |
366 | comp_ellint_1l(long double __k) |
367 | { return __detail::__comp_ellint_1<long double>(__k); } |
368 | |
369 | /** |
370 | * Return the complete elliptic integral of the first kind |
371 | * @f$ K(k) @f$ for real modulus @c k. |
372 | * |
373 | * The complete elliptic integral of the first kind is defined as |
374 | * @f[ |
375 | * K(k) = F(k,\pi/2) = \int_0^{\pi/2}\frac{d\theta} |
376 | * {\sqrt{1 - k^2 sin^2\theta}} |
377 | * @f] |
378 | * where @f$ F(k,\phi) @f$ is the incomplete elliptic integral of the |
379 | * first kind and the modulus @f$ |k| <= 1 @f$. |
380 | * @see ellint_1 for details of the incomplete elliptic function |
381 | * of the first kind. |
382 | * |
383 | * @tparam _Tp The floating-point type of the modulus @c __k. |
384 | * @param __k The modulus, <tt> abs(__k) <= 1 </tt> |
385 | * @throw std::domain_error if <tt> abs(__k) > 1 </tt>. |
386 | */ |
387 | template<typename _Tp> |
388 | inline typename __gnu_cxx::__promote<_Tp>::__type |
389 | comp_ellint_1(_Tp __k) |
390 | { |
391 | typedef typename __gnu_cxx::__promote<_Tp>::__type __type; |
392 | return __detail::__comp_ellint_1<__type>(__k); |
393 | } |
394 | |
395 | // Complete elliptic integrals of the second kind |
396 | |
397 | /** |
398 | * Return the complete elliptic integral of the second kind @f$ E(k) @f$ |
399 | * for @c float modulus @c k. |
400 | * |
401 | * @see comp_ellint_2 for details. |
402 | */ |
403 | inline float |
404 | comp_ellint_2f(float __k) |
405 | { return __detail::__comp_ellint_2<float>(__k); } |
406 | |
407 | /** |
408 | * Return the complete elliptic integral of the second kind @f$ E(k) @f$ |
409 | * for long double modulus @c k. |
410 | * |
411 | * @see comp_ellint_2 for details. |
412 | */ |
413 | inline long double |
414 | comp_ellint_2l(long double __k) |
415 | { return __detail::__comp_ellint_2<long double>(__k); } |
416 | |
417 | /** |
418 | * Return the complete elliptic integral of the second kind @f$ E(k) @f$ |
419 | * for real modulus @c k. |
420 | * |
421 | * The complete elliptic integral of the second kind is defined as |
422 | * @f[ |
423 | * E(k) = E(k,\pi/2) = \int_0^{\pi/2}\sqrt{1 - k^2 sin^2\theta} |
424 | * @f] |
425 | * where @f$ E(k,\phi) @f$ is the incomplete elliptic integral of the |
426 | * second kind and the modulus @f$ |k| <= 1 @f$. |
427 | * @see ellint_2 for details of the incomplete elliptic function |
428 | * of the second kind. |
429 | * |
430 | * @tparam _Tp The floating-point type of the modulus @c __k. |
431 | * @param __k The modulus, @c abs(__k) <= 1 |
432 | * @throw std::domain_error if @c abs(__k) > 1. |
433 | */ |
434 | template<typename _Tp> |
435 | inline typename __gnu_cxx::__promote<_Tp>::__type |
436 | comp_ellint_2(_Tp __k) |
437 | { |
438 | typedef typename __gnu_cxx::__promote<_Tp>::__type __type; |
439 | return __detail::__comp_ellint_2<__type>(__k); |
440 | } |
441 | |
442 | // Complete elliptic integrals of the third kind |
443 | |
444 | /** |
445 | * @brief Return the complete elliptic integral of the third kind |
446 | * @f$ \Pi(k,\nu) @f$ for @c float modulus @c k. |
447 | * |
448 | * @see comp_ellint_3 for details. |
449 | */ |
450 | inline float |
451 | comp_ellint_3f(float __k, float __nu) |
452 | { return __detail::__comp_ellint_3<float>(__k, __nu); } |
453 | |
454 | /** |
455 | * @brief Return the complete elliptic integral of the third kind |
456 | * @f$ \Pi(k,\nu) @f$ for <tt>long double</tt> modulus @c k. |
457 | * |
458 | * @see comp_ellint_3 for details. |
459 | */ |
460 | inline long double |
461 | comp_ellint_3l(long double __k, long double __nu) |
462 | { return __detail::__comp_ellint_3<long double>(__k, __nu); } |
463 | |
464 | /** |
465 | * Return the complete elliptic integral of the third kind |
466 | * @f$ \Pi(k,\nu) = \Pi(k,\nu,\pi/2) @f$ for real modulus @c k. |
467 | * |
468 | * The complete elliptic integral of the third kind is defined as |
469 | * @f[ |
470 | * \Pi(k,\nu) = \Pi(k,\nu,\pi/2) = \int_0^{\pi/2} |
471 | * \frac{d\theta} |
472 | * {(1 - \nu \sin^2\theta)\sqrt{1 - k^2 \sin^2\theta}} |
473 | * @f] |
474 | * where @f$ \Pi(k,\nu,\phi) @f$ is the incomplete elliptic integral of the |
475 | * second kind and the modulus @f$ |k| <= 1 @f$. |
476 | * @see ellint_3 for details of the incomplete elliptic function |
477 | * of the third kind. |
478 | * |
479 | * @tparam _Tp The floating-point type of the modulus @c __k. |
480 | * @tparam _Tpn The floating-point type of the argument @c __nu. |
481 | * @param __k The modulus, @c abs(__k) <= 1 |
482 | * @param __nu The argument |
483 | * @throw std::domain_error if @c abs(__k) > 1. |
484 | */ |
485 | template<typename _Tp, typename _Tpn> |
486 | inline typename __gnu_cxx::__promote_2<_Tp, _Tpn>::__type |
487 | comp_ellint_3(_Tp __k, _Tpn __nu) |
488 | { |
489 | typedef typename __gnu_cxx::__promote_2<_Tp, _Tpn>::__type __type; |
490 | return __detail::__comp_ellint_3<__type>(__k, __nu); |
491 | } |
492 | |
493 | // Regular modified cylindrical Bessel functions |
494 | |
495 | /** |
496 | * Return the regular modified Bessel function @f$ I_{\nu}(x) @f$ |
497 | * for @c float order @f$ \nu @f$ and argument @f$ x >= 0 @f$. |
498 | * |
499 | * @see cyl_bessel_i for setails. |
500 | */ |
501 | inline float |
502 | cyl_bessel_if(float __nu, float __x) |
503 | { return __detail::__cyl_bessel_i<float>(__nu, __x); } |
504 | |
505 | /** |
506 | * Return the regular modified Bessel function @f$ I_{\nu}(x) @f$ |
507 | * for <tt>long double</tt> order @f$ \nu @f$ and argument @f$ x >= 0 @f$. |
508 | * |
509 | * @see cyl_bessel_i for setails. |
510 | */ |
511 | inline long double |
512 | cyl_bessel_il(long double __nu, long double __x) |
513 | { return __detail::__cyl_bessel_i<long double>(__nu, __x); } |
514 | |
515 | /** |
516 | * Return the regular modified Bessel function @f$ I_{\nu}(x) @f$ |
517 | * for real order @f$ \nu @f$ and argument @f$ x >= 0 @f$. |
518 | * |
519 | * The regular modified cylindrical Bessel function is: |
520 | * @f[ |
521 | * I_{\nu}(x) = i^{-\nu}J_\nu(ix) = \sum_{k=0}^{\infty} |
522 | * \frac{(x/2)^{\nu + 2k}}{k!\Gamma(\nu+k+1)} |
523 | * @f] |
524 | * |
525 | * @tparam _Tpnu The floating-point type of the order @c __nu. |
526 | * @tparam _Tp The floating-point type of the argument @c __x. |
527 | * @param __nu The order |
528 | * @param __x The argument, <tt> __x >= 0 </tt> |
529 | * @throw std::domain_error if <tt> __x < 0 </tt>. |
530 | */ |
531 | template<typename _Tpnu, typename _Tp> |
532 | inline typename __gnu_cxx::__promote_2<_Tpnu, _Tp>::__type |
533 | cyl_bessel_i(_Tpnu __nu, _Tp __x) |
534 | { |
535 | typedef typename __gnu_cxx::__promote_2<_Tpnu, _Tp>::__type __type; |
536 | return __detail::__cyl_bessel_i<__type>(__nu, __x); |
537 | } |
538 | |
539 | // Cylindrical Bessel functions (of the first kind) |
540 | |
541 | /** |
542 | * Return the Bessel function of the first kind @f$ J_{\nu}(x) @f$ |
543 | * for @c float order @f$ \nu @f$ and argument @f$ x >= 0 @f$. |
544 | * |
545 | * @see cyl_bessel_j for setails. |
546 | */ |
547 | inline float |
548 | cyl_bessel_jf(float __nu, float __x) |
549 | { return __detail::__cyl_bessel_j<float>(__nu, __x); } |
550 | |
551 | /** |
552 | * Return the Bessel function of the first kind @f$ J_{\nu}(x) @f$ |
553 | * for <tt>long double</tt> order @f$ \nu @f$ and argument @f$ x >= 0 @f$. |
554 | * |
555 | * @see cyl_bessel_j for setails. |
556 | */ |
557 | inline long double |
558 | cyl_bessel_jl(long double __nu, long double __x) |
559 | { return __detail::__cyl_bessel_j<long double>(__nu, __x); } |
560 | |
561 | /** |
562 | * Return the Bessel function @f$ J_{\nu}(x) @f$ of real order @f$ \nu @f$ |
563 | * and argument @f$ x >= 0 @f$. |
564 | * |
565 | * The cylindrical Bessel function is: |
566 | * @f[ |
567 | * J_{\nu}(x) = \sum_{k=0}^{\infty} |
568 | * \frac{(-1)^k (x/2)^{\nu + 2k}}{k!\Gamma(\nu+k+1)} |
569 | * @f] |
570 | * |
571 | * @tparam _Tpnu The floating-point type of the order @c __nu. |
572 | * @tparam _Tp The floating-point type of the argument @c __x. |
573 | * @param __nu The order |
574 | * @param __x The argument, <tt> __x >= 0 </tt> |
575 | * @throw std::domain_error if <tt> __x < 0 </tt>. |
576 | */ |
577 | template<typename _Tpnu, typename _Tp> |
578 | inline typename __gnu_cxx::__promote_2<_Tpnu, _Tp>::__type |
579 | cyl_bessel_j(_Tpnu __nu, _Tp __x) |
580 | { |
581 | typedef typename __gnu_cxx::__promote_2<_Tpnu, _Tp>::__type __type; |
582 | return __detail::__cyl_bessel_j<__type>(__nu, __x); |
583 | } |
584 | |
585 | // Irregular modified cylindrical Bessel functions |
586 | |
587 | /** |
588 | * Return the irregular modified Bessel function @f$ K_{\nu}(x) @f$ |
589 | * for @c float order @f$ \nu @f$ and argument @f$ x >= 0 @f$. |
590 | * |
591 | * @see cyl_bessel_k for setails. |
592 | */ |
593 | inline float |
594 | cyl_bessel_kf(float __nu, float __x) |
595 | { return __detail::__cyl_bessel_k<float>(__nu, __x); } |
596 | |
597 | /** |
598 | * Return the irregular modified Bessel function @f$ K_{\nu}(x) @f$ |
599 | * for <tt>long double</tt> order @f$ \nu @f$ and argument @f$ x >= 0 @f$. |
600 | * |
601 | * @see cyl_bessel_k for setails. |
602 | */ |
603 | inline long double |
604 | cyl_bessel_kl(long double __nu, long double __x) |
605 | { return __detail::__cyl_bessel_k<long double>(__nu, __x); } |
606 | |
607 | /** |
608 | * Return the irregular modified Bessel function @f$ K_{\nu}(x) @f$ |
609 | * of real order @f$ \nu @f$ and argument @f$ x @f$. |
610 | * |
611 | * The irregular modified Bessel function is defined by: |
612 | * @f[ |
613 | * K_{\nu}(x) = \frac{\pi}{2} |
614 | * \frac{I_{-\nu}(x) - I_{\nu}(x)}{\sin \nu\pi} |
615 | * @f] |
616 | * where for integral @f$ \nu = n @f$ a limit is taken: |
617 | * @f$ lim_{\nu \to n} @f$. |
618 | * For negative argument we have simply: |
619 | * @f[ |
620 | * K_{-\nu}(x) = K_{\nu}(x) |
621 | * @f] |
622 | * |
623 | * @tparam _Tpnu The floating-point type of the order @c __nu. |
624 | * @tparam _Tp The floating-point type of the argument @c __x. |
625 | * @param __nu The order |
626 | * @param __x The argument, <tt> __x >= 0 </tt> |
627 | * @throw std::domain_error if <tt> __x < 0 </tt>. |
628 | */ |
629 | template<typename _Tpnu, typename _Tp> |
630 | inline typename __gnu_cxx::__promote_2<_Tpnu, _Tp>::__type |
631 | cyl_bessel_k(_Tpnu __nu, _Tp __x) |
632 | { |
633 | typedef typename __gnu_cxx::__promote_2<_Tpnu, _Tp>::__type __type; |
634 | return __detail::__cyl_bessel_k<__type>(__nu, __x); |
635 | } |
636 | |
637 | // Cylindrical Neumann functions |
638 | |
639 | /** |
640 | * Return the Neumann function @f$ N_{\nu}(x) @f$ |
641 | * of @c float order @f$ \nu @f$ and argument @f$ x @f$. |
642 | * |
643 | * @see cyl_neumann for setails. |
644 | */ |
645 | inline float |
646 | cyl_neumannf(float __nu, float __x) |
647 | { return __detail::__cyl_neumann_n<float>(__nu, __x); } |
648 | |
649 | /** |
650 | * Return the Neumann function @f$ N_{\nu}(x) @f$ |
651 | * of <tt>long double</tt> order @f$ \nu @f$ and argument @f$ x @f$. |
652 | * |
653 | * @see cyl_neumann for setails. |
654 | */ |
655 | inline long double |
656 | cyl_neumannl(long double __nu, long double __x) |
657 | { return __detail::__cyl_neumann_n<long double>(__nu, __x); } |
658 | |
659 | /** |
660 | * Return the Neumann function @f$ N_{\nu}(x) @f$ |
661 | * of real order @f$ \nu @f$ and argument @f$ x >= 0 @f$. |
662 | * |
663 | * The Neumann function is defined by: |
664 | * @f[ |
665 | * N_{\nu}(x) = \frac{J_{\nu}(x) \cos \nu\pi - J_{-\nu}(x)} |
666 | * {\sin \nu\pi} |
667 | * @f] |
668 | * where @f$ x >= 0 @f$ and for integral order @f$ \nu = n @f$ |
669 | * a limit is taken: @f$ lim_{\nu \to n} @f$. |
670 | * |
671 | * @tparam _Tpnu The floating-point type of the order @c __nu. |
672 | * @tparam _Tp The floating-point type of the argument @c __x. |
673 | * @param __nu The order |
674 | * @param __x The argument, <tt> __x >= 0 </tt> |
675 | * @throw std::domain_error if <tt> __x < 0 </tt>. |
676 | */ |
677 | template<typename _Tpnu, typename _Tp> |
678 | inline typename __gnu_cxx::__promote_2<_Tpnu, _Tp>::__type |
679 | cyl_neumann(_Tpnu __nu, _Tp __x) |
680 | { |
681 | typedef typename __gnu_cxx::__promote_2<_Tpnu, _Tp>::__type __type; |
682 | return __detail::__cyl_neumann_n<__type>(__nu, __x); |
683 | } |
684 | |
685 | // Incomplete elliptic integrals of the first kind |
686 | |
687 | /** |
688 | * Return the incomplete elliptic integral of the first kind @f$ E(k,\phi) @f$ |
689 | * for @c float modulus @f$ k @f$ and angle @f$ \phi @f$. |
690 | * |
691 | * @see ellint_1 for details. |
692 | */ |
693 | inline float |
694 | ellint_1f(float __k, float __phi) |
695 | { return __detail::__ellint_1<float>(__k, __phi); } |
696 | |
697 | /** |
698 | * Return the incomplete elliptic integral of the first kind @f$ E(k,\phi) @f$ |
699 | * for <tt>long double</tt> modulus @f$ k @f$ and angle @f$ \phi @f$. |
700 | * |
701 | * @see ellint_1 for details. |
702 | */ |
703 | inline long double |
704 | ellint_1l(long double __k, long double __phi) |
705 | { return __detail::__ellint_1<long double>(__k, __phi); } |
706 | |
707 | /** |
708 | * Return the incomplete elliptic integral of the first kind @f$ F(k,\phi) @f$ |
709 | * for @c real modulus @f$ k @f$ and angle @f$ \phi @f$. |
710 | * |
711 | * The incomplete elliptic integral of the first kind is defined as |
712 | * @f[ |
713 | * F(k,\phi) = \int_0^{\phi}\frac{d\theta} |
714 | * {\sqrt{1 - k^2 sin^2\theta}} |
715 | * @f] |
716 | * For @f$ \phi= \pi/2 @f$ this becomes the complete elliptic integral of |
717 | * the first kind, @f$ K(k) @f$. @see comp_ellint_1. |
718 | * |
719 | * @tparam _Tp The floating-point type of the modulus @c __k. |
720 | * @tparam _Tpp The floating-point type of the angle @c __phi. |
721 | * @param __k The modulus, <tt> abs(__k) <= 1 </tt> |
722 | * @param __phi The integral limit argument in radians |
723 | * @throw std::domain_error if <tt> abs(__k) > 1 </tt>. |
724 | */ |
725 | template<typename _Tp, typename _Tpp> |
726 | inline typename __gnu_cxx::__promote_2<_Tp, _Tpp>::__type |
727 | ellint_1(_Tp __k, _Tpp __phi) |
728 | { |
729 | typedef typename __gnu_cxx::__promote_2<_Tp, _Tpp>::__type __type; |
730 | return __detail::__ellint_1<__type>(__k, __phi); |
731 | } |
732 | |
733 | // Incomplete elliptic integrals of the second kind |
734 | |
735 | /** |
736 | * @brief Return the incomplete elliptic integral of the second kind |
737 | * @f$ E(k,\phi) @f$ for @c float argument. |
738 | * |
739 | * @see ellint_2 for details. |
740 | */ |
741 | inline float |
742 | ellint_2f(float __k, float __phi) |
743 | { return __detail::__ellint_2<float>(__k, __phi); } |
744 | |
745 | /** |
746 | * @brief Return the incomplete elliptic integral of the second kind |
747 | * @f$ E(k,\phi) @f$. |
748 | * |
749 | * @see ellint_2 for details. |
750 | */ |
751 | inline long double |
752 | ellint_2l(long double __k, long double __phi) |
753 | { return __detail::__ellint_2<long double>(__k, __phi); } |
754 | |
755 | /** |
756 | * Return the incomplete elliptic integral of the second kind |
757 | * @f$ E(k,\phi) @f$. |
758 | * |
759 | * The incomplete elliptic integral of the second kind is defined as |
760 | * @f[ |
761 | * E(k,\phi) = \int_0^{\phi} \sqrt{1 - k^2 sin^2\theta} |
762 | * @f] |
763 | * For @f$ \phi= \pi/2 @f$ this becomes the complete elliptic integral of |
764 | * the second kind, @f$ E(k) @f$. @see comp_ellint_2. |
765 | * |
766 | * @tparam _Tp The floating-point type of the modulus @c __k. |
767 | * @tparam _Tpp The floating-point type of the angle @c __phi. |
768 | * @param __k The modulus, <tt> abs(__k) <= 1 </tt> |
769 | * @param __phi The integral limit argument in radians |
770 | * @return The elliptic function of the second kind. |
771 | * @throw std::domain_error if <tt> abs(__k) > 1 </tt>. |
772 | */ |
773 | template<typename _Tp, typename _Tpp> |
774 | inline typename __gnu_cxx::__promote_2<_Tp, _Tpp>::__type |
775 | ellint_2(_Tp __k, _Tpp __phi) |
776 | { |
777 | typedef typename __gnu_cxx::__promote_2<_Tp, _Tpp>::__type __type; |
778 | return __detail::__ellint_2<__type>(__k, __phi); |
779 | } |
780 | |
781 | // Incomplete elliptic integrals of the third kind |
782 | |
783 | /** |
784 | * @brief Return the incomplete elliptic integral of the third kind |
785 | * @f$ \Pi(k,\nu,\phi) @f$ for @c float argument. |
786 | * |
787 | * @see ellint_3 for details. |
788 | */ |
789 | inline float |
790 | ellint_3f(float __k, float __nu, float __phi) |
791 | { return __detail::__ellint_3<float>(__k, __nu, __phi); } |
792 | |
793 | /** |
794 | * @brief Return the incomplete elliptic integral of the third kind |
795 | * @f$ \Pi(k,\nu,\phi) @f$. |
796 | * |
797 | * @see ellint_3 for details. |
798 | */ |
799 | inline long double |
800 | ellint_3l(long double __k, long double __nu, long double __phi) |
801 | { return __detail::__ellint_3<long double>(__k, __nu, __phi); } |
802 | |
803 | /** |
804 | * @brief Return the incomplete elliptic integral of the third kind |
805 | * @f$ \Pi(k,\nu,\phi) @f$. |
806 | * |
807 | * The incomplete elliptic integral of the third kind is defined by: |
808 | * @f[ |
809 | * \Pi(k,\nu,\phi) = \int_0^{\phi} |
810 | * \frac{d\theta} |
811 | * {(1 - \nu \sin^2\theta) |
812 | * \sqrt{1 - k^2 \sin^2\theta}} |
813 | * @f] |
814 | * For @f$ \phi= \pi/2 @f$ this becomes the complete elliptic integral of |
815 | * the third kind, @f$ \Pi(k,\nu) @f$. @see comp_ellint_3. |
816 | * |
817 | * @tparam _Tp The floating-point type of the modulus @c __k. |
818 | * @tparam _Tpn The floating-point type of the argument @c __nu. |
819 | * @tparam _Tpp The floating-point type of the angle @c __phi. |
820 | * @param __k The modulus, <tt> abs(__k) <= 1 </tt> |
821 | * @param __nu The second argument |
822 | * @param __phi The integral limit argument in radians |
823 | * @return The elliptic function of the third kind. |
824 | * @throw std::domain_error if <tt> abs(__k) > 1 </tt>. |
825 | */ |
826 | template<typename _Tp, typename _Tpn, typename _Tpp> |
827 | inline typename __gnu_cxx::__promote_3<_Tp, _Tpn, _Tpp>::__type |
828 | ellint_3(_Tp __k, _Tpn __nu, _Tpp __phi) |
829 | { |
830 | typedef typename __gnu_cxx::__promote_3<_Tp, _Tpn, _Tpp>::__type __type; |
831 | return __detail::__ellint_3<__type>(__k, __nu, __phi); |
832 | } |
833 | |
834 | // Exponential integrals |
835 | |
836 | /** |
837 | * Return the exponential integral @f$ Ei(x) @f$ for @c float argument @c x. |
838 | * |
839 | * @see expint for details. |
840 | */ |
841 | inline float |
842 | expintf(float __x) |
843 | { return __detail::__expint<float>(__x); } |
844 | |
845 | /** |
846 | * Return the exponential integral @f$ Ei(x) @f$ |
847 | * for <tt>long double</tt> argument @c x. |
848 | * |
849 | * @see expint for details. |
850 | */ |
851 | inline long double |
852 | expintl(long double __x) |
853 | { return __detail::__expint<long double>(__x); } |
854 | |
855 | /** |
856 | * Return the exponential integral @f$ Ei(x) @f$ for @c real argument @c x. |
857 | * |
858 | * The exponential integral is given by |
859 | * \f[ |
860 | * Ei(x) = -\int_{-x}^\infty \frac{e^t}{t} dt |
861 | * \f] |
862 | * |
863 | * @tparam _Tp The floating-point type of the argument @c __x. |
864 | * @param __x The argument of the exponential integral function. |
865 | */ |
866 | template<typename _Tp> |
867 | inline typename __gnu_cxx::__promote<_Tp>::__type |
868 | expint(_Tp __x) |
869 | { |
870 | typedef typename __gnu_cxx::__promote<_Tp>::__type __type; |
871 | return __detail::__expint<__type>(__x); |
872 | } |
873 | |
874 | // Hermite polynomials |
875 | |
876 | /** |
877 | * Return the Hermite polynomial @f$ H_n(x) @f$ of nonnegative order n |
878 | * and float argument @c x. |
879 | * |
880 | * @see hermite for details. |
881 | */ |
882 | inline float |
883 | hermitef(unsigned int __n, float __x) |
884 | { return __detail::__poly_hermite<float>(__n, __x); } |
885 | |
886 | /** |
887 | * Return the Hermite polynomial @f$ H_n(x) @f$ of nonnegative order n |
888 | * and <tt>long double</tt> argument @c x. |
889 | * |
890 | * @see hermite for details. |
891 | */ |
892 | inline long double |
893 | hermitel(unsigned int __n, long double __x) |
894 | { return __detail::__poly_hermite<long double>(__n, __x); } |
895 | |
896 | /** |
897 | * Return the Hermite polynomial @f$ H_n(x) @f$ of order n |
898 | * and @c real argument @c x. |
899 | * |
900 | * The Hermite polynomial is defined by: |
901 | * @f[ |
902 | * H_n(x) = (-1)^n e^{x^2} \frac{d^n}{dx^n} e^{-x^2} |
903 | * @f] |
904 | * |
905 | * The Hermite polynomial obeys a reflection formula: |
906 | * @f[ |
907 | * H_n(-x) = (-1)^n H_n(x) |
908 | * @f] |
909 | * |
910 | * @tparam _Tp The floating-point type of the argument @c __x. |
911 | * @param __n The order |
912 | * @param __x The argument |
913 | */ |
914 | template<typename _Tp> |
915 | inline typename __gnu_cxx::__promote<_Tp>::__type |
916 | hermite(unsigned int __n, _Tp __x) |
917 | { |
918 | typedef typename __gnu_cxx::__promote<_Tp>::__type __type; |
919 | return __detail::__poly_hermite<__type>(__n, __x); |
920 | } |
921 | |
922 | // Laguerre polynomials |
923 | |
924 | /** |
925 | * Returns the Laguerre polynomial @f$ L_n(x) @f$ of nonnegative degree @c n |
926 | * and @c float argument @f$ x >= 0 @f$. |
927 | * |
928 | * @see laguerre for more details. |
929 | */ |
930 | inline float |
931 | laguerref(unsigned int __n, float __x) |
932 | { return __detail::__laguerre<float>(__n, __x); } |
933 | |
934 | /** |
935 | * Returns the Laguerre polynomial @f$ L_n(x) @f$ of nonnegative degree @c n |
936 | * and <tt>long double</tt> argument @f$ x >= 0 @f$. |
937 | * |
938 | * @see laguerre for more details. |
939 | */ |
940 | inline long double |
941 | laguerrel(unsigned int __n, long double __x) |
942 | { return __detail::__laguerre<long double>(__n, __x); } |
943 | |
944 | /** |
945 | * Returns the Laguerre polynomial @f$ L_n(x) @f$ |
946 | * of nonnegative degree @c n and real argument @f$ x >= 0 @f$. |
947 | * |
948 | * The Laguerre polynomial is defined by: |
949 | * @f[ |
950 | * L_n(x) = \frac{e^x}{n!} \frac{d^n}{dx^n} (x^ne^{-x}) |
951 | * @f] |
952 | * |
953 | * @tparam _Tp The floating-point type of the argument @c __x. |
954 | * @param __n The nonnegative order |
955 | * @param __x The argument <tt> __x >= 0 </tt> |
956 | * @throw std::domain_error if <tt> __x < 0 </tt>. |
957 | */ |
958 | template<typename _Tp> |
959 | inline typename __gnu_cxx::__promote<_Tp>::__type |
960 | laguerre(unsigned int __n, _Tp __x) |
961 | { |
962 | typedef typename __gnu_cxx::__promote<_Tp>::__type __type; |
963 | return __detail::__laguerre<__type>(__n, __x); |
964 | } |
965 | |
966 | // Legendre polynomials |
967 | |
968 | /** |
969 | * Return the Legendre polynomial @f$ P_l(x) @f$ of nonnegative |
970 | * degree @f$ l @f$ and @c float argument @f$ |x| <= 0 @f$. |
971 | * |
972 | * @see legendre for more details. |
973 | */ |
974 | inline float |
975 | legendref(unsigned int __l, float __x) |
976 | { return __detail::__poly_legendre_p<float>(__l, __x); } |
977 | |
978 | /** |
979 | * Return the Legendre polynomial @f$ P_l(x) @f$ of nonnegative |
980 | * degree @f$ l @f$ and <tt>long double</tt> argument @f$ |x| <= 0 @f$. |
981 | * |
982 | * @see legendre for more details. |
983 | */ |
984 | inline long double |
985 | legendrel(unsigned int __l, long double __x) |
986 | { return __detail::__poly_legendre_p<long double>(__l, __x); } |
987 | |
988 | /** |
989 | * Return the Legendre polynomial @f$ P_l(x) @f$ of nonnegative |
990 | * degree @f$ l @f$ and real argument @f$ |x| <= 0 @f$. |
991 | * |
992 | * The Legendre function of order @f$ l @f$ and argument @f$ x @f$, |
993 | * @f$ P_l(x) @f$, is defined by: |
994 | * @f[ |
995 | * P_l(x) = \frac{1}{2^l l!}\frac{d^l}{dx^l}(x^2 - 1)^{l} |
996 | * @f] |
997 | * |
998 | * @tparam _Tp The floating-point type of the argument @c __x. |
999 | * @param __l The degree @f$ l >= 0 @f$ |
1000 | * @param __x The argument @c abs(__x) <= 1 |
1001 | * @throw std::domain_error if @c abs(__x) > 1 |
1002 | */ |
1003 | template<typename _Tp> |
1004 | inline typename __gnu_cxx::__promote<_Tp>::__type |
1005 | legendre(unsigned int __l, _Tp __x) |
1006 | { |
1007 | typedef typename __gnu_cxx::__promote<_Tp>::__type __type; |
1008 | return __detail::__poly_legendre_p<__type>(__l, __x); |
1009 | } |
1010 | |
1011 | // Riemann zeta functions |
1012 | |
1013 | /** |
1014 | * Return the Riemann zeta function @f$ \zeta(s) @f$ |
1015 | * for @c float argument @f$ s @f$. |
1016 | * |
1017 | * @see riemann_zeta for more details. |
1018 | */ |
1019 | inline float |
1020 | riemann_zetaf(float __s) |
1021 | { return __detail::__riemann_zeta<float>(__s); } |
1022 | |
1023 | /** |
1024 | * Return the Riemann zeta function @f$ \zeta(s) @f$ |
1025 | * for <tt>long double</tt> argument @f$ s @f$. |
1026 | * |
1027 | * @see riemann_zeta for more details. |
1028 | */ |
1029 | inline long double |
1030 | riemann_zetal(long double __s) |
1031 | { return __detail::__riemann_zeta<long double>(__s); } |
1032 | |
1033 | /** |
1034 | * Return the Riemann zeta function @f$ \zeta(s) @f$ |
1035 | * for real argument @f$ s @f$. |
1036 | * |
1037 | * The Riemann zeta function is defined by: |
1038 | * @f[ |
1039 | * \zeta(s) = \sum_{k=1}^{\infty} k^{-s} \hbox{ for } s > 1 |
1040 | * @f] |
1041 | * and |
1042 | * @f[ |
1043 | * \zeta(s) = \frac{1}{1-2^{1-s}}\sum_{k=1}^{\infty}(-1)^{k-1}k^{-s} |
1044 | * \hbox{ for } 0 <= s <= 1 |
1045 | * @f] |
1046 | * For s < 1 use the reflection formula: |
1047 | * @f[ |
1048 | * \zeta(s) = 2^s \pi^{s-1} \sin(\frac{\pi s}{2}) \Gamma(1-s) \zeta(1-s) |
1049 | * @f] |
1050 | * |
1051 | * @tparam _Tp The floating-point type of the argument @c __s. |
1052 | * @param __s The argument <tt> s != 1 </tt> |
1053 | */ |
1054 | template<typename _Tp> |
1055 | inline typename __gnu_cxx::__promote<_Tp>::__type |
1056 | riemann_zeta(_Tp __s) |
1057 | { |
1058 | typedef typename __gnu_cxx::__promote<_Tp>::__type __type; |
1059 | return __detail::__riemann_zeta<__type>(__s); |
1060 | } |
1061 | |
1062 | // Spherical Bessel functions |
1063 | |
1064 | /** |
1065 | * Return the spherical Bessel function @f$ j_n(x) @f$ of nonnegative order n |
1066 | * and @c float argument @f$ x >= 0 @f$. |
1067 | * |
1068 | * @see sph_bessel for more details. |
1069 | */ |
1070 | inline float |
1071 | sph_besself(unsigned int __n, float __x) |
1072 | { return __detail::__sph_bessel<float>(__n, __x); } |
1073 | |
1074 | /** |
1075 | * Return the spherical Bessel function @f$ j_n(x) @f$ of nonnegative order n |
1076 | * and <tt>long double</tt> argument @f$ x >= 0 @f$. |
1077 | * |
1078 | * @see sph_bessel for more details. |
1079 | */ |
1080 | inline long double |
1081 | sph_bessell(unsigned int __n, long double __x) |
1082 | { return __detail::__sph_bessel<long double>(__n, __x); } |
1083 | |
1084 | /** |
1085 | * Return the spherical Bessel function @f$ j_n(x) @f$ of nonnegative order n |
1086 | * and real argument @f$ x >= 0 @f$. |
1087 | * |
1088 | * The spherical Bessel function is defined by: |
1089 | * @f[ |
1090 | * j_n(x) = \left(\frac{\pi}{2x} \right) ^{1/2} J_{n+1/2}(x) |
1091 | * @f] |
1092 | * |
1093 | * @tparam _Tp The floating-point type of the argument @c __x. |
1094 | * @param __n The integral order <tt> n >= 0 </tt> |
1095 | * @param __x The real argument <tt> x >= 0 </tt> |
1096 | * @throw std::domain_error if <tt> __x < 0 </tt>. |
1097 | */ |
1098 | template<typename _Tp> |
1099 | inline typename __gnu_cxx::__promote<_Tp>::__type |
1100 | sph_bessel(unsigned int __n, _Tp __x) |
1101 | { |
1102 | typedef typename __gnu_cxx::__promote<_Tp>::__type __type; |
1103 | return __detail::__sph_bessel<__type>(__n, __x); |
1104 | } |
1105 | |
1106 | // Spherical associated Legendre functions |
1107 | |
1108 | /** |
1109 | * Return the spherical Legendre function of nonnegative integral |
1110 | * degree @c l and order @c m and float angle @f$ \theta @f$ in radians. |
1111 | * |
1112 | * @see sph_legendre for details. |
1113 | */ |
1114 | inline float |
1115 | sph_legendref(unsigned int __l, unsigned int __m, float __theta) |
1116 | { return __detail::__sph_legendre<float>(__l, __m, __theta); } |
1117 | |
1118 | /** |
1119 | * Return the spherical Legendre function of nonnegative integral |
1120 | * degree @c l and order @c m and <tt>long double</tt> angle @f$ \theta @f$ |
1121 | * in radians. |
1122 | * |
1123 | * @see sph_legendre for details. |
1124 | */ |
1125 | inline long double |
1126 | sph_legendrel(unsigned int __l, unsigned int __m, long double __theta) |
1127 | { return __detail::__sph_legendre<long double>(__l, __m, __theta); } |
1128 | |
1129 | /** |
1130 | * Return the spherical Legendre function of nonnegative integral |
1131 | * degree @c l and order @c m and real angle @f$ \theta @f$ in radians. |
1132 | * |
1133 | * The spherical Legendre function is defined by |
1134 | * @f[ |
1135 | * Y_l^m(\theta,\phi) = (-1)^m[\frac{(2l+1)}{4\pi} |
1136 | * \frac{(l-m)!}{(l+m)!}] |
1137 | * P_l^m(\cos\theta) \exp^{im\phi} |
1138 | * @f] |
1139 | * |
1140 | * @tparam _Tp The floating-point type of the angle @c __theta. |
1141 | * @param __l The order <tt> __l >= 0 </tt> |
1142 | * @param __m The degree <tt> __m >= 0 </tt> and <tt> __m <= __l </tt> |
1143 | * @param __theta The radian polar angle argument |
1144 | */ |
1145 | template<typename _Tp> |
1146 | inline typename __gnu_cxx::__promote<_Tp>::__type |
1147 | sph_legendre(unsigned int __l, unsigned int __m, _Tp __theta) |
1148 | { |
1149 | typedef typename __gnu_cxx::__promote<_Tp>::__type __type; |
1150 | return __detail::__sph_legendre<__type>(__l, __m, __theta); |
1151 | } |
1152 | |
1153 | // Spherical Neumann functions |
1154 | |
1155 | /** |
1156 | * Return the spherical Neumann function of integral order @f$ n >= 0 @f$ |
1157 | * and @c float argument @f$ x >= 0 @f$. |
1158 | * |
1159 | * @see sph_neumann for details. |
1160 | */ |
1161 | inline float |
1162 | sph_neumannf(unsigned int __n, float __x) |
1163 | { return __detail::__sph_neumann<float>(__n, __x); } |
1164 | |
1165 | /** |
1166 | * Return the spherical Neumann function of integral order @f$ n >= 0 @f$ |
1167 | * and <tt>long double</tt> @f$ x >= 0 @f$. |
1168 | * |
1169 | * @see sph_neumann for details. |
1170 | */ |
1171 | inline long double |
1172 | sph_neumannl(unsigned int __n, long double __x) |
1173 | { return __detail::__sph_neumann<long double>(__n, __x); } |
1174 | |
1175 | /** |
1176 | * Return the spherical Neumann function of integral order @f$ n >= 0 @f$ |
1177 | * and real argument @f$ x >= 0 @f$. |
1178 | * |
1179 | * The spherical Neumann function is defined by |
1180 | * @f[ |
1181 | * n_n(x) = \left(\frac{\pi}{2x} \right) ^{1/2} N_{n+1/2}(x) |
1182 | * @f] |
1183 | * |
1184 | * @tparam _Tp The floating-point type of the argument @c __x. |
1185 | * @param __n The integral order <tt> n >= 0 </tt> |
1186 | * @param __x The real argument <tt> __x >= 0 </tt> |
1187 | * @throw std::domain_error if <tt> __x < 0 </tt>. |
1188 | */ |
1189 | template<typename _Tp> |
1190 | inline typename __gnu_cxx::__promote<_Tp>::__type |
1191 | sph_neumann(unsigned int __n, _Tp __x) |
1192 | { |
1193 | typedef typename __gnu_cxx::__promote<_Tp>::__type __type; |
1194 | return __detail::__sph_neumann<__type>(__n, __x); |
1195 | } |
1196 | |
1197 | /// @} group mathsf |
1198 | |
1199 | _GLIBCXX_END_NAMESPACE_VERSION |
1200 | } // namespace std |
1201 | |
1202 | #ifndef __STRICT_ANSI__ |
1203 | namespace __gnu_cxx _GLIBCXX_VISIBILITY(default) |
1204 | { |
1205 | _GLIBCXX_BEGIN_NAMESPACE_VERSION |
1206 | |
1207 | /** @addtogroup mathsf |
1208 | * @{ |
1209 | */ |
1210 | |
1211 | // Airy functions |
1212 | |
1213 | /** |
1214 | * Return the Airy function @f$ Ai(x) @f$ of @c float argument x. |
1215 | */ |
1216 | inline float |
1217 | airy_aif(float __x) |
1218 | { |
1219 | float __Ai, __Bi, __Aip, __Bip; |
1220 | std::__detail::__airy<float>(__x, __Ai, __Bi, __Aip, __Bip); |
1221 | return __Ai; |
1222 | } |
1223 | |
1224 | /** |
1225 | * Return the Airy function @f$ Ai(x) @f$ of <tt>long double</tt> argument x. |
1226 | */ |
1227 | inline long double |
1228 | airy_ail(long double __x) |
1229 | { |
1230 | long double __Ai, __Bi, __Aip, __Bip; |
1231 | std::__detail::__airy<long double>(__x, __Ai, __Bi, __Aip, __Bip); |
1232 | return __Ai; |
1233 | } |
1234 | |
1235 | /** |
1236 | * Return the Airy function @f$ Ai(x) @f$ of real argument x. |
1237 | */ |
1238 | template<typename _Tp> |
1239 | inline typename __gnu_cxx::__promote<_Tp>::__type |
1240 | airy_ai(_Tp __x) |
1241 | { |
1242 | typedef typename __gnu_cxx::__promote<_Tp>::__type __type; |
1243 | __type __Ai, __Bi, __Aip, __Bip; |
1244 | std::__detail::__airy<__type>(__x, __Ai, __Bi, __Aip, __Bip); |
1245 | return __Ai; |
1246 | } |
1247 | |
1248 | /** |
1249 | * Return the Airy function @f$ Bi(x) @f$ of @c float argument x. |
1250 | */ |
1251 | inline float |
1252 | airy_bif(float __x) |
1253 | { |
1254 | float __Ai, __Bi, __Aip, __Bip; |
1255 | std::__detail::__airy<float>(__x, __Ai, __Bi, __Aip, __Bip); |
1256 | return __Bi; |
1257 | } |
1258 | |
1259 | /** |
1260 | * Return the Airy function @f$ Bi(x) @f$ of <tt>long double</tt> argument x. |
1261 | */ |
1262 | inline long double |
1263 | airy_bil(long double __x) |
1264 | { |
1265 | long double __Ai, __Bi, __Aip, __Bip; |
1266 | std::__detail::__airy<long double>(__x, __Ai, __Bi, __Aip, __Bip); |
1267 | return __Bi; |
1268 | } |
1269 | |
1270 | /** |
1271 | * Return the Airy function @f$ Bi(x) @f$ of real argument x. |
1272 | */ |
1273 | template<typename _Tp> |
1274 | inline typename __gnu_cxx::__promote<_Tp>::__type |
1275 | airy_bi(_Tp __x) |
1276 | { |
1277 | typedef typename __gnu_cxx::__promote<_Tp>::__type __type; |
1278 | __type __Ai, __Bi, __Aip, __Bip; |
1279 | std::__detail::__airy<__type>(__x, __Ai, __Bi, __Aip, __Bip); |
1280 | return __Bi; |
1281 | } |
1282 | |
1283 | // Confluent hypergeometric functions |
1284 | |
1285 | /** |
1286 | * Return the confluent hypergeometric function @f$ {}_1F_1(a;c;x) @f$ |
1287 | * of @c float numeratorial parameter @c a, denominatorial parameter @c c, |
1288 | * and argument @c x. |
1289 | * |
1290 | * @see conf_hyperg for details. |
1291 | */ |
1292 | inline float |
1293 | conf_hypergf(float __a, float __c, float __x) |
1294 | { return std::__detail::__conf_hyperg<float>(__a, __c, __x); } |
1295 | |
1296 | /** |
1297 | * Return the confluent hypergeometric function @f$ {}_1F_1(a;c;x) @f$ |
1298 | * of <tt>long double</tt> numeratorial parameter @c a, |
1299 | * denominatorial parameter @c c, and argument @c x. |
1300 | * |
1301 | * @see conf_hyperg for details. |
1302 | */ |
1303 | inline long double |
1304 | conf_hypergl(long double __a, long double __c, long double __x) |
1305 | { return std::__detail::__conf_hyperg<long double>(__a, __c, __x); } |
1306 | |
1307 | /** |
1308 | * Return the confluent hypergeometric function @f$ {}_1F_1(a;c;x) @f$ |
1309 | * of real numeratorial parameter @c a, denominatorial parameter @c c, |
1310 | * and argument @c x. |
1311 | * |
1312 | * The confluent hypergeometric function is defined by |
1313 | * @f[ |
1314 | * {}_1F_1(a;c;x) = \sum_{n=0}^{\infty} \frac{(a)_n x^n}{(c)_n n!} |
1315 | * @f] |
1316 | * where the Pochhammer symbol is @f$ (x)_k = (x)(x+1)...(x+k-1) @f$, |
1317 | * @f$ (x)_0 = 1 @f$ |
1318 | * |
1319 | * @param __a The numeratorial parameter |
1320 | * @param __c The denominatorial parameter |
1321 | * @param __x The argument |
1322 | */ |
1323 | template<typename _Tpa, typename _Tpc, typename _Tp> |
1324 | inline typename __gnu_cxx::__promote_3<_Tpa, _Tpc, _Tp>::__type |
1325 | conf_hyperg(_Tpa __a, _Tpc __c, _Tp __x) |
1326 | { |
1327 | typedef typename __gnu_cxx::__promote_3<_Tpa, _Tpc, _Tp>::__type __type; |
1328 | return std::__detail::__conf_hyperg<__type>(__a, __c, __x); |
1329 | } |
1330 | |
1331 | // Hypergeometric functions |
1332 | |
1333 | /** |
1334 | * Return the hypergeometric function @f$ {}_2F_1(a,b;c;x) @f$ |
1335 | * of @ float numeratorial parameters @c a and @c b, |
1336 | * denominatorial parameter @c c, and argument @c x. |
1337 | * |
1338 | * @see hyperg for details. |
1339 | */ |
1340 | inline float |
1341 | hypergf(float __a, float __b, float __c, float __x) |
1342 | { return std::__detail::__hyperg<float>(__a, __b, __c, __x); } |
1343 | |
1344 | /** |
1345 | * Return the hypergeometric function @f$ {}_2F_1(a,b;c;x) @f$ |
1346 | * of <tt>long double</tt> numeratorial parameters @c a and @c b, |
1347 | * denominatorial parameter @c c, and argument @c x. |
1348 | * |
1349 | * @see hyperg for details. |
1350 | */ |
1351 | inline long double |
1352 | hypergl(long double __a, long double __b, long double __c, long double __x) |
1353 | { return std::__detail::__hyperg<long double>(__a, __b, __c, __x); } |
1354 | |
1355 | /** |
1356 | * Return the hypergeometric function @f$ {}_2F_1(a,b;c;x) @f$ |
1357 | * of real numeratorial parameters @c a and @c b, |
1358 | * denominatorial parameter @c c, and argument @c x. |
1359 | * |
1360 | * The hypergeometric function is defined by |
1361 | * @f[ |
1362 | * {}_2F_1(a;c;x) = \sum_{n=0}^{\infty} \frac{(a)_n (b)_n x^n}{(c)_n n!} |
1363 | * @f] |
1364 | * where the Pochhammer symbol is @f$ (x)_k = (x)(x+1)...(x+k-1) @f$, |
1365 | * @f$ (x)_0 = 1 @f$ |
1366 | * |
1367 | * @param __a The first numeratorial parameter |
1368 | * @param __b The second numeratorial parameter |
1369 | * @param __c The denominatorial parameter |
1370 | * @param __x The argument |
1371 | */ |
1372 | template<typename _Tpa, typename _Tpb, typename _Tpc, typename _Tp> |
1373 | inline typename __gnu_cxx::__promote_4<_Tpa, _Tpb, _Tpc, _Tp>::__type |
1374 | hyperg(_Tpa __a, _Tpb __b, _Tpc __c, _Tp __x) |
1375 | { |
1376 | typedef typename __gnu_cxx::__promote_4<_Tpa, _Tpb, _Tpc, _Tp> |
1377 | ::__type __type; |
1378 | return std::__detail::__hyperg<__type>(__a, __b, __c, __x); |
1379 | } |
1380 | |
1381 | /// @} |
1382 | _GLIBCXX_END_NAMESPACE_VERSION |
1383 | } // namespace __gnu_cxx |
1384 | #endif // __STRICT_ANSI__ |
1385 | |
1386 | #endif // _GLIBCXX_BITS_SPECFUN_H |
1387 | |