7.9A. Nullstellensuche nach Newton

Wir haben früher bemerkt, daß zur Auffindung von Nullstellen einer gegebenen Funktion oft nur Näherungsverfahren helfen. Eine alte, aber wirkungsvolle Methode ist das Newton-Verfahren zur näherungsweisen Bestimmung von Nullstellen. Wie der Name sagt, stammt es von dem genialen Mathematiker und Physiker Isaak Newton, der unabhängig von Leibniz um 1670 die moderne Infinitesimalrechnung entwickelte. Die Grundidee basiert auf der Beobachtung, daß eine Tangente an die gegebene Kurve meist die x-Achse in der Nähe eines Schnittpunktes der Kurve mit der x-Achse schneidet, falls der Tangentialpunkt nicht allzu weit von der Nullstelle entfernt ist. Haben wir irgendwie nach k Schritten eine Näherung ak für die gesuchte Nullstelle konstruiert, so ergibt sich die x-Koordinate des Schnittpunktes der durch Ta

k

,1

f( x ) = f( ak ) + f ´( ak ) ( x − ak )

beschriebenen Tangente mit der x-Achse als Lösung der Gleichung 0 = f( ak ) + f ´( ak ) ( x − ak ) , und diese lautet an + 1 = ak −

f( ak ) f ´( ak )

(*)

sofern die Ableitung im Punkt ak nicht Null war. Es steht also zu hoffen, daß die durch (*) rekursiv definierte Folge (ak) gegen eine Nullstelle konvergiert, sofern der Startwert a0 nicht zu ungünstig gewählt wurde.

Dies bestätigt sich eindrucksvoll in Beispiel 1: Nullstellen des Polynoms f( x ) = x4 − x − 1. Es gibt, wie früher erwähnt, zwar monströse Formeln zur exakten Lösung von Gleichungen 4. Grades, aber auch das sehr versierte MAPLE streicht hier die Waffen. Und selbst wenn irgend ein anderes Programm die exakte Lösungsformel zustandebringen sollte - sie ist so unüberschaubar, daß niemand etwas damit anfangen kann. Versuchen wir es also näherungsweise - für die Praxis braucht man ohnehin Dezimaldarstellungen und nicht irgendwelche chaotischen Wurzelausdrücke. Stichprobenartig eingesetze Werte liefern zunächst  1 3 0 < f( −1 ) , f −  < 0 , f( 0 ) < 0 , f( 1 ) < 0 , 0 < f  ,  2 2 und der Zwischenwertsatz verrät uns, daß sowohl zwischen −1 und −

1

als auch zwischen 1 und

3

2 2 eine Nullstelle liegen muß. Es ist also vernünftig, diese Zahlen als Startwerte zu nehmen. Für den Startwert a0 = −1

berechnet MAPLE bei schrittweise verdoppelter Stellenzahl k = 0, ak = -1., f( ak ) = 1. k = 1, ak = -0.80, f( ak ) = 0.2096000000 k = 2, ak = -0.7312, f( ak ) = 0.01714043591 k = 3, ak = -0.72454848, f( ak ) = 0.0001425057630 k = 4, ak = -0.7244919629910424, f( ak ) = 0.1006055827 10-7 k = 5, ak = -0.72449195900051563148076470985776, f( ak ) = 0.5015091497 10-16 k = 6, ak = -0.7244919590005156115883722821870370600978455954423770171218710302, f( ak ) = 0.1246213427 10-32 Und für den Startwert b0 = 1 kommt heraus: k = 0, bk = 1., f( bk ) = -1. k = 1, bk = 1.3, f( bk ) = 0.8271604938 k = 2, bk = 1.236, f( bk ) = 0.09659632871 k = 3, bk = 1.2210590, f( bk ) = 0.001977477570 k = 4, bk = 1.220744225795051, f( bk ) = 0.8862015978 10-6 k = 5, bk = 1.2207440846057878724120649605478, f( bk ) = 0.1782394837 10-12 k = 6, bk = 1.220744084605759475361685350257557482895676462493199013748970373, f( bk ) = 0.7210194350 10-26

Mit Digits:=30 und fsolve liefert MAPLE die beiden Nullstellen auf 30 Stellen genau: -0.724491959000515611588372282187, 1.22074408460575947536168534911

Die Newton-Approximationen konvergieren so rasant gegen je eine Nullstelle, daß bereits nach drei Schritten die Tangenten im Bild sich nicht mehr von der Kurve unterscheiden lassen und die Zahl der korrekten Dezimalstellen nach dem Komma sich pro Schritt etwa verdoppelt. Das ist kein Zufall, wie die folgende raffinierte Anwendung des Mittelwertsatzes allgemein zeigt: Konvergenz des Newton-Verfahrens Ist f stetig differenzierbar und die Ableitung nirgends Null, so gibt es zu einer Nullstelle a von f und jeder Vorgabe ε > 0 ein δ > 0, so daß für Startwerte a0, die weniger als δ von a entfernt sind, alle Iterationen ak die folgende Abschätzung erfüllen: ak + 1 − a ≤ ak − a ε. Induktiv folgt dann ak − a ≤ a0 − a ε < δ ε . k

k

Falls f zweimal stetig differenzierbar ist, kann δ > 0 sogar so gewählt werden, daß die Iterationsfolgen ak "quadratisch" konvergieren, d.h. ak + 1 − a ≤ γ ak − a

2

mit einer geeigneten Konstante γ. Durch Logarithmieren zur Basis 10 wird daraus log( ak + 1 − a ) ≤ 2 log( ak − a ) + log( γ ) , und das ist im Wesentlichen (bis auf die additive Konstante log( γ )) die Aussage, daß sich die Anzahl der korrekten Nachkommastellen beim Schritt von ak zu ak + 1 ungefähr verdoppelt. So frappierend und effektiv das Newtonsche Verfahren ist - es klappt nicht immer.

Beispiel 2: Ein Halbkreis Es kann passieren, daß die Glieder der Iterationsfolge aus dem Definitionsbereich der Funktion hinausrutschen. Das ist zum Beispiel bei der Halbkreisfunktion 1 − x2 schon im ersten Schritt der Fall, wenn man nicht gerade mit einer der beiden Nullstellen startet.

Beispiel 3: Ein Flip-Flop Es kann auch vorkommen, daß man einen Startwert erwischt, von dem aus das Verfahren ständig zwischen zwei Werten hin und her springt, etwa bei der Funktion f( x ) = x3 − x mit der Ableitung f ´( x ) = 3 x2 − 1 und dem Startwert 1 a0 = . 5 Denn dann ergibt sich 2

2

2

2

2

5 a0 − 1 = 0 , 2 − 6 a0 = 1 − a0 , 2 a0 ( 1 − 3 a0 ) = a0 ( 1 − a0 ) , 3

a1 = a0 −

a0 − a0 2

3 a0 − 1

3

= −a0 und a2 = a1 −

a1 − a1 2

3 a1 − 1

= −a1 = a0 .

Das ist aber nicht schlimm - wir verschieben den Startwert einfach etwas und nehmen z.B.

a0 := 0.4 ... und schon konvergiert die Folge wie der Wirbelwind gegen 0. 0.4, -0.2461538462, 0.0364566877, -0.00009729639, 0.184 10-11, 0.

Beim Startwert a0 = .5 ergibt sich zufällig schon als erste Iteration die Nullstelle −1, aber das ist die von .5 am weitesten entfernte! Das Auffinden eines günstigen Startwertes für das Newton-Verfahren ist nicht immer einfach. Man wird zunächst versuchen, gewisse Bereiche einzugrenzen, innerhalb derer eine Nullstelle liegen muß. Hat man zum Beispiel zwei Stellen gefunden, an denen die Funktion verschiedenes Vorzeichen hat, so sichert der Zwischenwertsatz die Existenz einer dazwischen liegenden Nullstelle. Sollte die Funktion im gesamten Definitionsbereich keinen Vorzeichenwechsel haben, so muß bei einer Nullstelle ein Extremum liegen, und die Ableitung hat dort eine Nullstelle. Man kann dann mit f ´ ebenso verfahren wie mit f . Der folgende Satz ist deshalb nützlich, weil man unter den genannten Voraussetzungen sicher sein kann, daß das Verfahren konvergiert, ohne bereits nahe an einer Nullstelle zu starten. Konvergenzsatz Gegeben sei eine auf dem kompakten Intervall [a, b] zweimal stetig differenzierbare Funktion f , so daß zwar f , aber weder f ´ noch f ´´ in [a, b] eine Nullstelle hat. Liegen die Iterationswerte a1 für a0 = a und b1 für b0 = b ebenfalls in [a, b], so liegt für jeden Startwert x0 aus [a, b] auch die zugehörige Iterationsfolge ganz in diesem Intervall und konvergiert quadratisch gegen eine Nullstelle von f . Berechnung von Wurzeln Schließlich soll nicht unerwähnt bleiben, daß das uralte Babylonische Wurzelziehen (siehe 4.4) ein Spezialfall des Newton-Verfahrens ist: Für f( x ) = x2 − z lautet die Iteration 2

ak + 1 = ak −

ak − z 2 ak

=

ak 2

−

z 2 ak

.

Ebenso kann man mit dem Newton-Verfahren auch Wurzeln höheren Grades sehr schnell näherungsweise berechnen: Um die n-te Wurzel aus einer positiven Zahl z zu finden, sucht man die Nullstelle der Funktion f( x ) = xn − z. Daß es eine solche Nullstelle gibt, sichert der Zwischenwertsatz, denn f( 0 ) = −z ist kleiner als 0 und f( z + 1 )ist größer als 0. Daß es genau eine gibt, folgt aus der strengen Monotonie der Funktion im positiven Bereich. Wir leiten also ab: f ´( x ) = n x

(n − 1)

und erhalten die Rekursion n

ak + 1 = ak −

ak − z n ak

(n − 1)

 1 z =  1 −  ak − . (n − 1)  n n ak

Beispiel 4: Vierte Wurzel aus 5 Wir starten mit dem (schon recht guten) Wert a0 :=

3 2

und iterieren, wobei wir beim k-ten Schritt 2k + 1 Nachkommastellen von ak angeben: k = 0, ak = 1.50 f( ak ) = 0.06250000000 k = 1, ak = 1.495 f( ak ) = 0.0002887569371 k = 2, ak = 1.49535 f( ak ) = 0.6253121338 10-8 k = 3, ak = 1.495348781 f( ak ) = 0.2932614481 10-17 k = 4, ak = 1.49534878122122054 f( ak ) = 0.6450170771 10-36 k = 5, ak = 1.495348781221220541911898994140913 f( ak ) = 0.3120352723 10-73 k = 6, ak = 1.49534878122122054191189899414091339536345975761470634551659350005 f( ak ) = 0.7302450836 10-148 Das ist nach 6 Iterationen bereits auf 64 Stellen genau! MAPLE bestätigt: 5

(1 / 4)

= 1.495348781221220541911898994140913395363459757614706345516593500048

aber auch das ist natürlich nur ein Näherungswert, denn die vierte Wurzel aus 5 ist irrational.