[SciPy-User] numpy array root operation

[Friedrich Romstedt]

Am 22.03.2012 um 03:43 schrieb Odin Den <hate_pod at yahoo.com>:

> Hi,
> 5th root of -32 can be computed correctly as follows:
>>>> -32**(1/5)
>>>> -2.0

Warning, mathematician (physicist to be precise) speaking. 

Additional to the operator precedence issue pointed out by David, notice that powers to non-integer numbers are cumbersome to define. In fact, let q be a rational number q = n/d, and v a complex number v = |v| exp(i phi), e.g. –42 with |v| = 42 and phi = pi. Then the root v^(1/d) is d-fold, and can be defined as the set of all (complex) numbers whose d-th power is v, namely the set {r exp(2 pi i f/d + i phi/d), f = 0...(d – 1)}, with a real positive number r such that r^d = |v|. The d-th power of each of this numbers is v. v^q is then just (v^(1/d))^n. 

For real v phi is either 0 or pi, so for positive v phi/d in the set equation will be zero, so there's always a positive root, which we call just "root" in daily language and in Python. For negative v, phi = pi, and (2 pi)/d is just twice that large, so there's no longer a positive root. But for odd d, the term 2 pi f/d + phi/d will be, for f = (d – 1)/2, just pi, so there's then a negative root, amongst all this roots. For even d, there's neither a positive nor a negative root of v, but only d complex ones. 

Notice that taking the numbers in the set to the n-th power multiplies their angle with n. 

You can calculate the first root (f = 0) in Python and numpy by taking a complex number to the power, e.g. in Python (–42 + 0j) ** (0.2). But this will never be real, always complex for negative v and noninteger exponent. Notice that for numpy, the convention for phi is –pi < phi <= pi. Negative numbers have phi = pi. 

The "first" root can be defined for any kind of exponent, but for irrational numbers as exponent there will be an infinite set of roots, AISI. But don't bet on that. It cannot be finite because then the irrational number would have been rational. :-)

As I warned you, mathematically inclined people can speak an hour about the most obvious thing to CS people, but normally I found it worth doing it once and then just keeping in mind that there is no such thing as natural intuition. 

cu
Friedrich

> However, when I try to perform same operation on numpy arrays I get
> the following:
>>>> array([-32])**(1/5)
>>>> array([ nan])

So ya, it's Not A Number, but A Set Of Numbers. :-)

Maybe it could spit out the complex first root instead, as it will do for numpy.asarray(–42 + 0j) ** 0.2. But I'm not involved enough to be knowledgable about the design here. 

> Is there anyway to compute roots of numpy arrays? I have a huge matrix which
> contains both negative and positive values. What is the easiest way of making
> python compute the "nth" roots of each element of this matrix without knowing
> the value of "n" a priory?

The other proposed approach, by using just the modulus and keeping the sign, is, as I pointed out above, not in all cases mathematically valid. I would guess you screwed up your model if you ran across taking fractional powers of negative numbers.