std.mathspecial

pure nothrow @nogc @safe real gamma(real x);

The Gamma function, Γ(x)

Γ(x) is a generalisation of the factorial function to real and complex numbers. Like x!, Γ(x+1) = x * Γ(x).

Mathematically, if z.re > 0 then Γ(z) = ∫₀^∞ t^z-1e^-t dt

Special Values

x	Γ(x)
NAN	NAN
±0.0	±∞
integer > 0	(x-1)!
integer < 0	NAN
+∞	+∞
-∞	NAN

pure nothrow @nogc @safe real logGamma(real x);

Natural logarithm of the gamma function, Γ(x)

Returns the base e (2.718...) logarithm of the absolute value of the gamma function of the argument.

For reals, logGamma is equivalent to log(fabs(gamma(x))).

Special Values

x	logGamma(x)
NAN	NAN
integer <= 0	+∞
±∞	+∞

pure nothrow @nogc @safe real sgnGamma(real x);

The sign of Γ(x).

Returns -1 if Γ(x) < 0, +1 if Γ(x) > 0, NAN if sign is indeterminate.

Note that this function can be used in conjunction with logGamma(x) to evaluate gamma for very large values of x.

pure nothrow @nogc @safe real beta(real x, real y);

Beta function

The beta function is defined as

beta(x, y) = (Γ(x) * Γ(y)) / Γ(x + y)

pure nothrow @nogc @safe real digamma(real x);

Digamma function

The digamma function is the logarithmic derivative of the gamma function.

digamma(x) = d/dx logGamma(x)

pure nothrow @nogc @safe real betaIncomplete(real a, real b, real x);

Incomplete beta integral

Returns incomplete beta integral of the arguments, evaluated from zero to x. The regularized incomplete beta function is defined as

betaIncomplete(a, b, x) = Γ(a + b) / ( Γ(a) Γ(b) ) * ∫₀^x t^a-1(1-t)^b-1 dt

and is the same as the the cumulative distribution function.

The domain of definition is 0 <= x <= 1. In this implementation a and b are restricted to positive values. The integral from x to 1 may be obtained by the symmetry relation

betaIncompleteCompl(a, b, x ) = betaIncomplete( b, a, 1-x )

The integral is evaluated by a continued fraction expansion or, when b * x is small, by a power series.

pure nothrow @nogc @safe real betaIncompleteInverse(real a, real b, real y);

Inverse of incomplete beta integral

Given y, the function finds x such that

betaIncomplete(a, b, x) == y

Newton iterations or interval halving is used.

pure nothrow @nogc @safe real gammaIncomplete(real a, real x);

pure nothrow @nogc @safe real gammaIncompleteCompl(real a, real x);

Incomplete gamma integral and its complement

These functions are defined by

gammaIncomplete = ( ∫₀^x e^-t t^a-1 dt )/ Γ(a)

gammaIncompleteCompl(a,x) = 1 - gammaIncomplete(a,x) = (∫_x^∞ e^-t t^a-1 dt )/ Γ(a)

In this implementation both arguments must be positive. The integral is evaluated by either a power series or continued fraction expansion, depending on the relative values of a and x.

pure nothrow @nogc @safe real gammaIncompleteComplInverse(real a, real p);

Inverse of complemented incomplete gamma integral

Given a and p, the function finds x such that

gammaIncompleteCompl( a, x ) = p.

pure nothrow @nogc @safe real erf(real x);

Error function

The integral is

erf(x) = 2/ √(π) ∫₀^x exp( - t²) dt

The magnitude of x is limited to about 106.56 for IEEE 80-bit arithmetic; 1 or -1 is returned outside this range.

pure nothrow @nogc @safe real erfc(real x);

Complementary error function

erfc(x) = 1 - erf(x) = 2/ √(π) ∫_x^∞ exp( - t²) dt

This function has high relative accuracy for values of x far from zero. (For values near zero, use erf(x)).

pure nothrow @nogc @safe real normalDistribution(real x);

Normal distribution function.

The normal (or Gaussian, or bell-shaped) distribution is defined as:

normalDist(x) = 1/√ π ∫_-∞^x exp( - t²/2) dt = 0.5 + 0.5 * erf(x/sqrt(2)) = 0.5 * erfc(- x/sqrt(2))

To maintain accuracy at values of x near 1.0, use normalDistribution(x) = 1.0 - normalDistribution(-x).

References: http://www.netlib.org/cephes/ldoubdoc.html, G. Marsaglia, "Evaluating the Normal Distribution", Journal of Statistical Software 11, (July 2004).

pure nothrow @nogc @safe real normalDistributionInverse(real p);

Inverse of Normal distribution function

Returns the argument, x, for which the area under the Normal probability density function (integrated from minus infinity to x) is equal to p.