Natural logarithm of gamma function.

@param x Value.
@return Result.
423807
public static double Log(double x) {
        double p, q, w, z;

        double[] A =
                {
                        8.11614167470508450300E-4,
                        -5.95061904284301438324E-4,
                        7.93650340457716943945E-4,
                        -2.77777777730099687205E-3,
                        8.33333333333331927722E-2
                };

        double[] B =
                {
                        -1.37825152569120859100E3,
                        -3.88016315134637840924E4,
                        -3.31612992738871184744E5,
                        -1.16237097492762307383E6,
                        -1.72173700820839662146E6,
                        -8.53555664245765465627E5
                };

        double[] C =
                {
                        -3.51815701436523470549E2,
                        -1.70642106651881159223E4,
                        -2.20528590553854454839E5,
                        -1.13933444367982507207E6,
                        -2.53252307177582951285E6,
                        -2.01889141433532773231E6
                };

        if (x < -34.0) {
            q = -x;
            w = Log(q);
            p = Math.floor(q);

            if (p == q) {
                try {
                    throw new ArithmeticException("Overflow.");
                } catch (Exception e) {
                    e.printStackTrace();
                }
            }

            z = q - p;
            if (z > 0.5) {
                p += 1.0;
                z = p - q;
            }
            z = q * Math.sin(Math.PI * z);

            if (z == 0.0) {
                try {
                    throw new ArithmeticException("Overflow.");
                } catch (Exception e) {
                    e.printStackTrace();
                }
            }

            z = Constants.LogPI - Math.log(z) - w;
            return z;
        }

        if (x < 13.0) {
            z = 1.0;
            while (x >= 3.0) {
                x -= 1.0;
                z *= x;
            }
            while (x < 2.0) {
                if (x == 0.0) {
                    try {
                        throw new ArithmeticException("Overflow.");
                    } catch (Exception e) {
                        e.printStackTrace();
                    }
                }

                z /= x;
                x += 1.0;
            }
            if (z < 0.0) z = -z;
            if (x == 2.0) return Math.log(z);
            x -= 2.0;
            p = x * Special.Polevl(x, B, 5) / Special.P1evl(x, C, 6);
            return (Math.log(z) + p);
        }

        if (x > 2.556348e305) {
            try {
                throw new ArithmeticException("Overflow.");
            } catch (Exception e) {
                e.printStackTrace();
            }
        }

        q = (x - 0.5) * Math.log(x) - x + 0.91893853320467274178;
        if (x > 1.0e8) return (q);

        p = 1.0 / (x * x);
        if (x >= 1000.0) {
            q += ((7.9365079365079365079365e-4 * p
                    - 2.7777777777777777777778e-3) * p
                    + 0.0833333333333333333333) / x;
        } else {
            q += Special.Polevl(p, A, 4) / x;
        }

        return q;
    }
874309
public static double logpdf(double x, double sigma) {
    if(!(x > 0.) || x == Double.POSITIVE_INFINITY) {
      return x == x ? Double.NEGATIVE_INFINITY : Double.NaN;
    }
    final double xs = x / sigma, xs2 = xs * xs;
    return xs2 < Double.POSITIVE_INFINITY ? FastMath.log(xs / sigma) - .5 * xs2 : Double.NEGATIVE_INFINITY;
  }
1253026
@Override
    public double toReal(double x) {
        double unsignedReal = FastMath.exp(natlog(x));
        return (sign(x) == POSITIVE) ? unsignedReal : -unsignedReal;
    }
332469
public static double beta(double x, double y) {
        return Math.exp(Gamma.lgamma(x) + Gamma.lgamma(y) - Gamma.lgamma(x + y));
    }
765130
private function gamma($x)
    {
        $gamma = 0.577215664901532860606512090; # Euler's gamma constant
        if ($x < 0.001) {
            return 1.0/($x*(1.0 + $gamma*$x));
        }
        if ($x < 12.0) {
            # The algorithm directly approximates gamma over (1,2) and uses
            # reduction identities to reduce other arguments to this interval.
            $y = $x;
            $n = 0;
            $arg_was_less_than_one = ($y < 1.0);
            # Add or subtract integers as necessary to bring y into (1,2)
            # Will correct for this below
            if ($arg_was_less_than_one) {
                $y += 1.0;
            } else {
                $n = floor($y) - 1;  # will use n later
                $y -= $n;
            }
            # numerator coefficients for approximation over the interval (1,2)
            $p =
            array(
                -1.71618513886549492533811E+0,
                 2.47656508055759199108314E+1,
                -3.79804256470945635097577E+2,
                 6.29331155312818442661052E+2,
                 8.66966202790413211295064E+2,
                -3.14512729688483675254357E+4,
                -3.61444134186911729807069E+4,
                 6.64561438202405440627855E+4
            );

            # denominator coefficients for approximation over the interval (1,2)
            $q =
            array(
                -3.08402300119738975254353E+1,
                 3.15350626979604161529144E+2,
                -1.01515636749021914166146E+3,
                -3.10777167157231109440444E+3,
                 2.25381184209801510330112E+4,
                 4.75584627752788110767815E+3,
                -1.34659959864969306392456E+5,
                -1.15132259675553483497211E+5
            );

            $num = 0.0;
            $den = 1.0;

            $z = $y - 1;
            for ($i = 0; $i < 8; $i++) {
                $num = ($num + $p[$i])*$z;
                $den = $den*$z + $q[$i];
            }
            $result = $num/$den + 1.0;

            # Apply correction if argument was not initially in (1,2)
            if ($arg_was_less_than_one) {
                # Use identity gamma(z) = gamma(z+1)/z
                # The variable "result" now holds gamma of the original y + 1
                # Thus we use y-1 to get back the orginal y.
                $result /= ($y-1.0);
            } else {
                # Use the identity gamma(z+n) = z*(z+1)* ... *(z+n-1)*gamma(z)
                for ($i = 0; $i < $n; $i++) {
                    $result *= $y++;
                }
            }

            return $result;
        }

        ###########################################################################
        # Third interval: [12, infinity)

        if ($x > 171.624) {
            # Correct answer too large to display.

            return Double.POSITIVE_INFINITY;
        }

        return exp($this->log_gamma($x));
    }
1205324
public static double lgamma(double x) {
        double xcopy = x;
        double fg = 0.0;
        double first = x + LANCZOS_SMALL_GAMMA + 0.5;
        double second = LANCZOS_COEFF[0];

        if (x >= 0.0) {
            if (x >= 1.0 && x - (int) x == 0.0) {
                fg = Math.logFactorial((int) x - 1);
            } else {
                first -= (x + 0.5) * Math.log(first);
                for (int i = 1; i <= LANCZOS_N; i++) {
                    second += LANCZOS_COEFF[i] / ++xcopy;
                }
                fg = Math.log(Math.sqrt(2.0 * Math.PI) * second / x) - first;
            }
        } else {
            fg = Math.PI / (gamma(1.0 - x) * Math.sin(Math.PI * x));

            if (fg != 1.0 / 0.0 && fg != -1.0 / 0.0) {
                if (fg < 0) {
                    throw new IllegalArgumentException("The gamma function is negative: " + fg);
                } else {
                    fg = Math.log(fg);
                }
            }
        }
        return fg;
    }
790741
def log(x, base=None):
    
    _math = infer_math(x)
    if base is None:
        return _math.log(x)
    elif _math == math:
        return _math.log(x, base)
    else:
        # numpy has no option to set a base
        return _math.log(x) / _math.log(base)
768647
public static double logFactorial(int n) {
        if (n < 0) {
            throw new IllegalArgumentException(String.format("n has to be nonnegative: %d", n));
        }

        double f = 0.0;
        for (int i = 2; i <= n; i++) {
            f += Math.log(i);
        }

        return f;
    }
839219
protected static double Gser(double x, double A)
    {
        // Good for X<A+1.
        double T9 = 1 / A;
        double G = T9;
        double I = 1;
        while (T9 > G * 0.00001)
        {
            T9 = T9 * x / (A + I);
            G = G + T9;
            ++I;
        }
        G = G * Math.exp(A * Math.log(x) - x - LogGamma(A));

        return G;
    }
184228
public double Function1D(double x) {
        return Math.exp(x * x / (-2 * sqrSigma)) / (Math.sqrt(2 * Math.PI) * sigma);
    }
281096
function gamma(x) {
      var expx;

      expx = Math.exp(lgamma(x));

      return x > 0 && x === Math.round(x) ? Math.round(expx) : expx;
   }
887905
function gamma0_64(x)
	{
		var rotr1 = rotr_64(x, 1), rotr8 = rotr_64(x, 8), shr7 = shr_64(x, 7);

		return new Int_64(
				rotr1.highOrder ^ rotr8.highOrder ^ shr7.highOrder,
				rotr1.lowOrder ^ rotr8.lowOrder ^ shr7.lowOrder
			);
	}
82014
private static double log1pe(double x) {
        double y = 0.0;
        if (x > 15) {
            y = x;
        } else {
            y += Math.log1p(Math.exp(x));
        }

        return y;
    }
1094609
public static double sigmoid(double x) {
    if (x<0) {
      double num = Math.exp(x);
      return num / (1.0 + num);
    }
    else {
      double den = 1.0 + Math.exp(-x);
      return 1.0 / den;
    }
  }
531499
public static double logcdf(double x, double k, double theta, double shift) {
    final double e = FastMath.exp((x - shift) * theta);
    return e < Double.POSITIVE_INFINITY ? GammaDistribution.logregularizedGammaP(k, e) : 0.;
  }
1295989
function gamma( x ) {
	var sign;
	var q;
	var p;
	var z;
	if (
		(x%1 === 0 && x < 0) || // negative integer
		x === NINF ||
		x !== x // NaN
	) {
		return NaN;
	}
	if ( x === 0 ) {
		// Check for negative zero...
		if ( 1/x < 0 ) {
			return NINF;
		}
		return PINF;
	}
	if (
		x < -170.5674972726612 ||
		x > 171.61447887182298
	) {
		return PINF;
	}
	q = abs( x );
	if ( q > 33 ) {
		if ( x >= 0 ) {
			return stirlingApprox( x );
		}
		p = floor( q );

		// Check whether `x` is even...
		if ( (p&1) === 0 ) {
			sign = -1;
		} else {
			sign = 1;
		}
		z = q - p;
		if ( z > 0.5 ) {
			p += 1;
			z = q - p;
		}
		z = q * sin( PI * z );
		return sign * PI / ( abs(z)*stirlingApprox(q) );
	}
	// Reduce `x`...
	z = 1;
	while ( x >= 3 ) {
		x -= 1;
		z *= x;
	}
	while ( x < 0 ) {
		if ( x > -1e-9 ) {
			return smallApprox( x, z );
		}
		z /= x;
		x += 1;
	}
	while ( x < 2 ) {
		if ( x < 1e-9 ) {
			return smallApprox( x, z );
		}
		z /= x;
		x += 1;
	}
	if ( x === 2 ) {
		return z;
	}
	x -= 2;
	return z * rateval( x );
}
725069
public static double PowerSeries(double a, double b, double x) {
        double s, t, u, v, n, t1, z, ai;

        ai = 1.0 / a;
        u = (1.0 - b) * x;
        v = u / (a + 1.0);
        t1 = v;
        t = u;
        n = 2.0;
        s = 0.0;
        z = Constants.DoubleEpsilon * ai;
        while (Math.abs(v) > z) {
            u = (n - b) * x / n;
            t *= u;
            v = t / (a + n);
            s += v;
            n += 1.0;
        }
        s += t1;
        s += ai;

        u = a * Math.log(x);
        if ((a + b) < Gamma.GammaMax && Math.abs(u) < Constants.LogMax) {
            t = Gamma.Function(a + b) / (Gamma.Function(a) * Gamma.Function(b));
            s = s * t * Math.pow(x, a);
        } else {
            t = Gamma.Log(a + b) - Gamma.Log(a) - Gamma.Log(b) + u + Math.log(s);
            if (t < Constants.LogMin) s = 0.0;
            else s = Math.exp(t);
        }
        return s;
    }