[Numpy-discussion] roots and high-order polynomial
Mon Jul 6 10:13:10 CDT 2009
On Mon, 06 Jul 2009 16:53:42 +0200
Fabrice Silva <firstname.lastname@example.org> wrote:
> Le lundi 06 juillet 2009 à 08:16 -0600, Charles R Harris
>a écrit :
>> Double precision breaks down at about degree 25 if
>>things are well
>> scaled, so that is suspicious in itself. Also, the
>> isn't Hermitean so that property of the roots isn't
>>preserved by the
>> algorithm. If it were Hermitian then eigh would be used
>> eig. That said, there are other ways of computing
>> that might preserve the Hermitean property, but I don't
>>know that that
>> would solve your problem.
> I think there is a misunderstanding: I was referring to
>the fact the
> solution had to be real or complex conjugate, and not
>that the companion
> matrix would be a hermitian matrix (which it isn't due
IIRC, the coefficients of your polynomial are complex.
So, you cannot guarantee that the roots are complex
> Something I forgot to tell is that the roots are complex
> eigenfrequencies of a physical system, the real and
> expressing the damping and the frequency, respectively.
>If a positive
> frequency is solution then the opposite negative is
>solution too but
> with the «same-signed» damping.
>> The problem is floating point round off error in
>> coefficients. Even if you know the coefficients exactly
>> generally be represented exactly in double precision.
>> computational roundoff error just adds to that. If the
>> were all integers I would have more confidence in the no
>> Where do the coefficients come from? Perhaps there are
> Here is the construction:
> the coefficients are obtained from a modal analysis of a
>subsystem of a
> bigger system. A quantity called impedance depending of
>a variable X is
> the result of the combination of several terms:
> Z(X) = C1/(X-X1)+C1*/(X-X1*)+...+CN/(X-XN)+CN*/(X-XN*)
> where * denotes the complex conjugate.
> I have to get the solutions X of an equation
>Ze(X)=N(X)/D(X) with N and
> D are very-low order polynomial (orders 1 and 2). An
> using iterative root finding (for example with
> to expand the equation as a polynomial of order close to
> I do understand this approach might seem non-obvious but
>it is rather
> efficient for a low value of N (<10)...
> Laboratory of Mechanics and Acoustics - CNRS
> 31 chemin Joseph Aiguier, 13402 Marseille, France.
More information about the NumPy-Discussion