The C and C++ Include Header Files

/usr/include/c++/11/bits/specfun.h


$ cat -n /usr/include/c++/11/bits/specfun.h

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