Go to the first, previous, next, last section, table of contents.


Special Functions

This chapter describes the GSL special function library.

The gsl_sf_result struct

Almost all the special functions calculate an error estimate along with the value of the result. Therefore, structures are provided for amalgating a value and error estimate.

The gsl_sf_result struct contains value and error fields.

typedef struct
{
  double val;
  double err;
} gsl_sf_result;

In some cases, an overflow or underflow can be detected and handled by a function. In this case, it may be possible to return a scaling exponent as well as an error/value pair. The gsl_sf_result struct contains value and error fields as well as an exponent field such that the actual result is obtained as result * 10^(e10).

typedef struct
{
  double val;
  double err;
  int    e10;
} gsl_sf_result_e10;

Airy Functions

The Airy functions Ai(x) and Bi(x) are defined by the integral representations

Ai(x) = {1 / \pi} \int_0^\infty \cos({1/3} t^3 + xt ) dt,
Bi(x) = {1 / \pi} \int_0^\infty (e^{-t^3/3} + \sin({1/3} t^3 + xt) dt.

Airy Functions

Function: int gsl_sf_airy_Ai_impl (double x, gsl_mode_t mode, gsl_sf_result * result)
Function: int gsl_sf_airy_Ai_e (double x, gsl_mode_t mode, gsl_sf_result * result)

Function: int gsl_sf_airy_Bi_impl (double x, gsl_mode_t mode, gsl_sf_result * result)
Function: int gsl_sf_airy_Bi_e (double x, gsl_mode_t mode, gsl_sf_result * result)

Function: int gsl_sf_airy_Ai_scaled_impl (double x, gsl_mode_t mode, gsl_sf_result * result)
Function: int gsl_sf_airy_Ai_scaled_e (double x, gsl_mode_t mode, gsl_sf_result * result)

Function: int gsl_sf_airy_Bi_scaled_impl (double x, gsl_mode_t mode, gsl_sf_result * result)
Function: int gsl_sf_airy_Bi_scaled_e (double x, gsl_mode_t mode, gsl_sf_result * result)

Derivatives of Airy Functions

Function: int gsl_sf_airy_Ai_deriv_impl (double x, gsl_mode_t mode, gsl_sf_result * result)
Function: int gsl_sf_airy_Ai_deriv_e (double x, gsl_mode_t mode, gsl_sf_result * result)

Function: int gsl_sf_airy_Bi_deriv_impl (double x, gsl_mode_t mode, gsl_sf_result * result)
Function: int gsl_sf_airy_Bi_deriv_e (double x, gsl_mode_t mode, gsl_sf_result * result)

Function: int gsl_sf_airy_Ai_deriv_scaled_impl (double x, gsl_mode_t mode, gsl_sf_result * result)
Function: int gsl_sf_airy_Ai_deriv_scaled_e (double x, gsl_mode_t mode, gsl_sf_result * result)

Function: int gsl_sf_airy_Bi_deriv_scaled_e (double x, gsl_mode_t mode, gsl_sf_result * result)
Function: int gsl_sf_airy_Bi_deriv_scaled_e (double x, gsl_mode_t mode, gsl_sf_result * result)

Zeros of Airy Functions

Function: int gsl_sf_airy_zero_Ai_e (int s, gsl_sf_result * result)
Function: int gsl_sf_airy_zero_Ai_e (int s, gsl_sf_result * result)

Function: int gsl_sf_airy_zero_Bi_e (int s, gsl_sf_result * result)
Function: int gsl_sf_airy_zero_Bi_e (int s, gsl_sf_result * result)

Zeros of Derivatives of Airy Functions

Function: int gsl_sf_airy_zero_Ai_deriv_e (int s, gsl_sf_result * result)
Function: int gsl_sf_airy_zero_Ai_deriv_e (int s, gsl_sf_result * result)

Function: int gsl_sf_airy_zero_Bi_deriv_e (int s, gsl_sf_result * result)
Function: int gsl_sf_airy_zero_Bi_deriv_e (int s, gsl_sf_result * result)

Bessel Functions

Regular Cylindrical Bessel Functions

Function: int gsl_sf_bessel_J0_impl (double x, gsl_sf_result * result)
Function: int gsl_sf_bessel_J0_e (double x, gsl_sf_result * result)
Exceptional Return Values: none

Function: int gsl_sf_bessel_J1_impl (double x, gsl_sf_result * result)
Function: int gsl_sf_bessel_J1_e (double x, gsl_sf_result * result)
Exceptional Return Values: GSL_EUNDRFLW

Function: int gsl_sf_bessel_Jn_impl (int n, double x, gsl_sf_result * result)
Function: int gsl_sf_bessel_Jn_e (int n, double x, gsl_sf_result * result)
Exceptional Return Values: GSL_EUNDRFLW

Function: int gsl_sf_bessel_Jn_array_impl (int nmin, int nmax, double x, double * result_array)
Function: int gsl_sf_bessel_Jn_array_e (int nmin, int nmax, double x, double * result_array)
Conditions: nmin <= n <= nmax Exceptional Return Values: GSL_EDOM, GSL_EUNDRFLW

Irregular Cylindrical Bessel Functions

Function: int gsl_sf_bessel_Y0_impl (double x, gsl_sf_result * result)
Function: int gsl_sf_bessel_Y0_e (double x, gsl_sf_result * result)
Domain: x > 0.0 Exceptional Return Values: GSL_EDOM, GSL_EUNDRFLW

Function: int gsl_sf_bessel_Y1_impl (double x, gsl_sf_result * result)
Function: int gsl_sf_bessel_Y1_e (double x, gsl_sf_result * result)
Domain: x > 0.0 Exceptional Return Values: GSL_EDOM, GSL_EOVRFLW, GSL_EUNDRFLW

Function: int gsl_sf_bessel_Yn_impl (int n,double x, gsl_sf_result * result)
Function: int gsl_sf_bessel_Yn_e (int n,double x, gsl_sf_result * result)
Domain: x > 0.0 Exceptional Return Values: GSL_EDOM, GSL_EOVRFLW, GSL_EUNDRFLW

Function: int gsl_sf_bessel_Yn_array_impl (int nmin, int nmax, double x, double * result_array)
Function: int gsl_sf_bessel_Yn_array_e (int nmin, int nmax, double x, double * result_array)
Domain: x > 0.0 Conditions: nmin <= n <= nmax Exceptional Return Values: GSL_EDOM, GSL_EOVRFLW, GSL_EUNDRFLW

Regular Modified Cylindrical Bessel Functions

Function: int gsl_sf_bessel_I0_impl (double x, gsl_sf_result * result)
Function: int gsl_sf_bessel_I0_e (double x, gsl_sf_result * result)
Exceptional Return Values: GSL_EOVRFLW

Function: int gsl_sf_bessel_I1_impl (double x, gsl_sf_result * result)
Function: int gsl_sf_bessel_I1_e (double x, gsl_sf_result * result)
Exceptional Return Values: GSL_EOVRFLW, GSL_EUNDRFLW

Function: int gsl_sf_bessel_In_impl (int n, double x, gsl_sf_result * result)
Function: int gsl_sf_bessel_In_e (int n, double x, gsl_sf_result * result)
Exceptional Return Values: GSL_EOVRFLW, GSL_EUNDRFLW

Function: int gsl_sf_bessel_In_array_impl (int nmin, int nmax, double x, double * result_array)
Function: int gsl_sf_bessel_In_array_e (int nmin, int nmax, double x, double * result_array)
Domain: nmin >=0, nmax >= nmin Conditions: n=nmin,...,nmax, nmin >=0, nmax >= nmin Exceptional Return Values: GSL_EDOM, GSL_EOVRFLW, GSL_EUNDRFLW

Function: int gsl_sf_bessel_I0_scaled_impl (double x, gsl_sf_result * result)
Function: int gsl_sf_bessel_I0_scaled_e (double x, gsl_sf_result * result)
\exp(-|x|) I_0(x) Exceptional Return Values: none

Function: int gsl_sf_bessel_I1_scaled_impl (double x, gsl_sf_result * result)
Function: int gsl_sf_bessel_I1_scaled_e (double x, gsl_sf_result * result)
\exp(-|x|) I_1(x) Exceptional Return Values: GSL_EUNDRFLW

Function: int gsl_sf_bessel_In_scaled_impl (int n, double x, gsl_sf_result * result)
Function: int gsl_sf_bessel_In_scaled_e (int n, double x, gsl_sf_result * result)
\exp(-|x|) I_n(x) Exceptional Return Values: GSL_EUNDRFLW

Function: int gsl_sf_bessel_In_scaled_array_impl (int nmin, int nmax, double x, double * result_array)
Function: int gsl_sf_bessel_In_scaled_array_e (int nmin, int nmax, double x, double * result_array)
\exp(-|x|) I_n(x) Domain: nmin >=0, nmax >= nmin Conditions: n=nmin,...,nmax Exceptional Return Values: GSL_EUNDRFLW

Irregular Modified Cylindrical Bessel Functions

Function: int gsl_sf_bessel_K0_impl (double x, gsl_sf_result * result)
Function: int gsl_sf_bessel_K0_e (double x, gsl_sf_result * result)
Domain: x > 0.0 Exceptional Return Values: GSL_EDOM, GSL_EUNDRFLW

Function: int gsl_sf_bessel_K1_impl (double x, gsl_sf_result * result)
Function: int gsl_sf_bessel_K1_e (double x, gsl_sf_result * result)
Domain: x > 0.0 Exceptional Return Values: GSL_EDOM, GSL_EOVRFLW, GSL_EUNDRFLW

Function: int gsl_sf_bessel_Kn_impl (int n, double x, gsl_sf_result * result)
Function: int gsl_sf_bessel_Kn_e (int n, double x, gsl_sf_result * result)
Domain: x > 0.0 Exceptional Return Values: GSL_EDOM, GSL_EOVRFLW, GSL_EUNDRFLW

Function: int gsl_sf_bessel_Kn_array_impl (int nmin, int nmax, double x, double * result_array)
Function: int gsl_sf_bessel_Kn_array_e (int nmin, int nmax, double x, double * result_array)
Conditions: n=nmin,...,nmax Domain: x > 0.0, nmin >=0, nmax >= nmin Exceptional Return Values: GSL_EDOM, GSL_EOVRFLW, GSL_EUNDRFLW

Function: int gsl_sf_bessel_K0_scaled_impl (double x, gsl_sf_result * result)
Function: int gsl_sf_bessel_K0_scaled_e (double x, gsl_sf_result * result)
\exp(x) K_0(x) Domain: x > 0.0 Exceptional Return Values: GSL_EDOM

Function: int gsl_sf_bessel_K1_scaled_impl (double x, gsl_sf_result * result)
Function: int gsl_sf_bessel_K1_scaled_e (double x, gsl_sf_result * result)
\exp(x) K_1(x) Domain: x > 0.0 Exceptional Return Values: GSL_EDOM, GSL_EUNDRFLW

Function: int gsl_sf_bessel_Kn_scaled_impl (int n, double x, gsl_sf_result * result)
Function: int gsl_sf_bessel_Kn_scaled_e (int n, double x, gsl_sf_result * result)
\exp(x) K_n(x) Domain: x > 0.0 Exceptional Return Values: GSL_EDOM, GSL_EUNDRFLW

Function: int gsl_sf_bessel_Kn_scaled_array_impl (int nmin, int nmax, double x, double * result_array)
Function: int gsl_sf_bessel_Kn_scaled_array_e (int nmin, int nmax, double x, double * result_array)
\exp(x) K_n(x) Domain: x > 0.0, nmin >=0, nmax >= nmin Conditions: n=nmin,...,nmax Exceptional Return Values: GSL_EDOM, GSL_EUNDRFLW

Regular Spherical Bessel Functions

Function: int gsl_sf_bessel_j0_impl (double x, gsl_sf_result * result)
Function: int gsl_sf_bessel_j0_e (double x, gsl_sf_result * result)
j_0(x) = \sin(x)/x Exceptional Return Values: none

Function: int gsl_sf_bessel_j1_impl (double x, gsl_sf_result * result)
Function: int gsl_sf_bessel_j1_e (double x, gsl_sf_result * result)
j_1(x) = (\sin(x)/x - \cos(x))/x Exceptional Return Values: GSL_EUNDRFLW

Function: int gsl_sf_bessel_j2_impl (double x, gsl_sf_result * result)
Function: int gsl_sf_bessel_j2_e (double x, gsl_sf_result * result)
j_2(x) = ((3/x^2 - 1)\sin(x) - 3\cos(x)/x)/x Exceptional Return Values: GSL_EUNDRFLW

Function: int gsl_sf_bessel_jl_impl (int l, double x, gsl_sf_result * result)
Function: int gsl_sf_bessel_jl_e (int l, double x, gsl_sf_result * result)
Domain: l >= 0, x >= 0.0 Exceptional Return Values: GSL_EDOM, GSL_EUNDRFLW

Function: int gsl_sf_bessel_jl_array_impl (int lmax, double x, double * result_array)
Function: int gsl_sf_bessel_jl_array_e (int lmax, double x, double * result_array)
Domain: lmax >= 0 Conditions: l=0,1,...,lmax Exceptional Return Values: GSL_EDOM, GSL_EUNDRFLW

Function: int gsl_sf_bessel_jl_steed_array_impl (int lmax, double x, double * jl_x_array)
Uses Steed's method. Domain: lmax >= 0 Conditions: l=0,1,...,lmax Exceptional Return Values: GSL_EDOM, GSL_EUNDRFLW

Irregular Spherical Bessel Functions

Function: int gsl_sf_bessel_y0_impl (double x, gsl_sf_result * result)
Function: int gsl_sf_bessel_y0_e (double x, gsl_sf_result * result)
Exceptional Return Values: none

Function: int gsl_sf_bessel_y1_impl (double x, gsl_sf_result * result)
Function: int gsl_sf_bessel_y1_e (double x, gsl_sf_result * result)
Exceptional Return Values: GSL_EUNDRFLW

Function: int gsl_sf_bessel_y2_impl (double x, gsl_sf_result * result)
Function: int gsl_sf_bessel_y2_e (double x, gsl_sf_result * result)
Exceptional Return Values: GSL_EUNDRFLW

Function: int gsl_sf_bessel_yl_impl (int l, double x, gsl_sf_result * result)
Function: int gsl_sf_bessel_yl_e (int l, double x, gsl_sf_result * result)
Exceptional Return Values: GSL_EUNDRFLW

Function: int gsl_sf_bessel_yl_array_impl (int lmax, double x, double * result_array)
Function: int gsl_sf_bessel_yl_array_e (int lmax, double x, double * result_array)
Domain: lmax >= 0 Conditions: l=0,1,...,lmax Exceptional Return Values: GSL_EUNDRFLW

Regular Modified Spherical Bessel Functions

i_l(x) = \sqrt{\pi/(2x)} I_{l+1/2}(x)

Function: int gsl_sf_bessel_i0_scaled_impl (double x, gsl_sf_result * result)
Function: int gsl_sf_bessel_i0_scaled_e (double x, gsl_sf_result * result)
\exp(-|x|) i_0(x) Exceptional Return Values: none

Function: int gsl_sf_bessel_i1_scaled_impl (double x, gsl_sf_result * result)
Function: int gsl_sf_bessel_i1_scaled_e (double x, gsl_sf_result * result)
\exp(-|x|) i_1(x) Exceptional Return Values: GSL_EUNDRFLW

Function: int gsl_sf_bessel_i2_scaled_impl (double x, gsl_sf_result * result)
Function: int gsl_sf_bessel_i2_scaled_e (double x, gsl_sf_result * result)
\exp(-|x|) i_2(x) Exceptional Return Values: GSL_EUNDRFLW

Function: int gsl_sf_bessel_il_scaled_impl (int l, double x, gsl_sf_result * result)
Function: int gsl_sf_bessel_il_scaled_e (int l, double x, gsl_sf_result * result)
\exp(-|x|) i_l(x) Domain: l >= 0 Exceptional Return Values: GSL_EDOM, GSL_EUNDRFLW

Function: int gsl_sf_bessel_il_scaled_array_impl (int lmax, double x, double * result_array)
Function: int gsl_sf_bessel_il_scaled_array_e (int lmax, double x, double * result_array)
\exp(-|x|) i_l(x) Domain: lmax >= 0 Conditions: l=0,1,...,lmax Exceptional Return Values: GSL_EUNDRFLW

Irregular Modified Spherical Bessel Functions

k_l(x) = \sqrt{\pi/(2x)]} K_{l+1/2}(x)

Function: int gsl_sf_bessel_k0_scaled_impl (double x, gsl_sf_result * result)
Function: int gsl_sf_bessel_k0_scaled_e (double x, gsl_sf_result * result)
\exp(x) k_0(x) Domain: x > 0.0 Exceptional Return Values: GSL_EDOM, GSL_EUNDRFLW

Function: int gsl_sf_bessel_k1_scaled_impl (double x, gsl_sf_result * result)
Function: int gsl_sf_bessel_k1_scaled_e (double x, gsl_sf_result * result)
\exp(x) k_1(x) Domain: x > 0.0 Exceptional Return Values: GSL_EDOM, GSL_EUNDRFLW, GSL_EOVRFLW

Function: int gsl_sf_bessel_k2_scaled_impl (double x, gsl_sf_result * result)
Function: int gsl_sf_bessel_k2_scaled_e (double x, gsl_sf_result * result)
\exp(x) k_2(x) Domain: x > 0.0 Exceptional Return Values: GSL_EDOM, GSL_EUNDRFLW, GSL_EOVRFLW

Function: int gsl_sf_bessel_kl_scaled_impl (int l, double x, gsl_sf_result * result)
Function: int gsl_sf_bessel_kl_scaled_e (int l, double x, gsl_sf_result * result)
\exp(x) k_l(x) Domain: x > 0.0 Exceptional Return Values: GSL_EDOM, GSL_EUNDRFLW

Function: int gsl_sf_bessel_kl_scaled_array_impl (int lmax, double x, double * result_array)
Function: int gsl_sf_bessel_kl_scaled_array_e (int lmax, double x, double * result_array)
\exp(x) k_l(x) Domain: lmax >= 0 Conditions: l=0,1,...,lmax Exceptional Return Values: GSL_EDOM, GSL_EUNDRFLW

Regular Bessel Function, Fractional Order

Function: int gsl_sf_bessel_Jnu_impl (double nu, double x, gsl_sf_result * result)
Function: int gsl_sf_bessel_Jnu_e (double nu, double x, gsl_sf_result * result)
Exceptional Return Values: GSL_EDOM, GSL_EUNDRFLW

Function: int gsl_sf_bessel_sequence_Jnu_impl (double nu, gsl_mode_t mode, size_t size, double * v)
Function: int gsl_sf_bessel_sequence_Jnu_e (double nu, gsl_mode_t mode, size_t size, double * v)
Regular cylindrical Bessel function J_nu(x) evaluated at a series of x values. The array contains the x values. They are assumed to be strictly ordered and positive. The array is over-written with the values of J_nu(x_i). Exceptional Return Values: GSL_EDOM, GSL_EINVAL

Irregular Bessel Functions, Fractional Order

Function: int gsl_sf_bessel_Ynu_impl (double nu, double x, gsl_sf_result * result)
Function: int gsl_sf_bessel_Ynu_e (double nu, double x, gsl_sf_result * result)
Exceptional Return Values:

Regular Modified Bessel Functions, Fractional Order

Function: int gsl_sf_bessel_Inu_scaled_impl (double nu, double x, gsl_sf_result * result)
Function: int gsl_sf_bessel_Inu_scaled_e (double nu, double x, gsl_sf_result * result)
\exp(-|x|) I_{\nu}(x) Domain: x >= 0, nu >= 0 Exceptional Return Values: GSL_EDOM

Function: int gsl_sf_bessel_Inu_impl (double nu, double x, gsl_sf_result * result)
Function: int gsl_sf_bessel_Inu_e (double nu, double x, gsl_sf_result * result)
I_{\nu}(x) Domain: x >= 0, nu >= 0 Exceptional Return Values: GSL_EDOM, GSL_EOVRFLW

Irregular Modified Bessel Functions, Fractional Order

Function: int gsl_sf_bessel_Knu_scaled_impl (double nu, double x, gsl_sf_result * result)
Function: int gsl_sf_bessel_Knu_scaled_e (double nu, double x, gsl_sf_result * result)
\exp(+|x|) K_{\nu}(x) Domain: x > 0, nu >= 0 Exceptional Return Values: GSL_EDOM

Function: int gsl_sf_bessel_Knu_impl (double nu, double x, gsl_sf_result * result)
Function: int gsl_sf_bessel_Knu_e (double nu, double x, gsl_sf_result * result)
K_{\nu}(x) Domain: x > 0, nu >= 0 Exceptional Return Values: GSL_EDOM, GSL_EUNDRFLW

Function: int gsl_sf_bessel_lnKnu_impl (double nu, double x, gsl_sf_result * result)
Function: int gsl_sf_bessel_lnKnu_e (double nu, double x, gsl_sf_result * result)
\ln(K_{\nu}(x)) Domain: x > 0, nu >= 0 Exceptional Return Values: GSL_EDOM

Zeros of Regular Bessel Functions

Function: int gsl_sf_bessel_zero_J0_impl (int s, gsl_sf_result * result)
Function: int gsl_sf_bessel_zero_J0_e (int s, gsl_sf_result * result)
s'th positive zero of the Bessel function J_0(x). Exceptional Return Values:

Function: int gsl_sf_bessel_zero_J1_impl (int s, gsl_sf_result * result)
Function: int gsl_sf_bessel_zero_J1_e (int s, gsl_sf_result * result)
s'th positive zero of the Bessel function J_1(x). Exceptional Return Values:

Function: int gsl_sf_bessel_zero_Jnu_impl (double nu, int s, gsl_sf_result * result)
Function: int gsl_sf_bessel_zero_Jnu_e (double nu, int s, gsl_sf_result * result)
s'th positive zero of the Bessel function J_nu(x). Exceptional Return Values:

Chebyshev Polynomials

The Chebyshev polynomials T_n(x) = \cos(n \arccos x) provide an orthogonal basis of polynomials on the interval [-1,1], with the weight function 1 \over \sqrt{1-x^2}. The first few such polynomials are

 T_0(x) = 1,
 T_1(x) = x,
 T_2(x) = 2 x^2 - 1.

The gsl_sf_cheb_series struct

typedef struct
{
  double * c;   /* coefficients                */
  int order;    /* order of expansion          */
  double a;     /* lower interval point        */
  double b;     /* upper interval point        */
  double * cp;  /* coefficients of derivative  */
  double * ci;  /* coefficients of integral    */

  /* The following exists (mostly) for the benefit
   * of the implementation.  It is an effective single
   * precision order, for use in single precision
   * evaluation.  Users can use it if they like, but
   * only they know how to calculate it, since it is
   * specific to the approximated function.  By default,
   * order_sp = order.
   * It is used explicitly only by the gsl_sf_cheb_eval_mode
   * functions, which are not meant for casual use.

  int order_sp;
} gsl_sf_cheb_struct

Creation/Calculation of Chebyshev Series

Function: gsl_sf_cheb_series * gsl_sf_cheb_new (double (*func)(double),
double a, double b, int order) Calculate a Chebyshev series of specified order over a specified interval, for a given function. Return 0 on failure.

Function: int gsl_sf_cheb_calc_impl (gsl_sf_cheb_series * cs, double (*func)(double))
Function: int gsl_sf_cheb_calc_e (gsl_sf_cheb_series * cs, double (*func)(double))
Calculate a Chebyshev series, but do not allocate a new cheb_series struct. Instead use the one provided. Uses the interval (a,b) and the order with which it was initially created; if you want to change these, then use gsl_sf_cheb_new() instead. Exceptional Return Values: GSL_EFAULT, GSL_ENOMEM

Chebyshev Series Evaluation

Function: int gsl_sf_cheb_eval_impl (const gsl_sf_cheb_series * cs, double x, gsl_sf_result * result)
Function: int gsl_sf_cheb_eval_e (const gsl_sf_cheb_series * cs, double x, gsl_sf_result * result)
Evaluate a Chebyshev series at a given point. No errors can occur for a struct obtained from gsl_sf_cheb_new().

Function: int gsl_sf_cheb_eval_n_impl (const gsl_sf_cheb_series * cs, int order, double x, gsl_sf_result * result)
Function: int gsl_sf_cheb_eval_n_e (const gsl_sf_cheb_series * cs, int order, double x, gsl_sf_result * result)
Evaluate a Chebyshev series at a given point, to (at most) the given order. No errors can occur for a struct obtained from gsl_sf_cheb_new().

Function: int gsl_sf_cheb_eval_mode_impl (const gsl_sf_cheb_series * cs, double x, gsl_mode_t mode, gsl_sf_result * result)
Function: int gsl_sf_cheb_eval_mode_e (const gsl_sf_cheb_series * cs, double x, gsl_mode_t mode, gsl_sf_result * result)
Evaluate a Chebyshev series at a given point, using the default order for double precision mode(s) and the single precision order for other modes. No errors can occur for a struct obtained from gsl_sf_cheb_new().

Function: int gsl_sf_cheb_eval_deriv_impl (gsl_sf_cheb_series * cs, double x, gsl_sf_result * result)
Function: int gsl_sf_cheb_eval_deriv_e (gsl_sf_cheb_series * cs, double x, gsl_sf_result * result)
Evaluate derivative of a Chebyshev series at a given point.

Function: int gsl_sf_cheb_eval_integ_impl (gsl_sf_cheb_series * cs, double x, gsl_sf_result * result)
Function: int gsl_sf_cheb_eval_integ_e (gsl_sf_cheb_series * cs, double x, gsl_sf_result * result)
Evaluate integal of a Chebyshev series at a given point. The integral is fixed by the condition that it equals zero at the left end-point, i.e. it is precisely \int_a^x dt cs(t; a,b) .

Delete a Chebyshev Expansion

Function: void gsl_sf_cheb_free (gsl_sf_cheb_series * cs)
Free a Chebyshev series previously calculated with gsl_sf_cheb_new(). No error can occur.

Clausen Functions

The Clausen integral is defined to be Cl_2(x) = - \int_0^x dt \log(2 \sin(t/2)) .

It is related to the dilogarithm by Cl_2(theta) = Im Li_2(\exp(i theta)) .

Function: int gsl_sf_clausen_impl (double x, gsl_sf_result * result)
Function: int gsl_sf_clausen_e (double x, gsl_sf_result * result)

Coulomb Wave Functions

Normalized Hydrogenic Bound States

/* R_1 := 2Z sqrt(Z) exp(-Z r) */

Function: int gsl_sf_hydrogenicR_1_impl (double Z, double r, gsl_sf_result * result)
Function: int gsl_sf_hydrogenicR_1_e (double Z, double r, gsl_sf_result * result)

/* R_n := norm exp(-Z r/n) (2Z/n)^l Laguerre[n-l-1, 2l+1, 2Z/n r] * * normalization such that psi(n,l,r) = R_n Y_(lm) */

Function: int gsl_sf_hydrogenicR_impl (int n, int l, double Z, double r, gsl_sf_result * result)
Function: int gsl_sf_hydrogenicR_e (int n, int l, double Z, double r, gsl_sf_result * result)

Coulomb Wave Functions

/* Coulomb wave functions F_(lam_F)(eta,x), G_(lam_G)(eta,x) * and their derivatives; lam_G := lam_F - k_lam_G * * lam_F, lam_G > -0.5 * x > 0.0 * * Conventions of Abramowitz+Stegun. * * Because their can be a large dynamic range of values, * overflows are handled gracefully. If an overflow occurs, * GSL_EOVRFLW is signalled and exponent(s) are returned * through exp_F, exp_G. These are such that * * F_L(eta,x) = fc[k_L] * exp(exp_F) * G_L(eta,x) = gc[k_L] * exp(exp_G) * F_L'(eta,x) = fcp[k_L] * exp(exp_F) * G_L'(eta,x) = gcp[k_L] * exp(exp_G) */ int gsl_sf_coulomb_wave_FG_impl(const double eta, const double x, const double lam_F, const int k_lam_G, gsl_sf_result * F, gsl_sf_result * Fp, gsl_sf_result * G, gsl_sf_result * Gp, double * exp_F, double * exp_G)

/* F_L(eta,x) */

Function: int gsl_sf_coulomb_wave_F_array_impl (double lam_min, int kmax, double eta, double x, double * fc_array, double * F_exponent)
Function: int gsl_sf_coulomb_wave_F_array_e (double lam_min, int kmax, double eta, double x, double * fc_array, double * F_exponent)

/* F_L(eta,x), G_L(eta,x) */

Function: int gsl_sf_coulomb_wave_FG_array_impl (double lam_min, int kmax, double eta, double x, double * fc_array, double * gc_array, double * F_exponent, double * G_exponent)
Function: int gsl_sf_coulomb_wave_FG_array_e (double lam_min, int kmax, double eta, double x, double * fc_array, double * gc_array, double * F_exponent, double * G_exponent)

/* F_L(eta,x), G_L(eta,x), F'_L(eta,x), G'_L(eta,x) */

Function: int gsl_sf_coulomb_wave_FGp_impl (double lam_min, int kmax, double eta, double x, gsl_sf_result * fc, gsl_sf_result * fcp, gsl_sf_result * gc, gsl_sf_result * gcp, double * F_exponent, double * G_exponent)
Function: int gsl_sf_coulomb_wave_FGp_e (double lam_min, int kmax, double eta, double x, gsl_sf_result * fc, gsl_sf_result * fcp, gsl_sf_result * gc, gsl_sf_result * gcp, double * F_exponent, double * G_exponent)

/* Coulomb wave function divided by the argument, * F(xi, eta)/xi. This is the function which reduces to * spherical Bessel functions in the limit eta->0. */

Function: int gsl_sf_coulomb_wave_sphF_array_impl (double lam_min, int kmax, double eta, double x, double * fc_array, double * F_exponent)

Coulomb Wave Function Normalization Constant

[Abramowitz+Stegun 14.1.8, 14.1.9]

Function: int gsl_sf_coulomb_CL_impl (double L, double eta, gsl_sf_result * result)
Function: int gsl_sf_coulomb_CL_list (double Lmin, int kmax, double eta, double * cl)

Coupling Coefficients

Since the arguments of the standard coupling coefficient functions are integer or half-integer, the arguments of the following functions are by convention integers equal to twice the actual spin value.

3j Symbols

ja ma
jb jc
mb mc
Function: int gsl_sf_coupling_3j_impl (int two_ja, int two_jb, int two_jc, int two_ma, int two_mb, int two_mc, gsl_sf_result * result)
Function: int gsl_sf_coupling_3j_e (int two_ja, int two_jb, int two_jc, int two_ma, int two_mb, int two_mc, gsl_sf_result * result)
Exceptional Return Values: GSL_EDOM, GSL_EOVRFLW

6j Symbols

ja jd
jb jc
je jf
Function: int gsl_sf_coupling_6j_impl (int two_ja, int two_jb, int two_jc, int two_jd, int two_je, int two_jf, gsl_sf_result * result)
Function: int gsl_sf_coupling_6j_e (int two_ja, int two_jb, int two_jc, int two_jd, int two_je, int two_jf, gsl_sf_result * result)
Exceptional Return Values: GSL_EDOM, GSL_EOVRFLW

9j Symbols

ja jd jg
jb jc
je jf
jh ji
Function: int gsl_sf_coupling_9j_impl (int two_ja, int two_jb, int two_jc, int two_jd, int two_je, int two_jf, int two_jg, int two_jh, int two_ji, gsl_sf_result * result)
Function: int gsl_sf_coupling_9j_e (int two_ja, int two_jb, int two_jc, int two_jd, int two_je, int two_jf, int two_jg, int two_jh, int two_ji, gsl_sf_result * result)
Exceptional Return Values: GSL_EDOM, GSL_EOVRFLW

Dawson Function

The Dawson integral is defined by \exp(-x^2) \int_0^x dt \exp(t^2) .

Function: int gsl_sf_dawson_impl (double x, gsl_sf_result * result)
Function: int gsl_sf_dawson_e (double x, gsl_sf_result * result)
Exceptional Return Values: GSL_EUNDRFLW

Debye Functions

The Debye integrals are defined by D_n(x) = n/x^n \int_0^x dt t^n/(e^t - 1) .

Function: int gsl_sf_debye_1_impl (double x, gsl_sf_result * result)
Function: int gsl_sf_debye_1_e (double x, gsl_sf_result * result)
Exceptional Return Values: GSL_EDOM

Function: int gsl_sf_debye_2_impl (double x, gsl_sf_result * result)
Function: int gsl_sf_debye_2_e (double x, gsl_sf_result * result)
Exceptional Return Values: GSL_EDOM, GSL_EUNDRFLW

Function: int gsl_sf_debye_3_impl (double x, gsl_sf_result * result)
Function: int gsl_sf_debye_3_e (double x, gsl_sf_result * result)
Exceptional Return Values: GSL_EDOM, GSL_EUNDRFLW

Function: int gsl_sf_debye_4_impl (double x, gsl_sf_result * result)
Function: int gsl_sf_debye_4_e (double x, gsl_sf_result * result)
Exceptional Return Values: GSL_EDOM, GSL_EUNDRFLW

Dilogarithm

Real Argument

In Lewin's notation this is Li_2(x) , the real part of the dilogarithm of a real x. It is defined by the integral representation Li_2(x) = - Re \int_0^x ds \log(1-s) / s . Note that Im Li_2(x) = ( 0 for x <= 1, -Pi*log(x) for x > 1 ) .

Function: int gsl_sf_dilog_impl (double x, gsl_sf_result * result)
Function: int gsl_sf_dilog_e (double x, gsl_sf_result * result)

Complex Argument

This is the full complex-valued dilogarithm for complex argument. The argument is z = r \exp(i theta) .

Function: int gsl_sf_complex_dilog_impl (double r, double theta, gsl_sf_result * result_re, gsl_sf_result * result_im)
Function: int gsl_sf_complex_dilog_e (double r, double theta, gsl_sf_result * result_re, gsl_sf_result * result_im)

Elementary Operations

Function: int gsl_sf_multiply_impl (double x, double y, gsl_sf_result * result)
Function: int gsl_sf_multiply_e (double x, double y, gsl_sf_result * result)
Multiplication. Exceptional Return Values: GSL_EOVRFLW, GSL_EUNDRFLW

Function: int gsl_sf_multiply_err_impl (double x, double dx, double y, double dy, gsl_sf_result * result)
Function: int gsl_sf_multiply_err_e (double x, double dx, double y, double dy, gsl_sf_result * result)
Multiplication of quantities with associated errors. Exceptional Return Values: GSL_EOVRFLW, GSL_EUNDRFLW

Elliptic Integrals

Legendre Forms

The Legendre forms of elliptic integrals are defined by F(phi,k) = \int_0^\phi dt 1/\sqrt(1 - k^2 \sin(t)^2) ; E(phi,k) = \int_0^\phi dt \sqrt(1 - k^2 \sin(t)^2) ; P(phi,k,n) = \int_0^\phi dt (1 + n \sin(t)^2)^(-1)/\sqrt(1 - k^2 \sin(t)^2) .

The complete Legendre forms are denoted by K(k) = F(\pi/2, k) ; E(k) = E(\pi/2, k) .

Carlson Forms

The Carlson symmetric forms are defined by RC(x,y) = 1/2 \int_0^\infty dt (t+x)^(-1/2) (t+y)^(-1) ; RD(x,y,z) = 3/2 \int_0^\infty dt (t+x)^(-1/2) (t+y)^(-1/2) (t+z)^(-3/2) ; RF(x,y,z) = 1/2 \int_0^\infty dt (t+x)^(-1/2) (t+y)^(-1/2) (t+z)^(-1/2) ; RJ(x,y,z,p) = 3/2 \int_0^\infty dt (t+x)^(-1/2) (t+y)^(-1/2) (t+z)^(-1/2) (t+p)^(-1) .

Legendre Form of Complete Elliptic Integrals

Function: int gsl_sf_ellint_Kcomp_impl (double k, gsl_mode_t mode, gsl_sf_result * result);
Function: int gsl_sf_ellint_Kcomp_e (double k, gsl_mode_t mode, gsl_sf_result * result);
Exceptional Return Values: GSL_EDOM

Function: int gsl_sf_ellint_Ecomp_impl (double k, gsl_mode_t mode, gsl_sf_result * result);
Function: int gsl_sf_ellint_Ecomp_e (double k, gsl_mode_t mode, gsl_sf_result * result);
Exceptional Return Values: GSL_EDOM

Legendre Form of Incomplete Elliptic Integrals

Function: int gsl_sf_ellint_F_impl (double phi, double k, gsl_mode_t mode, gsl_sf_result * result);
Function: int gsl_sf_ellint_F_e (double phi, double k, gsl_mode_t mode, gsl_sf_result * result);
Exceptional Return Values: GSL_EDOM

Function: int gsl_sf_ellint_E_impl (double phi, double k, gsl_mode_t mode, gsl_sf_result * result);
Function: int gsl_sf_ellint_E_e (double phi, double k, gsl_mode_t mode, gsl_sf_result * result);
Exceptional Return Values: GSL_EDOM

Function: int gsl_sf_ellint_P_impl (double phi, double k, double n, gsl_mode_t mode, gsl_sf_result * result);
Function: int gsl_sf_ellint_P_e (double phi, double k, double n, gsl_mode_t mode, gsl_sf_result * result);
Exceptional Return Values: GSL_EDOM

Function: int gsl_sf_ellint_D_impl (double phi, double k, double n, gsl_mode_t mode, gsl_sf_result * result);
Function: int gsl_sf_ellint_D_e (double phi, double k, double n, gsl_mode_t mode, gsl_sf_result * result);
Exceptional Return Values: GSL_EDOM

Carlson Forms

Function: int gsl_sf_ellint_RC_impl (double x, double y, gsl_mode_t mode, gsl_sf_result * result);
Function: int gsl_sf_ellint_RC_e (double x, double y, gsl_mode_t mode, gsl_sf_result * result);
Exceptional Return Values: GSL_EDOM

Function: int gsl_sf_ellint_RD_impl (double x, double y, double z, gsl_mode_t mode, gsl_sf_result * result);
Function: int gsl_sf_ellint_RD_e (double x, double y, double z, gsl_mode_t mode, gsl_sf_result * result);
Exceptional Return Values: GSL_EDOM

Function: int gsl_sf_ellint_RF_impl (double x, double y, double z, gsl_mode_t mode, gsl_sf_result * result);
Function: int gsl_sf_ellint_RF_e (double x, double y, double z, gsl_mode_t mode, gsl_sf_result * result);
Exceptional Return Values: GSL_EDOM

Function: int gsl_sf_ellint_RJ_impl (double x, double y, double z, double p, gsl_mode_t mode, gsl_sf_result * result);
Function: int gsl_sf_ellint_RJ_e (double x, double y, double z, double p, gsl_mode_t mode, gsl_sf_result * result);
Exceptional Return Values: GSL_EDOM

Elliptic Functions (Jacobi)

Function: int gsl_sf_elljac_impl (double u, double m, double * sn, double * cn, double * dn)
Function: int gsl_sf_elljac_e (double u, double m, double * sn, double * cn, double * dn)
Jacobian elliptic functions sn, dn, cn, calculated by descending Landen transformations. Exceptional Return Values: GSL_EDOM

Error Function

Complementary Error Function

We have erfc(x) = 2/\sqrt(Pi) \int_x^\infty \exp(-t^2) .

Function: int gsl_sf_erfc_impl (double x, gsl_sf_result * result)
Function: int gsl_sf_erfc_e (double x, gsl_sf_result * result)
Exceptional Return Values: none

Log Complementary Error Function

Function: int gsl_sf_log_erfc_impl (double x, gsl_sf_result * result)
Function: int gsl_sf_log_erfc_e (double x, gsl_sf_result * result)
Exceptional Return Values: none

Error Function

We have erf(x) = 2/\sqrt(Pi) \int_0^x dt \exp(-t^2) .

Function: int gsl_sf_erf_impl (double x, gsl_sf_result * result)
Function: int gsl_sf_erf_e (double x, gsl_sf_result * result)
Exceptional Return Values: none

Probability functions

Z(x) : Abramowitz+Stegun 26.2.1 Q(x) : Abramowitz+Stegun 26.2.3

Function: int gsl_sf_erf_Z_impl (double x, gsl_sf_result * result)
Function: int gsl_sf_erf_Q_impl (double x, gsl_sf_result * result)
Function: int gsl_sf_erf_Z_e (double x, gsl_sf_result * result)
Function: int gsl_sf_erf_Q_e (double x, gsl_sf_result * result)
Exceptional Return Values: none

Exponential Function

Exponential Function

/* Provide an exp() function with GSL semantics, * i.e. with proper error checking, etc. *

Function: int gsl_sf_exp_impl (double x, gsl_sf_result * result);
Function: int gsl_sf_exp_e (double x, gsl_sf_result * result);
Exceptional Return Values: GSL_EOVRFLW, GSL_EUNDRFLW

/* Exp(x) *

Function: int gsl_sf_exp_e10_impl (double x, gsl_sf_result_e10 * result);
Function: int gsl_sf_exp_e10_e (double x, gsl_sf_result_e10 * result);
Exceptional Return Values: GSL_EOVRFLW, GSL_EUNDRFLW

/* Exponentiate and multiply by a given factor: y * Exp(x) *

Function: int gsl_sf_exp_mult_impl (double x, double y, gsl_sf_result * result);
Function: int gsl_sf_exp_mult_e (double x, double y, gsl_sf_result * result);
Exceptional Return Values: GSL_EOVRFLW, GSL_EUNDRFLW

/* Exponentiate and multiply by a given factor: y * Exp(x) *

Function: int gsl_sf_exp_mult_e10_impl (const double x, const double y, gsl_sf_result_e10 * result);
Function: int gsl_sf_exp_mult_e10_e (const double x, const double y, gsl_sf_result_e10 * result);
Exceptional Return Values: GSL_EOVRFLW, GSL_EUNDRFLW

Relative Exponential Functions

/* exp(x)-1 *

Function: int gsl_sf_expm1_impl (double x, gsl_sf_result * result);
Function: int gsl_sf_expm1_e (double x, gsl_sf_result * result);
Exceptional Return Values: GSL_EOVRFLW

/* (exp(x)-1)/x = 1 + x/2 + x^2/(2*3) + x^3/(2*3*4) + ... *

Function: int gsl_sf_exprel_impl (double x, gsl_sf_result * result);
Function: int gsl_sf_exprel_e (double x, gsl_sf_result * result);
Exceptional Return Values: GSL_EOVRFLW

/* 2(exp(x)-1-x)/x^2 = 1 + x/3 + x^2/(3*4) + x^3/(3*4*5) + ... *

Function: int gsl_sf_exprel_2_impl (double x, gsl_sf_result * result);
Function: int gsl_sf_exprel_2_e (double x, gsl_sf_result * result);
Exceptional Return Values: GSL_EOVRFLW

/* Similarly for the N-th generalization of * the above. The so-called N-relative exponential * * exprel_N(x) = N!/x^N (exp(x) - Sum[x^k/k!, (k,0,N-1)]) * = 1 + x/(N+1) + x^2/((N+1)(N+2)) + ... * = 1F1(1,1+N,x) */

Function: int gsl_sf_exprel_n_impl (int n, double x, gsl_sf_result * result);
Function: int gsl_sf_exprel_n_e (int n, double x, gsl_sf_result * result);
Exceptional Return Values:

Exponentiation With Error Estimate

/* Exponentiate a quantity with an associated error. */

Function: int gsl_sf_exp_err_impl (double x, double dx, gsl_sf_result * result);
Function: int gsl_sf_exp_err_e (double x, double dx, gsl_sf_result * result);
Exceptional Return Values:

/* Exponentiate a quantity with an associated error. */

Function: int gsl_sf_exp_err_e10_impl (double x, double dx, gsl_sf_result_e10 * result);
Function: int gsl_sf_exp_err_e10_e (double x, double dx, gsl_sf_result_e10 * result);
Exceptional Return Values:

/* Exponentiate and multiply by a given factor: y * Exp(x), * for quantities with associated errors. *

Function: int gsl_sf_exp_mult_err_impl (double x, double dx, double y, double dy, gsl_sf_result * result);
Function: int gsl_sf_exp_mult_err_e (double x, double dx, double y, double dy, gsl_sf_result * result);
Exceptional Return Values: GSL_EOVRFLW, GSL_EUNDRFLW

/* Exponentiate and multiply by a given factor: y * Exp(x), * for quantities with associated errors. *

Function: int gsl_sf_exp_mult_err_e10_impl (double x, double dx, double y, double dy, gsl_sf_result_e10 * result);
Function: int gsl_sf_exp_mult_err_e10_e (double x, double dx, double y, double dy, gsl_sf_result_e10 * result);
Exceptional Return Values: GSL_EOVRFLW, GSL_EUNDRFLW

Exponential Integrals

Exponential Integrals

We have E_1(x) := Re \int_1^\infty dt \exp(-xt)/t ; E_2(x) := Re \int_1^\infty dt \exp(-xt)/t^2 .

Function: int gsl_sf_expint_E1_impl (double x, gsl_sf_result * result)
Function: int gsl_sf_expint_E1_e (double x, gsl_sf_result * result)
Domain: x != 0.0 Exceptional Return Values: GSL_EDOM, GSL_EOVRFLW, GSL_EUNDRFLW

Function: int gsl_sf_expint_E2_impl (double x, gsl_sf_result * result)
Function: int gsl_sf_expint_E2_e (double x, gsl_sf_result * result)
Domain: x != 0.0 Exceptional Return Values: GSL_EDOM, GSL_EOVRFLW, GSL_EUNDRFLW

Ei(x)

We have Ei(x) := PV \int_(-x)^\infty dt \exp(-t)/t .

Function: int gsl_sf_expint_Ei_impl (double x, gsl_sf_result * result)
Function: int gsl_sf_expint_Ei_e (double x, gsl_sf_result * result)
Domain: x != 0.0 Exceptional Return Values: GSL_EDOM, GSL_EOVRFLW, GSL_EUNDRFLW

Hyperbolic Integrals

We have Shi(x) = \int_0^x dt \sinh(t)/t ; Chi(x) := Re[ M_EULER + log(x) + \int_0^x dt (Cosh[t]-1)/t .

Function: int gsl_sf_Shi_impl (double x, gsl_sf_result * result)
Function: int gsl_sf_Shi_e (double x, gsl_sf_result * result)
Exceptional Return Values: GSL_EOVRFLW, GSL_EUNDRFLW

Function: int gsl_sf_Chi_impl (double x, gsl_sf_result * result)
Function: int gsl_sf_Chi_e (double x, gsl_sf_result * result)
Domain: x != 0.0 Exceptional Return Values: GSL_EDOM, GSL_EOVRFLW, GSL_EUNDRFLW

Ei_3(x)

We have Ei_3(x) = \int_0^x dt \exp(-t^3) .

Function: int gsl_sf_expint_3_impl (double x, gsl_sf_result * result)
Function: int gsl_sf_expint_3_e (double x, gsl_sf_result * result)
Domain: x >= 0.0 Exceptional Return Values: GSL_EDOM

Trigonometric Integrals

We have Si(x) = \int_0^x dt \sin(t)/t ; Ci(x) = -\int_x^\infty \cos(t)/t .

Function: int gsl_sf_Si_impl (const double x, gsl_sf_result * result)
Function: int gsl_sf_Si_e (double x, gsl_sf_result * result)
Exceptional Return Values: none

Function: int gsl_sf_Ci_impl (const double x, gsl_sf_result * result)
Function: int gsl_sf_Ci_e (double x, gsl_sf_result * result)
Domain: x > 0.0 Exceptional Return Values: GSL_EDOM

Arctangent Integral

We have AtanInt(x) = \int_0^x dt \arctan(t)/t .

Function: int gsl_sf_atanint_impl (double x, gsl_sf_result * result);
Function: int gsl_sf_atanint_e (double x, gsl_sf_result * result);
Domain: Exceptional Return Values:

Fermi-Dirac Function

Complete Fermi-Dirac Integrals

F_j(x) := 1/Gamma[j+1] Integral[ t^j /(Exp[t-x] + 1), (t,0,Infinity)]

/* Complete integral F_(-1)(x) = e^x / (1 + e^x) *

Function: int gsl_sf_fermi_dirac_m1_impl (double x, gsl_sf_result * result)
Function: int gsl_sf_fermi_dirac_m1_e (double x, gsl_sf_result * result)
Exceptional Return Values: GSL_EUNDRFLW

/* Complete integral F_0(x) = ln(1 + e^x) *

Function: int gsl_sf_fermi_dirac_0_impl (double x, gsl_sf_result * result)
Function: int gsl_sf_fermi_dirac_0_e (double x, gsl_sf_result * result)
Exceptional Return Values: GSL_EUNDRFLW

/* Complete integral F_1(x) *

Function: int gsl_sf_fermi_dirac_1_impl (double x, gsl_sf_result * result)
Function: int gsl_sf_fermi_dirac_1_e (double x, gsl_sf_result * result)
Exceptional Return Values: GSL_EUNDRFLW, GSL_EOVRFLW

/* Complete integral F_2(x) *

Function: int gsl_sf_fermi_dirac_2_impl (double x, gsl_sf_result * result)
Function: int gsl_sf_fermi_dirac_2_e (double x, gsl_sf_result * result)
Exceptional Return Values: GSL_EUNDRFLW, GSL_EOVRFLW

/* Complete integral F_j(x) * for integer j *

Function: int gsl_sf_fermi_dirac_int_impl (int j, double x, gsl_sf_result * result)
Function: int gsl_sf_fermi_dirac_int_e (int j, double x, gsl_sf_result * result)
Exceptional Return Values: GSL_EUNDRFLW, GSL_EOVRFLW

/* Complete integral F_(-1/2)(x) *

Function: int gsl_sf_fermi_dirac_mhalf_impl (double x, gsl_sf_result * result)
Function: int gsl_sf_fermi_dirac_mhalf_e (double x, gsl_sf_result * result)
Exceptional Return Values: GSL_EUNDRFLW, GSL_EOVRFLW

/* Complete integral F_(1/2)(x) *

Function: int gsl_sf_fermi_dirac_half_impl (double x, gsl_sf_result * result)
Function: int gsl_sf_fermi_dirac_half_e (double x, gsl_sf_result * result)
Exceptional Return Values: GSL_EUNDRFLW, GSL_EOVRFLW

/* Complete integral F_(3/2)(x) *

Function: int gsl_sf_fermi_dirac_3half_impl (double x, gsl_sf_result * result)
Function: int gsl_sf_fermi_dirac_3half_e (double x, gsl_sf_result * result)
Exceptional Return Values: GSL_EUNDRFLW, GSL_EOVRFLW

Incomplete Fermi-Dirac Integrals:

F_j(x,b) := 1/Gamma[j+1] Integral[ t^j /(Exp[t-x] + 1), (t,b,Infinity)]

/* Incomplete integral F_0(x,b) = ln(1 + e^(b-x)) - (b-x) *

Function: int gsl_sf_fermi_dirac_inc_0_impl (double x, double b, gsl_sf_result * result)
Function: int gsl_sf_fermi_dirac_inc_0_e (double x, double b, gsl_sf_result * result)
Exceptional Return Values: GSL_EUNDRFLW, GSL_EDOM

Gamma Function

/* Log[Gamma(x)], x not a negative integer * Uses real Lanczos method. * Returns the real part of Log[Gamma[x]] when x < 0, * i.e. Log[|Gamma[x]|]. * * exceptions: GSL_EDOM, GSL_EROUND */ int gsl_sf_lngamma_impl(double x, gsl_sf_result * result); int gsl_sf_lngamma_e(double x, gsl_sf_result * result);

/* Log[Gamma(x)], x not a negative integer * Uses real Lanczos method. Determines * the sign of Gamma[x] as well as Log[|Gamma[x]|] for x < 0. * So Gamma[x] = sgn * Exp[result_lg]. * * exceptions: GSL_EDOM, GSL_EROUND */ int gsl_sf_lngamma_sgn_impl(double x, gsl_sf_result * result_lg, double *sgn); int gsl_sf_lngamma_sgn_e(double x, gsl_sf_result * result_lg, double * sgn);

/* Gamma(x), x not a negative integer * Uses real Lanczos method. * * exceptions: GSL_EDOM, GSL_EOVRFLW, GSL_EROUND */ int gsl_sf_gamma_impl(double x, gsl_sf_result * result); int gsl_sf_gamma_e(double x, gsl_sf_result * result);

/* Regulated Gamma Function, x > 0 * Gamma^*(x) = Gamma(x)/(Sqrt[2Pi] x^(x-1/2) exp(-x)) * = (1 + 1/(12x) + ...), x->Inf * A useful suggestion of Temme. * * exceptions: GSL_EDOM */ int gsl_sf_gammastar_impl(double x, gsl_sf_result * result); int gsl_sf_gammastar_e(double x, gsl_sf_result * result);

/* 1/Gamma(x) * Uses real Lanczos method. * * exceptions: GSL_EUNDRFLW, GSL_EROUND */ int gsl_sf_gammainv_impl(double x, gsl_sf_result * result); int gsl_sf_gammainv_e(double x, gsl_sf_result * result);

/* Log[Gamma(z)] for z complex, z not a negative integer * Uses complex Lanczos method. Note that the phase part (arg) * is not well-determined when |z| is very large, due * to inevitable roundoff in restricting to (-Pi,Pi]. * This will raise the GSL_ELOSS exception when it occurs. * The absolute value part (lnr), however, never suffers. * * Calculates: * lnr = log|Gamma(z)| * arg = arg(Gamma(z)) in (-Pi, Pi] * * exceptions: GSL_EDOM, GSL_ELOSS */ int gsl_sf_lngamma_complex_impl(double zr, double zi, gsl_sf_result * lnr, gsl_sf_result * arg); int gsl_sf_lngamma_complex_e(double zr, double zi, gsl_sf_result * lnr, gsl_sf_result * arg);

/* x^n / n! * * x >= 0.0, n >= 0 * exceptions: GSL_EDOM, GSL_EOVRFLW, GSL_EUNDRFLW */ int gsl_sf_taylorcoeff_impl(int n, double x, gsl_sf_result * result); int gsl_sf_taylorcoeff_e(int n, double x, gsl_sf_result * result);

/* n! * * exceptions: GSL_EDOM, GSL_OVRFLW */ int gsl_sf_fact_impl(unsigned int n, gsl_sf_result * result); int gsl_sf_fact_e(unsigned int n, gsl_sf_result * result);

/* n!! = n(n-2)(n-4) ... * * exceptions: GSL_EDOM, GSL_OVRFLW */ int gsl_sf_doublefact_impl(unsigned int n, gsl_sf_result * result); int gsl_sf_doublefact_e(unsigned int n, gsl_sf_result * result);

/* log(n!) * Faster than ln(Gamma(n+1)) for n < 170; defers for larger n. * * exceptions: none */ int gsl_sf_lnfact_impl(unsigned int n, gsl_sf_result * result); int gsl_sf_lnfact_e(unsigned int n, gsl_sf_result * result);

/* log(n!!) * * exceptions: none */ int gsl_sf_lndoublefact_impl(unsigned int n, gsl_sf_result * result); int gsl_sf_lndoublefact_e(unsigned int n, gsl_sf_result * result);

/* log(n choose m) * * exceptions: GSL_EDOM */ int gsl_sf_lnchoose_impl(unsigned int n, unsigned int m, gsl_sf_result * result); int gsl_sf_lnchoose_e(unsigned int n, unsigned int m, gsl_sf_result * result);

/* n choose m * * exceptions: GSL_EDOM, GSL_EOVRFLW */ int gsl_sf_choose_impl(unsigned int n, unsigned int m, gsl_sf_result * result); int gsl_sf_choose_e(unsigned int n, unsigned int m, gsl_sf_result * result);

/* Logarithm of Pochhammer (Apell) symbol * log( (a)_x ) * where (a)_x := Gamma[a + x]/Gamma[a] * * a > 0, a+x > 0 * * exceptions: GSL_EDOM */ int gsl_sf_lnpoch_impl(double a, double x, gsl_sf_result * result); int gsl_sf_lnpoch_e(double a, double x, gsl_sf_result * result);

/* Logarithm of Pochhammer (Apell) symbol, with sign information. * result = log( |(a)_x| ) * sgn = sgn( (a)_x ) * where (a)_x := Gamma[a + x]/Gamma[a] * * a != neg integer, a+x != neg integer * * exceptions: GSL_EDOM */ int gsl_sf_lnpoch_sgn_impl(double a, double x, gsl_sf_result * result, double * sgn); int gsl_sf_lnpoch_sgn_e(double a, double x, gsl_sf_result * result, double * sgn);

/* Pochhammer (Apell) symbol * (a)_x := Gamma[a + x]/Gamma[x] * * a != neg integer, a+x != neg integer * * exceptions: GSL_EDOM, GSL_EOVRFLW */ int gsl_sf_poch_impl(double a, double x, gsl_sf_result * result); int gsl_sf_poch_e(double a, double x, gsl_sf_result * result);

/* Relative Pochhammer (Apell) symbol * ((a,x) - 1)/x * where (a,x) = (a)_x := Gamma[a + x]/Gamma[a] * * exceptions: GSL_EDOM */ int gsl_sf_pochrel_impl(double a, double x, gsl_sf_result * result); int gsl_sf_pochrel_e(double a, double x, gsl_sf_result * result);

/* Normalized Incomplete Gamma Function * * Q(a,x) = 1/Gamma(a) Integral[ t^(a-1) e^(-t), (t,x,Infinity) ] * * a > 0, x >= 0 * * exceptions: GSL_EDOM */ int gsl_sf_gamma_inc_Q_impl(double a, double x, gsl_sf_result * result); int gsl_sf_gamma_inc_Q_e(double a, double x, gsl_sf_result * result);

/* Complementary Normalized Incomplete Gamma Function * * P(a,x) = 1/Gamma(a) Integral[ t^(a-1) e^(-t), (t,0,x) ] * * a > 0, x >= 0 * * exceptions: GSL_EDOM */ int gsl_sf_gamma_inc_P_impl(double a, double x, gsl_sf_result * result); int gsl_sf_gamma_inc_P_e(double a, double x, gsl_sf_result * result);

/* Logarithm of Beta Function * Log[B(a,b)] * * a > 0, b > 0 * exceptions: GSL_EDOM */ int gsl_sf_lnbeta_impl(double a, double b, gsl_sf_result * result); int gsl_sf_lnbeta_e(double a, double b, gsl_sf_result * result);

/* Beta Function * B(a,b) * * a > 0, b > 0 * exceptions: GSL_EDOM, GSL_EOVRFLW, GSL_EUNDRFLW */ int gsl_sf_beta_impl(double a, double b, gsl_sf_result * result); int gsl_sf_beta_e(double a, double b, gsl_sf_result * result);

/* The maximum x such that gamma(x) is not * considered an overflow. */ #define GSL_SF_GAMMA_XMAX 171.0

Gegenbauer Functions

Function: int gsl_sf_gegenpoly_1_impl (double lambda, double x, gsl_sf_result * result)
Function: int gsl_sf_gegenpoly_2_impl (double lambda, double x, gsl_sf_result * result)
Function: int gsl_sf_gegenpoly_3_impl (double lambda, double x, gsl_sf_result * result)
Function: int gsl_sf_gegenpoly_1_e (double lambda, double x, gsl_sf_result * result)
Function: int gsl_sf_gegenpoly_2_e (double lambda, double x, gsl_sf_result * result)
Function: int gsl_sf_gegenpoly_3_e (double lambda, double x, gsl_sf_result * result)
Evaluate Gegenbauer polynomials using explicit representations. Exceptional Return Values: none

Function: int gsl_sf_gegenpoly_n_impl (int n, double lambda, double x, gsl_sf_result * result)
Function: int gsl_sf_gegenpoly_n_e (int n, double lambda, double x, gsl_sf_result * result)
Evaluate Gegenbauer polynomials. Domain: lambda > -1/2, n >= 0 Exceptional Return Values: GSL_EDOM

Function: int gsl_sf_gegenpoly_array_impl (int nmax, double lambda, double x, double * result_array);
Function: int gsl_sf_gegenpoly_array_e (int nmax, double lambda, double x, double * result_array);
Calculate array of Gegenbauer polynomials. Conditions: n = 0, 1, 2, ... nmax Domain: lambda > -1/2, nmax >= 0 Exceptional Return Values: GSL_EDOM

Hypergeometric Functions

/* Hypergeometric function related to Bessel functions * 0F1[c,x] = * Gamma[c] x^(1/2(1-c)) I_(c-1)(2 Sqrt[x]) * Gamma[c] (-x)^(1/2(1-c)) J_(c-1)(2 Sqrt[-x]) * * exceptions: GSL_EOVRFLW, GSL_EUNDRFLW */ int gsl_sf_hyperg_0F1_impl(double c, double x, gsl_sf_result * result); int gsl_sf_hyperg_0F1_e(double c, double x, gsl_sf_result * result);

/* Confluent hypergeometric function for integer parameters. * 1F1[m,n,x] = M(m,n,x) * * exceptions: */ int gsl_sf_hyperg_1F1_int_impl(int m, int n, double x, gsl_sf_result * result); int gsl_sf_hyperg_1F1_int_e(int m, int n, double x, gsl_sf_result * result);

/* Confluent hypergeometric function. * 1F1[a,b,x] = M(a,b,x) * * exceptions: */ int gsl_sf_hyperg_1F1_impl(double a, double b, double x, gsl_sf_result * result); int gsl_sf_hyperg_1F1_e(double a, double b, double x, gsl_sf_result * result);

/* Confluent hypergeometric function for integer parameters. * U(m,n,x) * * exceptions: */ int gsl_sf_hyperg_U_int_impl(int m, int n, double x, gsl_sf_result * result); int gsl_sf_hyperg_U_int_e(int m, int n, double x, gsl_sf_result * result);

/* Confluent hypergeometric function for integer parameters. * U(m,n,x) * * exceptions: */ int gsl_sf_hyperg_U_int_e10_impl(int m, int n, double x, gsl_sf_result_e10 * result); int gsl_sf_hyperg_U_int_e10_e(int m, int n, double x, gsl_sf_result_e10 * result);

/* Confluent hypergeometric function. * U(a,b,x) * * exceptions: */ int gsl_sf_hyperg_U_impl(double a, double b, double x, gsl_sf_result * result); int gsl_sf_hyperg_U_e(double a, double b, double x, gsl_sf_result * result);

/* Confluent hypergeometric function. * U(a,b,x) * * exceptions: */ int gsl_sf_hyperg_U_e10_impl(double a, double b, double x, gsl_sf_result_e10 * result); int gsl_sf_hyperg_U_e10_e(double a, double b, double x, gsl_sf_result_e10 * result);

/* Gauss hypergeometric function 2F1[a,b,c,x] * |x| < 1 * * exceptions: */ int gsl_sf_hyperg_2F1_impl(double a, double b, double c, double x, gsl_sf_result * result); int gsl_sf_hyperg_2F1_e(double a, double b, double c, double x, gsl_sf_result * result);

/* Gauss hypergeometric function * 2F1[aR + I aI, aR - I aI, c, x] * |x| < 1 * * exceptions: */ int gsl_sf_hyperg_2F1_conj_impl(double aR, double aI, double c, double x, gsl_sf_result * result); int gsl_sf_hyperg_2F1_conj_e(double aR, double aI, double c, double x, gsl_sf_result * result);

/* Renormalized Gauss hypergeometric function * 2F1[a,b,c,x] / Gamma[c] * |x| < 1 * * exceptions: */ int gsl_sf_hyperg_2F1_renorm_impl(double a, double b, double c, double x, gsl_sf_result * result); int gsl_sf_hyperg_2F1_renorm_e(double a, double b, double c, double x, gsl_sf_result * result);

/* Renormalized Gauss hypergeometric function * 2F1[aR + I aI, aR - I aI, c, x] / Gamma[c] * |x| < 1 * * exceptions: */ int gsl_sf_hyperg_2F1_conj_renorm_impl(double aR, double aI, double c, double x, gsl_sf_result * result); int gsl_sf_hyperg_2F1_conj_renorm_e(double aR, double aI, double c, double x, gsl_sf_result * result);

/* Mysterious hypergeometric function. The series representation * is a divergent hypergeometric series. However, for x < 0 we * have 2F0(a,b,x) = (-1/x)^a U(a,1+a-b,-1/x) * * exceptions: GSL_EDOM */ int gsl_sf_hyperg_2F0_impl(double a, double b, double x, gsl_sf_result * result); int gsl_sf_hyperg_2F0_e(double a, double b, double x, gsl_sf_result * result);

Laguerre Functions

The Laguerre polynomials are defined in terms of confluent hypergeometric functions as L^a_n(x) = (a+1)_n / n! 1F1(-n,a+1,x) .

Function: int gsl_sf_laguerre_1_impl (double a, double x, gsl_sf_result * result)
Function: int gsl_sf_laguerre_2_impl (double a, double x, gsl_sf_result * result)
Function: int gsl_sf_laguerre_3_impl (double a, double x, gsl_sf_result * result)
Function: int gsl_sf_laguerre_1_e (double a, double x, gsl_sf_result * result)
Function: int gsl_sf_laguerre_2_e (double a, double x, gsl_sf_result * result)
Function: int gsl_sf_laguerre_3_e (double a, double x, gsl_sf_result * result)
Evaluate generalized Laguerre polynomials using explicit representations. Exceptional Return Values: none

Function: int gsl_sf_laguerre_n_impl (const int n, const double a, const double x, gsl_sf_result * result)
Function: int gsl_sf_laguerre_n_e (int n, double a, double x, gsl_sf_result * result)
Domain: a > -1.0, n >= 0 Evaluate generalized Laguerre polynomials. Exceptional Return Values: GSL_EDOM

Legendre Functions and Spherical Harmonics

Legendre Polynomials

/* P_l(x), l=1,2,3 *

Function: int gsl_sf_legendre_P1_impl (double x, gsl_sf_result * result)
Function: int gsl_sf_legendre_P2_impl (double x, gsl_sf_result * result)
Function: int gsl_sf_legendre_P3_impl (double x, gsl_sf_result * result)
Function: int gsl_sf_legendre_P1_e (double x, gsl_sf_result * result)
Function: int gsl_sf_legendre_P2_e (double x, gsl_sf_result * result)
Function: int gsl_sf_legendre_P3_e (double x, gsl_sf_result * result)
Exceptional Return Values: none

/* P_l(x) l >= 0; |x| <= 1 *

Function: int gsl_sf_legendre_Pl_impl (int l, double x, gsl_sf_result * result)
Function: int gsl_sf_legendre_Pl_e (int l, double x, gsl_sf_result * result)
Exceptional Return Values: GSL_EDOM

/* P_l(x) for l=0,...,lmax; |x| <= 1 *

Function: int gsl_sf_legendre_Pl_array_impl (int lmax, double x, double * result_array)
Function: int gsl_sf_legendre_Pl_array_e (int lmax, double x, double * result_array)
Exceptional Return Values: GSL_EDOM

/* Q_0(x), x > -1, x != 1 *

Function: int gsl_sf_legendre_Q0_impl (double x, gsl_sf_result * result)
Function: int gsl_sf_legendre_Q0_e (double x, gsl_sf_result * result)
Exceptional Return Values: GSL_EDOM

/* Q_1(x), x > -1, x != 1 *

Function: int gsl_sf_legendre_Q1_impl (double x, gsl_sf_result * result)
Function: int gsl_sf_legendre_Q1_e (double x, gsl_sf_result * result)
Exceptional Return Values: GSL_EDOM

/* Q_l(x), x > -1, x != 1, l >= 0 *

Function: int gsl_sf_legendre_Ql_impl (int l, double x, gsl_sf_result * result)
Function: int gsl_sf_legendre_Ql_e (int l, double x, gsl_sf_result * result)
Exceptional Return Values: GSL_EDOM

Associated Legendre Polynomials and Spherical Harmonics

/* P_l^m(x) m >= 0; l >= m; |x| <= 1.0 * * Note that this function grows combinatorially with l. * Therefore we can easily generate an overflow for l larger * than about 150. * * There is no trouble for small m, but when m and l are both large, * then there will be trouble. Rather than allow overflows, these * functions refuse to calculate when they can sense that l and m are * too big. * * If you really want to calculate a spherical harmonic, then DO NOT * use this. Instead use legendre_sphPlm() below, which uses a similar * recursion, but with the normalized functions. *

Function: int gsl_sf_legendre_Plm_impl (int l, int m, double x, gsl_sf_result * result)
Function: int gsl_sf_legendre_Plm_e (int l, int m, double x, gsl_sf_result * result)
Exceptional Return Values: GSL_EDOM, GSL_EOVRFLW

/* P_l^m(x) m >= 0; l >= m; |x| <= 1.0 * l=|m|,...,lmax *

Function: int gsl_sf_legendre_Plm_array_impl (int lmax, int m, double x, double * result_array)
Function: int gsl_sf_legendre_Plm_array_e (int lmax, int m, double x, double * result_array)
Exceptional Return Values: GSL_EDOM, GSL_EOVRFLW

/* P_l^m(x), normalized properly for use in spherical harmonics * m >= 0; l >= m; |x| <= 1.0 * * There is no overflow problem, as there is for the * standard normalization of P_l^m(x). * * Specifically, it returns: * * sqrt((2l+1)/(4pi)) sqrt((l-m)!/(l+m)!) P_l^m(x) *

Function: int gsl_sf_legendre_sphPlm_impl (int l, int m, double x, gsl_sf_result * result)
Function: int gsl_sf_legendre_sphPlm_e (int l, int m, double x, gsl_sf_result * result)
Exceptional Return Values: GSL_EDOM

/* P_l^m(x), normalized properly for use in spherical harmonics * m >= 0; l >= m; |x| <= 1.0 * l=|m|,...,lmax *

Function: int gsl_sf_legendre_sphPlm_array_impl (int lmax, int m, double x, double * result_array)
Function: int gsl_sf_legendre_sphPlm_array_e (int lmax, int m, double x, double * result_array)
Exceptional Return Values: GSL_EDOM

/* size of result_array[] needed for the array versions of Plm * (lmax - m + 1) */

Function: int gsl_sf_legendre_array_size (const int lmax, const int m)
Exceptional Return Values: none

Conical Functions

/* Irregular Spherical Conical Function * P^(1/2)_(-1/2 + I lambda)(x) * * x > -1.0

Function: int gsl_sf_conicalP_half_impl (double lambda, double x, gsl_sf_result * result)
Function: int gsl_sf_conicalP_half_e (double lambda, double x, gsl_sf_result * result)
Exceptional Return Values: GSL_EDOM

/* Regular Spherical Conical Function * P^(-1/2)_(-1/2 + I lambda)(x) * * x > -1.0

Function: int gsl_sf_conicalP_mhalf_impl (double lambda, double x, gsl_sf_result * result)
Function: int gsl_sf_conicalP_mhalf_e (double lambda, double x, gsl_sf_result * result)
Exceptional Return Values: GSL_EDOM

/* Conical Function * P^(0)_(-1/2 + I lambda)(x) * * x > -1.0

Function: int gsl_sf_conicalP_0_impl (double lambda, double x, gsl_sf_result * result)
Function: int gsl_sf_conicalP_0_e (double lambda, double x, gsl_sf_result * result)
Exceptional Return Values: GSL_EDOM

/* Conical Function * P^(1)_(-1/2 + I lambda)(x) * * x > -1.0

Function: int gsl_sf_conicalP_1_impl (double lambda, double x, gsl_sf_result * result)
Function: int gsl_sf_conicalP_1_e (double lambda, double x, gsl_sf_result * result)
Exceptional Return Values: GSL_EDOM

/* Regular Spherical Conical Function * P^(-1/2-l)_(-1/2 + I lambda)(x) * * x > -1.0, l >= -1

Function: int gsl_sf_conicalP_sph_reg_impl (int l, double lambda, double x, gsl_sf_result * result)
Function: int gsl_sf_conicalP_sph_reg_e (int l, double lambda, double x, gsl_sf_result * result)
Exceptional Return Values: GSL_EDOM

/* Regular Cylindrical Conical Function * P^(-m)_(-1/2 + I lambda)(x) * * x > -1.0, m >= -1

Function: int gsl_sf_conicalP_cyl_reg_impl (int m, double lambda, double x, gsl_sf_result * result)
Function: int gsl_sf_conicalP_cyl_reg_e (int m, double lambda, double x, gsl_sf_result * result)
Exceptional Return Values: GSL_EDOM

Radial Functions for Hyperbolic Space

/* The following spherical functions are specializations * of Legendre functions which give the regular eigenfunctions * of the Laplacian on a 3-dimensional hyperbolic space. * Of particular interest is the flat limit, which is * Flat-Lim := (lambda->Inf, eta->0, lambda*eta fixed). */ /* Zeroth radial eigenfunction of the Laplacian on the * 3-dimensional hyperbolic space. * * legendre_H3d_0(lambda,eta) := sin(lambda*eta)/(lambda*sinh(eta)) * * Normalization: * Flat-Lim legendre_H3d_0(lambda,eta) = j_0(lambda*eta) * * eta >= 0.0

Function: int gsl_sf_legendre_H3d_0_impl (double lambda, double eta, gsl_sf_result * result)
Function: int gsl_sf_legendre_H3d_0_e (double lambda, double eta, gsl_sf_result * result)
Exceptional Return Values: GSL_EDOM

/* First radial eigenfunction of the Laplacian on the * 3-dimensional hyperbolic space. * * legendre_H3d_1(lambda,eta) := * 1/sqrt(lambda^2 + 1) sin(lam eta)/(lam sinh(eta)) * (coth(eta) - lambda cot(lambda*eta)) * * Normalization: * Flat-Lim legendre_H3d_1(lambda,eta) = j_1(lambda*eta) * * eta >= 0.0

Function: int gsl_sf_legendre_H3d_1_impl (double lambda, double eta, gsl_sf_result * result)
Function: int gsl_sf_legendre_H3d_1_e (double lambda, double eta, gsl_sf_result * result)
Exceptional Return Values: GSL_EDOM

/* l'th radial eigenfunction of the Laplacian on the * 3-dimensional hyperbolic space. * * Normalization: * Flat-Lim legendre_H3d_l(l,lambda,eta) = j_l(lambda*eta) * * eta >= 0.0, l >= 0

Function: int gsl_sf_legendre_H3d_impl (int l, double lambda, double eta, gsl_sf_result * result)
Function: int gsl_sf_legendre_H3d_e (int l, double lambda, double eta, gsl_sf_result * result)
Exceptional Return Values: GSL_EDOM

/* Array of H3d(ell), 0 <= ell <= lmax */

Function: int gsl_sf_legendre_H3d_array_impl (int lmax, double lambda, double eta, double * result_array)
Function: int gsl_sf_legendre_H3d_array_e (int lmax, double lambda, double eta, double * result_array)
Exceptional Return Values:

Logarithm and Related Functions

Function: int gsl_sf_log_impl (double x, gsl_sf_result * result)
Function: int gsl_sf_log_e (double x, gsl_sf_result * result)
Domain: x > 0.0 Exceptional Return Values: GSL_EDOM

Function: int gsl_sf_log_abs_impl (double x, gsl_sf_result * result)
Function: int gsl_sf_log_abs_e (double x, gsl_sf_result * result)
\log(|x|) Domain: x != 0.0 Exceptional Return Values: GSL_EDOM

/* Complex Logarithm * exp(lnr + I theta) = zr + I zi * Returns argument in [-pi,pi]. *

Function: int gsl_sf_complex_log_impl (double zr, double zi, gsl_sf_result * lnr, gsl_sf_result * theta)
Function: int gsl_sf_complex_log_e (double zr, double zi, gsl_sf_result * lnr, gsl_sf_result * theta)
Exceptional Return Values: GSL_EDOM

Function: int gsl_sf_log_1plusx_impl (double x, gsl_sf_result * result)
Function: int gsl_sf_log_1plusx_e (double x, gsl_sf_result * result)
\log(1 + x) Domain: x > -1.0 Exceptional Return Values: GSL_EDOM

Function: int gsl_sf_log_1plusx_mx_impl (double x, gsl_sf_result * result)
Function: int gsl_sf_log_1plusx_mx_e (double x, gsl_sf_result * result)
\log(1 + x) - x Domain: x > -1.0 Exceptional Return Values: GSL_EDOM

Polynomial Manipulation

Function: double gsl_sf_poly_eval (const double c[], const int len, const double x)
Evaluate a polynomial, c[0] + c[1] x + c[2] x^2 + ... + c[len-1] x^(len-1) . Exceptional Return Values: none

Power Function

A common complaint about the standard C library is its lack of a function for calculating (small) integer powers. GSL provides a simple function to fill this gap. For reasons of efficiency, these functions do not check for overflow or underflow conditions.

int gsl_sf_pow_int_impl(double x, int n, gsl_sf_result * result); int gsl_sf_pow_int_e(double x, int n, gsl_sf_result * result);

Function: double gsl_sf_pow_int (double x, int n)

double gsl_sf_pow_2(const double x); double gsl_sf_pow_3(const double x); double gsl_sf_pow_4(const double x); double gsl_sf_pow_5(const double x); double gsl_sf_pow_6(const double x); double gsl_sf_pow_7(const double x); double gsl_sf_pow_8(const double x); double gsl_sf_pow_9(const double x);

#include <gsl/gsl_sf_pow_int.h>
double y = gsl_sf_pow_int(3.0, 12)

Psi (Digamma) Function

The poly-gamma functions are defined by psi(m,x) = (d/dx)^m psi(0,x) = (d/dx)^(m+1) log(gamma(x)) , where \psi(0,x) is called the di-gamma function.

Digamma Function

Function: int gsl_sf_psi_int_impl (int n, gsl_sf_result * result)
Function: int gsl_sf_psi_int_e (int n, gsl_sf_result * result)
Domain: n integer, n > 0 Exceptional Return Values: GSL_EDOM

Function: int gsl_sf_psi_impl (double x, gsl_sf_result * result)
Function: int gsl_sf_psi_e (double x, gsl_sf_result * result)
Domain: x != 0.0 Exceptional Return Values: GSL_EDOM, GSL_ELOSS

/* Di-Gamma Function Re[psi(1 + I y)] * * exceptions: none */

Function: int gsl_sf_psi_1piy_impl (double y, gsl_sf_result * result)
Function: int gsl_sf_psi_1piy_e (double y, gsl_sf_result * result)
Computes Re psi(1 + I y) . Exceptional Return Values: none

Trigamma Function

Function: int gsl_sf_psi_1_int_impl (int n, gsl_sf_result * result)
Function: int gsl_sf_psi_1_int_e (int n, gsl_sf_result * result)
Domain: n integer, n > 0 Exceptional Return Values: GSL_EDOM

Polygamma Function

Function: int gsl_sf_psi_n_impl (int m, double x, gsl_sf_result * result)
Function: int gsl_sf_psi_n_e (int m, double x, gsl_sf_result * result)
Domain: m >= 0, x > 0.0 Exceptional Return Values: GSL_EDOM

Synchrotron Functions

The first synchrotron function is defined by the integral representation synchrotron_1(x) = x \int_x^\infty dt K_(5/3)(t) .

Function: int gsl_sf_synchrotron_1_impl (double x, gsl_sf_result * result)
Function: int gsl_sf_synchrotron_1_e (double x, gsl_sf_result * result)
Domain: x >= 0.0 Exceptional Return Values: GSL_EDOM, GSL_EUNDRFLW

The second synchroton function is defined by synchrotron_2(x) = x * K_(2/3)(x) .

Function: int gsl_sf_synchrotron_2_impl (double x, gsl_sf_result * result)
Function: int gsl_sf_synchrotron_2_e (double x, gsl_sf_result * result)
Domain: x >= 0.0 Exceptional Return Values: GSL_EDOM, GSL_EUNDRFLW

Transport Functions

The transport functions are defined by the integral representations J(n,x) := \int_0^x dt t^n e^t /(e^t - 1)^2 .

Function: int gsl_sf_transport_2_impl (double x, gsl_sf_result * result)
Function: int gsl_sf_transport_2_e (double x, gsl_sf_result * result)
Calculates J(2,x). Exceptional Return Values: GSL_EDOM

Function: int gsl_sf_transport_3_impl (double x, gsl_sf_result * result)
Function: int gsl_sf_transport_3_e (double x, gsl_sf_result * result)
Calculates J(3,x). Exceptional Return Values: GSL_EDOM, GSL_EUNDRFLW

Function: int gsl_sf_transport_4_impl (double x, gsl_sf_result * result)
Function: int gsl_sf_transport_4_e (double x, gsl_sf_result * result)
Calculates J(4,x). Exceptional Return Values: GSL_EDOM, GSL_EUNDRFLW

Function: int gsl_sf_transport_5_impl (double x, gsl_sf_result * result)
Function: int gsl_sf_transport_5_e (double x, gsl_sf_result * result)
Calculates J(5,x). Exceptional Return Values: GSL_EDOM, GSL_EUNDRFLW

Trigonometric Functions

Trigonometric Functions

These simple trogonometric functions are important because we want to control the error estimate, and trying to guess the error for the standard library implementation every time it is used would be a little goofy.

Function: int gsl_sf_sin_impl (double x, gsl_sf_result * result);
Function: int gsl_sf_sin_e (double x, gsl_sf_result * result);
Exceptional Return Values:

Function: int gsl_sf_cos_impl (double x, gsl_sf_result * result);
Function: int gsl_sf_cos_e (double x, gsl_sf_result * result);
Exceptional Return Values:

Function: int gsl_sf_hypot_impl (double x, double y, gsl_sf_result * result);
Function: int gsl_sf_hypot_e (double x, double y, gsl_sf_result * result);
Exceptional Return Values:

Function: int gsl_sf_complex_sin_impl (double zr, double zi, gsl_sf_result * szr, gsl_sf_result * szi);
Function: int gsl_sf_complex_sin_e (double zr, double zi, gsl_sf_result * szr, gsl_sf_result * szi);
Sin(z) Exceptional Return Values: GSL_EOVRFLW

Function: int gsl_sf_complex_cos_impl (double zr, double zi, gsl_sf_result * czr, gsl_sf_result * czi);
Function: int gsl_sf_complex_cos_e (double zr, double zi, gsl_sf_result * czr, gsl_sf_result * czi);
Cos(z) Exceptional Return Values: GSL_EOVRFLW

Function: int gsl_sf_complex_logsin_impl (double zr, double zi, gsl_sf_result * lszr, gsl_sf_result * lszi);
Function: int gsl_sf_complex_logsin_e (double zr, double zi, gsl_sf_result * lszr, gsl_sf_result * lszi);
Log(Sin(z)) Exceptional Return Values: GSL_EDOM, GSL_ELOSS

Function: int gsl_sf_sinc_impl (double x, gsl_sf_result * result);
Function: int gsl_sf_sinc_e (double x, gsl_sf_result * result);
Sinc(x) = sin(pi x) / (pi x) Exceptional Return Values: none

Function: int gsl_sf_lnsinh_impl (double x, gsl_sf_result * result);
Function: int gsl_sf_lnsinh_e (double x, gsl_sf_result * result);
Log(Sinh(x)) Domain: x > 0 Exceptional Return Values: GSL_EDOM

Function: int gsl_sf_lncosh_impl (double x, gsl_sf_result * result);
Function: int gsl_sf_lncosh_e (double x, gsl_sf_result * result);
Log(Cosh(x)) Exceptional Return Values: none

Conversion Functions

Function: int gsl_sf_polar_to_rect_impl (double r, double theta, gsl_sf_result * x, gsl_sf_result * y);
Function: int gsl_sf_polar_to_rect_e (double r, double theta, gsl_sf_result * x, gsl_sf_result * y);
Convert polar to rectlinear coordinates. Exceptional Return Values: GSL_ELOSS

Function: int gsl_sf_rect_to_polar_impl (double x, double y, gsl_sf_result * r, gsl_sf_result * theta)
Function: int gsl_sf_rect_to_polar_e (double x, double y, gsl_sf_result * r, gsl_sf_result * theta)
Convert rectilinear to polar coordinates. Return argument in range [-pi, pi]. Exceptional Return Values: GSL_EDOM

Restriction Functions

Function: int gsl_sf_angle_restrict_symm_impl (double * theta);
Function: int gsl_sf_angle_restrict_symm_e (double * theta);
Force an angle to lie in the range (-pi,pi]. Exceptional Return Values: GSL_ELOSS

Function: int gsl_sf_angle_restrict_pos_impl (double * theta);
Function: int gsl_sf_angle_restrict_pos_e (double * theta);
Force an angle to lie in the range [0, 2pi). Exceptional Return Values: GSL_ELOSS

 Trigonometric Functions With Error Estimate

Function: int gsl_sf_sin_err_impl (double x, double dx, gsl_sf_result * result);
Function: int gsl_sf_sin_err_e (double x, double dx, gsl_sf_result * result);

Function: int gsl_sf_cos_err_impl (double x, double dx, gsl_sf_result * result);
Function: int gsl_sf_cos_err_e (double x, double dx, gsl_sf_result * result);

Zeta Functions

Riemann Zeta Function

The Riemann zeta function is defined by

Function: int gsl_sf_zeta_int_impl (int n, gsl_sf_result * result)
Function: int gsl_sf_zeta_int_e (int n, gsl_sf_result * result)
Domain: n integer, n != 1 Exceptional Return Values: GSL_EDOM, GSL_EOVRFLW

Function: int gsl_sf_zeta_impl (double s, gsl_sf_result * result)
Function: int gsl_sf_zeta_e (double s, gsl_sf_result * result)
Domain: s != 1.0 Exceptional Return Values: GSL_EDOM, GSL_EOVRFLW

Hurwitz Zeta Function

The Hurwitz zeta function is defined by zeta(s,q) = \sum_0^\infty (k+q)^(-s) .

Function: int gsl_sf_hzeta_impl (double s, double q, gsl_sf_result * result)
Function: int gsl_sf_hzeta_e (double s, double q, gsl_sf_result * result)
Domain: s > 1.0, q > 0.0 Exceptional Return Values: GSL_EDOM, GSL_EUNDRFLW, GSL_EOVRFLW

Eta Function

The eta function is defined by \eta(s) = (1-2^(1-s)) zeta(s) .

Function: int gsl_sf_eta_int_impl (int n, gsl_sf_result * result)
Function: int gsl_sf_eta_int_e (int n, gsl_sf_result * result)
Exceptional Return Values: GSL_EUNDRFLW, GSL_EOVRFLW

Function: int gsl_sf_eta_impl (double s, gsl_sf_result * result)
Function: int gsl_sf_eta_e (double s, gsl_sf_result * result)
Exceptional Return Values: GSL_EUNDRFLW, GSL_EOVRFLW


Go to the first, previous, next, last section, table of contents.