Iterating over np.arange gives incorrect results starting from i = 10

import numpy as np
from math import pi, sqrt, exp, factorial

def pyas(x, sred):
    return sred**x*exp(-sred)/factorial(x)

for i, j in zip(np.arange(20), range(20)):
    print(i, j, pyas(i,10), pyas(j,10))

The question is, why do the values of the pyas() function start to differ? The output:

0 0 4.5399929762484854e-05 4.5399929762484854e-05
1 1 0.00045399929762484856 0.00045399929762484856
2 2 0.0022699964881242427 0.0022699964881242427
3 3 0.007566654960414142 0.007566654960414142
4 4 0.018916637401035354 0.018916637401035354
5 5 0.03783327480207071 0.03783327480207071
6 6 0.06305545800345118 0.06305545800345118
7 7 0.09007922571921599 0.09007922571921599
8 8 0.11259903214901998 0.11259903214901998
9 9 0.1251100357211333 0.1251100357211333
10 10 0.01764133335640144 0.1251100357211333
11 11 0.001382752728810601 0.11373639611012118
12 12 -6.89413134691794e-05 0.09478033009176766
13 13 9.595669338819968e-06 0.07290794622443666
14 14 1.4396571374678832e-07 0.05207710444602619
15 15 -5.313583114618023e-08 0.03471806963068413
16 16 4.0683489446471356e-09 0.021698793519177577
17 17 2.0030859032126935e-10 0.012763996187751515
18 18 -1.0541774694098002e-11 0.007091108993195286
19 19 -7.394475413970681e-13 0.0037321626279975192

CodePudding user response：

The difference is that np.arange is iterating fixed size np.int64 values while range is iterating python's unbounded int. We can see this by choosing the highest value in your range: 19.

>>> foo=10**np.int64(19)
>>> bar=10**19
>>> type(foo), type(bar)
(<class 'numpy.int64'>, <class 'int'>)
>>> foo, bar
(-8446744073709551616, 10000000000000000000)

The numpy version overflowed and wrapped to a negative number with a smaller absolute value. The python int is not limited to 64 bits and can expand indefinately.