Evaluating a polynomial on one of its roots produces "a not so small" number

37 Views Asked by At

I am dealing with the roots of the following polynomial:

$p(x) = 882x^{12} + 6867x^{11} + 23566x^{10} + 53265x^9 + 84904x^8 + 101863x^7 + 92002x^6 + 63175x^5 + 31576x^4 + 11089x^3 + 2302x^2 + 243x + 10$

The matlab function roots() outputs:

 -3.333333333333349e+00 + 0.000000000000000e+00i
 -5.002524795857115e-01 + 8.661787610509997e-01i
 -5.002524795857115e-01 - 8.661787610509997e-01i
 -4.998466746602097e-01 + 8.662777735465724e-01i
 -4.998466746602097e-01 - 8.662777735465724e-01i
 -5.001532474054582e-01 + 8.657728922922465e-01i
 -5.001532474054582e-01 - 8.657728922922465e-01i
 -4.997475983486135e-01 + 8.658721882479348e-01i
 -4.997475983486135e-01 - 8.658721882479348e-01i
 -1.666666666666945e-01 + 0.000000000000000e+00i
 -1.428571428571326e-01 + 9.525799160446818e-09i
 -1.428571428571326e-01 - 9.525799160446818e-09i

So we have six (3 pairs of complex conjugate roots) which are very close to each other, two other complex roots and 2 real roots. Now as far as I now problems can arise when the roots are very close to each other, what bothers me is the fact that evaluating the polynomial on all the roots using polyval(p, roots(p)) outputs:

  2.407688793226726e-06 + 0.000000000000000e+00i
 -1.117346215551152e-10 - 1.041549069213943e-10i
 -1.117346215551152e-10 + 1.041549069213943e-10i
 -9.761613739556196e-11 - 1.169686569824080e-10i
 -9.761613739556196e-11 + 1.169686569824080e-10i
 -1.019273554447864e-10 - 1.033395591321096e-10i
 -1.019273554447864e-10 + 1.033395591321096e-10i
 -1.032542940038184e-10 - 1.166284846476628e-10i
 -1.032542940038184e-10 + 1.166284846476628e-10i
 -2.664535259100376e-14 + 0.000000000000000e+00i
 -3.552713678800501e-15 + 6.035109749186890e-21i
 -3.552713678800501e-15 - 6.035109749186890e-21i

Is this result normal? The polynomial evaluated on the first root results in ~$2.41 · 10^{-06}$

Could someone explain me this result? Should I trust the numerical aproximation of this first root?

1

There are 1 best solutions below

2
On

You should read the IEEE floating point specification. THat's what Matlab probably works with. Floats are not real numbers --- they're a very discrete sampling of the reals. When you get an answer like this one, there's a good chance that the number you get has the property that the floats on either side of it, $x_L$ and $x_R$, have the property that $p(x_L)$ and $p(x_R)$ have opposite signs, so there HAS to be a zero of the polynomial between them. (At least to the degree that you trust the floating point computation of $p$ at those points.) And that means that your particular "root" is one of the two floats that best approximate the root that must be in that interval.