Faster Factorial (Python recipe) by Ivo Danihelka

ActiveState Code (http://code.activestate.com/recipes/577241/)

All you're testing is how fast python can do the first calculation and then return an element in a list 999 times.

By that same logic instead of the complicated process you're using you may as well combine both ideas to make the far simpler:

def mfactorial(n, _factcache={}):
  if n not in _factcache:
    _factcache[n] = math.factorial(n)
  return _factcache[n]

factorial:      0.00080 s
mfactorial:      0.00059 s
math.factorial: 0.25890 s

You are right that the test does not measure all usages. Your mfactorial() has a disadvantage. It does not store the computed factorial(1) .. factorial(n-1).

For what it's worth, Python 3.2 will include a much faster factorial function than earlier versions of Python (see Issue8692).

To compute factorials faster for large arguments, use Stirling's Approximation (http://en.wikipedia.org/wiki/Stirling%27s_approximation):

fN = float(N)
Log(N!) = N*math.log(fN) + math.log(2*math.pi*fN)/2 - fN + \
      (fN**-1)/12 - (fN**-3)/360 + (fN**-5)/1260 - (fN**-7)/1680 + (fN**-9)/1188

As N grows, the negative powers of N improve the approximation.
This approximation is so good that:

• For N > 11, the (fN^-9) term hits float's limit of precision.
• For N > 17, the (fN^-7) term hits float's limit of precision.
• For N > 32, the (fN^-5) term hits float's limit of precision.
• For N around 94, the (fN^-3) term hits float's limit of precision.
• For N around 670, the (fN^-1) term hits float's limit of precision.

When a term hits float's limit of precision, you can omit it from the approximation.

This program shows the accuracy of Stirling's Approximation:

import math

# Create output file.
sFile = 'C:\Documents\Desktop\Factorial_Float.tsv'
try:
    oFile = open(sFile, 'w')
except:
    print "Problem opening file " + sFile + ".\nGoodbye"
    quit()

# Column headers.
sOutput = "N\tExact\t0\t1\t2\t3\t4\t5\n"
oFile.write( sOutput )

# Calculate exact & approx values of Ln( N! ) for N = 1 to 1000.
# Long variables have prefix "i", float variables have prefix "f".

# Exact value of Ln( N! ).
fExact = 0
for iN in range(1,1000):
    # Float version of iN.
    fN = float(iN)
    fExact = fExact + math.log(fN)
    sOutput = str(iN) + "\t" + str(fExact)

    # f1 = 1st term of Stirling's Approx.
    f1 = math.log( 2*math.pi*fN )/2 + fN * math.log( fN ) - fN
    # fDiff always = difference between exact & approx values of Ln( N! )
    fDiff = f1-fExact
    sOutput = sOutput + "\t" + str(abs(fDiff))

    # f2 = 2nd term of Stirling's Approx.
    # dDelta always = next term in Stirling's Approx.
    fDelta = (fN**-1)/12
    f2 = f1 + fDelta
    fDiff = fDiff + fDelta
    sOutput = sOutput + "\t" + str(abs(fDiff))

    # f3 = 3rd term of Stirling's Approx.
    fDelta = -(fN**-3)/360
    f3 = f2 + fDelta
    fDiff = fDiff + fDelta
    sOutput = sOutput + "\t" + str(abs(fDiff))

    # f4 = 4th term of Stirling's Approx.
    fDelta = (fN**-5)/1260
    f4 = f3 + fDelta
    fDiff = fDiff + fDelta
    sOutput = sOutput + "\t" + str(abs(fDiff))

    # f5 = 5th term of Stirling's Approx.
    fDelta = -(fN**-7)/1680
    f5 = f4 + fDelta
    fDiff = fDiff + fDelta
    sOutput = sOutput + "\t" + str(abs(fDiff))

    # f6 = 6th term of Stirling's Approx.
    fDelta = (fN**-9)/1188
    f6 = f5 + fDelta
    fDiff = fDiff + fDelta
    sOutput = sOutput + "\t" + str(abs(fDiff))

    oFile.write( sOutput + "\n" )

oFile.close()

print "Goodbye"

David Lambert this method seems pretty accurate ,...but when i give value to N greater than 170 i get an overflowError

Example code:

import math

for N in xrange(1,10):
    fN=float(N)
    FactorialOfN=N*math.log(fN) + math.log(2*math.pi*fN)/2 - fN + (fN**-1)/12 \
    - (fN**-3)/360 + (fN**-5)/1260 - (fN**-7)/1680 + (fN**-9)/1188
    Result = math.e **FactorialOfN
    print "Factorial of %d" %N ," = " , Result


for N in xrange(10,190,10):
    fN=float(N)
    FactorialOfN=N*math.log(fN) + math.log(2*math.pi*fN)/2 - fN + (fN**-1)/12 \
    - (fN**-3)/360 + (fN**-5)/1260 - (fN**-7)/1680 + (fN**-9)/1188
    Result = math.e **FactorialOfN
    print "Factorial of %d" %N ," = " , Result

Example Error:

~$ ./Stirlingfactorial.py 
    Factorial of 1  =  1.00053439504
    Factorial of 2  =  2.00000108765
    Factorial of 3  =  6.00000004849
    Factorial of 4  =  24.0000000092
    Factorial of 5  =  120.000000004
    Factorial of 6  =  720.000000003
    Factorial of 7  =  5040.0
    Factorial of 8  =  40320.0
    Factorial of 9  =  362880.0
    Factorial of 10  =  3628800.0
    Factorial of 20  =  2.43290200818e+18
    Factorial of 30  =  2.65252859812e+32
    Factorial of 40  =  8.15915283248e+47
    Factorial of 50  =  3.04140932017e+64
    Factorial of 60  =  8.32098711274e+81
    Factorial of 70  =  1.197857167e+100
    Factorial of 80  =  7.15694570463e+118
    Factorial of 90  =  1.48571596448e+138
    Factorial of 100  =  9.33262154439e+157
    Factorial of 110  =  1.58824554152e+178
    Factorial of 120  =  6.68950291345e+198
    Factorial of 130  =  6.46685548922e+219
    Factorial of 140  =  1.34620124757e+241
    Factorial of 150  =  5.71338395644e+262
    Factorial of 160  =  4.71472363599e+284
    Factorial of 170  =  7.25741561531e+306
    Traceback (most recent call last):
      File "./Stirlingfactorial.py", line 16, in <module>
        Result = math.e **FactorialOfN
    OverflowError: (34, 'Numerical result out of range')

Can you suggest a way to get rid of this error :) ?

Thanks in advance :-)

i understand it's a really big number for the Computer... to represent it ... so i assume we can keep it in it's "log" representation.... for "mathematistic" manipulations and "stuff" ...