any Python equivalent of Math::Polynomial::Solve?

Alex Renelt renelt at web.de
Sun Feb 27 11:21:09 EST 2005


Just wrote:
> In article <1109485184.587035.176540 at l41g2000cwc.googlegroups.com>,
>  "Carl Banks" <invalidemail at aerojockey.com> wrote:

>>It should be pretty easy to set up a Numeric matrix and call
>>LinearAlgebra.eigenvalues.  For example, here is a simple quintic
>>solver:
>>
>>. from Numeric import *
>>. from LinearAlgebra import *
>>.
>>. def quinticroots(p):
>>.     cm = zeros((5,5),Float32)
>>.     cm[0,1] = cm[1,2] = cm[2,3] = cm[3,4] = 1.0
>>.     cm[4,0] = -p[0]
>>.     cm[4,1] = -p[1]
>>.     cm[4,2] = -p[2]
>>.     cm[4,3] = -p[3]
>>.     cm[4,4] = -p[4]
>>.     return eigenvalues(cm)
>>
>>
>>now-you-can-find-all-five-Lagrange-points-ly yr's,
> 
> 
> Wow, THANKS. This was the answer I was secretly hoping for... "Great 
> need for speed", no, not really, but this Numeric-based version is about 
> 9 times faster than what I translated from Perl code yesterday, so from 
> where I'm standing your version is blazingly fast...
> 
> Thanks again,
> 
> Just

in addition:
I'm writing a class for polynomial manipulation. The generalization of 
the above code is:

definitions:
1.) p = array([a_0, a_i, ..., a_n]) represents your polynomial
P(x) = \sum _{i=0} ^n a_i x^i

2.) deg(p) is its degree

3.) monic(p) makes P monic, i.e. monic(p) = p / p[:-1]

then you get:
from numarray import *
import numarray.linear_algebra as la

def roots(p):
	p = monic(p); n = deg(p)
	M = asarray(zeros((n,n)), typecode = 'f8')
# or 'c16' if you need complex coefficients
	M[:-1,1:] = identity(n-1)
	M[-1,:] = -p[:-1]
	return la.eigenvalues(M)

Alex



More information about the Python-list mailing list