gamma function

double gamma(num x)

Returns an approximation of the gamma function, for details see https://en.wikipedia.org/wiki/Gamma_function.

This uses a Lanczos approximation from https://en.wikipedia.org/wiki/Lanczos_approximation.

Implementation

double gamma(num x) {
  const g = 7;
  const p = [
    0.99999999999980993,
    676.5203681218851,
    -1259.1392167224028,
    771.32342877765313,
    -176.61502916214059,
    12.507343278686905,
    -0.13857109526572012,
    9.9843695780195716e-6,
    1.5056327351493116e-7,
  ];
  if (x < 0.5) {
    if (x.roundToDouble() == x) {
      return double.nan;
    } else {
      return pi / (sin(pi * x) * gamma(1 - x));
    }
  } else if (x > 100.0) {
    return exp(gammaLn(x));
  } else {
    x -= 1.0;
    var y = p[0];
    for (var i = 1; i < g + 2; i++) {
      y += p[i] / (x + i);
    }
    final t = x + g + 0.5;
    return sqrt(2.0 * pi) * pow(t, x + 0.5) * exp(-t) * y;
  }
}