Numpy failing to properly square array

Question

I'm trying to map a simple quadratic function, where zs is a numpy array and R is a constant

Ns = -np.square(zs) + 2*zs*R+ 3*R**2

It works fine most of the time, but for some reason whenever I have the evaluation set up as following the code breaks:

>>>zs = np.array(range(80262,80268)
>>>R = 26756
>>>Ns = -np.square(zs) + 2*zs*R+ 3*R**2
>>>print Ns
array([    642108,    535095,    428080,    321063,    214044
       4295074319], dtype=int64)

That last value in the array should be 107023. Whenever I go above 80267, the squaring function breaks completely and starts giving me absolutely ridiculous answers. Is this just a data type error, or is something else going on here that I don't know about?

Answer 1

The trouble is that zs = np.array(range(80262,80268)) creates an array of int32 values.

np.square(zs) returns an array of the same datatype as zs and the final squared value in the array overflows the four bytes of memory it's been allocated.

You see that Ns = -np.square(zs) + 2*zs*R+ 3*R**2 has a datatype of int64 because NumPy has given this array more memory in order to accommodate the larger numbers. However, it's too late: you already have an overflowed value in np.square(zs).

To solve the issue, create zs using the np.int64 datatype:

zs = np.arange(80262, 80268, dtype=np.int64)

Be aware that the same problem will occur again if the numbers in zs get large enough!