[SciPy-user] help with precision for big numbers

[Johann Cohen-Tanugi]

Hello,
I am computing :

In [22]: for i in range(6):
    s=(1.+Toff/Ton)**i*sp.factorial(Non+Noff-i)/sp.factorial(Non-i)
    print "%.14g"%s
   ....:    
   ....:    
4.3585218122217e+42
9.7493251062853e+42
1.7917678573714e+43
2.5383377979428e+43
2.4658138608587e+43
1.2329069304293e+43

A colleague using GSL and C code with double precision and long double ( 
I am not sure whether he has a 64bit machine) obtained the following 
values :
4.3585218122216e+42
1.0131651581042e+43
1.9350541758386e+43
2.8488297588735e+43
2.8759614708627e+43
1.4943721368208e+43

Close but not identical...... I was wondering if there is a way to 
increase numerical accuracy within scipy, assuming the standard behavior 
is not optimal with this respect. Or any other thoughts about these 
discrepancies? Or some nifty tricks to recover lost precision by 
organizing the computation differently?

thanks in advance,
Johann