Durand–Kerner method: Difference between revisions

Revision as of 21:07, 6 February 2025

Latest revision as of 05:14, 4 August 2025

(4 intermediate revisions by 4 users not shown)

Line 1:

{{Short description|Root-finding algorithm for polynomials}}

{{refstyle|date=November 2020}}

In [[numerical analysis]], the '''Weierstrass method''' or '''Durand–Kerner method''', discovered by [[Karl Weierstrass]] in 1891 and rediscovered independently by Durand in 1960 and Kerner in 1966, is a [[root-finding algorithm]] for solving [[polynomial]] [[equation (mathematics)|equations]].<ref name="Petković">{{cite book |last1=Petković |first1=Miodrag |title=Iterative methods for simultaneous inclusion of polynomial zeros |date=1989 |publisher=Springer |___location=Berlin [u.a.] |isbn=978-3-540-51485-5 |pages=31–32}}</ref> In other words, the method can be used to solve numerically the equation ''f''(''x''){{=}}0, where ''f'' is a given polynomial, which can be taken to be scaled so that the leading coefficient is&nbsp;1.

: ''f''(''x'') = 0,

where ''f'' is a given polynomial, which can be taken to be scaled so that the leading coefficient is&nbsp;1.

==Explanation==

Line 14 ⟶ 10:

: <math>f(x) = x^4 + ax^3 + bx^2 + cx + d</math>

for all ''x''. The known numbers ''a'', ''b'', ''c'', ''d'' are the [[coefficient]]s.

for all ''x''.

TheLet knownthe (potentially complex) numbers ''aP'', ''bQ'', ''cR'', ''dS'' arebe the [[coefficient]]sroots of this polynomial ''f''. Then

Let the (potentially complex) numbers ''P'', ''Q'', ''R'', ''S'' be the roots of this polynomial ''f''.

Then

: <math>f(x) = (x - P)(x - Q)(x - R)(x - S)</math>

Line 30 ⟶ 22:

So if used as a [[fixed point (mathematics)|fixed-point]] [[iteration]]

: <math>x_1 := x_0 - \frac{f(x_0)}{(x_0 - Q)(x_0 - R)(x_0 - S)},</math>

it is strongly stable in that every initial point ''x''₀ ≠ ''Q'', ''R'', ''S'' delivers after one iteration the root ''P'' = ''x''₁. Furthermore, if one replaces the zeros ''Q'', ''R'' and ''S'' by approximations ''q'' ≈ ''Q'', ''r'' ≈ ''R'', ''s'' ≈ ''S'', such that ''q'', ''r'', ''s'' are not equal to ''P'', then ''P'' is still a fixed point of the perturbed fixed-point iteration

delivers after one iteration the root ''P'' = ''x''₁.

Furthermore, if one replaces the zeros ''Q'', ''R'' and ''S''

by approximations ''q'' ≈ ''Q'', ''r'' ≈ ''R'', ''s'' ≈ ''S'',

such that ''q'', ''r'', ''s'' are not equal to ''P'', then ''P''

is still a fixed point of the perturbed fixed-point iteration

: <math>x_{k+1} := x_k - \frac{f(x_k)}{(x_k - q)(x_k - r)(x_k - s)},</math>

Line 43 ⟶ 29:

: <math>P - \frac{f(P)}{(P - q)(P - r)(P - s)} = P - 0 = P.</math>

Note that the denominator is still different from zero. This fixed-point iteration is a [[contraction mapping]] for ''x'' around ''P''.

This fixed-point iteration is a [[contraction mapping]]

for ''x'' around ''P''.

The clue to the method now is to combine the fixed-point iteration for ''P'' with similar iterations for ''Q'', ''R'', ''S'' into a simultaneous iteration for all roots.

the fixed-point iteration for ''P'' with similar iterations

for ''Q'', ''R'', ''S'' into a simultaneous iteration for all roots.

Initialize ''p'', ''q'', ''r'', ''s'':

Line 66 ⟶ 48:

: <math>s_n = s_{n-1} - \frac{f(s_{n-1})}{(s_{n-1} - p_n)(s_{n-1} - q_n)(s_{n-1} - r_n)}.</math>

Re-iterate until the numbers ''p'', ''q'', ''r'', ''s'' essentially stop changing relative to the desired precision. They then have the values ''P'', ''Q'', ''R'', ''S'' in some order and in the chosen precision. So the problem is solved.

essentially stop changing relative to the desired precision.

They then have the values ''P'', ''Q'', ''R'', ''S'' in some order

and in the chosen precision. So the problem is solved.

Note that [[complex number]] arithmetic must be used, and that the roots are found simultaneously rather than one at a time.

and that the roots are found simultaneously rather than one at a time.

== Variations ==

This iteration procedure, like the [[Gauss–Seidel method]] for linear equations, computes one number at a time based on the already computed numbers. A variant of this procedure, like the [[Jacobi method]], computes a vector of root approximations at a time. Both variants are effective root-finding algorithms.

computes one number at a time based on the already computed numbers.

A variant of this procedure, like the [[Jacobi method]],

computes a vector of root approximations at a time.

Both variants are effective root-finding algorithms.

One could also choose the initial values for ''p'', ''q'', ''r'', ''s'' by some other procedure, even randomly, but in a way that

by some other procedure, even randomly, but in a way that

* they are inside some not-too-large circle containing also the roots of ''f''(''x''), e.g. the circle around the origin with radius <math>1 + \max\big(|a|, |b|, |c|, |d|\big)</math>, (where 1, ''a'', ''b'', ''c'', ''d'' are the coefficients of ''f''(''x''))

and that

* they are not too close to each other,

which may increasingly become a concern as the degree of the polynomial increases.

as the degree of the polynomial increases.

If the coefficients are real and the polynomial has odd degree, then it must have at least one real root. To find this, use a real value of ''p''₀ as the initial guess and make ''q''₀ and ''r''₀, etc., [[complex conjugate]] pairs. Then the iteration will preserve these properties; that is, ''p''_''n'' will always be real, and ''q''_''n'' and ''r''_''n'', etc., will always be conjugate. In this way, the ''p''_''n'' will converge to a real root ''P''. Alternatively, make all of the initial guesses real; they will remain so.

Line 183 ⟶ 155:

In the [[quotient ring]] (algebra) of [[residue class]]es modulo &fnof;(''X''), the multiplication by ''X'' defines an [[endomorphism]] that has the zeros of &fnof;(''X'') as [[eigenvalue]]s with the corresponding multiplicities. Choosing a basis, the multiplication operator is represented by its coefficient matrix ''A'', the [[companion matrix]] of &fnof;(''X'') for this basis.

Since every polynomial can be reduced modulo &fnof;(''X'') to a polynomial of degree ''n''&nbsp;&minus;&nbsp;1 or lower, the space of residue classes can be identified with the space of polynomials of degree bounded by ''n''&nbsp;&minus;&nbsp;1. A problem-specific basis can be taken from [[Lagrange interpolation]] as the set of ''n'' polynomials

A problem specific basis can be taken from [[Lagrange interpolation]] as the set of ''n'' polynomials

:<math>b_k(X)=\prod_{1\le j\le n,\;j\ne k}(X-z_j),\quad k=1,\dots,n,</math>

Line 220 ⟶ 191:

: <math>\max_{1 \le k \le n} |w_k| \le \frac{1}{5n} \min_{1 \le j < k \le n} |z_k - z_j|,</math>

then this inequality also holds for all iterates, all inclusion disks <math>D\big(z_k + w_k, (n - 1) |w_k|\big)</math> are disjoint, and linear convergence with a contraction factor of 1/2 holds. Further, the inclusion disks can in this case be chosen as

are disjoint, and linear convergence with a contraction factor of 1/2 holds. Further, the inclusion disks can in this case be chosen as

: <math>D\left(z_k + w_k, \tfrac14 |w_k|\right),\quad k = 1, \dots, n,</math>

Line 245 ⟶ 215:

* [[Victor Pan]] (May 2002): [https://web.archive.org/web/20060907205721/http://www.cs.gc.cuny.edu/tr/techreport.php?id=26 ''Univariate Polynomial Root-Finding with Lower Computational Precision and Higher Convergence Rates'']. Tech-Report, City University of New York

* {{cite journal|first= Arnold|last= Neumaier|title= Enclosing clusters of zeros of polynomials|journal= Journal of Computational and Applied Mathematics|volume= 156 |year=2003|issue= 2|url=https://www.mat.univie.ac.at/~neum/papers.html#polzer|doi= 10.1016/S0377-0427(03)00380-7|pages= 389–401|bibcode= 2003JCoAM.156..389N|doi-access= free}}

* Jan Verschelde, ''[httphttps://www2.math.uic.edu/~jan/mcs471f03/Project_Two/proj2/node2.html The method of Weierstrass (also known as the Durand–Kerner method)]'', 2003.

* Bernhard Reinke, Dierk Schleicher, and Michael Stoll, ``"[https://arxiv.org/abs/2004.04777 The Weierstrass root finder is not generally convergent]''", 2020

** Bernhard Reinke, Dierk Schleicher and Michael Stoll: "The Weierstrass-Durand-Kerner root finder is not generally convergent", Math. Comp. vol.92 (2023), pp.839-866. DOI: https://doi.org/10.1090/mcom/3783 .

==External links==