Methode van Bairstow

De methode van Bairstow is een numerieke methode om de nulpunten van een reële veelterm van willekeurige graad te vinden. De methode werd voor het eerst beschreven door Leonard Bairstow in de appendix van diens boek Applied Aerodynamics (1920), en is naar hem genoemd. Een reële veelterm heeft wortels die reëel zijn of anders voorkomen in paren van complex toegevoegde getallen. De methode zoekt naar kwadratische factoren, die ofwel reducibel zijn en twee reële wortels hebben, ofwel irreducibel zijn en twee complex toegevoegde wortels hebben. Op deze manier maakt de methode alleen gebruik van reële getallen. Een alternatief is de methode van Müller.

Beschrijving van de methode

Om de wortels van de veelterm

$p(x)=\sum _{i=0}^{n}a_{i}x^{i}$

te vinden zoekt deze methode de twee coëfficiënten $u$ en $v$ van de kwadratische veelterm

$k(x)=x^{2}+ux+v$

waarvan de wortels ook wortels zijn van de veelterm $p(x)$. Dan kunnen enerzijds twee wortels bepaald worden en kan anderzijds de oorspronkelijke veelterm $p(x)$ in graad worden verlaagd door hem te delen door de kwadratische veelterm $k(x)$. Dit proces wordt herhaald tot een veelterm van graad 1 of 2 overblijft.

Deling van $p(x)$ door $k(x)$ geeft als quotiënt

$q(x)=\sum _{i=0}^{n-2}b_{i}x^{i}$

en als rest

$r(x)=cx+d$,

zodat

$p(x)=k(x)q(x)+r(x)=(x^{2}+ux+v)\left(\sum _{i=0}^{n-2}b_{i}x^{i}\right)+(cx+d)$

Het gaat er nu om waarden van $u$ en $v$ te vinden waarvoor $r(x)=0$. Daartoe wordt in een soortgelijke procedure als de methode van Newton-Raphson een verbetering $u+\Delta u$ en $v+\Delta v$ van $u$ en $v$ gezocht, waarvoor de rest bij deling door $k(x)+\Delta k(x)=x^{2}+(u+\Delta u)x+v+\Delta v$ gelijk is aan 0, of althans kleiner is dan $r(x)$. Dus:

$p(x)=(k(x)+\Delta k(x))(q(x)+\Delta q(x))=$

$=k(x)(q(x)+\Delta q(x))+\Delta k(x)(q(x)+\Delta q(x))$

$=k(x)q(x)+k(x)\Delta q(x)+\Delta k(x)q(x)+\Delta k(x)\Delta q(x)$,

zodat

$r(x)=p(x)-k(x)q(x)=k(x)\Delta q(x)+\Delta k(x)q(x)+\Delta k(x)\Delta q(x)$

$\approx k(x)\Delta q(x)+\Delta k(x)q(x)$

Door de veelterm $q(x)$ een tweede keer te delen door $k(x)$ kan de uitdrukking voor de restterm vereenvoudigd worden:

$r(x)\approx k(x)\Delta q(x)+\Delta k(x)(k(x)f(x)+g(x))=k(x)m(x)+(\Delta ux+\Delta v)(g_{1}x+g_{0})$

Omdat in de restterm $r(x)$ geen kwadratische term voorkomt, moet

$cx+d\approx (x^{2}+ux+v)(-\Delta u\,g_{1})+(\Delta u\,x+\Delta v)(g_{1}x+g_{0})$

Daaruit volgen de twee lineaire vergelijkingen:

$c=-u\Delta u\,g_{1}+\Delta u\,g_{0}+\Delta v\,g_{1}$

$d=-v\Delta u\,g_{1}+\Delta v\,g_{0}$

In matrixvorm:

${\begin{bmatrix}d\\c\end{bmatrix}}={\begin{bmatrix}g_{0}&-g_{1}v\\g_{1}&g_{0}-g_{1}u\end{bmatrix}}{\begin{bmatrix}\Delta v\\\Delta u\end{bmatrix}}$

Een verbeterde benadering wordt dus verkregen door de uitgangswaarden $u$ en $v$ te corrigeren met:

${\begin{bmatrix}\Delta v\\\Delta u\end{bmatrix}}={\begin{bmatrix}g_{0}&-g_{1}v\\g_{1}&g_{0}-g_{1}u\end{bmatrix}}^{-1}{\begin{bmatrix}c\\d\end{bmatrix}}$

De onbekenden $c,\,d,\,g_{1},\,g_{0},(b_{i})$ en $(f_{i})$ zijn functies van $u$ en $v$. Ze worden recursief bepaald via:

${\begin{aligned}b_{n-2}&=a_{n}&f_{n-4}&=b_{n-2}\\b_{n-3}&=a_{n-1}-ub_{n-2}&f_{n-5}&=b_{n-3}-uf_{n-4}\\b_{i}&=a_{i+2}-ub_{i+1}-vb_{i+2}&&\quad (i=n-4,\ldots ,0)\\&&f_{i}&=b_{i+2}-uf_{i+1}-vf_{i+2}\quad (i=n-6,\ldots ,0)\\c&=a_{1}-ub_{0}-vb_{1}&g_{1}&=b_{1}-uf_{0}-vf_{1}\\d&=a_{0}-vb_{0}&g_{0}&=b_{0}-vf_{0}\end{aligned}}$

Deze iteraties worden herhaald tot wanneer convergentie bereikt is. De methode kan op eenvoudige manier geprogrammeerd worden in de standaard programmeertalen.

Voorbeeld

Bepaal de wortels van de veelterm

$p(x)=6\,x^{5}+11\,x^{4}-33\,x^{3}-33\,x^{2}+11\,x+6$

Neem als eerste kwadratische veelterm $k(x)$ een veelterm met als coëfficiënten de drie hoogste macht-coëfficiënten van $p(x)$, genormaliseerd op de hoogste graad, dus:

$k(x)=x^{2}+{\tfrac {11}{6}}\,x-{\tfrac {33}{6}}$,

zodat:

$u={\frac {a_{n-1}}{a_{n}}}={\tfrac {11}{6}};\quad v={\frac {a_{n-2}}{a_{n}}}=-{\tfrac {33}{6}}$

De eerste stap geeft als betere benadering:

$v+\Delta v=-5{,}500000000000+5{,}460103215562=-0{,}039896784438$

$u+\Delta u=\ \ \,1{,}833333333333+1{,}145692735213=\quad 2{,}979026068546$

Berekening 1e stap 

Berekening

$b_{3}=6$

$b_{2}=11-{\tfrac {11}{6}}\times 6=0$

$b_{1}=-33+{\tfrac {11}{6}}\times 0+{\tfrac {33}{6}}\times 6=0$

$b_{0}=-33-{\tfrac {11}{6}}\times 0+{\tfrac {33}{6}}\times 0=-33$

$f_{1}=6$

$f_{0}=0-11\times 6=-11$

$c=11+{\tfrac {11}{6}}\times 33+{\tfrac {33}{6}}\times 0=71{,}5$

$d=6-{\tfrac {33}{6}}\times 33=-175{,}5$

$g_{1}=0+{\tfrac {11}{6}}\times 11+{\tfrac {33}{6}}\times 6=53{,}166666666667$

$g_{0}=-33-{\tfrac {33}{6}}\times 11=-93{,}5$

$-g_{1}v=-53{,}166666666667\times (-{\tfrac {33}{6}})=292{,}416666666667$

$g_{0}-g_{1}u=-93{,}5-53{,}166666666667\times {\tfrac {11}{6}}=-190{,}972222222222$

De vergelijkingen zijn dus:

${\begin{bmatrix}d\\c\end{bmatrix}}={\begin{bmatrix}-93{,}5&292{,}416666666667\\53{,}16666666667&-190{,}972222222222\end{bmatrix}}{\begin{bmatrix}\Delta v\\\Delta u\end{bmatrix}}$

De determinant van de matrix is:

$\det =2309{,}083333333333$,

zodat

$\Delta v=(190{,}972222222222\times 175{,}5-292{,}416666666667\times 71{,}5)/2309{,}083333333333=5{,}460103215562$

$\Delta u=(53{,}1667\times 175{,}5-93{,}5\times 71{,}5)/2309{,}083333333333=1{,}145692735213$

De resultaten van de opeenvolgende iteraties staan in de volgende tabel. De daarin genoemde stapgrootte is gelijk aan:

$s={\sqrt {\Delta v^{2}+\Delta u^{2}}}$

Voor de eerste stap is dus:

$s={\sqrt {5{,}460103215562^{2}+1{,}145692735213^{2}}}=5{,}579008780071$

De kolom "wortels" geeft de wortels van de kwadratische veelterm $k(x)=x^{2}+ux+v$ in de vorm $w_{12}=m\pm h$. Omdat

$u=-(w_{1}+w_{2})=-2m$ en $v=w_{1}w_{2}=m^{2}-h^{2}$

geldt:

$m=-{\tfrac {1}{2}}u$ en $h={\sqrt {{\tfrac {1}{4}}u^{2}-v}}$

Iteratiestappen van Bairstow

Nr	$u$	$v$	stapgrootte	wortels
0	1,833333333333	−5,500000000000	5,579008780071	−0,916666666667±2.517990821623
1	2,979026068546	−0,039896784438	2,048558558641	−1,489513034273±1,502845921479
2	3,635306053091	1,900693009946	1,799922838287	−1,817653026545±1,184554563945
3	3,064938039761	0,193530875538	1,256481376254	−1,532469019881±1,467968126819
4	3,461834191232	1,385679731101	0,428931413521	−1,730917095616±1,269013105052
5	3,326244386565	0,978742927192	0,022431883898	−1,663122193282±1,336874153612
6	3,333340909351	1,000022701147	0,000023931927	−1,666670454676±1,333329555414
7	3,333333333340	1,000000000020	0,000000000021	−1,666666666670±1,333333333330
8	3,333333333333	1,000000000000	0,000000000000	−1,666666666667±1,333333333333

Na acht iteraties wordt een kwadratische veelterm $k(x)$ gevonden met wortels −1/3 en −3, en dit binnen de gewenste nauwkeurigheid. Dit zijn dus ook oplossingen van de oorspronkelijke veelterm $p(x)$. De superlineaire snelheid van convergentie is merkbaar vanaf de vierde stap. Na deling van de oorspronkelijke veelterm $p(x)$ door de kwadratische veelterm

$k(x)=(x+3)(x+{\tfrac {1}{3}})$

vindt men als quotiënt een nieuwe veelterm:

$6x^{3}-9x^{2}-9x+6$

waarop de methode opnieuw kan worden toegepast.

Zwakkere punten van de methode

De methode van Bairstow erft de lokale kwadratische convergentiesnelheid van de Newton-Raphson methode, behalve indien wortels optreden met een multipliciteit groter dan 1. Dan daalt de snelheid tot lineair. Een tweede probleem kan optreden bij veeltermen met een oneven graad, die slechts één reële wortel hebben. Kwadratische factoren die in deze wortel een kleine functiewaarde hebben, hebben de neiging te divergeren naar oneindig. Deze neiging kan worden onderdrukt door in dat geval de veelterm uit te breiden met een extra reële wortel, waardoor het convergentiegedrag bevorderd wordt, zij het ten koste van de snelheid van convergentie.

Externe links

• Bairstow's Algorithm on Mathworld
• Numerical Recipes in Fortran 77 Online
• Example polynomial root solver (deg(P)⇐10) using Bairstow's Method
• LinBairstowSolve, an open-source C++ implementation of the Lin-Bairstow method available as a method of the VTK library