Why can't I raise the negative force in numpy?
I simulate the riemann theta function:
import numpy as np
def theta(s, n=100):
a_range = np.arange(2, n + 1)
return 1 + sum(1/(a_range ** s))
It doesn't work for the negative s
; for example theta(-2)
results in this error:
1 def theta(s, n=100):
2 a_range = np.arange(1)
----> 3 return 1 + sum(1/(a_range ** s))
4
5 theta(-2)
ValueError: Integers to negative integer powers are not allowed.
Why? x^-1
should be 1/x
if I name my math correctly.
source to share
In NumPy, the logic used to select the output dtype of a type operation a_range ** s
is based on dtypes, not values. This means it a_range ** -2
must have the same output dtype as a_range ** 2
.
It's important for something like to numpy.array([2]) ** 2
give integer output, so means it numpy.array([2]) ** -2
should give integers or nothing. They didn't choose anything; integers with negative integer powers is a bug in NumPy.
If you want floating point output, enter floating point input:
a_range = np.arange(2, n + 1, dtype=float)
or
a_range = np.arange(2, n + 1).astype(float)
There are a few strange aspects of NumPy type rules that you cannot expect from the above description. One is that, for operations using both scalars and arrays, the dtype scalar can actually be "demoted" based on its value, before the input data types are used to select the dtype result:
>>> (numpy.array([1], dtype='int8') + numpy.int32(1)).dtype
dtype('int8')
>>> (numpy.array([1], dtype='int8') + numpy.array([1], dtype='int32')).dtype
dtype('int32')
Here the scalar numpy.int32(1)
is being "demoted" to int8, but the arrays are not being promoted. (This is actually a bit more complicated than just downgrading to int8, especially for signed / unsigned handling; see more details.)
Second, when uint64s are involved, NumPy can suddenly end up okay with negative values:
>>> numpy.arange(5, dtype='uint64') ** -2
__main__:1: RuntimeWarning: divide by zero encountered in power
array([ inf, 1. , 0.25 , 0.11111111, 0.0625 ])
This is because NumPy cannot find an integer dtype large enough for both uint64 and negative values, so it refuses and forces the inputs to float. The same can be seen with a positive exponent of the signed dtype, if you avoid "downgrading" the scalar type:
>>> numpy.arange(5, dtype='uint64') ** numpy.array([2], dtype='int32')
array([ 0., 1., 4., 9., 16.])
source to share