Numerieke nulpunten - Deel 2
[ooO]
In het vorige artikel hebben we twee numerieke methoden behandeld voor het benaderen van nulpunten van functies, de bisectiemethode en de regula-falsimethode. In dit artikel bespreken we twee andere methoden, de secantmethode en de methode van Newton. Deze zijn in het algemeen sneller, maar komen niet voor alle startpunten op een nulpuntsbenadering uit.
De methoden uit Pythagoras 65–3 zijn veilig in de zin dat de benaderingen van het nulpunt in ieder geval in de buurt van het nulpunt blijven, en wel binnen het interval dat aan het begin is gekozen. In dit artikel bekijken we de secantmethode en de methode van Newton (ook wel de methode van Newton-Raphson genoemd). Deze lijken op de regula-falsimethode, want elke volgende benadering wordt bepaald door een rechte lijn met de $x$-as te snijden. Het verschil is dat het niet altijd zo hoeft te zijn dat het gezochte nulpunt daadwerkelijk steeds beter wordt benaderd. De methoden kunnen ontsporen. Hun kracht is echter dat ze, als de functie een beetje netjes is en je de startwaarden met beleid kiest, heel snel naar een nulpunt kunnen convergeren. We passen ze toe op dezelfde derdegraadsveelterm $f(x) = x^3 - 3x + 1$ als in het vorige artikel.
De seCantmethode
De secantmethode gebruikt een secant of snijlijn ('secare' betekent 'snijden'), een rechte lijn die de grafiek van de gegeven functie $f$ snijdt in (minstens) twee punten. De methode start met twee punten $(a, f(a))$ en $(b, f(b))$ in de buurt van het gezochte nulpunt. Het punt waar de secant door deze twee punten de $x$-as snijdt wordt genomen als een benadering van het nulpunt van $f$. Deze secant heeft richtingscoëfficiënt $\frac{f(b) – f(a)}{b – a}$ en als vergelijking:
| $y=f(a)+\frac{f(b)-f(a)}{b-a}(x-a)$. | $(\star)$ |
Het snijpunt $(c, 0)$ van de secant en de $x$-as volgt voor $y = 0$, met $f(a) \neq f(b)$:
| $c=a-f(a)\left(\frac{b-a}{f(b)-f(a)}\right)=\frac{af(b)-bf(a)}{f(b)-f(a)}$. | $(\star\star)$ |
Zie figuur 1. Als $f(c) = 0$, dan is het nulpunt gevonden. Anders wordt $(\star\star)$ opnieuw toegepast, en wel voor de twee punten $(b, f(b))$ en $(c, f(c))$. Enzovoort.
De formule $(\star\star)$ is dezelfde als in de regula-falsimethode, maar wordt anders gebruikt. De startpunten hoeven niet aan weerszijden van het nulpunt te liggen, en de volgende benadering is steeds gebaseerd op de laatste twee benaderingen ongeacht hun ligging ($a$ en $b$ krijgen de waarden van $b$ en $c$, in die volgorde). De gedachte erachter is dat de verbindingslijn tussen $(a, f(a))$ en $(b, f(b))$ op de grafiek van $f$ lijkt en doorgetrokken tot op de $x$-as in de buurt van een nulpunt zou moeten uitkomen. De secantmethode is dus vooral geschikt op intervallen waar de afgeleide van de functie weinig varieert en de secant dus veel op de grafiek van $f$ lijkt.
In het kader staat de secantmethode in pseudocode. Hierin is steeds $c$ de volgende benadering op basis van de op een na laatste benadering $a$ en de laatste benadering $b$ tot dan toe. De berekening kan worden gestopt als $|f(c)|$ en/of het verschil $|c - b|$ tussen de twee laatste benaderingen niet groter is dan een vooraf gekozen klein getal. De Python-programma's (te vinden via de knop [Bekijk oplossing] onder deze pagina) voeren een vooraf gekozen aantal iteraties uit. Ook moet je bij het programmeren rekening houden met de mogelijkheid dat de secant door $(a, f(a))$ en $(b, f(b))$ horizontaal kan lopen en de $x$-as niet snijdt, en de berekening dan laten stoppen (als $f(a) = f(b)$ wordt de noemer in $(\star\star)$ nul).
{Gegeven: functie f, startwaardenn a en b, stopcriterium}
fa = f(a)
fb = f(b)
c = (a*fb - b*fa)/(fb – fa)
fc = f(c)
Zolang het stopcriterium nog niet geldt
a krijgt de waarde van b
fa krijgt de waarde van fb
b krijgt de waarde van c
fb krijgt de waarde van fc
c = (a*fb - b*fa)/(fb – fa), voor fa ≠ fb
fc = f(c)
{Uitkomst: c, voldaan aan het stopcriterium}
Voor $f(x) = x^3 - 3x + 1$ en startwaarden $a = -0{,}9$ en $b = -0{,}4$ geeft figuur 2 een aantal iteraties van de secantmethode gericht op het middelste nulpunt.
Opgave 1In de secantmethode doet de volgorde van de startwaarden $a$ en $b$ ertoe. In de tweede keer dat $(\star\star)$ wordt toegepast, telt de startwaarde $b$, en niet $a$, als de op een na laatste benadering. Omkering van de waarden van $a$ en $b$ geeft een ander verloop van de iteraties. Ga dit na met de Python-programma's hieronder voor de secantmethode voor $f(x) = x^3 - 3x + 1$ met eerst de startwaarden $a = 0$ en $b = 1{,}5$ en daarna de startwaarden $a = 1{,}5$ en $b = 0$. Je zult zien dat in beide gevallen het rechter nulpunt wordt gevonden, maar dat er in het tweede geval grote uitschieters optreden en er meer iteraties nodig zijn om tot een goede benadering te komen. Je kiest de startwaarden $a$ en $b$ meestal zodanig dat $b$ dichter bij het nulpunt ligt dan $a$ – dus altijd eerst de grafiek schetsen of laten plotten om zicht te krijgen op het nulpunt voordat je de startpunten kiest. En wat gebeurt er voor de startwaarden $a = 0$ en $b = 2$? Voor $a = 2$ en $b = 0$ gaat het wel goed! |
De methode van Newton
Worden in de secantmethode de startpunten $a$ en $b$ heel dicht bij elkaar gekozen, dan valt de secant door de punten $(a, f(a))$ en $(b, f(b))$ bijna samen met een raaklijn aan de grafiek van de functie $f$. In de limiet voor $b$ gaat naar $a$ gaat de secant over in de raaklijn in het punt $(a, f(a))$. De methode van Newton neemt het nulpunt van deze raaklijn als benadering van het nulpunt van $f$.
De richtingscoëfficiënt in $(\star)$ wordt voor $b$ gaat naar $a$ gelijk aan $f'(a)$, de waarde van de afgeleide van $f$ in het punt $a$. De vergelijking van de raaklijn is $y = f(a) + f'(a)(x - a)$. Voor de methode van Newton is dus de afgeleide $f'$ van $f$ nodig. Het snijpunt $(c, 0)$ van de raaklijn en de $x$-as volgt voor $y = 0$:
| $c=N(a)=a-\frac{f(a)}{f'(a)}$. | $(\star\star\star)$ |
Zie figuur 3. Als $f(c) = 0$, dan is het nulpunt gevonden. Anders wordt $(\star\star\star)$ opnieuw toegepast in het punt $(c, f(c))$. Enzovoort.
In de andere drie methoden is de nulpuntsbenadering gebaseerd op de functiewaarden in twee verschillende punten, in de methode van Newton op de functiewaarde en de waarde van de afgeleide in een en hetzelfde punt.
De afgeleide van een veelterm
$p(x)=a_nx^n+a_{n-1}x^{n-1}+\cdots+a_2x^2+a_1x+a0$
is zoals bekend
$p'(x)=na_nx^{n-1}+(n-1)x^{n-2}+\cdots+2a_2x+a_1$.
In het kader hieronder staat de methode van Newton in pseudocode. Hierin is steeds $c$ de volgende benadering op basis van de vorige benadering $a$. De berekening kan worden gestopt als $|f(c)|$ en/of het verschil $|c - a|$ met de vorige benadering niet groter is dan een vooraf gekozen klein getal, of na een vast aantal iteraties. Ook moet je rekening houden met de mogelijkheid dat de raaklijn in $(a, f(a))$ horizontaal kan lopen en de $x$-as niet snijdt, en de berekening dan laten stoppen (als $f'(a)=0$ is de noemer in $(\star\star\star)$ nul).
{Gegeven: functie f en afgeleide f’, startwaarde a, stopcriterium}
fa = f(a)
dfdxa = f’(a)
c = a – fa/dfdxa, voor dfdxa ≠ 0
fc = f(c)
dfdxc = f’(c)
Zolang het stopcriterium nog niet geldt
a krijgt de waarde van c
fa krijgt de waarde van fc
dfdxa krijgt de waarde van dfdxc
c = a – fa/dfdxa, voor dfdxa ≠ 0
fc = f(c)
dfdxc = f’(c)
{Uitkomst: c, voldaan aan het stopcriterium}
Voor $f(x) = x^3 - 3x + 1$, $f'(x) = 3x^2 - 3$ en startpunt $a = -0{,}6$ geeft figuur 4 een aantal iteraties van de methode van Newton gericht op de middelste nulwaarde $r = 2\cos\left(\frac{4\pi}{9}\right)$.
In de buurt van een nulwaarde $r$ neemt de fout snel af. Als $f'$ niet veel verandert nabij $r$, dan kan worden afgeleid dat $c - a$ ongeveer gelijk is aan $r - a$. Elke volgende iteratiestap is dan van dezelfde orde van grootte als de fout tot dan toe. Tabel 1 laat dit zien voor de berekening van figuur 4.
Opgave 2Experimenteer met de Pythonprogramma's hieronder via de knop [Bekijk oplossing] voor de methode van Newton. Bereken een benadering van de $n$de-machtswortel $\sqrt[n]{a}$ van een positief getal $a$ door het nulpunt te bepalen van $f(x) = xn - a, f'(x) = nx^{n-1}$. Benader $r$ als nulwaarde van $f(x) = \sin(x$) met $f'(x) = \cos(x)$. Neem een waarde in de buurt van $3$; merk op dat je met een startwaarde dichtbij $\tfrac{\pi}{2}$ uit kunt komen op een veelvoud van $\pi$). |
In het algemeen is de methode van Newton sneller dan de secantmethode, die weer sneller is dan de regula-falsimethode. Je moet wel altijd met beleid te werk gaan. In opgave 2 zie je dat je soms niet het nulpunt vindt dat je dacht te vinden. Probeer de startwaarde zo te kiezen dat de raaklijn aldaar al min of meer in de richting van het nulpunt wijst. In sommige gevallen leidt de methode van Newton niet naar een nulpunt. We noemden al de horizontale raaklijn. Ook kan de berekening tussen punten heen en weer springen. Of de berekening kan ontsporen en uitkomsten geven die steeds verder van een nulpunt afwijken. Beide gedragingen treden op bij de functie in opgave 3.
Opgave 3Beschouw de functie $f(x) = \tan^{-1}(x)$, die afgeleide $f'(x)=\frac{1}{x^2+1}$ heeft. Gebruik de Python-programma's onder de knop [Bekijk oplossing] voor de methode van Newton met verschillende startwaarden. Er is een startwaarde $s \approx 1{,}39$ waarvoor de berekening heen en weer blijft springen tussen $s$ en $-s$; er geldt $N(\pm s) = \mp s$ ofwel $N(N(\pm s)) = \pm s$. Ga na dat $s$ voldoet aan $2s - (s^2 + 1)tan^{-1}(s) = 0$ en bepaal $s$ met de bisectiemethode. Voor startpunten op het interval $(-s, s)$ wordt de nulwaarde $r = 0$ gevonden. Voor startwaarden op $(-\infty, -s)$ en $(s,\infty)$ ontspoort de berekening. |
De besproken methoden zijn de eenvoudigste nulpuntsberekeningen. Er bestaan verfijndere methoden, maar we laten het hierbij. Ook komen we niet toe aan het diepgaand analyseren van het complexe, interessante gedrag van de methode van Newton voor veeltermen. We geven een doorkijkje in opgave 4.
Opgave 4Een nulpuntsberekening voor $f(x) = x^3 - 3x$ levert een benadering van $\sqrt{3}$. Gebruik de Pythonprogramma's onder knop [Bekijk oplossing] voor de methode van Newton met verschillende startwaarden. Pas op voor $f'(\pm 1) = 0$. Ga na dat voor de startpunten $x =\pm\sqrt{\tfrac{3}{5}}$ geldt dat $N(N(\pm\sqrt{\tfrac{3}{5}})) = \pm\sqrt{\tfrac{3}{5}}$. Startwaarden op $(-\infty, -1)$, $(-\sqrt{\tfrac{3}{5}},\sqrt{\tfrac{3}{5}})$ en $(1,\infty)$ leiden naar resp. de nulwaarden $-\sqrt{3}$, $0$ en $\sqrt{3}$. Op $(-1, -\sqrt{\tfrac{3}{5}})$ en $(\sqrt{\tfrac{3}{5}},1)$ heeft het gedrag zoals dat heet 'gevoelige afhankelijkheid van de beginvoorwaarde'. Daar wisselen intervallen waarin de iteraties naar $\sqrt{3}$ of naar $-\sqrt{3}$ gaan elkaar af. Vergelijk bijvoorbeeld de uitkomsten voor de startwaarden $x = \sqrt{\tfrac{3}{5}} + 0{,}0005$ en $\sqrt{\tfrac{3}{5}} + 0{,}001$, die heel dicht bij elkaar liggen. |
| $\color{white}{iteratie}$ | $\color{white}{c}$ | $\color{white}{f(c)}$ | $\color{white}{c-a}$ | $\color{white}{r-a}$ |
| $0$ | $-0{,}6$ | $2{,}583999999999$ | ||
| $1$ | $0{,}745833333333$ | $-0{,}822617259838$ | $1{,}345833333333$ | $0{,}947296355334$ |
| $2$ | $0{,}127880698688$ | $0{,}618449197500$ | $-0{,}617952634645$ | $-0{,}398536977999$ |
| $3$ | $0{,}337457743757$ | $0{,}026055690162$ | $0{,}209577045068$ | $0{,}219415656645$ |
| $4$ | $0{,}347259133854$ | $9{,}8197502691\times10^{-5\ }$ | $0{,}009801390098$ | $0{,}009838611577$ |
| $5$ | $0{,}347296354787$ | $1{,}4433275686\times10^{-9\ }$ | $3{,}7220932419\times10^{-5\ }$ | $3{,}7221479517\times10^{-5\ }$ |
| $6$ | $0{,}347296355334$ | $-2{,}2204460493\times10^{-16}$ | $5,4709720088\times10^{-10}$ | $5{,}4709731190\times10^{-10}$ |
Tabel 1 Afname van de fout in figuur 4
Bekijk oplossing