Wilfried Rohm, HTL Saalfelden

Numerische Stabilität

Mathematische Inhalte :

- Fehlerrechnung
- Vorbereitung für die Betrachtung numerischer Aspekte bei komplizierteren numerischen Algorithmen (Reihenentwicklungen, Lineare Gleichungssysteme) Kurzzusammenfassung : Das Thema "numerische Stabilität" kommt im herkömmlichen Mathematikunterricht zu kurz. Es werden im Hinblick auf später zu besprechende, kompliziertere numerische Algorithmen und als Schnittstelle zum EDV-Unterricht die Probleme und Lösungsmöglichkeiten an Hand einfacher Formeln und Iterationen aufgezeigt. Zentrales Thema ist das Problem der "numerischen Auslöschung". Lehrplanbezug: Die verwendeten Beispiele passen in den 2.Jahrgang, sind als Einstieg in das Thema "Numerische Methoden" aber auch im 3.Jahrgang denkbar. Zeitaufwand: Je nach der zur Verfügung stehenden Zeit und der gewünschten Ausführlichkeit kann das Thema in 1-3 Unterrichtseinheiten behandelt werden. Mediales Umfeld: verwendete Medien: Demonstrationsrechner(PC) und Taschenrechner (in den angeführten Beispielen ist es ein SHARP-1403)

verwendete Software: Turbo-Pascal, Version 6.0

Dateinamen zum Herunterladen: RO-SCHLF.EXE (Demoprogramm zum Thema Rechner-genauigkeit am PC); RO-PINST.EXE (instabile Berechnung von PI); RO-PSTAB.EXE (stabile Berechnung von PI). Alle Programme liegen auch als Source-Code (*.PAS) vor.

 Anmerkungen: Es könnte beim Durchlesen dieses Artikels das Mißverständnis entstehen, im Mathematikunterricht solle EDV ("Programmieren") betrieben werden. Es ist mir aber vielmehr darum gegangen, numerisch und damit auch mathematisch interessante Probleme, die auch beim Programmieren auftreten können, in einem für den Mathematikunterricht relevanten Zusammenhang darzustellen. Die angegebenen Pascal-Programme können nach Besprechung des mathematischen Hintergrundes (bzw. parallel dazu) mittels Demonstrationsrechner in der Klasse vorgeführt und besprochen, oder in Abstimmung mit dem EDV-Lehrer (insbesondere im 2.Jahrgang!) von den Schülern erstellt werden. Hingegen kann auch im Mathematikunterricht das Gelernte beispielsweise in ein numerisch stabiles Programm zum Lösen von quadratischen Gleichungen umgesetzt werden (siehe Anhang).
  1. Rechnergenauigkeit:

Das Problem der "numerischen Stabiltät" von Algorithmen bzw. Formeln ist letztlich auf die begrenzte Rechengenauigkeit jeder Rechenmaschine (Taschenrechner,PC) zurückzuführen. Daher müssen zunächst die Auswirkungen dieser begrenzten Genauigkeit vor Augen geführt werden. Dafür gibt es viele Möglichkeiten, von denen hier einige einfach nachvollziehende dargestellt werden:

Rechnergenauigkeit am Taschenrechner

Bei neueren Taschenrechner ist es bemerkenswerter Weise gar nicht mehr so einfach, die internen Rechenungenauigkeiten zu demonstrieren, da bei Verwendung eingebauter Funktionen Algorithmen ablaufen, die die internen Rechenungenauigkeiten "verschleiern" :

Beispiele:

- (Ö 2)2 ergibt bei den meisten Taschenrechnern tatsächlich "2". Drückt man hingegen vor dem Quadrieren die "Enter" bzw. "=" Taste, so weicht das Ergebnis (meistens) vom exakten Wert ab!

- Mehrmaliges Wurzelziehen und anschließendes Quadrieren zeigt aber die Kumulation der Rechenungenauigkeit auch bei ausschließlicher Verwendung der Funktionstasten auf:: 10-maliges Wurzelziehen und anschließendes 10-maliges Quadrieren liefert am Taschenrechner (PC-1403) das Ergebnis 1.999999963.

- "199*199*199 - 199^3" liefert am PC-1403 das bemerkenswerte Ergebnis "0.00024". Dieses Beispiel stammt aus der Bedienungsanleitung und demonstriert, daß numerische Algorithmen (hier das Potenzieren) "fehlerhaft" arbeiten.
 

 Rechnergenauigkeit am PC :

Im EDV-Unterricht kommen die Schüler meist viel direkter mit dem Problem der internen Rechenungenauigkeiten in Kontakt. Das folgende Beispiel verwende ich daher seit Jahren als Einstieg in die Diskussion über die Folgen der internen Rechenungenauigkeit (Anzumerken ist, daß Taschenrechnerprogramme sich anders verhalten können !) :

Aufgabenstellung:

Es soll eine Tabelle von Funktionswerten mit dem PC in einer höheren Programmiersprache (zB.PASCAL) erstellt werden. Als Eingabe sind Anfangswert, Endwert und Schrittweite der unabhängigen Variablen (hier x) vorzusehen. In beinahe allen Fällen werden Schüler ohne weitere Vorkenntnisse ungefähr folgendermaßen vorgehen :


Dieses Programm (-> RO-SCHLF.EXE bzw. RO-SCHLF.PAS) verhält sich überraschender Weise bei verschiedenen Werten von "step" unterschiedlich :
 

- xa:=-1 , xb:=1 , step := 0.5/0.25/0.125/... führt zu "ordnungsgemäßem" Abbruch.

- xa:=-1 , xb:=1 , step := 0.1 führt zu Endlosschleife, da 0.1 "intern" als Dualzahl nicht exakt dargestellt werden kann, während 0.5/0.25/0.125/... exakt umgewandelt werden.

 Begründung: Es gilt beispielsweise: (0.1)10 = (0.0001100110011001100110011....)2

d.h. die interne Darstellung der (exakten) Dezimalzahl 0.1 kann gar nicht exakt sein.

 Bemerkung : Mit diesem PASCAL-Beispiel läßt sich auch ganz gut die Fehler-"Kumulierung" bei vielen Operationen zeigen.Das Programm liefert zum Beispiel bei "step:=0.1" nach vielen Iterationen den Wert "x = 479.90000032" ,d.h. bereits einen deutlich "sichtbaren" Fehler. Aus diesem Beispiel ergibt sich eine allgemeingültige Grundregel für numerische Berechnungen in Programmen (REAL-Zahlen sind Zahlen in Gleitkommadarstellung):
 


Ein weiteres Beispiel:

Ein zum obigen Beispiel ähnliches Problem tritt bei den Fallunterscheidungen der Lösungen einer quadratischen Gleichung auf :

Ist D=(B2 - 4AC) gleich 0, so liegt bekanntlich eine "reelle Doppellösung" vor. Eine Abfrage wie

"IF D=0 THEN ...."

wird aber nicht immer das gewünschte Ergebnis liefern !
 

 2. Numerische Auslöschung
 

Damit wird ein spezieller Rundungsfehler bezeichnet,der auftritt, wenn zwei nahezu gleichgrosse Zahlen voneinander subtrahiert werden .

Beispiel:

Die Nullen nach dem Dezimalpunkt im Resultat sind nur richtig, wenn die beiden Operanden 1.2345 und 1.2344 "exakte" Zahlen waren. Waren diese auch schon durch Rundungsfehler verfälscht, so ist das Resultat 1.xxxx * 10-4
  und die mit "xxxx" bezeichneten Stellen sind nicht bekannt. Das Resultat (nur einer einzigen Rechenoperation !!) hat also einen großen relativen Fehler. Daher :
 


Beispiel 1: (aus [2])  soll auf 3 Dezimalstellen genau berechnet werden.

a) Nenner wird zuerst berechnet:  ergibt 120.67 » 121 (aufgerundet)

bzw.» 120 (abgerundet)

Das Ergebnis ist daher entweder ¥ oder 380.

b) Der umgeformte Ausdruck  ergibt den wahren Wert 1140.
 

 Beispiel 2: Der angeführte Taschenrechner (PC-1403) verarbeitet Zahlen mit 10-stelliger Genauigkeit. Daher ist logisch, daß 1E11 + 1 - 1E11 zur Ausgabe "0" führt. Beispiel 3: (Iteration am PC) Berechnung von PI durch Approximation des Kreises durch eingeschriebene n-Ecke nach einer Iterationsformel (frei nach Archimedes).
 
  An ....Fläche des eingeschriebenen n-Ecks.

Es gilt : 

Für r=1 und n=6 ist : A6 = 

allgemein für r=1: An =  und  an

 Wir verwenden die folgende Identität zur Berechnung von A2n aus An

(Rekursionsformel):


 

Das folgende Struktogramm ist das Grundgerüst des Programms "RO-PINST.EXE", das die nachfolgenden Ergebnisse liefert:

(Anmerkung: sqrt( ) ... Quadratwurzelfunktion)
 
 


 
 

 Der Grund für die fehlerhafte Berechnung ab einem gewissen n ist "numerische Auslöschung" in der verwendeten Rekursionsformel. Wir können jedoch daraus eine stabile Berechnung machen, indem in der Formel die Auslöschung durch Umformung behoben wird : Gemäß 

erhält man für a=1 und  die numerisch stabile Formel:


  Ferner läßt sich das Abbruchkriterium verbessern: Die Folge der Flächen An ist monoton wachsend und beschränkt - in einer endlichen Arithmetik kann jedoch An nicht für alle n monoton wachsend sein, da es ja nur endlich viele Zahlen gibt! Daher sollte die Iteration abgebrochen werden, sobald A2n < An gilt! Dadurch läuft das Programm maschinen-unabhängig bis zur bestmöglichen Näherung für p .(-> Programm RO-PSTAB.EXE)

(Anmerkung: sqrt( ) ... Quadratwurzelfunktion)
 
 


 


  Bemerkung : Ein entsprechendes Taschenrechnerprogramm liefert auch in der unstabilen Form richtige Ergebnisse, was mit einer anderen internen Verarbeitung zusammenhängen dürfte. Auch ein PASCAL-Programm "RO-PINST" liefert bis zu den angegebenen n-Werten richtige Ergebnisse, wenn bei Vorhandensein eines numerischen Coprozessors mit der weitaus größeren Genauigkeit des Typs "Extended" anstatt "Real" gearbeitet wird.
 
  Im zweiten Jahrgang werden folgende Formeln besprochen, die (bei bestimmten Werten) unstabil sind und durch stabile Varianten ersetzt werden können :
 
 

Beispiel 4: (Stabile Formel zur Lösung einer quadratischen Gleichung)

Die Berechnung der Lösungen einer quadratischen Gleichung gemäß der Formel (mit p=B/A und q= C/A) führt bei ½ p½ >>½ q½ zu Ungenauigkeiten wegen numerischer Auslöschung.

Zum Beispiel liefert am Taschenrechner (PC-1403) die Gleichung x² - 108x +1 = 0

die Lösungen x1=108 und die "falsche" Lösung x2=0.

Nach einem numerisch stabilen Algorithmus berechne man zuerst die betragsmäßig größere Lösung und anschließend die zweite Lösung mit Hilfe des Satz von Vieta:

Ein weiteres Problem könnte bei großen Werten von p die Überschreitung des Rechenbereiches durch p² sein ( Der Datentyp real hat zB. in PASCAL eine Obergrenze von 1.7*1039).

Hier hilft die Umwandlung:

Allerdings muß man dabei  über  berechnen ! Bei Berücksichtigung der angeführten Punkte kann man ein numerisch stabiles Programm zur Lösung quadratischer Gleichungen schreiben. (siehe Taschenrechnerprogramm im ANHANG)
  Beispiel 5 : (Kosinussatz) ist instabil für kleine Werte von g und a» b wegen numerischer Auslöschung.
 

Eine algebraisch äquivalente, aber stabile Formel erhält man durch Anwendung der Summensätze:


  3.Ausblick

In diesem Artikel wurde das Thema "Numerische Stabilität" auf einer sehr elementaren Ebene mit einfachen Beispielen aus den beiden ersten Jahrgängen behandelt. Wirklich bedeutend wird dieses Thema bei der Behandlung komplizierterer numerischer Verfahren im 3. und 4.Jahrgang. Dazu folgende Beispiele :

- Potenzreihenentwicklungen, speziell von alternierenden Reihen, können zum bereits be-schriebenen Problem der numerischen Auslöschung führen. Geplant ist für die nächste Ausgabe dieser Zeitschrift ein Artikel "Taylorreihen und Computereinsatz", in dem unter anderem darauf näher eingegangen werden soll.

- Das numerische Lösen von Gleichungssystemen ist ein Paradebeispiel für das Suchen nach numerisch stabilen Algorithmen, weil sehr viele Berechnungsschritte durchzuführen sind. In vielen Fällen (zB. Gauß-Algorithmus) wird diese Stabilität durch die sogenannte "PIVOT"-Strategie erreicht. Es erscheint aber sinnvoll, den Schüler schon VOR Besprechung derartiger Algorithmen für das Problem der "numerischen Stabilität" sensibilisiert zu haben.

- Hinweis auf die Regeln der Fehlerfortpflanzung:

Addition/Subtraktion : Es addieren sich im ungünstigsten Fall die absoluten Fehler der Operanden.

Multiplikation/Division : Es addieren sich im ungünstigsten Fall die relativen Fehler der Operanden.

Numerisch stabile Algorithmen sind gegen das Akkumulieren von Rundungsfehlern weitgehend unempfindlich.

- Schlecht konditionierte Problemstellungen können auch mit numerisch stabilen Algorithmen NICHT befriedigend gelöst werden: Darunter versteht man Aufgabenstellungen, die bei nur leicht veränderten Ausgangswerten zu sehr unterschiedlichen Ergebniswerten führen.

Das Resultat der numerischen Rechnung wird durch Rundungsfehler in allen Fällen so weit verfälscht, daß es das exakte Resultat leicht geänderter Anfangsdaten ist. Daher erhält man bei schlechter Kondition der Aufgabenstellung Lösungen, die sehr stark von der exakten Lösung abweichen. Beispiele :

- Gleichungssysteme mit "beinahe" singulären Koeffizientenmatrizen

- Berechnung mehrfacher Nullstellen von Polynomen, zum Beispiel mit dem Newton-verfahren (wodurch sich eine Schnittstelle zur Chaostheorie ergibt)
 

 4. LITERATURHINWEISE :
 
  [1] Walter Gander: Computermathematik,Birkhäuser Verlag,1985. (viele numerische Algorithmen mit Stabilitätsbetrachtungen und ihre Umsetzung in PASCAL-Programme)

[2] Walter Riemer: EDV mit Pascal,Österr.Gerwerbeverlag,1990. (Anhang über die Problematik numerischer Methoden)
 

 5. ANHANG :

Struktogramm für ein numerisch stabiles Programm zum Lösen einer quadratischen Gleichung :


 
  (Anmerkung : squ( )...Quadratfunktion ; sqr( )...Wurzelfunktion; abs()...Absolutbetrag)
 
  Umsetzung des obigen Algorithmus auf ein Taschenrechnerprogramm am Beispiel PC-1403 (etwas altertümliches BASIC ! ; das Programm ist ohne formatierte Ausgabe und ohne Fallunterscheidung für spezielle Fälle, etwa A=0, ausgestattet.)

Source- Text der PASCAL-Programme (in etwas vereinfachter Form):

1) Demo-Programm für Schleifenende (->  RO-SCHLF.PAS):


 
 

2) Berechnung von PI -unstabile Version (-> RO-PINST.PAS):


 

3) Berechnung von PI - stabile Version (->  RO-PSTAB.PAS)