2

I've implemented the Jenkins Traub algorithm in c++ (Github repo). While the majority of the solutions work well, it seems that a small portion of the roots are unstable. Here is a figure showing polynomials evaluated at the computed roots (termed "residual"). This value should always be close to zero since we are finding polynomial roots, however, sometimes the value is unstable.

My intuition is that these unstable roots occur in polynomials with poor conditioning, and that it could be fixed by rescaling the polynomial somehow. Methods exist for balancing and scaling companion matrices, but does anything exist for the polynomial itself?

Thanks very much for the feedback, I realized a flaw in my experiments. For evaluating the Jenkins-Traub method I was taking the real value of all roots whereas for the companion matrix I was only consider roots where the imaginary component was 0. This has a large effect on the residuals. After correcting this error here are my updated results that I am now very happy with:

Considering the real values of roots (ignoring the imaginary component):

Considering only real roots (only roots where imaginary part == 0):

When the histogram starts on the y-axis it means that portion of roots evaluated to 0 within machine precision.

grg
  • 1,017
kip622
  • 131
  • 4
  • Doesn't Jenkins–Traub already scale the polynomial's coefficients by rescaling the variable, in an attempt to better condition the polynomial? – wltrup Aug 05 '15 at 00:09
  • The original RPOLY alg does some rescaling but it's not clear to me why the scaling is chosen that way. The original Jenkins Traub paper does not suggest a way to properly scale the polynomial.

    I'm looking for either a) a principled way to rescale coefficients or b) a reasoning as to why rescaling them in the way the original RPOLY scales is correct.

     – kip622 Aug 05 '15 at 03:34
  • Sorry, I'm afraid I don't have an answer to either (a) or (b). For (a), I tried to play with the idea of rescaling by attempting to set the sum of the roots to $N$ and/or the product of the roots to 1, where $N$ is the order of the polynomial, in order to see if I could have the scaled roots, in magnitude, be near 1 but I wasn't successful in getting anything substantial. As for (b), I haven't tried to understand the implementation you linked to. Perhaps someone else may be able to assist you. – wltrup Aug 05 '15 at 06:40
  • Please see this question. Maybe you wouldn't mind sharing the specifics of the polynomials you are having trouble with. – user5108_Dan Aug 05 '15 at 15:48

0 Answers0