Methods of polynomial scaling

Question

This is a followup on this question Polynomial scaling for numeric root finding (warning, very large question; btw, half-band is $\frac12\textrm{sinc}{\frac12n}$). Among others, I have also seen this.

Dealing with the polynomial scaling problem (while knowing nothing about it), I thought maybe I should reverse the problem, so I imposed some known roots (r), constructed the polynomial (p), then scaled the roots (r/10), and then constructed another polynomial from these scaled roots (q). Here's how it looks in Octave:

r=[1 2 3 4];
p=poly(r)
q=poly(r/10)
s=roots(q)
s=s*10

And this is the result:

p = 1  -10   35  -50   24
q = 1.0000000  -1.0000000   0.3500000  -0.0500000   0.0024000
s = 0.40000   0.30000   0.20000   0.10000
s = 4.0000   3.0000   2.0000   1.0000

Despite q having very different coefficients, the roots s can be scaled back to form the originals, r. I suppose it's no surprise, since I am dividing by a factor with the power of x of every coefficient. Then I tried with logarithm:

q=poly(log(r))
s=roots(q)
s=exp(s)

And this also works:

q = 1.00000  -3.17805   3.24541  -1.05566   0.00000
s = 1.38629   1.09861   0.69315   0.00000
s = 4.0000   3.0000   2.0000   1.0000

Which brings me to my question: is there a relation between the ratio, or the coefficients of the two polynomials, p and q, and the method of scaling, which allows pre-scaling a polynomial, finding the scaled polynomial's roots, and then scaling back these roots to find the original roots? As mentioned in the linked question, this is for a numerical root finder which fails when very large differences between coeeficients cause ill-formed polymonials. The roots are not known, they need to be found.

Here's an example which is closer to what I have. That's two complex conjugate roots, and two real roots which are light years apart from themselves, and from the rest:

r=[-1-i -1+i 1e16 1e-16]
p=poly(r)
q=poly(log(r))
s=roots(q)
s=exp(s)

Notice how p can get pretty scary in terms of machine precision, while q looks inoffensive, yet the roots come out just fine:

p = 1.0000e+00  -1.0000e+16  -2.0000e+16  -2.0000e+16   2.0000e+00
q = 1.00000     -0.69315  -1351.61415    940.79891  -7698.20774
s = -36.84136 -  0.00000i   36.84136 -  0.00000i    0.34657 -  2.35619i    0.34657 +  2.35619i
s = 1.0000e-16 - 0.0000e+00i   1.0000e+16 - 0.0000e+00i  -1.0000e+00 - 1.0000e+00i  -1.0000e+00 + 1.0000e+00i

I realize this question (and the linked one) are not quite so mathematically-friendly, since my aim is numerical root finding, so if moderators decide this should be moved somewhere else, so be it (and I'm sorry for taking up your time).

score 1 · Answer 1 · answered Mar 30 '20 at 09:25

It turns out that scaling the polynomial is not practical:

a) the roots can be all scaled (down), which results in the first case in the OP, where the polynomial's coefficients drop proportionally with $x^n$, so they end up less than machine precision for $x^n$.

b) the roots can be all scaled by some function, say logarithm, which results in the 2nd case where for low orders the polynomial does look tamed, but for high orders, the coefficients "blow up" to the point where Octave says it cannot plot the coefficients because they exceed the allowed range. In fact, as early as 60th order, some go beyond 1e16.

I also tried gradually scaling the roots (r), to try to avoid overflowing the resulting coefficients, but that results in not knowing how to scale the found roots (s). Even so, in all the cases, the resulting polynomials seem to have wildly different relations with the original.

Methods of polynomial scaling

1 Answers1

Linked