10

How can I make a function that calculates the factorial (or the gamma function) of decimal numbers in JavaScript? For example, how could I calculate 2.33!?

Peter Olson
  • 139,199
  • 49
  • 202
  • 242
Mich'
  • 683
  • 1
  • 8
  • 17
  • What have you tried? At least there's a possible duplicated: http://stackoverflow.com/questions/3959211/fast-factorial-function-in-javascript – Rubens Mariuzzo Mar 16 '13 at 20:11
  • 1
    function fact(n) { var x = 1 for (var i = 2; i <= n; i++) { x = x * i } return x} But this can't compute factorials with decimals and i don't know how i can make a gamma function for it – Mich' Mar 16 '13 at 20:14
  • 1
    I found that before and an answer but it's doesn't really goes exacly correct – Mich' Mar 16 '13 at 20:16
  • You are right @Mich, voting up. – Rubens Mariuzzo Mar 16 '13 at 20:16
  • the array with constants doesn't seem to be correct – Mich' Mar 16 '13 at 20:17
  • It's not trivial, so what is your background (did you read and understand http://en.wikipedia.org/wiki/Gamma_function)? What do you need this for in javascript? – Bergi Mar 16 '13 at 20:59

6 Answers6

12

I might have found an existing solution... It's an implementation of Lanczos method, I found it at the swedish wikipedia (http://sv.wikipedia.org/wiki/Gammafunktionen). It was written in python and says to be correct up to 15 decimals. I ported it to js, cross checked some random values against (http://www.efunda.com/math/gamma/findgamma.cfm).

http://jsfiddle.net/Fzy9C/

var g = 7;
var C = [0.99999999999980993, 676.5203681218851, -1259.1392167224028,771.32342877765313, -176.61502916214059, 12.507343278686905, -0.13857109526572012, 9.9843695780195716e-6, 1.5056327351493116e-7];

function gamma(z) {

    if (z < 0.5) return Math.PI / (Math.sin(Math.PI * z) * gamma(1 - z));
    else {
        z -= 1;

        var x = C[0];
        for (var i = 1; i < g + 2; i++)
        x += C[i] / (z + i);

        var t = z + g + 0.5;
        return Math.sqrt(2 * Math.PI) * Math.pow(t, (z + 0.5)) * Math.exp(-t) * x;
    }
}

(and ofcourse it does not support imaginary numbers, since js does not)

apelsinapa
  • 575
  • 3
  • 9
7

As an alternative to the other answers here, here's a much simpler approximation for the gamma function, proposed in 2007 by Gergő Nemes. (See the wikipedia page on Stirling's approximation).

enter image description here

This can be implemented directly in JavaScript in a single line:

function gamma(z) {
  return Math.sqrt(2 * Math.PI / z) * Math.pow((1 / Math.E) * (z + 1 / (12 * z - 1 / (10 * z))), z);
}

You can see this in action on jsFiddle.

This is accurate to 8 digits for z > 8, but it is still accurate to a handful of digits for smaller z. It is not quite as accurate as Lanczos approximation, but it is simpler and also slightly faster.

Note that the gamma function and the factorial function are slightly different. The factorial function can be defined in terms of the gamma function thus:

function factorial(n) {
  return gamma(n + 1);
}
Peter Olson
  • 139,199
  • 49
  • 202
  • 242
  • Is there something incorrect here? I'd like to know why this was downvoted. – Peter Olson Mar 16 '13 at 21:59
  • @PeterOhlson I'm not sure, but I get 3.11... instead of 1.27... when i do factorial(2.44) in jsfiddle. And i don't get other number to return a good value either... sorry about the downvote, i voted up but ran into this problem so ill vote up again if i see it functional. It't looks like a good solution. – apelsinapa Mar 16 '13 at 22:00
  • @apelsinapa [3.11 is the correct answer](http://www.wolframalpha.com/input/?i=2.44%21). 1.27 is factorial(1.44), not factorial(2.44). – Peter Olson Mar 16 '13 at 22:02
  • @PeterOhlson not accordning to: http://www.danielsoper.com/statcalc3/calc.aspx?id=30 and http://www.efunda.com/math/gamma/findgamma.cfm and http://www.miniwebtool.com/gamma-function-calculator/?number=2.44 ? – apelsinapa Mar 16 '13 at 22:08
  • @apelsinapa That's because that's the *gamma* function, and the jsFiddle is the *factorial* function. They are slightly different. On the jsFiddle page you'll see the gamma function definition inside of the factorial function. If you like, [here's the gamma function by itself](http://jsfiddle.net/5EGXG/1/) – Peter Olson Mar 16 '13 at 22:09
  • @PeterOhlson Haha, sorry I think it might be this beer =) By the way, its an far more elegant solution – apelsinapa Mar 16 '13 at 22:10
4

This is not a trivial problem. There is not a simple closed-form formula for the gamma function. That said, there are some numerical approximations that should suit your needs.

The following answer will be using a technique called Lanczos approximation. The formula is as follows:

enter image description here

where g is an arbitrarily chosen constant that controls how accurate the approximation will be. For larger g, the approximation will be more accurate. Ag(z) is defined thus:

enter image description here

The hardest part is finding Ag(z), since pn is also defined with a complicated formula dependent on g.

I can't take too much credit for the following code, since I am just writing a port of the Python program on the wikipedia page.

function gamma(n) {  // accurate to about 15 decimal places
  //some magic constants 
  var g = 7, // g represents the precision desired, p is the values of p[i] to plug into Lanczos' formula
      p = [0.99999999999980993, 676.5203681218851, -1259.1392167224028, 771.32342877765313, -176.61502916214059, 12.507343278686905, -0.13857109526572012, 9.9843695780195716e-6, 1.5056327351493116e-7];
  if(n < 0.5) {
    return Math.PI / Math.sin(n * Math.PI) / gamma(1 - n);
  }
  else {
    n--;
    var x = p[0];
    for(var i = 1; i < g + 2; i++) {
      x += p[i] / (n + i);
    }
    var t = n + g + 0.5;
    return Math.sqrt(2 * Math.PI) * Math.pow(t, (n + 0.5)) * Math.exp(-t) * x;
  }
}

and of course, by definition of the gamma function:

function factorial(n) {
  return gamma(n + 1);
}

You can see this in action on jsFiddle.

Peter Olson
  • 139,199
  • 49
  • 202
  • 242
2

Just to complete @apelsinapa answer to correct the calculation for an integer (we didn't get an integer solution when inputing an integer number).

@apelsinapa's great solution:

var g = 7;
var C = [0.99999999999980993, 676.5203681218851, -1259.1392167224028,771.32342877765313, -176.61502916214059, 12.507343278686905, -0.13857109526572012, 9.9843695780195716e-6, 1.5056327351493116e-7];

function gamma(z) {
    if (z < 0.5) return Math.PI / (Math.sin(Math.PI * z) * gamma(1 - z));
    else {
        z -= 1;

        var x = C[0];
        for (var i = 1; i < g + 2; i++)
        x += C[i] / (z + i);

        var t = z + g + 0.5;
        return Math.sqrt(2 * Math.PI) * Math.pow(t, (z + 0.5)) * Math.exp(-t) * x;
    }
}

And to get a correct answer for integer:

 function factorialOfNumber(number) {
    if (number % 1 != 0 || number<0){
        return gamma(number + 1);
    }
    else {
        if(number == 0) {
           return 1;
        }
        for(var i = number; --i; ) {
           number *= i;
        }
        return number;
    }
}     
vntrp
  • 21
  • 3
0

Here's a version I wrote a few years ago ... a bit messy but tested :)

var
    M_GAMMA = [76.18009172947146, -86.50532032941677, 24.01409824083091, -1.231739572450155, 0.1208650973866179e-2, -0.5395239384953e-5],
    M_GAMMAS = 6;

function gammaFn(x) // Modified to JavaScript from "Numerical Recipies in C" by P. Mainwaring
{
    var i = 0, n = ++x, tmp = x + 5.5, ser = 1.000000000190015;
    for (tmp -= (x + 0.5) * Math.log(tmp); i < M_GAMMAS; ++i) ser += M_GAMMA [i] / ++n;
    return Math.log(2.5066282746310005 * ser / x) - tmp;
}

function fact(x) { return x > 1 ? Math.exp(gammaFn(x)) : 1 }

function combin(n, k) { return (Math.exp(gammaFn(n) - gammaFn(n - k) - gammaFn(k)) + 0.5) | 0 } // Ms Excel COMBIN() n! / k!(n - k)!

n = 49; k = 6; alert(fact(n) + ' ' + fact(k) + ' ' + combin(n, k)); // Lottery odds! (13,983,816)

The gamma and gammaLn functions are then:

function gammaLn(x) { return gammaFn(--x) }

function gamma(x) { return Math.exp(gammaLn(x)) }

:-)

0

If you're just looking for the function to compute factorials of real numbers then you only need this bit of code from Lanczos approximation:

function = factorial(z) {

  var g = 7;
  var C = [0.99999999999980993, 676.5203681218851, -1259.1392167224028, 771.32342877765313, -176.61502916214059, 12.507343278686905, -0.13857109526572012, 9.9843695780195716e-6, 1.5056327351493116e-7];
  var x = C[0];

  for (var i = 1; i < g + 2; i++)
  x += C[i] / (z + i);

  var t = z + g + 0.5;
  return Math.sqrt(2 * Math.PI) * Math.pow(t, (z + 0.5)) * Math.exp(-t) * x;
}

Works great for negative numbers in addition to decimals.