[Numpy-discussion] optimizing power() for complex and real cases

[Tim Hochberg]

Gary Ruben wrote:

> Tim Hochberg wrote:
> <big snip>
>
>> As I've been thinking about this some more, I think the correct thing 
>> to do is not to mess with  the power ufuncs at all. Rather in 
>> x.__pow__ (since I don't know that there's anywhere else to do it), 
>> after the above checks check the types of the array and in the cases 
>> where the first argument is a float or complex and the second 
>> argument is some sort of integer array. This would be dispatched to 
>> some other helper function instead of the normal pow_ufunc. In other 
>> words, optimize:
>>
>> A**2, A**2.0, A**(2.0+0j), etc
>>
>> and
>>
>> A**array([1,2,3])
>>
>> but not
>>
>> A**array[1.0, 2.0, 3.0]
>>
>> I think that this takes care of the optimization slowing down power 
>> for general floats and optimizes the only array-array case that 
>> really matter.
>
>
> I think this might still be a tiny bit dangerous despite not observing 
> monotonicity problems and would be a bit more conservative and change 
> it to:
>
> optimize:
>
> A**2, A**(2+0j), etc


I'm guessing here that you did not mean to include (2+0j) on both lists 
and that, in fact you wanted not to optimize on complex exponents:
.
So, optimize:

A**-1, A**0, A**1, A**2, etc.

>
> and
>
> A**array([1,2,3])
>
> but not
>
> A**array[1.0, 2.0, 3.0], A**2.0, A**(2.0+0j)


That makes sense. It's safer and easier to explain: "numpy optimizes 
raising matrices (and possibly scalars) to integer powers)". The only 
sticking point that I see is if David is still interested in optimizing 
A**0.5, that's not going to mesh with this. On the other hand, perhaps 
he can be persuaded that sqrt(A) is just as good. After all, it's only 
one more character long ;)

-tim