What's up with this code?

[Gregory Petrosyan]

Hello, it's me again.
I am trying to optimise small module for working with polynomials, and
I have encountered problem: new, optimised, code, doesn't work in some
specific case. Here's the old version of the code:

(3x^2 + 2x + 1 <=> poly([3, 2, 1]), btw)

def __add__(self, other):
"Return sum of two polynomials"
rcoefs1 = list(reversed(s elf.coefs))
rcoefs2 = list(reversed(o ther.coefs))
coefs = [c1+c2 for c1,c2 in zip(rcoefs1, rcoefs2)]
minlen = min(len(rcoefs1 ), len(rcoefs2))
coefs.extend(rc oefs1[minlen:] + rcoefs2[minlen:])
return Polynomial(reve rsed(coefs))

Here's the new version:

def __add__(self, other):
rcoefs1 = reversed(self.c oefs)
rcoefs2 = reversed(other. coefs)
coefs = [c1+c2 for c1,c2 in it.izip(rcoefs1 , rcoefs2)]
coefs.extend(it .chain(rcoefs1, rcoefs2)) #? -- here is magic
if coefs[-1] != 0:
return Polynomial(reve rsed(coefs), False, False)
else:
return Polynomial(reve rsed(coefs), False)

Last four lines are ok (tested), "magic" seems to occure in #? line.
What's the magic? Here it is:

poly([1,2,1]) + poly([1,-2,1]) is ok,
poly([1,-2,1]) + poly([5,4,3,2,1]) is ok,
poly([1]) + poly([2]) is ok,
poly([1,2,1]) + poly([5,4,3,2,1]) is ok, but

poly([1,0]) + poly([-1]) gives poly([-1]),
poly([1,-3]) + poly([-2]) gives poly([-5]) etc.

So, "magic" happens when processing something like poly([a,b]) +
poly([c]) and it gives poly([b+c]) as result. Can anybody help me with
it?

Regards, Gregory.

[bearophileHUGS@lycos.com]

Re: What's up with this code?

Gregory Petrosyan:[color=blue]
> coefs.extend(it .chain(rcoefs1, rcoefs2)) #? -- here is magic[/color]

Can't you just do a couple of extend? Something like:

coefs.extend(rc oefs1)
coefs.extend(rc oefs2)

This looks simpler and probably faster too.

Bye,
bearophile

[Gerard Flanagan]

Re: What's up with this code?

Gregory Petrosyan wrote:[color=blue]
> Hello, it's me again.
> I am trying to optimise small module for working with polynomials, and
> I have encountered problem: new, optimised, code, doesn't work in some
> specific case. Here's the old version of the code:
>
> (3x^2 + 2x + 1 <=> poly([3, 2, 1]), btw)
>
> def __add__(self, other):
> "Return sum of two polynomials"
> rcoefs1 = list(reversed(s elf.coefs))
> rcoefs2 = list(reversed(o ther.coefs))
> coefs = [c1+c2 for c1,c2 in zip(rcoefs1, rcoefs2)]
> minlen = min(len(rcoefs1 ), len(rcoefs2))
> coefs.extend(rc oefs1[minlen:] + rcoefs2[minlen:])
> return Polynomial(reve rsed(coefs))
>[/color]

Gregory

It seems that zip/izip will only generate values while all of the
iterable parameters have values, eg. if rcoefs1 has length 3 and
rcoefs2 has length 5, zip will have length 3. This is the point of the
lines

minlen = min(len(rcoefs1 ), len(rcoefs2))
coefs.extend(rc oefs1[minlen:] + rcoefs2[minlen:])

in the old code - to pad out 'coefs' with the unused values in the
longer 'rcoef'. I think that's the source of the results you're
getting.

Here's an interactive session:
[color=blue][color=green][color=darkred]
>>> import itertools as it
>>> a = [1,2,3]
>>> b = [0,0,0,4,5]
>>> c = [ alpha+beta for alpha,beta in it.izip(a,b) ]
>>> c[/color][/color][/color]
[1, 2, 3][color=blue][color=green][color=darkred]
>>> inf = min( len(a), len(b) )
>>> inf[/color][/color][/color]
3[color=blue][color=green][color=darkred]
>>> a[inf:][/color][/color][/color]
[][color=blue][color=green][color=darkred]
>>> b[inf:][/color][/color][/color]
[4, 5][color=blue][color=green][color=darkred]
>>> c.extend( a[inf:] + b[inf:] )
>>> c[/color][/color][/color]
[1, 2, 3, 4, 5]

[...]
[color=blue]
> Here's the new version:
>
> def __add__(self, other):
> rcoefs1 = reversed(self.c oefs)
> rcoefs2 = reversed(other. coefs)
> coefs = [c1+c2 for c1,c2 in it.izip(rcoefs1 , rcoefs2)]
> coefs.extend(it .chain(rcoefs1, rcoefs2)) #? -- here is magic
> if coefs[-1] != 0:
> return Polynomial(reve rsed(coefs), False, False)
> else:
> return Polynomial(reve rsed(coefs), False)
>
> Last four lines are ok (tested), "magic" seems to occure in #? line.
> What's the magic? Here it is:
>
> poly([1,2,1]) + poly([1,-2,1]) is ok,
> poly([1,-2,1]) + poly([5,4,3,2,1]) is ok,
> poly([1]) + poly([2]) is ok,
> poly([1,2,1]) + poly([5,4,3,2,1]) is ok, but
>
> poly([1,0]) + poly([-1]) gives poly([-1]),
> poly([1,-3]) + poly([-2]) gives poly([-5]) etc.
>
> So, "magic" happens when processing something like poly([a,b]) +
> poly([c]) and it gives poly([b+c]) as result. Can anybody help me with
> it?
>
> Regards, Gregory.[/color]


I don't really understand your "magic" line.

Here's some code I was playing about with, it pads the shorter
coefficient list *before* applying izip - nothing distinguished about
it, but maybe instructive as a comparison to what you've got.

import itertools as it

class Polynomial(obje ct):
def __init__(self, *coeffs):
self.coefs = list(coeffs)
#print self.coefs

def __add__(self, other):
rcoefs1 = self.coefs[:]
rcoefs2 = other.coefs[:]
d1 = len(rcoefs1)
d2 = len(rcoefs2)
sup = max( d1, d2 )
#pad the shorter of the two lists with zeros
#so that it is the same length as the other
rcoefs1 += [0] * (sup - d1)
rcoefs2 += [0] * (sup - d2)
coefs = [c1+c2 for c1,c2 in it.izip(rcoefs1 , rcoefs2 )]
return Polynomial(*coe fs)

def __str__(self):
return str(tuple(self. coefs))

print
print Polynomial(1,2, 1) + Polynomial(1,-2,1)
a = Polynomial(1,-2,1)
b = Polynomial(5,4, 3,2,1)
c = a+b
print a
print b
print c
print Polynomial(1) + Polynomial(2)
print Polynomial(1,0) + Polynomial(-1)
print Polynomial(1,-3) + Polynomial(-2)

Gerard

[Scott David Daniels]

Re: What's up with this code?

Gregory Petrosyan wrote:[color=blue]
> Hello, it's me again.
> I am trying to optimise small module for working with polynomials, and
> I have encountered problem: new, optimised, code, doesn't work in some
> specific case. Here's the old version of the code:
>
> (3x^2 + 2x + 1 <=> poly([3, 2, 1]), btw)
>
> def __add__(self, other):
> "Return sum of two polynomials"
> rcoefs1 = list(reversed(s elf.coefs))
> rcoefs2 = list(reversed(o ther.coefs))
> coefs = [c1+c2 for c1,c2 in zip(rcoefs1, rcoefs2)]
> minlen = min(len(rcoefs1 ), len(rcoefs2))
> coefs.extend(rc oefs1[minlen:] + rcoefs2[minlen:])
> return Polynomial(reve rsed(coefs))
>
> Here's the new version:
>
> def __add__(self, other):
> rcoefs1 = reversed(self.c oefs)
> rcoefs2 = reversed(other. coefs)
> coefs = [c1+c2 for c1,c2 in it.izip(rcoefs1 , rcoefs2)]
> coefs.extend(it .chain(rcoefs1, rcoefs2)) #? -- here is magic
> if coefs[-1] != 0:
> return Polynomial(reve rsed(coefs), False, False)
> else:
> return Polynomial(reve rsed(coefs), False)
>
> Last four lines are ok (tested), "magic" seems to occure in #? line.
> What's the magic? Here it is:
>
> poly([1,2,1]) + poly([1,-2,1]) is ok,
> poly([1,-2,1]) + poly([5,4,3,2,1]) is ok,
> poly([1]) + poly([2]) is ok,
> poly([1,2,1]) + poly([5,4,3,2,1]) is ok, but
>
> poly([1,0]) + poly([-1]) gives poly([-1]),
> poly([1,-3]) + poly([-2]) gives poly([-5]) etc.
>
> "magic" happens when processing something like poly([a,b]) +
> poly([c]) and it gives poly([b+c]) as result. Can anybody help me with
> it?[/color]

OK, the trick is to think about how zip/izip _must_ work.
When the first list runs out first it stops. When the
second list runs out, it has already read ahead one element
of the first list.

def __add__(self, other):
shorter = reversed(self.c oefs)
longer = reversed(other. coefs)
if len(self.coefs) > len(other.coefs ):
shorter, longer = longer, shorter # insure shorter first
coefs = [c1+c2 for c1,c2 in it.izip(shorter , longer)]
if len(self.coefs) == len(other.coefs ):
while coefs and coefs[-1] == 0:
coefs.pop()
else:
coefs.extend(lo nger)
...

--Scott David Daniels
scott.daniels@a cm.org

[plahey@alumni.caltech.edu]

Re: What's up with this code?

Scott,

this was a really clever catch, but I don't agree with the solution.
The problem is your assumption of how zip/izip _must_ work. I don't
think there is any contract that states that they will always advance
the first argument first and the second argument second... that is an
implementation detail (the end result, ignoring iterator state, is the
same regardless of the order in which you advance things). Given this,
I think using izip is simply a bad idea.

How about something like this:

def __add__(self, other):
rcoefs1 = reversed(self.c oefs)
rcoefs2 = reversed(other. coefs)
minlen = min(len(rcoefs1 ), len(rcoefs2))
coefs = [rcoefs1.next()+ rcoefs2.next() for i in xrange(minlen)]
coefs.extend(rc oefs1)
coefs.extend(rc oefs2)
if coefs[-1] != 0:
return Polynomial(reve rsed(coefs), False, False)
else:
return Polynomial(reve rsed(coefs), False)

[Gregory Petrosyan]

Re: What's up with this code?

That's it! Thanks a lot!

[Gregory Petrosyan]

Re: What's up with this code?

Oh, by the way, is there any documentation about time of execution of
standart functions (is len(list) O(1) or O(n), etc) and is there any C
version of decimal?

[Scott David Daniels]

Re: What's up with this code?

Gregory Petrosyan wrote:[color=blue]
> Oh, by the way, is there any documentation about time of execution of
> standard functions?[/color]
Nope. The reason is (A) that's a lot of work, (B) there are too many
caveats in practice to make such information useful, (C) any such spec
might corner the possible implementations in a later version.
[color=blue]
> is len(list) O(1) or O(n), etc) and is there any C version of decimal?[/color]
len(list) is O(1), list[n] is O(1), and if you (and others) continue
using decimal and it gets popular, rest assured that it will get faster.
Also, dictionary access is \Theta(1) (but only because you run on
machines with quite finite limitations).

The general Python rule when writing is "only worry about fast enough"
and optimize only after the code works and it needs to be faster.
Further, to determine speed, measure (use the timeit module) rather than
guess. Nobody's intuition is great WRT performance, and you'd be
shocked at the number of hours people spend speeding up a chunk of code
that cannot possibly substantially affect an applications performance.

--Scott David Daniels
scott.daniels@a cm.org

[Gregory Petrosyan]

Re: What's up with this code?

Again, thanks a lot!
I'm optimising my polynomial module as hard as possible, because when
it used native floats/integers it was very fast, but in some specific
cases -- too inaccurate. So, I migrated to decimals, and module began
to be very slow. Current success is nearly 4x speedup, and now it seems
to me that decimals are the bottleneck. Is anybody working on C version
of decimals? Is it a really hard work? (I am thinking about
re-implementing decimals in C by myself)

[Kent Johnson]

Re: What's up with this code?

Gregory Petrosyan wrote:[color=blue]
> Again, thanks a lot!
> I'm optimising my polynomial module as hard as possible, because when
> it used native floats/integers it was very fast, but in some specific
> cases -- too inaccurate. So, I migrated to decimals, and module began
> to be very slow. Current success is nearly 4x speedup, and now it seems
> to me that decimals are the bottleneck. Is anybody working on C version
> of decimals? Is it a really hard work? (I am thinking about
> re-implementing decimals in C by myself)[/color]

Maybe gmpy would help?

http://gmpy.sourceforge.net/

Kent