[Tutor] strange behavior of matrix**matrix

[Steven D'Aprano]

On Wed, Jul 20, 2016 at 10:27:50AM +0200, AB wrote:
> Hello
> 
> W dniu 2016-07-17 o 17:23, Steven D'Aprano pisze:
> >[...]
> >What result did you expect? 2**-1 as an int32 cannot be 0.5, as that's a
> >float.
> 
> I expected 0.5: as 2^(-1) is in fact 1/2, and as in Python 3 division of 
> two integers 1/2 produces float 0.5, I naturally expected the case of 
> arrays to be consistent with this behavior.

Ah, but the exponentiation operator ** is not the division operator / 
and is not guaranteed to give the same results.

It looks to me that numpy has array[ints]/int coerce to array[floats], 
but array[ints]**int remains an array[ints]. In that case, they have to 
choose between 2**-1 rounds down to 0 or rounds up to 1, and they chose 
rounding down.

Would I make the same decision? Probably not.


> >I'm not really sure about the rules that numpy uses for coercing from
> >one type to another, but I'm not surprised by this result. I don't know
> >if it is documented anywhere, but it seems like the sort of thing numpy
> >would do.
> 
> I'm only learning Python, so the behavior of version 3 is natural to me.

That's why it was changed :-)



> >Here's another similar example:
> >
> >py> np.array([0])**-1
> >__main__:1: RuntimeWarning: divide by zero encountered in power
> >__main__:1: RuntimeWarning: invalid value encountered in power
> >array([-2147483648])
> 
> In my eye it's not similar - 1/0 should always produce an error, so the 
> behavior is exactly as expected.

The point is that the numpy functions seem to me to be designed to be as 
fast as possible, not as correct as possible. There's no integer value 
to represent INFINITY, like for floats, and for some reason the numpy 
people didn't want to halt the calculation with an error, so they have 
to return *something*, and it has to be a 32-bit signed integer:

-2147483648 is 32 "1" bits (including the sign), so that makes a good 
error value, at least if you think like a C programmer.

[...]
> I thought that this may be a 'sin of omission', to be corrected in some 
> future versions of numpy/scipy.

If you report it as a bug on the numpy bug tracker, they will hopefully 
either explain why they think its not a bug, or fix it.



-- 
Steve