[Numpy-discussion] Unpleasant behavior with poly1d and numpy scalar multiplication

[Charles R Harris]

On Sat, Feb 13, 2010 at 1:11 AM, Fernando Perez <fperez.net at gmail.com>wrote:

> Mmh, today I got bitten by this again.  It took me a while to figure
> out what was going on while trying to construct a pedagogical example
> manipulating numpy poly1d objects, and after searching for 'poly1d
> multiplication float' in my gmail inbox, the *only* post I found was
> this old one of mine, so I guess I'll just resuscitate it:
>
> On Tue, Jul 31, 2007 at 2:54 PM, Fernando Perez <fperez.net at gmail.com>
> wrote:
> > Hi all,
> >
> > consider this little script:
> >
> > from numpy import poly1d, float, float32
> > p=poly1d([1.,2.])
> > three=float(3)
> > three32=float32(3)
> >
> > print 'three*p:',three*p
> > print 'three32*p:',three32*p
> > print 'p*three32:',p*three32
> >
> >
> > which produces when run:
> >
> > In [3]: run pol1d.py
> > three*p:
> > 3 x + 6
> > three32*p: [ 3.  6.]
> > p*three32:
> > 3 x + 6
> >
> >
> > The fact that multiplication between poly1d objects and numbers is:
> >
> > - non-commutative when the numbers are numpy scalars
> > - different for the same number if it is a python float vs a numpy scalar
> >
> > is rather unpleasant, and I can see this causing hard to find bugs,
> > depending on whether your code gets a parameter that came as a python
> > float or a numpy one.
> >
> > This was found today by a colleague on numpy 1.0.4.dev3937.   It feels
> > like a bug to me, do others agree? Or is it consistent with a part of
> > the zen of numpy I've missed thus far?
>
> Tim H. mentioned how it might be tricky to fix. I'm wondering if there
> are any new ideas since on this front, because it's really awkward to
> explain to new students that poly1d objects have this kind of odd
> behavior regarding operations with scalars.
>
> The same underlying problem happens for addition, but in this case the
> answer (depending on the order of operations) changes even more:
>
> In [560]: p
> Out[560]: poly1d([ 1.,  2.])
>
> In [561]: print(p)
>
> 1 x + 2
>
> In [562]: p+3
> Out[562]: poly1d([ 1.,  5.])
>
> In [563]: p+three32
> Out[563]: poly1d([ 1.,  5.])
>
> In [564]: three32+p
> Out[564]: array([ 4.,  5.])  # !!!
>
> I'm ok with teaching students that in floating point, basic algebraic
> operations may not be exactly associative and that ignoring this fact
> can lead to nasty surprises.  But explaining that a+b and b+a give
> completely different *types* of answer is kind of defeating my 'python
> is the simple language you want to learn' :)
>
> Is this really unfixable, or does one of our resident gurus have some
> ideas on how to approach the problem?
>
>
The new polynomials don't have that problem.

In [1]: from numpy.polynomial import Polynomial as Poly

In [2]: p = Poly([1,2])

In [3]: 3*p
Out[3]: Polynomial([ 3.,  6.], [-1.,  1.])

In [4]: p*3
Out[4]: Polynomial([ 3.,  6.], [-1.,  1.])

In [5]: float32(3)*p
Out[5]: Polynomial([ 3.,  6.], [-1.,  1.])

In [6]: p*float32(3)
Out[6]: Polynomial([ 3.,  6.], [-1.,  1.])

In [7]: 3.*p
Out[7]: Polynomial([ 3.,  6.], [-1.,  1.])

In [8]: p*3.
Out[8]: Polynomial([ 3.,  6.], [-1.,  1.])

In [9]: p + float32(3)
Out[9]: Polynomial([ 4.,  2.], [-1.,  1.])

In [10]: float32(3) + p
Out[10]: Polynomial([ 4.,  2.], [-1.,  1.])

They are only in the removed 1.4 release, unfortunately. You could just pull
that folder and run them as a separate module. They do have a problem with
ndarrays behaving differently on the left and right, but __array_priority__
can be use to fix that. I haven't made that last fix because I'm not quite
sure how I want them to behave.

Chuck
-------------- next part --------------
An HTML attachment was scrubbed...
URL: <http://mail.python.org/pipermail/numpy-discussion/attachments/20100213/1de126ad/attachment.html>