libstdc++
specfun.h
Go to the documentation of this file.
00001 // Mathematical Special Functions for -*- C++ -*-
00002 
00003 // Copyright (C) 2006-2017 Free Software Foundation, Inc.
00004 //
00005 // This file is part of the GNU ISO C++ Library.  This library is free
00006 // software; you can redistribute it and/or modify it under the
00007 // terms of the GNU General Public License as published by the
00008 // Free Software Foundation; either version 3, or (at your option)
00009 // any later version.
00010 
00011 // This library is distributed in the hope that it will be useful,
00012 // but WITHOUT ANY WARRANTY; without even the implied warranty of
00013 // MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
00014 // GNU General Public License for more details.
00015 
00016 // Under Section 7 of GPL version 3, you are granted additional
00017 // permissions described in the GCC Runtime Library Exception, version
00018 // 3.1, as published by the Free Software Foundation.
00019 
00020 // You should have received a copy of the GNU General Public License and
00021 // a copy of the GCC Runtime Library Exception along with this program;
00022 // see the files COPYING3 and COPYING.RUNTIME respectively.  If not, see
00023 // <http://www.gnu.org/licenses/>.
00024 
00025 /** @file bits/specfun.h
00026  *  This is an internal header file, included by other library headers.
00027  *  Do not attempt to use it directly. @headername{cmath}
00028  */
00029 
00030 #ifndef _GLIBCXX_BITS_SPECFUN_H
00031 #define _GLIBCXX_BITS_SPECFUN_H 1
00032 
00033 #pragma GCC visibility push(default)
00034 
00035 #include <bits/c++config.h>
00036 
00037 #define __STDCPP_MATH_SPEC_FUNCS__ 201003L
00038 
00039 #define __cpp_lib_math_special_functions 201603L
00040 
00041 #if __cplusplus <= 201403L && __STDCPP_WANT_MATH_SPEC_FUNCS__ == 0
00042 # error include <cmath> and define __STDCPP_WANT_MATH_SPEC_FUNCS__
00043 #endif
00044 
00045 #include <bits/stl_algobase.h>
00046 #include <limits>
00047 #include <type_traits>
00048 
00049 #include <tr1/gamma.tcc>
00050 #include <tr1/bessel_function.tcc>
00051 #include <tr1/beta_function.tcc>
00052 #include <tr1/ell_integral.tcc>
00053 #include <tr1/exp_integral.tcc>
00054 #include <tr1/hypergeometric.tcc>
00055 #include <tr1/legendre_function.tcc>
00056 #include <tr1/modified_bessel_func.tcc>
00057 #include <tr1/poly_hermite.tcc>
00058 #include <tr1/poly_laguerre.tcc>
00059 #include <tr1/riemann_zeta.tcc>
00060 
00061 namespace std _GLIBCXX_VISIBILITY(default)
00062 {
00063 _GLIBCXX_BEGIN_NAMESPACE_VERSION
00064 
00065   /**
00066    * @defgroup mathsf Mathematical Special Functions
00067    * @ingroup numerics
00068    *
00069    * A collection of advanced mathematical special functions,
00070    * defined by ISO/IEC IS 29124.
00071    * @{
00072    */
00073 
00074   /**
00075    * @mainpage Mathematical Special Functions
00076    *
00077    * @section intro Introduction and History
00078    * The first significant library upgrade on the road to C++2011,
00079    * <a href="http://www.open-std.org/JTC1/SC22/WG21/docs/papers/2005/n1836.pdf">
00080    * TR1</a>, included a set of 23 mathematical functions that significantly
00081    * extended the standard transcendental functions inherited from C and declared
00082    * in @<cmath@>.
00083    *
00084    * Although most components from TR1 were eventually adopted for C++11 these
00085    * math functions were left behind out of concern for implementability.
00086    * The math functions were published as a separate international standard
00087    * <a href="http://www.open-std.org/JTC1/SC22/WG21/docs/papers/2010/n3060.pdf">
00088    * IS 29124 - Extensions to the C++ Library to Support Mathematical Special
00089    * Functions</a>.
00090    *
00091    * For C++17 these functions were incorporated into the main standard.
00092    *
00093    * @section contents Contents
00094    * The following functions are implemented in namespace @c std:
00095    * - @ref assoc_laguerre "assoc_laguerre - Associated Laguerre functions"
00096    * - @ref assoc_legendre "assoc_legendre - Associated Legendre functions"
00097    * - @ref beta "beta - Beta functions"
00098    * - @ref comp_ellint_1 "comp_ellint_1 - Complete elliptic functions of the first kind"
00099    * - @ref comp_ellint_2 "comp_ellint_2 - Complete elliptic functions of the second kind"
00100    * - @ref comp_ellint_3 "comp_ellint_3 - Complete elliptic functions of the third kind"
00101    * - @ref cyl_bessel_i "cyl_bessel_i - Regular modified cylindrical Bessel functions"
00102    * - @ref cyl_bessel_j "cyl_bessel_j - Cylindrical Bessel functions of the first kind"
00103    * - @ref cyl_bessel_k "cyl_bessel_k - Irregular modified cylindrical Bessel functions"
00104    * - @ref cyl_neumann "cyl_neumann - Cylindrical Neumann functions or Cylindrical Bessel functions of the second kind"
00105    * - @ref ellint_1 "ellint_1 - Incomplete elliptic functions of the first kind"
00106    * - @ref ellint_2 "ellint_2 - Incomplete elliptic functions of the second kind"
00107    * - @ref ellint_3 "ellint_3 - Incomplete elliptic functions of the third kind"
00108    * - @ref expint "expint - The exponential integral"
00109    * - @ref hermite "hermite - Hermite polynomials"
00110    * - @ref laguerre "laguerre - Laguerre functions"
00111    * - @ref legendre "legendre - Legendre polynomials"
00112    * - @ref riemann_zeta "riemann_zeta - The Riemann zeta function"
00113    * - @ref sph_bessel "sph_bessel - Spherical Bessel functions"
00114    * - @ref sph_legendre "sph_legendre - Spherical Legendre functions"
00115    * - @ref sph_neumann "sph_neumann - Spherical Neumann functions"
00116    *
00117    * The hypergeometric functions were stricken from the TR29124 and C++17
00118    * versions of this math library because of implementation concerns.
00119    * However, since they were in the TR1 version and since they are popular
00120    * we kept them as an extension in namespace @c __gnu_cxx:
00121    * - @ref conf_hyperg "conf_hyperg - Confluent hypergeometric functions"
00122    * - @ref hyperg "hyperg - Hypergeometric functions"
00123    *
00124    * @section general General Features
00125    *
00126    * @subsection promotion Argument Promotion
00127    * The arguments suppled to the non-suffixed functions will be promoted
00128    * according to the following rules:
00129    * 1. If any argument intended to be floating point is given an integral value
00130    * That integral value is promoted to double.
00131    * 2. All floating point arguments are promoted up to the largest floating
00132    *    point precision among them.
00133    *
00134    * @subsection NaN NaN Arguments
00135    * If any of the floating point arguments supplied to these functions is
00136    * invalid or NaN (std::numeric_limits<Tp>::quiet_NaN),
00137    * the value NaN is returned.
00138    *
00139    * @section impl Implementation
00140    *
00141    * We strive to implement the underlying math with type generic algorithms
00142    * to the greatest extent possible.  In practice, the functions are thin
00143    * wrappers that dispatch to function templates. Type dependence is
00144    * controlled with std::numeric_limits and functions thereof.
00145    *
00146    * We don't promote @c float to @c double or @c double to <tt>long double</tt>
00147    * reflexively.  The goal is for @c float functions to operate more quickly,
00148    * at the cost of @c float accuracy and possibly a smaller domain of validity.
00149    * Similaryly, <tt>long double</tt> should give you more dynamic range
00150    * and slightly more pecision than @c double on many systems.
00151    *
00152    * @section testing Testing
00153    *
00154    * These functions have been tested against equivalent implementations
00155    * from the <a href="http://www.gnu.org/software/gsl">
00156    * Gnu Scientific Library, GSL</a> and
00157    * <a href="http://www.boost.org/doc/libs/1_60_0/libs/math/doc/html/index.html>Boost</a>
00158    * and the ratio
00159    * @f[
00160    *   \frac{|f - f_{test}|}{|f_{test}|}
00161    * @f]
00162    * is generally found to be within 10^-15 for 64-bit double on linux-x86_64 systems
00163    * over most of the ranges of validity.
00164    * 
00165    * @todo Provide accuracy comparisons on a per-function basis for a small
00166    *       number of targets.
00167    *
00168    * @section bibliography General Bibliography
00169    *
00170    * @see Abramowitz and Stegun: Handbook of Mathematical Functions,
00171    * with Formulas, Graphs, and Mathematical Tables
00172    * Edited by Milton Abramowitz and Irene A. Stegun,
00173    * National Bureau of Standards  Applied Mathematics Series - 55
00174    * Issued June 1964, Tenth Printing, December 1972, with corrections
00175    * Electronic versions of A&S abound including both pdf and navigable html.
00176    * @see for example  http://people.math.sfu.ca/~cbm/aands/
00177    *
00178    * @see The old A&S has been redone as the
00179    * NIST Digital Library of Mathematical Functions: http://dlmf.nist.gov/
00180    * This version is far more navigable and includes more recent work.
00181    *
00182    * @see An Atlas of Functions: with Equator, the Atlas Function Calculator
00183    * 2nd Edition, by Oldham, Keith B., Myland, Jan, Spanier, Jerome
00184    *
00185    * @see Asymptotics and Special Functions by Frank W. J. Olver,
00186    * Academic Press, 1974
00187    *
00188    * @see Numerical Recipes in C, The Art of Scientific Computing,
00189    * by William H. Press, Second Ed., Saul A. Teukolsky,
00190    * William T. Vetterling, and Brian P. Flannery,
00191    * Cambridge University Press, 1992
00192    *
00193    * @see The Special Functions and Their Approximations: Volumes 1 and 2,
00194    * by Yudell L. Luke, Academic Press, 1969
00195    */
00196 
00197   // Associated Laguerre polynomials
00198 
00199   /**
00200    * Return the associated Laguerre polynomial of order @c n,
00201    * degree @c m: @f$ L_n^m(x) @f$ for @c float argument.
00202    *
00203    * @see assoc_laguerre for more details.
00204    */
00205   inline float
00206   assoc_laguerref(unsigned int __n, unsigned int __m, float __x)
00207   { return __detail::__assoc_laguerre<float>(__n, __m, __x); }
00208 
00209   /**
00210    * Return the associated Laguerre polynomial of order @c n,
00211    * degree @c m: @f$ L_n^m(x) @f$.
00212    *
00213    * @see assoc_laguerre for more details.
00214    */
00215   inline long double
00216   assoc_laguerrel(unsigned int __n, unsigned int __m, long double __x)
00217   { return __detail::__assoc_laguerre<long double>(__n, __m, __x); }
00218 
00219   /**
00220    * Return the associated Laguerre polynomial of nonnegative order @c n,
00221    * nonnegative degree @c m and real argument @c x: @f$ L_n^m(x) @f$.
00222    *
00223    * The associated Laguerre function of real degree @f$ \alpha @f$,
00224    * @f$ L_n^\alpha(x) @f$, is defined by
00225    * @f[
00226    *     L_n^\alpha(x) = \frac{(\alpha + 1)_n}{n!}
00227    *                     {}_1F_1(-n; \alpha + 1; x)
00228    * @f]
00229    * where @f$ (\alpha)_n @f$ is the Pochhammer symbol and
00230    * @f$ {}_1F_1(a; c; x) @f$ is the confluent hypergeometric function.
00231    *
00232    * The associated Laguerre polynomial is defined for integral
00233    * degree @f$ \alpha = m @f$ by:
00234    * @f[
00235    *     L_n^m(x) = (-1)^m \frac{d^m}{dx^m} L_{n + m}(x)
00236    * @f]
00237    * where the Laguerre polynomial is defined by:
00238    * @f[
00239    *     L_n(x) = \frac{e^x}{n!} \frac{d^n}{dx^n} (x^ne^{-x})
00240    * @f]
00241    * and @f$ x >= 0 @f$.
00242    * @see laguerre for details of the Laguerre function of degree @c n
00243    *
00244    * @tparam _Tp The floating-point type of the argument @c __x.
00245    * @param __n The order of the Laguerre function, <tt>__n >= 0</tt>.
00246    * @param __m The degree of the Laguerre function, <tt>__m >= 0</tt>.
00247    * @param __x The argument of the Laguerre function, <tt>__x >= 0</tt>.
00248    * @throw std::domain_error if <tt>__x < 0</tt>.
00249    */
00250   template<typename _Tp>
00251     inline typename __gnu_cxx::__promote<_Tp>::__type
00252     assoc_laguerre(unsigned int __n, unsigned int __m, _Tp __x)
00253     {
00254       typedef typename __gnu_cxx::__promote<_Tp>::__type __type;
00255       return __detail::__assoc_laguerre<__type>(__n, __m, __x);
00256     }
00257 
00258   // Associated Legendre functions
00259 
00260   /**
00261    * Return the associated Legendre function of degree @c l and order @c m
00262    * for @c float argument.
00263    *
00264    * @see assoc_legendre for more details.
00265    */
00266   inline float
00267   assoc_legendref(unsigned int __l, unsigned int __m, float __x)
00268   { return __detail::__assoc_legendre_p<float>(__l, __m, __x); }
00269 
00270   /**
00271    * Return the associated Legendre function of degree @c l and order @c m.
00272    *
00273    * @see assoc_legendre for more details.
00274    */
00275   inline long double
00276   assoc_legendrel(unsigned int __l, unsigned int __m, long double __x)
00277   { return __detail::__assoc_legendre_p<long double>(__l, __m, __x); }
00278 
00279 
00280   /**
00281    * Return the associated Legendre function of degree @c l and order @c m.
00282    *
00283    * The associated Legendre function is derived from the Legendre function
00284    * @f$ P_l(x) @f$ by the Rodrigues formula:
00285    * @f[
00286    *   P_l^m(x) = (1 - x^2)^{m/2}\frac{d^m}{dx^m}P_l(x)
00287    * @f]
00288    * @see legendre for details of the Legendre function of degree @c l
00289    *
00290    * @tparam _Tp The floating-point type of the argument @c __x.
00291    * @param  __l  The degree <tt>__l >= 0</tt>.
00292    * @param  __m  The order <tt>__m <= l</tt>.
00293    * @param  __x  The argument, <tt>abs(__x) <= 1</tt>.
00294    * @throw std::domain_error if <tt>abs(__x) > 1</tt>.
00295    */
00296   template<typename _Tp>
00297     inline typename __gnu_cxx::__promote<_Tp>::__type
00298     assoc_legendre(unsigned int __l, unsigned int __m, _Tp __x)
00299     {
00300       typedef typename __gnu_cxx::__promote<_Tp>::__type __type;
00301       return __detail::__assoc_legendre_p<__type>(__l, __m, __x);
00302     }
00303 
00304   // Beta functions
00305 
00306   /**
00307    * Return the beta function, @f$ B(a,b) @f$, for @c float parameters @c a, @c b.
00308    *
00309    * @see beta for more details.
00310    */
00311   inline float
00312   betaf(float __a, float __b)
00313   { return __detail::__beta<float>(__a, __b); }
00314 
00315   /**
00316    * Return the beta function, @f$B(a,b)@f$, for long double
00317    * parameters @c a, @c b.
00318    *
00319    * @see beta for more details.
00320    */
00321   inline long double
00322   betal(long double __a, long double __b)
00323   { return __detail::__beta<long double>(__a, __b); }
00324 
00325   /**
00326    * Return the beta function, @f$B(a,b)@f$, for real parameters @c a, @c b.
00327    *
00328    * The beta function is defined by
00329    * @f[
00330    *   B(a,b) = \int_0^1 t^{a - 1} (1 - t)^{b - 1} dt
00331    *          = \frac{\Gamma(a)\Gamma(b)}{\Gamma(a+b)}
00332    * @f]
00333    * where @f$ a > 0 @f$ and @f$ b > 0 @f$
00334    *
00335    * @tparam _Tpa The floating-point type of the parameter @c __a.
00336    * @tparam _Tpb The floating-point type of the parameter @c __b.
00337    * @param __a The first argument of the beta function, <tt> __a > 0 </tt>.
00338    * @param __b The second argument of the beta function, <tt> __b > 0 </tt>.
00339    * @throw std::domain_error if <tt> __a < 0 </tt> or <tt> __b < 0 </tt>.
00340    */
00341   template<typename _Tpa, typename _Tpb>
00342     inline typename __gnu_cxx::__promote_2<_Tpa, _Tpb>::__type
00343     beta(_Tpa __a, _Tpb __b)
00344     {
00345       typedef typename __gnu_cxx::__promote_2<_Tpa, _Tpb>::__type __type;
00346       return __detail::__beta<__type>(__a, __b);
00347     }
00348 
00349   // Complete elliptic integrals of the first kind
00350 
00351   /**
00352    * Return the complete elliptic integral of the first kind @f$ E(k) @f$
00353    * for @c float modulus @c k.
00354    *
00355    * @see comp_ellint_1 for details.
00356    */
00357   inline float
00358   comp_ellint_1f(float __k)
00359   { return __detail::__comp_ellint_1<float>(__k); }
00360 
00361   /**
00362    * Return the complete elliptic integral of the first kind @f$ E(k) @f$
00363    * for long double modulus @c k.
00364    *
00365    * @see comp_ellint_1 for details.
00366    */
00367   inline long double
00368   comp_ellint_1l(long double __k)
00369   { return __detail::__comp_ellint_1<long double>(__k); }
00370 
00371   /**
00372    * Return the complete elliptic integral of the first kind
00373    * @f$ K(k) @f$ for real modulus @c k.
00374    *
00375    * The complete elliptic integral of the first kind is defined as
00376    * @f[
00377    *   K(k) = F(k,\pi/2) = \int_0^{\pi/2}\frac{d\theta}
00378    *                                         {\sqrt{1 - k^2 sin^2\theta}}
00379    * @f]
00380    * where @f$ F(k,\phi) @f$ is the incomplete elliptic integral of the
00381    * first kind and the modulus @f$ |k| <= 1 @f$.
00382    * @see ellint_1 for details of the incomplete elliptic function
00383    * of the first kind.
00384    *
00385    * @tparam _Tp The floating-point type of the modulus @c __k.
00386    * @param  __k  The modulus, <tt> abs(__k) <= 1 </tt>
00387    * @throw std::domain_error if <tt> abs(__k) > 1 </tt>.
00388    */
00389   template<typename _Tp>
00390     inline typename __gnu_cxx::__promote<_Tp>::__type
00391     comp_ellint_1(_Tp __k)
00392     {
00393       typedef typename __gnu_cxx::__promote<_Tp>::__type __type;
00394       return __detail::__comp_ellint_1<__type>(__k);
00395     }
00396 
00397   // Complete elliptic integrals of the second kind
00398 
00399   /**
00400    * Return the complete elliptic integral of the second kind @f$ E(k) @f$
00401    * for @c float modulus @c k.
00402    *
00403    * @see comp_ellint_2 for details.
00404    */
00405   inline float
00406   comp_ellint_2f(float __k)
00407   { return __detail::__comp_ellint_2<float>(__k); }
00408 
00409   /**
00410    * Return the complete elliptic integral of the second kind @f$ E(k) @f$
00411    * for long double modulus @c k.
00412    *
00413    * @see comp_ellint_2 for details.
00414    */
00415   inline long double
00416   comp_ellint_2l(long double __k)
00417   { return __detail::__comp_ellint_2<long double>(__k); }
00418 
00419   /**
00420    * Return the complete elliptic integral of the second kind @f$ E(k) @f$
00421    * for real modulus @c k.
00422    *
00423    * The complete elliptic integral of the second kind is defined as
00424    * @f[
00425    *   E(k) = E(k,\pi/2) = \int_0^{\pi/2}\sqrt{1 - k^2 sin^2\theta}
00426    * @f]
00427    * where @f$ E(k,\phi) @f$ is the incomplete elliptic integral of the
00428    * second kind and the modulus @f$ |k| <= 1 @f$.
00429    * @see ellint_2 for details of the incomplete elliptic function
00430    * of the second kind.
00431    *
00432    * @tparam _Tp The floating-point type of the modulus @c __k.
00433    * @param  __k  The modulus, @c abs(__k) <= 1
00434    * @throw std::domain_error if @c abs(__k) > 1.
00435    */
00436   template<typename _Tp>
00437     inline typename __gnu_cxx::__promote<_Tp>::__type
00438     comp_ellint_2(_Tp __k)
00439     {
00440       typedef typename __gnu_cxx::__promote<_Tp>::__type __type;
00441       return __detail::__comp_ellint_2<__type>(__k);
00442     }
00443 
00444   // Complete elliptic integrals of the third kind
00445 
00446   /**
00447    * @brief Return the complete elliptic integral of the third kind
00448    * @f$ \Pi(k,\nu) @f$ for @c float modulus @c k.
00449    *
00450    * @see comp_ellint_3 for details.
00451    */
00452   inline float
00453   comp_ellint_3f(float __k, float __nu)
00454   { return __detail::__comp_ellint_3<float>(__k, __nu); }
00455 
00456   /**
00457    * @brief Return the complete elliptic integral of the third kind
00458    * @f$ \Pi(k,\nu) @f$ for <tt>long double</tt> modulus @c k.
00459    *
00460    * @see comp_ellint_3 for details.
00461    */
00462   inline long double
00463   comp_ellint_3l(long double __k, long double __nu)
00464   { return __detail::__comp_ellint_3<long double>(__k, __nu); }
00465 
00466   /**
00467    * Return the complete elliptic integral of the third kind
00468    * @f$ \Pi(k,\nu) = \Pi(k,\nu,\pi/2) @f$ for real modulus @c k.
00469    *
00470    * The complete elliptic integral of the third kind is defined as
00471    * @f[
00472    *   \Pi(k,\nu) = \Pi(k,\nu,\pi/2) = \int_0^{\pi/2}
00473    *                 \frac{d\theta}
00474    *               {(1 - \nu \sin^2\theta)\sqrt{1 - k^2 \sin^2\theta}}
00475    * @f]
00476    * where @f$ \Pi(k,\nu,\phi) @f$ is the incomplete elliptic integral of the
00477    * second kind and the modulus @f$ |k| <= 1 @f$.
00478    * @see ellint_3 for details of the incomplete elliptic function
00479    * of the third kind.
00480    *
00481    * @tparam _Tp The floating-point type of the modulus @c __k.
00482    * @tparam _Tpn The floating-point type of the argument @c __nu.
00483    * @param  __k  The modulus, @c abs(__k) <= 1
00484    * @param  __nu  The argument
00485    * @throw std::domain_error if @c abs(__k) > 1.
00486    */
00487   template<typename _Tp, typename _Tpn>
00488     inline typename __gnu_cxx::__promote_2<_Tp, _Tpn>::__type
00489     comp_ellint_3(_Tp __k, _Tpn __nu)
00490     {
00491       typedef typename __gnu_cxx::__promote_2<_Tp, _Tpn>::__type __type;
00492       return __detail::__comp_ellint_3<__type>(__k, __nu);
00493     }
00494 
00495   // Regular modified cylindrical Bessel functions
00496 
00497   /**
00498    * Return the regular modified Bessel function @f$ I_{\nu}(x) @f$
00499    * for @c float order @f$ \nu @f$ and argument @f$ x >= 0 @f$.
00500    *
00501    * @see cyl_bessel_i for setails.
00502    */
00503   inline float
00504   cyl_bessel_if(float __nu, float __x)
00505   { return __detail::__cyl_bessel_i<float>(__nu, __x); }
00506 
00507   /**
00508    * Return the regular modified Bessel function @f$ I_{\nu}(x) @f$
00509    * for <tt>long double</tt> order @f$ \nu @f$ and argument @f$ x >= 0 @f$.
00510    *
00511    * @see cyl_bessel_i for setails.
00512    */
00513   inline long double
00514   cyl_bessel_il(long double __nu, long double __x)
00515   { return __detail::__cyl_bessel_i<long double>(__nu, __x); }
00516 
00517   /**
00518    * Return the regular modified Bessel function @f$ I_{\nu}(x) @f$
00519    * for real order @f$ \nu @f$ and argument @f$ x >= 0 @f$.
00520    *
00521    * The regular modified cylindrical Bessel function is:
00522    * @f[
00523    *  I_{\nu}(x) = i^{-\nu}J_\nu(ix) = \sum_{k=0}^{\infty}
00524    *            \frac{(x/2)^{\nu + 2k}}{k!\Gamma(\nu+k+1)}
00525    * @f]
00526    *
00527    * @tparam _Tpnu The floating-point type of the order @c __nu.
00528    * @tparam _Tp The floating-point type of the argument @c __x.
00529    * @param  __nu  The order
00530    * @param  __x   The argument, <tt> __x >= 0 </tt>
00531    * @throw std::domain_error if <tt> __x < 0 </tt>.
00532    */
00533   template<typename _Tpnu, typename _Tp>
00534     inline typename __gnu_cxx::__promote_2<_Tpnu, _Tp>::__type
00535     cyl_bessel_i(_Tpnu __nu, _Tp __x)
00536     {
00537       typedef typename __gnu_cxx::__promote_2<_Tpnu, _Tp>::__type __type;
00538       return __detail::__cyl_bessel_i<__type>(__nu, __x);
00539     }
00540 
00541   // Cylindrical Bessel functions (of the first kind)
00542 
00543   /**
00544    * Return the Bessel function of the first kind @f$ J_{\nu}(x) @f$
00545    * for @c float order @f$ \nu @f$ and argument @f$ x >= 0 @f$.
00546    *
00547    * @see cyl_bessel_j for setails.
00548    */
00549   inline float
00550   cyl_bessel_jf(float __nu, float __x)
00551   { return __detail::__cyl_bessel_j<float>(__nu, __x); }
00552 
00553   /**
00554    * Return the Bessel function of the first kind @f$ J_{\nu}(x) @f$
00555    * for <tt>long double</tt> order @f$ \nu @f$ and argument @f$ x >= 0 @f$.
00556    *
00557    * @see cyl_bessel_j for setails.
00558    */
00559   inline long double
00560   cyl_bessel_jl(long double __nu, long double __x)
00561   { return __detail::__cyl_bessel_j<long double>(__nu, __x); }
00562 
00563   /**
00564    * Return the Bessel function @f$ J_{\nu}(x) @f$ of real order @f$ \nu @f$
00565    * and argument @f$ x >= 0 @f$.
00566    *
00567    * The cylindrical Bessel function is:
00568    * @f[
00569    *    J_{\nu}(x) = \sum_{k=0}^{\infty}
00570    *              \frac{(-1)^k (x/2)^{\nu + 2k}}{k!\Gamma(\nu+k+1)}
00571    * @f]
00572    *
00573    * @tparam _Tpnu The floating-point type of the order @c __nu.
00574    * @tparam _Tp The floating-point type of the argument @c __x.
00575    * @param  __nu  The order
00576    * @param  __x   The argument, <tt> __x >= 0 </tt>
00577    * @throw std::domain_error if <tt> __x < 0 </tt>.
00578    */
00579   template<typename _Tpnu, typename _Tp>
00580     inline typename __gnu_cxx::__promote_2<_Tpnu, _Tp>::__type
00581     cyl_bessel_j(_Tpnu __nu, _Tp __x)
00582     {
00583       typedef typename __gnu_cxx::__promote_2<_Tpnu, _Tp>::__type __type;
00584       return __detail::__cyl_bessel_j<__type>(__nu, __x);
00585     }
00586 
00587   // Irregular modified cylindrical Bessel functions
00588 
00589   /**
00590    * Return the irregular modified Bessel function @f$ K_{\nu}(x) @f$
00591    * for @c float order @f$ \nu @f$ and argument @f$ x >= 0 @f$.
00592    *
00593    * @see cyl_bessel_k for setails.
00594    */
00595   inline float
00596   cyl_bessel_kf(float __nu, float __x)
00597   { return __detail::__cyl_bessel_k<float>(__nu, __x); }
00598 
00599   /**
00600    * Return the irregular modified Bessel function @f$ K_{\nu}(x) @f$
00601    * for <tt>long double</tt> order @f$ \nu @f$ and argument @f$ x >= 0 @f$.
00602    *
00603    * @see cyl_bessel_k for setails.
00604    */
00605   inline long double
00606   cyl_bessel_kl(long double __nu, long double __x)
00607   { return __detail::__cyl_bessel_k<long double>(__nu, __x); }
00608 
00609   /**
00610    * Return the irregular modified Bessel function @f$ K_{\nu}(x) @f$
00611    * of real order @f$ \nu @f$ and argument @f$ x @f$.
00612    *
00613    * The irregular modified Bessel function is defined by:
00614    * @f[
00615    *    K_{\nu}(x) = \frac{\pi}{2}
00616    *                 \frac{I_{-\nu}(x) - I_{\nu}(x)}{\sin \nu\pi}
00617    * @f]
00618    * where for integral @f$ \nu = n @f$ a limit is taken:
00619    * @f$ lim_{\nu \to n} @f$.
00620    * For negative argument we have simply:
00621    * @f[
00622    *    K_{-\nu}(x) = K_{\nu}(x)
00623    * @f]
00624    *
00625    * @tparam _Tpnu The floating-point type of the order @c __nu.
00626    * @tparam _Tp The floating-point type of the argument @c __x.
00627    * @param  __nu  The order
00628    * @param  __x   The argument, <tt> __x >= 0 </tt>
00629    * @throw std::domain_error if <tt> __x < 0 </tt>.
00630    */
00631   template<typename _Tpnu, typename _Tp>
00632     inline typename __gnu_cxx::__promote_2<_Tpnu, _Tp>::__type
00633     cyl_bessel_k(_Tpnu __nu, _Tp __x)
00634     {
00635       typedef typename __gnu_cxx::__promote_2<_Tpnu, _Tp>::__type __type;
00636       return __detail::__cyl_bessel_k<__type>(__nu, __x);
00637     }
00638 
00639   // Cylindrical Neumann functions
00640 
00641   /**
00642    * Return the Neumann function @f$ N_{\nu}(x) @f$
00643    * of @c float order @f$ \nu @f$ and argument @f$ x @f$.
00644    *
00645    * @see cyl_neumann for setails.
00646    */
00647   inline float
00648   cyl_neumannf(float __nu, float __x)
00649   { return __detail::__cyl_neumann_n<float>(__nu, __x); }
00650 
00651   /**
00652    * Return the Neumann function @f$ N_{\nu}(x) @f$
00653    * of <tt>long double</tt> order @f$ \nu @f$ and argument @f$ x @f$.
00654    *
00655    * @see cyl_neumann for setails.
00656    */
00657   inline long double
00658   cyl_neumannl(long double __nu, long double __x)
00659   { return __detail::__cyl_neumann_n<long double>(__nu, __x); }
00660 
00661   /**
00662    * Return the Neumann function @f$ N_{\nu}(x) @f$
00663    * of real order @f$ \nu @f$ and argument @f$ x >= 0 @f$.
00664    *
00665    * The Neumann function is defined by:
00666    * @f[
00667    *    N_{\nu}(x) = \frac{J_{\nu}(x) \cos \nu\pi - J_{-\nu}(x)}
00668    *                      {\sin \nu\pi}
00669    * @f]
00670    * where @f$ x >= 0 @f$ and for integral order @f$ \nu = n @f$
00671    * a limit is taken: @f$ lim_{\nu \to n} @f$.
00672    *
00673    * @tparam _Tpnu The floating-point type of the order @c __nu.
00674    * @tparam _Tp The floating-point type of the argument @c __x.
00675    * @param  __nu  The order
00676    * @param  __x   The argument, <tt> __x >= 0 </tt>
00677    * @throw std::domain_error if <tt> __x < 0 </tt>.
00678    */
00679   template<typename _Tpnu, typename _Tp>
00680     inline typename __gnu_cxx::__promote_2<_Tpnu, _Tp>::__type
00681     cyl_neumann(_Tpnu __nu, _Tp __x)
00682     {
00683       typedef typename __gnu_cxx::__promote_2<_Tpnu, _Tp>::__type __type;
00684       return __detail::__cyl_neumann_n<__type>(__nu, __x);
00685     }
00686 
00687   // Incomplete elliptic integrals of the first kind
00688 
00689   /**
00690    * Return the incomplete elliptic integral of the first kind @f$ E(k,\phi) @f$
00691    * for @c float modulus @f$ k @f$ and angle @f$ \phi @f$.
00692    *
00693    * @see ellint_1 for details.
00694    */
00695   inline float
00696   ellint_1f(float __k, float __phi)
00697   { return __detail::__ellint_1<float>(__k, __phi); }
00698 
00699   /**
00700    * Return the incomplete elliptic integral of the first kind @f$ E(k,\phi) @f$
00701    * for <tt>long double</tt> modulus @f$ k @f$ and angle @f$ \phi @f$.
00702    *
00703    * @see ellint_1 for details.
00704    */
00705   inline long double
00706   ellint_1l(long double __k, long double __phi)
00707   { return __detail::__ellint_1<long double>(__k, __phi); }
00708 
00709   /**
00710    * Return the incomplete elliptic integral of the first kind @f$ F(k,\phi) @f$
00711    * for @c real modulus @f$ k @f$ and angle @f$ \phi @f$.
00712    *
00713    * The incomplete elliptic integral of the first kind is defined as
00714    * @f[
00715    *   F(k,\phi) = \int_0^{\phi}\frac{d\theta}
00716    *                                 {\sqrt{1 - k^2 sin^2\theta}}
00717    * @f]
00718    * For  @f$ \phi= \pi/2 @f$ this becomes the complete elliptic integral of
00719    * the first kind, @f$ K(k) @f$.  @see comp_ellint_1.
00720    *
00721    * @tparam _Tp The floating-point type of the modulus @c __k.
00722    * @tparam _Tpp The floating-point type of the angle @c __phi.
00723    * @param  __k  The modulus, <tt> abs(__k) <= 1 </tt>
00724    * @param  __phi  The integral limit argument in radians
00725    * @throw std::domain_error if <tt> abs(__k) > 1 </tt>.
00726    */
00727   template<typename _Tp, typename _Tpp>
00728     inline typename __gnu_cxx::__promote_2<_Tp, _Tpp>::__type
00729     ellint_1(_Tp __k, _Tpp __phi)
00730     {
00731       typedef typename __gnu_cxx::__promote_2<_Tp, _Tpp>::__type __type;
00732       return __detail::__ellint_1<__type>(__k, __phi);
00733     }
00734 
00735   // Incomplete elliptic integrals of the second kind
00736 
00737   /**
00738    * @brief Return the incomplete elliptic integral of the second kind
00739    * @f$ E(k,\phi) @f$ for @c float argument.
00740    *
00741    * @see ellint_2 for details.
00742    */
00743   inline float
00744   ellint_2f(float __k, float __phi)
00745   { return __detail::__ellint_2<float>(__k, __phi); }
00746 
00747   /**
00748    * @brief Return the incomplete elliptic integral of the second kind
00749    * @f$ E(k,\phi) @f$.
00750    *
00751    * @see ellint_2 for details.
00752    */
00753   inline long double
00754   ellint_2l(long double __k, long double __phi)
00755   { return __detail::__ellint_2<long double>(__k, __phi); }
00756 
00757   /**
00758    * Return the incomplete elliptic integral of the second kind
00759    * @f$ E(k,\phi) @f$.
00760    *
00761    * The incomplete elliptic integral of the second kind is defined as
00762    * @f[
00763    *   E(k,\phi) = \int_0^{\phi} \sqrt{1 - k^2 sin^2\theta}
00764    * @f]
00765    * For  @f$ \phi= \pi/2 @f$ this becomes the complete elliptic integral of
00766    * the second kind, @f$ E(k) @f$.  @see comp_ellint_2.
00767    *
00768    * @tparam _Tp The floating-point type of the modulus @c __k.
00769    * @tparam _Tpp The floating-point type of the angle @c __phi.
00770    * @param  __k  The modulus, <tt> abs(__k) <= 1 </tt>
00771    * @param  __phi  The integral limit argument in radians
00772    * @return  The elliptic function of the second kind.
00773    * @throw std::domain_error if <tt> abs(__k) > 1 </tt>.
00774    */
00775   template<typename _Tp, typename _Tpp>
00776     inline typename __gnu_cxx::__promote_2<_Tp, _Tpp>::__type
00777     ellint_2(_Tp __k, _Tpp __phi)
00778     {
00779       typedef typename __gnu_cxx::__promote_2<_Tp, _Tpp>::__type __type;
00780       return __detail::__ellint_2<__type>(__k, __phi);
00781     }
00782 
00783   // Incomplete elliptic integrals of the third kind
00784 
00785   /**
00786    * @brief Return the incomplete elliptic integral of the third kind
00787    * @f$ \Pi(k,\nu,\phi) @f$ for @c float argument.
00788    *
00789    * @see ellint_3 for details.
00790    */
00791   inline float
00792   ellint_3f(float __k, float __nu, float __phi)
00793   { return __detail::__ellint_3<float>(__k, __nu, __phi); }
00794 
00795   /**
00796    * @brief Return the incomplete elliptic integral of the third kind
00797    * @f$ \Pi(k,\nu,\phi) @f$.
00798    *
00799    * @see ellint_3 for details.
00800    */
00801   inline long double
00802   ellint_3l(long double __k, long double __nu, long double __phi)
00803   { return __detail::__ellint_3<long double>(__k, __nu, __phi); }
00804 
00805   /**
00806    * @brief Return the incomplete elliptic integral of the third kind
00807    * @f$ \Pi(k,\nu,\phi) @f$.
00808    *
00809    * The incomplete elliptic integral of the third kind is defined by:
00810    * @f[
00811    *   \Pi(k,\nu,\phi) = \int_0^{\phi}
00812    *                     \frac{d\theta}
00813    *                     {(1 - \nu \sin^2\theta)
00814    *                      \sqrt{1 - k^2 \sin^2\theta}}
00815    * @f]
00816    * For  @f$ \phi= \pi/2 @f$ this becomes the complete elliptic integral of
00817    * the third kind, @f$ \Pi(k,\nu) @f$.  @see comp_ellint_3.
00818    *
00819    * @tparam _Tp The floating-point type of the modulus @c __k.
00820    * @tparam _Tpn The floating-point type of the argument @c __nu.
00821    * @tparam _Tpp The floating-point type of the angle @c __phi.
00822    * @param  __k  The modulus, <tt> abs(__k) <= 1 </tt>
00823    * @param  __nu  The second argument
00824    * @param  __phi  The integral limit argument in radians
00825    * @return  The elliptic function of the third kind.
00826    * @throw std::domain_error if <tt> abs(__k) > 1 </tt>.
00827    */
00828   template<typename _Tp, typename _Tpn, typename _Tpp>
00829     inline typename __gnu_cxx::__promote_3<_Tp, _Tpn, _Tpp>::__type
00830     ellint_3(_Tp __k, _Tpn __nu, _Tpp __phi)
00831     {
00832       typedef typename __gnu_cxx::__promote_3<_Tp, _Tpn, _Tpp>::__type __type;
00833       return __detail::__ellint_3<__type>(__k, __nu, __phi);
00834     }
00835 
00836   // Exponential integrals
00837 
00838   /**
00839    * Return the exponential integral @f$ Ei(x) @f$ for @c float argument @c x.
00840    *
00841    * @see expint for details.
00842    */
00843   inline float
00844   expintf(float __x)
00845   { return __detail::__expint<float>(__x); }
00846 
00847   /**
00848    * Return the exponential integral @f$ Ei(x) @f$
00849    * for <tt>long double</tt> argument @c x.
00850    *
00851    * @see expint for details.
00852    */
00853   inline long double
00854   expintl(long double __x)
00855   { return __detail::__expint<long double>(__x); }
00856 
00857   /**
00858    * Return the exponential integral @f$ Ei(x) @f$ for @c real argument @c x.
00859    *
00860    * The exponential integral is given by
00861    * \f[
00862    *   Ei(x) = -\int_{-x}^\infty \frac{e^t}{t} dt
00863    * \f]
00864    *
00865    * @tparam _Tp The floating-point type of the argument @c __x.
00866    * @param  __x  The argument of the exponential integral function.
00867    */
00868   template<typename _Tp>
00869     inline typename __gnu_cxx::__promote<_Tp>::__type
00870     expint(_Tp __x)
00871     {
00872       typedef typename __gnu_cxx::__promote<_Tp>::__type __type;
00873       return __detail::__expint<__type>(__x);
00874     }
00875 
00876   // Hermite polynomials
00877 
00878   /**
00879    * Return the Hermite polynomial @f$ H_n(x) @f$ of nonnegative order n
00880    * and float argument @c x.
00881    *
00882    * @see hermite for details.
00883    */
00884   inline float
00885   hermitef(unsigned int __n, float __x)
00886   { return __detail::__poly_hermite<float>(__n, __x); }
00887 
00888   /**
00889    * Return the Hermite polynomial @f$ H_n(x) @f$ of nonnegative order n
00890    * and <tt>long double</tt> argument @c x.
00891    *
00892    * @see hermite for details.
00893    */
00894   inline long double
00895   hermitel(unsigned int __n, long double __x)
00896   { return __detail::__poly_hermite<long double>(__n, __x); }
00897 
00898   /**
00899    * Return the Hermite polynomial @f$ H_n(x) @f$ of order n
00900    * and @c real argument @c x.
00901    *
00902    * The Hermite polynomial is defined by:
00903    * @f[
00904    *   H_n(x) = (-1)^n e^{x^2} \frac{d^n}{dx^n} e^{-x^2}
00905    * @f]
00906    *
00907    * The Hermite polynomial obeys a reflection formula:
00908    * @f[
00909    *   H_n(-x) = (-1)^n H_n(x)
00910    * @f]
00911    *
00912    * @tparam _Tp The floating-point type of the argument @c __x.
00913    * @param __n The order
00914    * @param __x The argument
00915    */
00916   template<typename _Tp>
00917     inline typename __gnu_cxx::__promote<_Tp>::__type
00918     hermite(unsigned int __n, _Tp __x)
00919     {
00920       typedef typename __gnu_cxx::__promote<_Tp>::__type __type;
00921       return __detail::__poly_hermite<__type>(__n, __x);
00922     }
00923 
00924   // Laguerre polynomials
00925 
00926   /**
00927    * Returns the Laguerre polynomial @f$ L_n(x) @f$ of nonnegative degree @c n
00928    * and @c float argument  @f$ x >= 0 @f$.
00929    *
00930    * @see laguerre for more details.
00931    */
00932   inline float
00933   laguerref(unsigned int __n, float __x)
00934   { return __detail::__laguerre<float>(__n, __x); }
00935 
00936   /**
00937    * Returns the Laguerre polynomial @f$ L_n(x) @f$ of nonnegative degree @c n
00938    * and <tt>long double</tt> argument @f$ x >= 0 @f$.
00939    *
00940    * @see laguerre for more details.
00941    */
00942   inline long double
00943   laguerrel(unsigned int __n, long double __x)
00944   { return __detail::__laguerre<long double>(__n, __x); }
00945 
00946   /**
00947    * Returns the Laguerre polynomial @f$ L_n(x) @f$
00948    * of nonnegative degree @c n and real argument @f$ x >= 0 @f$.
00949    *
00950    * The Laguerre polynomial is defined by:
00951    * @f[
00952    *     L_n(x) = \frac{e^x}{n!} \frac{d^n}{dx^n} (x^ne^{-x})
00953    * @f]
00954    *
00955    * @tparam _Tp The floating-point type of the argument @c __x.
00956    * @param __n The nonnegative order
00957    * @param __x The argument <tt> __x >= 0 </tt>
00958    * @throw std::domain_error if <tt> __x < 0 </tt>.
00959    */
00960   template<typename _Tp>
00961     inline typename __gnu_cxx::__promote<_Tp>::__type
00962     laguerre(unsigned int __n, _Tp __x)
00963     {
00964       typedef typename __gnu_cxx::__promote<_Tp>::__type __type;
00965       return __detail::__laguerre<__type>(__n, __x);
00966     }
00967 
00968   // Legendre polynomials
00969 
00970   /**
00971    * Return the Legendre polynomial @f$ P_l(x) @f$ of nonnegative
00972    * degree @f$ l @f$ and @c float argument @f$ |x| <= 0 @f$.
00973    *
00974    * @see legendre for more details.
00975    */
00976   inline float
00977   legendref(unsigned int __l, float __x)
00978   { return __detail::__poly_legendre_p<float>(__l, __x); }
00979 
00980   /**
00981    * Return the Legendre polynomial @f$ P_l(x) @f$ of nonnegative
00982    * degree @f$ l @f$ and <tt>long double</tt> argument @f$ |x| <= 0 @f$.
00983    *
00984    * @see legendre for more details.
00985    */
00986   inline long double
00987   legendrel(unsigned int __l, long double __x)
00988   { return __detail::__poly_legendre_p<long double>(__l, __x); }
00989 
00990   /**
00991    * Return the Legendre polynomial @f$ P_l(x) @f$ of nonnegative
00992    * degree @f$ l @f$ and real argument @f$ |x| <= 0 @f$.
00993    *
00994    * The Legendre function of order @f$ l @f$ and argument @f$ x @f$,
00995    * @f$ P_l(x) @f$, is defined by:
00996    * @f[
00997    *   P_l(x) = \frac{1}{2^l l!}\frac{d^l}{dx^l}(x^2 - 1)^{l}
00998    * @f]
00999    *
01000    * @tparam _Tp The floating-point type of the argument @c __x.
01001    * @param __l The degree @f$ l >= 0 @f$
01002    * @param __x The argument @c abs(__x) <= 1
01003    * @throw std::domain_error if @c abs(__x) > 1
01004    */
01005   template<typename _Tp>
01006     inline typename __gnu_cxx::__promote<_Tp>::__type
01007     legendre(unsigned int __l, _Tp __x)
01008     {
01009       typedef typename __gnu_cxx::__promote<_Tp>::__type __type;
01010       return __detail::__poly_legendre_p<__type>(__l, __x);
01011     }
01012 
01013   // Riemann zeta functions
01014 
01015   /**
01016    * Return the Riemann zeta function @f$ \zeta(s) @f$
01017    * for @c float argument @f$ s @f$.
01018    *
01019    * @see riemann_zeta for more details.
01020    */
01021   inline float
01022   riemann_zetaf(float __s)
01023   { return __detail::__riemann_zeta<float>(__s); }
01024 
01025   /**
01026    * Return the Riemann zeta function @f$ \zeta(s) @f$
01027    * for <tt>long double</tt> argument @f$ s @f$.
01028    *
01029    * @see riemann_zeta for more details.
01030    */
01031   inline long double
01032   riemann_zetal(long double __s)
01033   { return __detail::__riemann_zeta<long double>(__s); }
01034 
01035   /**
01036    * Return the Riemann zeta function @f$ \zeta(s) @f$
01037    * for real argument @f$ s @f$.
01038    *
01039    * The Riemann zeta function is defined by:
01040    * @f[
01041    *    \zeta(s) = \sum_{k=1}^{\infty} k^{-s} \hbox{ for } s > 1
01042    * @f]
01043    * and
01044    * @f[
01045    *    \zeta(s) = \frac{1}{1-2^{1-s}}\sum_{k=1}^{\infty}(-1)^{k-1}k^{-s}
01046    *              \hbox{ for } 0 <= s <= 1
01047    * @f]
01048    * For s < 1 use the reflection formula:
01049    * @f[
01050    *    \zeta(s) = 2^s \pi^{s-1} \sin(\frac{\pi s}{2}) \Gamma(1-s) \zeta(1-s)
01051    * @f]
01052    *
01053    * @tparam _Tp The floating-point type of the argument @c __s.
01054    * @param __s The argument <tt> s != 1 </tt>
01055    */
01056   template<typename _Tp>
01057     inline typename __gnu_cxx::__promote<_Tp>::__type
01058     riemann_zeta(_Tp __s)
01059     {
01060       typedef typename __gnu_cxx::__promote<_Tp>::__type __type;
01061       return __detail::__riemann_zeta<__type>(__s);
01062     }
01063 
01064   // Spherical Bessel functions
01065 
01066   /**
01067    * Return the spherical Bessel function @f$ j_n(x) @f$ of nonnegative order n
01068    * and @c float argument @f$ x >= 0 @f$.
01069    *
01070    * @see sph_bessel for more details.
01071    */
01072   inline float
01073   sph_besself(unsigned int __n, float __x)
01074   { return __detail::__sph_bessel<float>(__n, __x); }
01075 
01076   /**
01077    * Return the spherical Bessel function @f$ j_n(x) @f$ of nonnegative order n
01078    * and <tt>long double</tt> argument @f$ x >= 0 @f$.
01079    *
01080    * @see sph_bessel for more details.
01081    */
01082   inline long double
01083   sph_bessell(unsigned int __n, long double __x)
01084   { return __detail::__sph_bessel<long double>(__n, __x); }
01085 
01086   /**
01087    * Return the spherical Bessel function @f$ j_n(x) @f$ of nonnegative order n
01088    * and real argument @f$ x >= 0 @f$.
01089    *
01090    * The spherical Bessel function is defined by:
01091    * @f[
01092    *  j_n(x) = \left(\frac{\pi}{2x} \right) ^{1/2} J_{n+1/2}(x)
01093    * @f]
01094    *
01095    * @tparam _Tp The floating-point type of the argument @c __x.
01096    * @param  __n  The integral order <tt> n >= 0 </tt>
01097    * @param  __x  The real argument <tt> x >= 0 </tt>
01098    * @throw std::domain_error if <tt> __x < 0 </tt>.
01099    */
01100   template<typename _Tp>
01101     inline typename __gnu_cxx::__promote<_Tp>::__type
01102     sph_bessel(unsigned int __n, _Tp __x)
01103     {
01104       typedef typename __gnu_cxx::__promote<_Tp>::__type __type;
01105       return __detail::__sph_bessel<__type>(__n, __x);
01106     }
01107 
01108   // Spherical associated Legendre functions
01109 
01110   /**
01111    * Return the spherical Legendre function of nonnegative integral
01112    * degree @c l and order @c m and float angle @f$ \theta @f$ in radians.
01113    *
01114    * @see sph_legendre for details.
01115    */
01116   inline float
01117   sph_legendref(unsigned int __l, unsigned int __m, float __theta)
01118   { return __detail::__sph_legendre<float>(__l, __m, __theta); }
01119 
01120   /**
01121    * Return the spherical Legendre function of nonnegative integral
01122    * degree @c l and order @c m and <tt>long double</tt> angle @f$ \theta @f$
01123    * in radians.
01124    *
01125    * @see sph_legendre for details.
01126    */
01127   inline long double
01128   sph_legendrel(unsigned int __l, unsigned int __m, long double __theta)
01129   { return __detail::__sph_legendre<long double>(__l, __m, __theta); }
01130 
01131   /**
01132    * Return the spherical Legendre function of nonnegative integral
01133    * degree @c l and order @c m and real angle @f$ \theta @f$ in radians.
01134    *
01135    * The spherical Legendre function is defined by
01136    * @f[
01137    *  Y_l^m(\theta,\phi) = (-1)^m[\frac{(2l+1)}{4\pi}
01138    *                              \frac{(l-m)!}{(l+m)!}]
01139    *                   P_l^m(\cos\theta) \exp^{im\phi}
01140    * @f]
01141    *
01142    * @tparam _Tp The floating-point type of the angle @c __theta.
01143    * @param __l The order <tt> __l >= 0 </tt>
01144    * @param __m The degree <tt> __m >= 0 </tt> and <tt> __m <= __l </tt>
01145    * @param __theta The radian polar angle argument
01146    */
01147   template<typename _Tp>
01148     inline typename __gnu_cxx::__promote<_Tp>::__type
01149     sph_legendre(unsigned int __l, unsigned int __m, _Tp __theta)
01150     {
01151       typedef typename __gnu_cxx::__promote<_Tp>::__type __type;
01152       return __detail::__sph_legendre<__type>(__l, __m, __theta);
01153     }
01154 
01155   // Spherical Neumann functions
01156 
01157   /**
01158    * Return the spherical Neumann function of integral order @f$ n >= 0 @f$
01159    * and @c float argument @f$ x >= 0 @f$.
01160    *
01161    * @see sph_neumann for details.
01162    */
01163   inline float
01164   sph_neumannf(unsigned int __n, float __x)
01165   { return __detail::__sph_neumann<float>(__n, __x); }
01166 
01167   /**
01168    * Return the spherical Neumann function of integral order @f$ n >= 0 @f$
01169    * and <tt>long double</tt> @f$ x >= 0 @f$.
01170    *
01171    * @see sph_neumann for details.
01172    */
01173   inline long double
01174   sph_neumannl(unsigned int __n, long double __x)
01175   { return __detail::__sph_neumann<long double>(__n, __x); }
01176 
01177   /**
01178    * Return the spherical Neumann function of integral order @f$ n >= 0 @f$
01179    * and real argument @f$ x >= 0 @f$.
01180    *
01181    * The spherical Neumann function is defined by
01182    * @f[
01183    *    n_n(x) = \left(\frac{\pi}{2x} \right) ^{1/2} N_{n+1/2}(x)
01184    * @f]
01185    *
01186    * @tparam _Tp The floating-point type of the argument @c __x.
01187    * @param  __n  The integral order <tt> n >= 0 </tt>
01188    * @param  __x  The real argument <tt> __x >= 0 </tt>
01189    * @throw std::domain_error if <tt> __x < 0 </tt>.
01190    */
01191   template<typename _Tp>
01192     inline typename __gnu_cxx::__promote<_Tp>::__type
01193     sph_neumann(unsigned int __n, _Tp __x)
01194     {
01195       typedef typename __gnu_cxx::__promote<_Tp>::__type __type;
01196       return __detail::__sph_neumann<__type>(__n, __x);
01197     }
01198 
01199   // @} group mathsf
01200 
01201 _GLIBCXX_END_NAMESPACE_VERSION
01202 } // namespace std
01203 
01204 namespace __gnu_cxx _GLIBCXX_VISIBILITY(default)
01205 {
01206 
01207   // Confluent hypergeometric functions
01208 
01209   /**
01210    * Return the confluent hypergeometric function @f$ {}_1F_1(a;c;x) @f$
01211    * of @c float numeratorial parameter @c a, denominatorial parameter @c c,
01212    * and argument @c x.
01213    *
01214    * @see conf_hyperg for details.
01215    */
01216   inline float
01217   conf_hypergf(float __a, float __c, float __x)
01218   { return std::__detail::__conf_hyperg<float>(__a, __c, __x); }
01219 
01220   /**
01221    * Return the confluent hypergeometric function @f$ {}_1F_1(a;c;x) @f$
01222    * of <tt>long double</tt> numeratorial parameter @c a,
01223    * denominatorial parameter @c c, and argument @c x.
01224    *
01225    * @see conf_hyperg for details.
01226    */
01227   inline long double
01228   conf_hypergl(long double __a, long double __c, long double __x)
01229   { return std::__detail::__conf_hyperg<long double>(__a, __c, __x); }
01230 
01231   /**
01232    * Return the confluent hypergeometric function @f$ {}_1F_1(a;c;x) @f$
01233    * of real numeratorial parameter @c a, denominatorial parameter @c c,
01234    * and argument @c x.
01235    *
01236    * The confluent hypergeometric function is defined by
01237    * @f[
01238    *    {}_1F_1(a;c;x) = \sum_{n=0}^{\infty} \frac{(a)_n x^n}{(c)_n n!}
01239    * @f]
01240    * where the Pochhammer symbol is @f$ (x)_k = (x)(x+1)...(x+k-1) @f$,
01241    * @f$ (x)_0 = 1 @f$
01242    *
01243    * @param __a The numeratorial parameter
01244    * @param __c The denominatorial parameter
01245    * @param __x The argument
01246    */
01247   template<typename _Tpa, typename _Tpc, typename _Tp>
01248     inline typename __gnu_cxx::__promote_3<_Tpa, _Tpc, _Tp>::__type
01249     conf_hyperg(_Tpa __a, _Tpc __c, _Tp __x)
01250     {
01251       typedef typename __gnu_cxx::__promote_3<_Tpa, _Tpc, _Tp>::__type __type;
01252       return std::__detail::__conf_hyperg<__type>(__a, __c, __x);
01253     }
01254 
01255   // Hypergeometric functions
01256 
01257   /**
01258    * Return the hypergeometric function @f$ {}_2F_1(a,b;c;x) @f$
01259    * of @ float numeratorial parameters @c a and @c b,
01260    * denominatorial parameter @c c, and argument @c x.
01261    *
01262    * @see hyperg for details.
01263    */
01264   inline float
01265   hypergf(float __a, float __b, float __c, float __x)
01266   { return std::__detail::__hyperg<float>(__a, __b, __c, __x); }
01267 
01268   /**
01269    * Return the hypergeometric function @f$ {}_2F_1(a,b;c;x) @f$
01270    * of <tt>long double</tt> numeratorial parameters @c a and @c b,
01271    * denominatorial parameter @c c, and argument @c x.
01272    *
01273    * @see hyperg for details.
01274    */
01275   inline long double
01276   hypergl(long double __a, long double __b, long double __c, long double __x)
01277   { return std::__detail::__hyperg<long double>(__a, __b, __c, __x); }
01278 
01279   /**
01280    * Return the hypergeometric function @f$ {}_2F_1(a,b;c;x) @f$
01281    * of real numeratorial parameters @c a and @c b,
01282    * denominatorial parameter @c c, and argument @c x.
01283    *
01284    * The hypergeometric function is defined by
01285    * @f[
01286    *    {}_2F_1(a;c;x) = \sum_{n=0}^{\infty} \frac{(a)_n (b)_n x^n}{(c)_n n!}
01287    * @f]
01288    * where the Pochhammer symbol is @f$ (x)_k = (x)(x+1)...(x+k-1) @f$,
01289    * @f$ (x)_0 = 1 @f$
01290    *
01291    * @param __a The first numeratorial parameter
01292    * @param __b The second numeratorial parameter
01293    * @param __c The denominatorial parameter
01294    * @param __x The argument
01295    */
01296   template<typename _Tpa, typename _Tpb, typename _Tpc, typename _Tp>
01297     inline typename __gnu_cxx::__promote_4<_Tpa, _Tpb, _Tpc, _Tp>::__type
01298     hyperg(_Tpa __a, _Tpb __b, _Tpc __c, _Tp __x)
01299     {
01300       typedef typename __gnu_cxx::__promote_4<_Tpa, _Tpb, _Tpc, _Tp>
01301                 ::__type __type;
01302       return std::__detail::__hyperg<__type>(__a, __b, __c, __x);
01303     }
01304 
01305 } // namespace __gnu_cxx
01306 
01307 #pragma GCC visibility pop
01308 
01309 #endif // _GLIBCXX_BITS_SPECFUN_H