Precision issue in python

Mark Dickinson dickinsm at gmail.com
Sat Feb 20 12:17:12 EST 2010


On Feb 20, 3:37 pm, mukesh tiwari <mukeshtiwari.ii... at gmail.com>
wrote:
> I don't know if is possible to import this decimal module but kindly
> tell me.Also a bit about log implementation

The decimal module is part of the standard library;  I don't know what
the rules are for SPOJ, but you're already importing the math module,
so at least *some* imports are obviously permitted.  Here's a quick
demonstration of how you might use the decimal module:  you can
probably find ways to tweak it for speed and accuracy.

from decimal import Decimal as D
from decimal import getcontext

getcontext().prec = 100   # working precision = 100 digits

pi =
D('3.14159265358979323846264338327950288419716939937510582097494459'
 
'2307816406286208998628034825342117067982148086513282306647093844')
half = D('0.5')
log = D.ln

def logfac(x):
    """Approximation to log(x!), using first few terms of Stirling's
series."""

    x = D(x)
    return log(2*pi)/2 + (x + half)*log(x) - x + \
        1/(12*x) - 1/(360*x**3) + 1/(1260*x**5)

def fac_first_digits(n, k):
    """Get first k decimal digits of n!."""

    log10_nfac = logfac(n)/log(D(10))
    frac = log10_nfac - int(log10_nfac)
    return int(10**(frac - 1 + k))


With the above code I get, for example (with Python 2.6):

>>> fac_first_digits(12345, 15)
344364246918678
>>> from math import factorial
>>> int(str(factorial(12345))[:15])
344364246918678

For small n, you'll need more terms of Stirling's series than I've
used above.  And there's always a small chance that the intermediate
errors involved in the computation might change those first 15
digits;  with a working precision of 100 and lots of terms of
Stirling's series it's a *very* small chance, but it's there.  If you
want something 100% watertight, you'd need to  compute strict upper
and lower bounds for the true result, and increase precision until
those bounds match sufficiently far that you can be sure of the first
k digits being correct.

--
Mark



More information about the Python-list mailing list