Namespaces

Types in MathNet.Numerics

Type SpecialFunctions

Namespace MathNet.Numerics

This partial implementation of the SpecialFunctions class contains all methods related to the Airy functions.

Static Functions

Public Static Functions

Complex AiryAi(Complex z)

Returns the Airy function Ai.

AiryAi(z) is a solution to the Airy equation, y'' - y * z = 0.

Parameters
Complex z

The value to compute the Airy function of.

Return
Complex

The Airy function Ai.

double AiryAi(double z)

Returns the Airy function Ai.

AiryAi(z) is a solution to the Airy equation, y'' - y * z = 0.

Parameters
double z

The value to compute the Airy function of.

Return
double

The Airy function Ai.

Complex AiryAiPrime(Complex z)

Returns the derivative of the Airy function Ai.

AiryAiPrime(z) is defined as d/dz AiryAi(z).

Parameters
Complex z

The value to compute the derivative of the Airy function of.

Return
Complex

The derivative of the Airy function Ai.

double AiryAiPrime(double z)

Returns the derivative of the Airy function Ai.

AiryAiPrime(z) is defined as d/dz AiryAi(z).

Parameters
double z

The value to compute the derivative of the Airy function of.

Return
double

The derivative of the Airy function Ai.

Complex AiryAiPrimeScaled(Complex z)

Returns the exponentially scaled derivative of Airy function Ai

ScaledAiryAiPrime(z) is given by Exp(zta) * AiryAiPrime(z), where zta = (2/3) * z * Sqrt(z).

Parameters
Complex z

The value to compute the derivative of the Airy function of.

Return
Complex

The exponentially scaled derivative of Airy function Ai.

double AiryAiPrimeScaled(double z)

Returns the exponentially scaled derivative of the Airy function Ai.

ScaledAiryAiPrime(z) is given by Exp(zta) * AiryAiPrime(z), where zta = (2/3) * z * Sqrt(z).

Parameters
double z

The value to compute the derivative of the Airy function of.

Return
double

The exponentially scaled derivative of the Airy function Ai.

Complex AiryAiScaled(Complex z)

Returns the exponentially scaled Airy function Ai.

ScaledAiryAi(z) is given by Exp(zta) * AiryAi(z), where zta = (2/3) * z * Sqrt(z).

Parameters
Complex z

The value to compute the Airy function of.

Return
Complex

The exponentially scaled Airy function Ai.

double AiryAiScaled(double z)

Returns the exponentially scaled Airy function Ai.

ScaledAiryAi(z) is given by Exp(zta) * AiryAi(z), where zta = (2/3) * z * Sqrt(z).

Parameters
double z

The value to compute the Airy function of.

Return
double

The exponentially scaled Airy function Ai.

double AiryBi(double z)

Returns the Airy function Bi.

AiryBi(z) is a solution to the Airy equation, y'' - y * z = 0.

Parameters
double z

The value to compute the Airy function of.

Return
double

The Airy function Bi.

Complex AiryBi(Complex z)

Returns the Airy function Bi.

AiryBi(z) is a solution to the Airy equation, y'' - y * z = 0.

Parameters
Complex z

The value to compute the Airy function of.

Return
Complex

The Airy function Bi.

Complex AiryBiPrime(Complex z)

Returns the derivative of the Airy function Bi.

AiryBiPrime(z) is defined as d/dz AiryBi(z).

Parameters
Complex z

The value to compute the derivative of the Airy function of.

Return
Complex

The derivative of the Airy function Bi.

double AiryBiPrime(double z)

Returns the derivative of the Airy function Bi.

AiryBiPrime(z) is defined as d/dz AiryBi(z).

Parameters
double z

The value to compute the derivative of the Airy function of.

Return
double

The derivative of the Airy function Bi.

Complex AiryBiPrimeScaled(Complex z)

Returns the exponentially scaled derivative of the Airy function Bi.

ScaledAiryBiPrime(z) is given by Exp(-Abs(zta.Real)) * AiryBiPrime(z) where zta = (2 / 3) * z * Sqrt(z).

Parameters
Complex z

The value to compute the derivative of the Airy function of.

Return
Complex

The exponentially scaled derivative of the Airy function Bi.

double AiryBiPrimeScaled(double z)

Returns the exponentially scaled derivative of the Airy function Bi.

ScaledAiryBiPrime(z) is given by Exp(-Abs(zta.Real)) * AiryBiPrime(z) where zta = (2 / 3) * z * Sqrt(z).

Parameters
double z

The value to compute the derivative of the Airy function of.

Return
double

The exponentially scaled derivative of the Airy function Bi.

double AiryBiScaled(double z)

Returns the exponentially scaled Airy function Bi.

ScaledAiryBi(z) is given by Exp(-Abs(zta.Real)) * AiryBi(z) where zta = (2 / 3) * z * Sqrt(z).

Parameters
double z

The value to compute the Airy function of.

Return
double

The exponentially scaled Airy function Bi.

Complex AiryBiScaled(Complex z)

Returns the exponentially scaled Airy function Bi.

ScaledAiryBi(z) is given by Exp(-Abs(zta.Real)) * AiryBi(z) where zta = (2 / 3) * z * Sqrt(z).

Parameters
Complex z

The value to compute the Airy function of.

Return
Complex

The exponentially scaled Airy function Bi(z).

double BesselI(double n, double z)

Returns the modified Bessel function of the first kind.

BesselI(n, z) is a solution to the modified Bessel differential equation.

Parameters
double n

The order of the modified Bessel function.

double z

The value to compute the modified Bessel function of.

Return
double

The modified Bessel function of the first kind.

Complex BesselI(double n, Complex z)

Returns the modified Bessel function of the first kind.

BesselI(n, z) is a solution to the modified Bessel differential equation.

Parameters
double n

The order of the modified Bessel function.

Complex z

The value to compute the modified Bessel function of.

Return
Complex

The modified Bessel function of the first kind.

double BesselI0(double x)

Returns the modified Bessel function of first kind, order 0 of the argument. The function is defined as . The range is partitioned into the two intervals [0, 8] and (8, infinity). Chebyshev polynomial expansions are employed in each interval.
Parameters
double x

The value to compute the Bessel function of.

double BesselI0MStruveL0(double x)

Returns the difference between the Bessel I0 and Struve L0 functions.
Parameters
double x

The value to compute the function of.

double BesselI1(double x)

Returns the modified Bessel function of first kind, order 1 of the argument. The function is defined as . The range is partitioned into the two intervals [0, 8] and (8, infinity). Chebyshev polynomial expansions are employed in each interval.
Parameters
double x

The value to compute the Bessel function of.

double BesselI1MStruveL1(double x)

Returns the difference between the Bessel I1 and Struve L1 functions.
Parameters
double x

The value to compute the function of.

Complex BesselIScaled(double n, Complex z)

Returns the exponentially scaled modified Bessel function of the first kind.

ScaledBesselI(n, z) is given by Exp(-Abs(z.Real)) * BesselI(n, z).

Parameters
double n

The order of the modified Bessel function.

Complex z

The value to compute the modified Bessel function of.

Return
Complex

The exponentially scaled modified Bessel function of the first kind.

double BesselIScaled(double n, double z)

Returns the exponentially scaled modified Bessel function of the first kind.

ScaledBesselI(n, z) is given by Exp(-Abs(z.Real)) * BesselI(n, z).

Parameters
double n

The order of the modified Bessel function.

double z

The value to compute the modified Bessel function of.

Return
double

The exponentially scaled modified Bessel function of the first kind.

Complex BesselJ(double n, Complex z)

Returns the Bessel function of the first kind.

BesselJ(n, z) is a solution to the Bessel differential equation.

Parameters
double n

The order of the Bessel function.

Complex z

The value to compute the Bessel function of.

Return
Complex

The Bessel function of the first kind.

double BesselJ(double n, double z)

Returns the Bessel function of the first kind.

BesselJ(n, z) is a solution to the Bessel differential equation.

Parameters
double n

The order of the Bessel function.

double z

The value to compute the Bessel function of.

Return
double

The Bessel function of the first kind.

double BesselJScaled(double n, double z)

Returns the exponentially scaled Bessel function of the first kind.

ScaledBesselJ(n, z) is given by Exp(-Abs(z.Imaginary)) * BesselJ(n, z).

Parameters
double n

The order of the Bessel function.

double z

The value to compute the Bessel function of.

Return
double

The exponentially scaled Bessel function of the first kind.

Complex BesselJScaled(double n, Complex z)

Returns the exponentially scaled Bessel function of the first kind.

ScaledBesselJ(n, z) is given by Exp(-Abs(z.Imaginary)) * BesselJ(n, z).

Parameters
double n

The order of the Bessel function.

Complex z

The value to compute the Bessel function of.

Return
Complex

The exponentially scaled Bessel function of the first kind.

Complex BesselK(double n, Complex z)

Returns the modified Bessel function of the second kind.

BesselK(n, z) is a solution to the modified Bessel differential equation.

Parameters
double n

The order of the modified Bessel function.

Complex z

The value to compute the modified Bessel function of.

Return
Complex

The modified Bessel function of the second kind.

double BesselK(double n, double z)

Returns the modified Bessel function of the second kind.

BesselK(n, z) is a solution to the modified Bessel differential equation.

Parameters
double n

The order of the modified Bessel function.

double z

The value to compute the modified Bessel function of.

Return
double

The modified Bessel function of the second kind.

double BesselK0(double x)

Returns the modified Bessel function of the second kind of order 0 of the argument. The range is partitioned into the two intervals [0, 8] and (8, infinity). Chebyshev polynomial expansions are employed in each interval.
Parameters
double x

The value to compute the Bessel function of.

double BesselK0e(double x)

Returns the exponentially scaled modified Bessel function of the second kind of order 0 of the argument.
Parameters
double x

The value to compute the Bessel function of.

double BesselK1(double x)

Returns the modified Bessel function of the second kind of order 1 of the argument. The range is partitioned into the two intervals [0, 2] and (2, infinity). Chebyshev polynomial expansions are employed in each interval.
Parameters
double x

The value to compute the Bessel function of.

double BesselK1e(double x)

Returns the exponentially scaled modified Bessel function of the second kind of order 1 of the argument..
Parameters
double x

The value to compute the Bessel function of.

double BesselKScaled(double n, double z)

Returns the exponentially scaled modified Bessel function of the second kind.

ScaledBesselK(n, z) is given by Exp(z) * BesselK(n, z).

Parameters
double n

The order of the modified Bessel function.

double z

The value to compute the modified Bessel function of.

Return
double

The exponentially scaled modified Bessel function of the second kind.

Complex BesselKScaled(double n, Complex z)

Returns the exponentially scaled modified Bessel function of the second kind.

ScaledBesselK(n, z) is given by Exp(z) * BesselK(n, z).

Parameters
double n

The order of the modified Bessel function.

Complex z

The value to compute the modified Bessel function of.

Return
Complex

The exponentially scaled modified Bessel function of the second kind.

double BesselY(double n, double z)

Returns the Bessel function of the second kind.

BesselY(n, z) is a solution to the Bessel differential equation.

Parameters
double n

The order of the Bessel function.

double z

The value to compute the Bessel function of.

Return
double

The Bessel function of the second kind.

Complex BesselY(double n, Complex z)

Returns the Bessel function of the second kind.

BesselY(n, z) is a solution to the Bessel differential equation.

Parameters
double n

The order of the Bessel function.

Complex z

The value to compute the Bessel function of.

Return
Complex

The Bessel function of the second kind.

Complex BesselYScaled(double n, Complex z)

Returns the exponentially scaled Bessel function of the second kind.

ScaledBesselY(n, z) is given by Exp(-Abs(z.Imaginary)) * Y(n, z).

Parameters
double n

The order of the Bessel function.

Complex z

The value to compute the Bessel function of.

Return
Complex

The exponentially scaled Bessel function of the second kind.

double BesselYScaled(double n, double z)

Returns the exponentially scaled Bessel function of the second kind.

ScaledBesselY(n, z) is given by Exp(-Abs(z.Imaginary)) * BesselY(n, z).

Parameters
double n

The order of the Bessel function.

double z

The value to compute the Bessel function of.

Return
double

The exponentially scaled Bessel function of the second kind.

double Beta(double z, double w)

Computes the Euler Beta function.
Parameters
double z

The first Beta parameter, a positive real number.

double w

The second Beta parameter, a positive real number.

Return
double

The Euler Beta function evaluated at z,w.

double BetaIncomplete(double a, double b, double x)

Returns the lower incomplete (unregularized) beta function B(a,b,x) = int(t^(a-1)*(1-t)^(b-1),t=0..x) for real a > 0, b > 0, 1 >= x >= 0.
Parameters
double a

The first Beta parameter, a positive real number.

double b

The second Beta parameter, a positive real number.

double x

The upper limit of the integral.

Return
double

The lower incomplete (unregularized) beta function.

double BetaLn(double z, double w)

Computes the logarithm of the Euler Beta function.
Parameters
double z

The first Beta parameter, a positive real number.

double w

The second Beta parameter, a positive real number.

Return
double

The logarithm of the Euler Beta function evaluated at z,w.

double BetaRegularized(double a, double b, double x)

Returns the regularized lower incomplete beta function I_x(a,b) = 1/Beta(a,b) * int(t^(a-1)*(1-t)^(b-1),t=0..x) for real a > 0, b > 0, 1 >= x >= 0.
Parameters
double a

The first Beta parameter, a positive real number.

double b

The second Beta parameter, a positive real number.

double x

The upper limit of the integral.

Return
double

The regularized lower incomplete beta function.

double Binomial(int n, int k)

Computes the binomial coefficient: n choose k.
Parameters
int n

A nonnegative value n.

int k

A nonnegative value h.

Return
double

The binomial coefficient: n choose k.

double BinomialLn(int n, int k)

Computes the natural logarithm of the binomial coefficient: ln(n choose k).
Parameters
int n

A nonnegative value n.

int k

A nonnegative value h.

Return
double

The logarithmic binomial coefficient: ln(n choose k).

double DiGamma(double x)

Computes the Digamma function which is mathematically defined as the derivative of the logarithm of the gamma function. This implementation is based on Jose Bernardo Algorithm AS 103: Psi ( Digamma ) Function, Applied Statistics, Volume 25, Number 3, 1976, pages 315-317. Using the modifications as in Tom Minka's lightspeed toolbox.
Parameters
double x

The argument of the digamma function.

Return
double

The value of the DiGamma function at x.

double DiGammaInv(double p)

Computes the inverse Digamma function: this is the inverse of the logarithm of the gamma function. This function will only return solutions that are positive.

This implementation is based on the bisection method.

Parameters
double p

The argument of the inverse digamma function.

Return
double

The positive solution to the inverse DiGamma function at p.

double Erf(double x)

Calculates the error function.
Parameters
double x

The value to evaluate.

Return
double

the error function evaluated at given value.

double Erfc(double x)

Calculates the complementary error function.
Parameters
double x

The value to evaluate.

Return
double

the complementary error function evaluated at given value.

double ErfcInv(double z)

Calculates the complementary inverse error function evaluated at z.
We have tested this implementation against the arbitrary precision mpmath library and found cases where we can only guarantee 9 significant figures correct.
Parameters
double z

value to evaluate.

Return
double

The complementary inverse error function evaluated at given value.

double ErfInv(double z)

Calculates the inverse error function evaluated at z.
Parameters
double z

value to evaluate.

Return
double

The inverse error function evaluated at given value.

double Expm1(double power)

Numerically stable exponential minus one, i.e. x -> exp(x)-1
Parameters
double power

A number specifying a power.

Return
double

Returns exp(power)-1.

double ExponentialIntegral(double x, int n)

Computes the generalized Exponential Integral function (En).

This implementation of the computation of the Exponential Integral function follows the derivation in "Handbook of Mathematical Functions, Applied Mathematics Series, Volume 55", Abramowitz, M., and Stegun, I.A. 1964, reprinted 1968 by Dover Publications, New York), Chapters 6, 7, and 26. AND "Advanced mathematical methods for scientists and engineers", Bender, Carl M.; Steven A. Orszag (1978). page 253

for x > 1 uses continued fraction approach that is often used to compute incomplete gamma. for 0 < x <= 1 uses Taylor series expansion

Our unit tests suggest that the accuracy of the Exponential Integral function is correct up to 13 floating point digits.

Parameters
double x

The argument of the Exponential Integral function.

int n

Integer power of the denominator term. Generalization index.

Return
double

The value of the Exponential Integral function.

double ExponentialMinusOne(double power)

Numerically stable exponential minus one, i.e. x -> exp(x)-1
Obsolete: Use Expm1 instead
Parameters
double power

A number specifying a power.

Return
double

Returns exp(power)-1.

BigInteger Factorial(BigInteger x)

Computes the factorial of an integer.

double Factorial(int x)

Computes the factorial function x -> x! of an integer number > 0. The function can represent all number up to 22! exactly, all numbers up to 170! using a double representation. All larger values will overflow.
If you need to multiply or divide various such factorials, consider using the logarithmic version FactorialLn instead so you can add instead of multiply and subtract instead of divide, and then exponentiate the result using Exp. This will also circumvent the problem that factorials become very large even for small parameters.
Return
double

A value value! for value > 0

double FactorialLn(int x)

Computes the logarithmic factorial function x -> ln(x!) of an integer number > 0.
Return
double

A value value! for value > 0

double FallingFactorial(double x, int n)

Computes the Falling Factorial (Pochhammer function) x -> x(n), n>= 0. see: https://en.wikipedia.org/wiki/Falling_and_rising_factorials
Return
double

The real value of the Falling Factorial for x and n

double Gamma(double z)

Computes the Gamma function.

This implementation of the computation of the gamma and logarithm of the gamma function follows the derivation in "An Analysis Of The Lanczos Gamma Approximation", Glendon Ralph Pugh, 2004. We use the implementation listed on p. 116 which should achieve an accuracy of 16 floating point digits. Although 16 digit accuracy should be sufficient for double values, improving accuracy is possible (see p. 126 in Pugh).

Our unit tests suggest that the accuracy of the Gamma function is correct up to 13 floating point digits.

Parameters
double z

The argument of the gamma function.

Return
double

The logarithm of the gamma function.

double GammaLn(double z)

Computes the logarithm of the Gamma function.

This implementation of the computation of the gamma and logarithm of the gamma function follows the derivation in "An Analysis Of The Lanczos Gamma Approximation", Glendon Ralph Pugh, 2004. We use the implementation listed on p. 116 which achieves an accuracy of 16 floating point digits. Although 16 digit accuracy should be sufficient for double values, improving accuracy is possible (see p. 126 in Pugh).

Our unit tests suggest that the accuracy of the Gamma function is correct up to 14 floating point digits.

Parameters
double z

The argument of the gamma function.

Return
double

The logarithm of the gamma function.

double GammaLowerIncomplete(double a, double x)

Returns the lower incomplete gamma function gamma(a,x) = int(exp(-t)t^(a-1),t=0..x) for real a > 0, x > 0.
Parameters
double a

The argument for the gamma function.

double x

The upper integral limit.

Return
double

The lower incomplete gamma function.

double GammaLowerRegularized(double a, double x)

Returns the lower incomplete regularized gamma function P(a,x) = 1/Gamma(a) * int(exp(-t)t^(a-1),t=0..x) for real a > 0, x > 0.
Parameters
double a

The argument for the gamma function.

double x

The upper integral limit.

Return
double

The lower incomplete gamma function.

double GammaLowerRegularizedInv(double a, double y0)

Returns the inverse P^(-1) of the regularized lower incomplete gamma function P(a,x) = 1/Gamma(a) * int(exp(-t)t^(a-1),t=0..x) for real a > 0, x > 0, such that P^(-1)(a,P(a,x)) == x.

double GammaUpperIncomplete(double a, double x)

Returns the upper incomplete gamma function Gamma(a,x) = int(exp(-t)t^(a-1),t=0..x) for real a > 0, x > 0.
Parameters
double a

The argument for the gamma function.

double x

The lower integral limit.

Return
double

The upper incomplete gamma function.

double GammaUpperRegularized(double a, double x)

Returns the upper incomplete regularized gamma function Q(a,x) = 1/Gamma(a) * int(exp(-t)t^(a-1),t=0..x) for real a > 0, x > 0.
Parameters
double a

The argument for the gamma function.

double x

The lower integral limit.

Return
double

The upper incomplete regularized gamma function.

double GeneralHarmonic(int n, double m)

Compute the generalized harmonic number of order n of m. (1 + 1/2^m + 1/3^m +... + 1/n^m)
Parameters
int n

The order parameter.

double m

The power parameter.

Return
double

General Harmonic number.

double GeneralizedHypergeometric(Double[] a, Double[] b, int z)

A generalized hypergeometric series is a power series in which the ratio of successive coefficients indexed by n is a rational function of n. This is the most common pFq(a1,..., ap; b1,...,bq; z) representation see: https://en.wikipedia.org/wiki/Generalized_hypergeometric_function
Parameters
Double[] a

The list of coefficients in the numerator

Double[] b

The list of coefficients in the denominator

int z

The variable in the power series

Return
double

The value of the Generalized HyperGeometric Function.

Complex HankelH1(double n, Complex z)

Returns the Hankel function of the first kind.

HankelH1(n, z) is defined as BesselJ(n, z) + j * BesselY(n, z).

Parameters
double n

The order of the Hankel function.

Complex z

The value to compute the Hankel function of.

Return
Complex

The Hankel function of the first kind.

Complex HankelH1Scaled(double n, Complex z)

Returns the exponentially scaled Hankel function of the first kind.

ScaledHankelH1(n, z) is given by Exp(-z * j) * HankelH1(n, z) where j = Sqrt(-1).

Parameters
double n

The order of the Hankel function.

Complex z

The value to compute the Hankel function of.

Return
Complex

The exponentially scaled Hankel function of the first kind.

Complex HankelH2(double n, Complex z)

Returns the Hankel function of the second kind.

HankelH2(n, z) is defined as BesselJ(n, z) - j * BesselY(n, z).

Parameters
double n

The order of the Hankel function.

Complex z

The value to compute the Hankel function of.

Return
Complex

The Hankel function of the second kind.

Complex HankelH2Scaled(double n, Complex z)

Returns the exponentially scaled Hankel function of the second kind.

ScaledHankelH2(n, z) is given by Exp(z * j) * HankelH2(n, z) where j = Sqrt(-1).

Parameters
double n

The order of the Hankel function.

Complex z

The value to compute the Hankel function of.

Return
Complex

The exponentially scaled Hankel function of the second kind.

double Harmonic(int t)

Computes the t 'th Harmonic number.
Parameters
int t

The Harmonic number which needs to be computed.

Return
double

The t'th Harmonic number.

Complex32 Hypotenuse(Complex32 a, Complex32 b)

Numerically stable hypotenuse of a right angle triangle, i.e. (a,b) -> sqrt(a^2 + b^2)
Parameters
Complex32 a

The length of side a of the triangle.

Complex32 b

The length of side b of the triangle.

Return
Complex32

Returns sqrt(a2 + b2) without underflow/overflow.

double Hypotenuse(double a, double b)

Numerically stable hypotenuse of a right angle triangle, i.e. (a,b) -> sqrt(a^2 + b^2)
Parameters
double a

The length of side a of the triangle.

double b

The length of side b of the triangle.

Return
double

Returns sqrt(a2 + b2) without underflow/overflow.

float Hypotenuse(float a, float b)

Numerically stable hypotenuse of a right angle triangle, i.e. (a,b) -> sqrt(a^2 + b^2)
Parameters
float a

The length of side a of the triangle.

float b

The length of side b of the triangle.

Return
float

Returns sqrt(a2 + b2) without underflow/overflow.

Complex Hypotenuse(Complex a, Complex b)

Numerically stable hypotenuse of a right angle triangle, i.e. (a,b) -> sqrt(a^2 + b^2)
Parameters
Complex a

The length of side a of the triangle.

Complex b

The length of side b of the triangle.

Return
Complex

Returns sqrt(a2 + b2) without underflow/overflow.

Complex KelvinBe(double nu, double x)

Returns the Kelvin function of the first kind.

KelvinBe(nu, x) is given by BesselJ(0, j * sqrt(j) * x) where j = sqrt(-1).

KelvinBer(nu, x) and KelvinBei(nu, x) are the real and imaginary parts of the KelvinBe(nu, x)

Parameters
double nu

the order of the the Kelvin function.

double x

The value to compute the Kelvin function of.

Return
Complex

The Kelvin function of the first kind.

double KelvinBei(double nu, double x)

Returns the Kelvin function bei.

KelvinBei(nu, x) is given by the imaginary part of BesselJ(nu, j * sqrt(j) * x) where j = sqrt(-1).

Parameters
double nu

the order of the the Kelvin function.

double x

The value to compute the Kelvin function of.

Return
double

The Kelvin function bei.

double KelvinBei(double x)

Returns the Kelvin function bei.

KelvinBei(x) is given by the imaginary part of BesselJ(0, j * sqrt(j) * x) where j = sqrt(-1).

KelvinBei(x) is equivalent to KelvinBei(0, x).

Parameters
double x

The value to compute the Kelvin function of.

Return
double

The Kelvin function bei.

double KelvinBeiPrime(double x)

Returns the derivative of the Kelvin function bei.
Parameters
double x

The value to compute the derivative of the Kelvin function of.

Return
double

The derivative of the Kelvin function bei.

double KelvinBeiPrime(double nu, double x)

Returns the derivative of the Kelvin function bei.
Parameters
double nu

The order of the Kelvin function.

double x

The value to compute the derivative of the Kelvin function of.

Return
double

the derivative of the Kelvin function bei.

double KelvinBer(double nu, double x)

Returns the Kelvin function ber.

KelvinBer(nu, x) is given by the real part of BesselJ(nu, j * sqrt(j) * x) where j = sqrt(-1).

Parameters
double nu

the order of the the Kelvin function.

double x

The value to compute the Kelvin function of.

Return
double

The Kelvin function ber.

double KelvinBer(double x)

Returns the Kelvin function ber.

KelvinBer(x) is given by the real part of BesselJ(0, j * sqrt(j) * x) where j = sqrt(-1).

KelvinBer(x) is equivalent to KelvinBer(0, x).

Parameters
double x

The value to compute the Kelvin function of.

Return
double

The Kelvin function ber.

double KelvinBerPrime(double nu, double x)

Returns the derivative of the Kelvin function ber.
Parameters
double nu

The order of the Kelvin function.

double x

The value to compute the derivative of the Kelvin function of.

Return
double

the derivative of the Kelvin function ber

double KelvinBerPrime(double x)

Returns the derivative of the Kelvin function ber.
Parameters
double x

The value to compute the derivative of the Kelvin function of.

Return
double

The derivative of the Kelvin function ber.

Complex KelvinKe(double nu, double x)

Returns the Kelvin function of the second kind

KelvinKe(nu, x) is given by Exp(-nu * pi * j / 2) * BesselK(nu, x * sqrt(j)) where j = sqrt(-1).

KelvinKer(nu, x) and KelvinKei(nu, x) are the real and imaginary parts of the KelvinBe(nu, x)

Parameters
double nu

The order of the Kelvin function.

double x

The value to calculate the kelvin function of,

double KelvinKei(double x)

Returns the Kelvin function kei.

KelvinKei(x) is given by the imaginary part of Exp(-nu * pi * j / 2) * BesselK(0, sqrt(j) * x) where j = sqrt(-1).

KelvinKei(x) is equivalent to KelvinKei(0, x).

Parameters
double x

The non-negative real value to compute the Kelvin function of.

Return
double

The Kelvin function kei.

double KelvinKei(double nu, double x)

Returns the Kelvin function kei.

KelvinKei(nu, x) is given by the imaginary part of Exp(-nu * pi * j / 2) * BesselK(nu, sqrt(j) * x) where j = sqrt(-1).

Parameters
double nu

the order of the the Kelvin function.

double x

The non-negative real value to compute the Kelvin function of.

Return
double

The Kelvin function kei.

double KelvinKeiPrime(double x)

Returns the derivative of the Kelvin function kei.
Parameters
double x

The value to compute the derivative of the Kelvin function of.

Return
double

The derivative of the Kelvin function kei.

double KelvinKeiPrime(double nu, double x)

Returns the derivative of the Kelvin function kei.
Parameters
double nu

The order of the Kelvin function.

double x

The value to compute the derivative of the Kelvin function of.

Return
double

The derivative of the Kelvin function kei.

double KelvinKer(double nu, double x)

Returns the Kelvin function ker.

KelvinKer(nu, x) is given by the real part of Exp(-nu * pi * j / 2) * BesselK(nu, sqrt(j) * x) where j = sqrt(-1).

Parameters
double nu

the order of the the Kelvin function.

double x

The non-negative real value to compute the Kelvin function of.

Return
double

The Kelvin function ker.

double KelvinKer(double x)

Returns the Kelvin function ker.

KelvinKer(x) is given by the real part of Exp(-nu * pi * j / 2) * BesselK(0, sqrt(j) * x) where j = sqrt(-1).

KelvinKer(x) is equivalent to KelvinKer(0, x).

Parameters
double x

The non-negative real value to compute the Kelvin function of.

Return
double

The Kelvin function ker.

double KelvinKerPrime(double nu, double x)

Returns the derivative of the Kelvin function ker.
Parameters
double nu

The order of the Kelvin function.

double x

The non-negative real value to compute the derivative of the Kelvin function of.

Return
double

The derivative of the Kelvin function ker.

double KelvinKerPrime(double x)

Returns the derivative of the Kelvin function ker.
Parameters
double x

The value to compute the derivative of the Kelvin function of.

Return
double

The derivative of the Kelvin function ker.

double Log1p(double x)

Computes ln(1+x) with good relative precision when |x| is small
Parameters
double x

The parameter for which to compute the log1p function. Range: x > 0.

double Logistic(double p)

Computes the logistic function. see: http://en.wikipedia.org/wiki/Logistic
Parameters
double p

The parameter for which to compute the logistic function.

Return
double

The logistic function of p.

double Logit(double p)

Computes the logit function, the inverse of the sigmoid logistic function. see: http://en.wikipedia.org/wiki/Logit
Parameters
double p

The parameter for which to compute the logit function. This number should be between 0 and 1.

Return
double

The logarithm of p divided by 1.0 - p.

double MarcumQ(double nu, double a, double b)

Returns the Marcum Q-function Q[ν](a,b).

References: A. Gil, J. Segura and N.M. Temme. Efficient and accurate algorithms for the computation and inversion of the incomplete gamma function ratios. SIAM J Sci Comput. (2012) 34(6), A2965-A2981

Parameters
double nu

The order of generalized Marcum Q-function. Range: 1≦ν≦10000

double a

The value to compute the Marcum Q-function of. Range: 0≦a≦10000

double b

The value to compute the Marcum Q-function of. Range: 0≦b≦10000

Return
double

The Marcum Q-function Q[ν](a,b)

double MarcumQ(double nu, double a, double b, Int32& err)

double Multinomial(int n, Int32[] ni)

Computes the multinomial coefficient: n choose n1, n2, n3,...
Parameters
int n

A nonnegative value n.

Int32[] ni

An array of nonnegative values that sum to n.

Return
double

The multinomial coefficient.

double RisingFactorial(double x, int n)

Computes the Rising Factorial (Pochhammer function) x -> (x)n, n>= 0. see: https://en.wikipedia.org/wiki/Falling_and_rising_factorials
Return
double

The real value of the Rising Factorial for x and n

Complex SphericalBesselJ(double n, Complex z)

Returns the spherical Bessel function of the first kind.

SphericalBesselJ(n, z) is given by Sqrt(pi/2) / Sqrt(z) * BesselJ(n + 1/2, z).

Parameters
double n

The order of the spherical Bessel function.

Complex z

The value to compute the spherical Bessel function of.

Return
Complex

The spherical Bessel function of the first kind.

double SphericalBesselJ(double n, double z)

Returns the spherical Bessel function of the first kind.

SphericalBesselJ(n, z) is given by Sqrt(pi/2) / Sqrt(z) * BesselJ(n + 1/2, z).

Parameters
double n

The order of the spherical Bessel function.

double z

The value to compute the spherical Bessel function of.

Return
double

The spherical Bessel function of the first kind.

Complex SphericalBesselY(double n, Complex z)

Returns the spherical Bessel function of the second kind.

SphericalBesselY(n, z) is given by Sqrt(pi/2) / Sqrt(z) * BesselY(n + 1/2, z).

Parameters
double n

The order of the spherical Bessel function.

Complex z

The value to compute the spherical Bessel function of.

Return
Complex

The spherical Bessel function of the second kind.

double SphericalBesselY(double n, double z)

Returns the spherical Bessel function of the second kind.

SphericalBesselY(n, z) is given by Sqrt(pi/2) / Sqrt(z) * BesselY(n + 1/2, z).

Parameters
double n

The order of the spherical Bessel function.

double z

The value to compute the spherical Bessel function of.

Return
double

The spherical Bessel function of the second kind.

double StruveL0(double x)

Returns the modified Struve function of order 0.
Parameters
double x

The value to compute the function of.

double StruveL1(double x)

Returns the modified Struve function of order 1.
Parameters
double x

The value to compute the function of.