[6.61] Find more accurate polynomial roots
solve() may return false solutions to ill-conditioned polynomials. An ill-conditioned polynomial is one in
which small changes in the coefficients cause large changes in the roots. Since zeros() uses solve(),
zeros() can return the same incorrect solutions. This tip gives an example, shows how to check the
results and gives alternative methods for finding the better solutions. In general (and as usual), finding
polynomial roots is not trivial, and the theoretic best attainable accuracy might be worse than you
would expect, even for polynomials of relatively low degree.
For some examples in this tip, I define a 'best' solution for a root as x such that f(x+e) and f(x-e) are of
opposite sign, or that one or both are zero, and e is the least significant digit of the floating-point
mantissa. I will call this the sign test. There are better tests for root solutions, and I will discuss some of
those as well.
The 'best' value for a polynomial root depends on what you want to do with it. Some common criteria
for a root x are
1. x such that |f(x)| is as 'small' as possible. Ideally, |f(x)| would be zero, but this is unlikely with
the finite resolution of floating-point arithmetic and its round-off errors.
2. x such that |f(x)| is as small as required. For problems based on physical measurements, it
may be a waste of time to find the true minimum, since measurement error prevents such an
accurate solution, anyway.
3. The polynomial coefficients can be reconstructed as accurately as possible from the roots.
This is the same as saying that the errors in coefficients of the reconstructed polynomial are
minimized. The polynomial is reconstructed with the roots z
0
, z
1
, ... z
n
as f(x) =
(x-z
0
)(x-z
1
)...(x-z
n
)
4. Some other various conditions for polynomial roots are met.
The conditions in criteria 4 may include these properties of polynomials:
✟
i=1
n
z
i
=
a
n−1
a
n
✟
i>j
z
i
z
j
=
a
n−2
a
n
z
1
z
2
z
3
...z
n
=
(
−1
)
n
a
0
a
n
for this polynomial, with roots z
1
, z
2
, ... z
n
f
(
x
)
= a
n
x
n
+ a
n−1
x
n−1
+ ... + a
1
x + a
0
Or, we may want to use the roots to evaluate the polynomial in this form:
f
(
x
)
= a
n
x − z
1
x − z
2
...
(
x − z
n
)
The point is that different root-finding algorithms can return slightly different roots (or fail completely to
find all or any roots), so, as usual, you really need to know why you want the roots to get useful results.
Getting back to the failure of solve(), I'll use this polynomial as an example:
[1]
y1
(
x
)
= x
3
+ 4.217E17 $ x
2
− 3.981E20 $ x − 6.494E22
Define this polynomial in the Y= editor, then
zeros(y1(x),x)
returns these candidate solutions:
{-4.217E17, -141.81969643465, 0}
6 - 116
To Next Page
Loading...
Need help?
General
