Newton's Method for Roots of Polynomials

The standard way to use Newton's Method for finding a root of a polynomial $p(x)$ is to use the iteration formula $$x_{n+1}=x_n-{p(x)\over p'(x)}$$ I recently thought of a new way of finding the roots. Letting $x_1,...,x_n$ denote the $n$ roots and setting $p(x)=x^n+\cdots +a_1x+a_0$, we have the formulas $$x_1+\cdots +x_n=-a_{n-1}$$ $$x_1x_2+\cdots +x_{n-1}x_n=a_{n-2}$$ and so on to $$x_1x_2\cdots x_n=(-1)^na_0$$ This gives a non-linear system of $n$ variables and $n$ equations, which can be solved using the multivariable version of Newton's method. We would use the iteration $$\bf x_{n+1}=x_n-[Df(x_n)]^{-1}f(x_n)$$ where $$f(x_1,...,x_n)=\begin{bmatrix}x_1+\cdots +x_n+a_{n-1} \\ x_1x_2+\cdots +x_{n-1}x_n-a_{n-2} \\ \vdots \\ x_1\cdots x_n-(-1)^na_0 \end{bmatrix}$$ and $Df$ is the derivative matrix of $f$.

My question is: is there any advantage to using this method for finding roots? The obvious plus is that using the multivariable Newton's method provides all $n$ roots simultaneously while the regular method only gives them one at a time. However, I'm also interested in how the rate of convergence of this method would compare with the other, and if one would tend to be more numerically stable than the other.

• I'd guess the computatin of all the $n^2$ derivatives is prohibitive, let alone the solving of the linear system each step. There are lots of algorithms for finding roots of polynomials, it is a well-studied area. – vonbrand Apr 7 '14 at 0:39
• Offhand I would expect the numerical stability to be quite bad. Also note that the recursion would be a little more stable if you implemented it as $$Df(x_n)(x_n-x_{n+1})=f(x_n)$$ and solved for $x_{n+1}$. – Alex Becker Apr 7 '14 at 0:40
• Laguerre Method for Polynomia Roots is a 'safe' alternative. – Felix Marin Apr 7 '14 at 18:19
• Nice write-up. This is a better description of Durand-Kerner than the Wikipedia page. – Ryan Budney Feb 18 '16 at 7:25
• Possible duplicate of Vieta's theorem – J. M. is a poor mathematician Jun 4 '16 at 6:24

Congratulation, you have reinvented the Durand-Kerner method (1960) that was developed by Weierstraß in 1891.

It has quite satisfying convergence, but may show a rather slow initial phase until the root approximations are all close enough to the exact roots. (Once, there were Java applets, but since they are now virtually unusable, and even sometimes no longer to be found, that demonstrate that nicely, especially for symmetric configurations of roots and initial values.)

• Haha thanks for the answer... lol – user140943 Apr 7 '14 at 17:52
• A fine answer, for sure. And the prevalence of Java applets for mathematics demonstration is an ongoing annoyance for me, especially with the availability of HTML5 and Canvas to handle it all in-browser without any plugins. – Steven Stadnicki Apr 7 '14 at 18:47
• This method is implemented in the mpmath library, in Python. I find it quite effective for polynomials with distinct roots, of order less than about 800 or so. After that convergence appears to become very slow. – Ryan Budney Feb 18 '16 at 7:19

Newton or Quasi-Newton methods are a good method for solving systems of nonlinear equations such as polynomial rootfinding. These methods can be made robust by using line-search and trust-region methods to ensure convergence, and continuation techniques to prevent getting stuck in local minima. The full matrix of second derivatives does not need to be formed, only the ability to apply it to vectors, or approximations to it's inverse are needed. A good reference is chapter 11 of Nocedal and Wright's book Numerical Optimization.

For the special case of numerical root finding you can convert the problem to finding the eigenvalues of the companion matrix, to which fast stable methods for solving sparse eigenvalue problems can be used.

• Companion matrix methods for polynomial roots usually reduce to the Bernoulli iteration for the forwards power method and a variant of Jenkins-Traub for the shifted inverse power method. Also, Durand-Kerner has an interpretation in terms of the companion matrix using an adapted basis in the polynomial vector space. – LutzL Apr 7 '14 at 21:39