I/T5 - WS 04/05 CR – Übung 2 Termin: 12.-16.11.04
Behandelt wird eine rundefehlerbehaftete Rekursion.
Das Integral
1 1
0
( )
tx dt
z x t a
−
=
∫
+ ist eine Funktion des reellen Parameters x und einer Konstanten a>0. Für ganzzahlige x= ∈n N schreibt man auch zn =z n( ). Da der Integrand1
( ) tx
f t t a
= −
+ für alle t∈
[
0,1]
positiv ist, sind alle Integrale zn >0.Die Integrale zn kann mit verschiedenen Techniken für positive ganzzahlige Werte von n berechnen:
(a) durch direkte Integration, (b) durch Vorwärtsrekursion, (c) durch Rückwärtsrekursion.
Geeignete Rekursionsformeln findet man durch Umformung des Integranden:
(
1 1)
1 1 1 1 1
1 1
0 0 0 0
( ) n n 1
n n
n
n n
t a t at dt
t dt t dt
z t dt
t a t a t a n
− − −
− +
+ −
= = = − = −
+ + +
∫ ∫ ∫
a∫
az .Aus dieser Gleichung liest man für n∈N die Vorwärtsrekursionsformel ab: 1 1
n n
z a
+ = −n z . Als exakt berechenbaren Anfangswert kann man für die Vorwärtsrekursion
[ ]
1
1
1 0
0
1 ln( ) t ln( 1) ln( )
z dt t a t a a
t a
=
= = + = = + −
∫
+ einsetzen und damit der Reihe nach z z2, 3,...u.s.w. berechnen.
Nach Umstellen der Terme gilt auch die Rückwärtsrekursionsformel 1 1 1
n n
z z
a n +
⎛ ⎞
= ⎜ − ⎟
⎝ ⎠. Da tn⎯⎯ →n→⎯∞ 0 für t∈
[
0,1]
gilt für große N 01 1
0 N 0
z + ≈
∫
dt= . Deshalb eignet sich zN+1=0 als (etwas ungenauer) Anfangswert für die Rückwärtsrekursion und man findet damit der Reihe nach zN−1,zN−2,...,z1.Sie sollen an diesem Beispiel die Fortpflanzung von Gleitpunktfehlern untersuchen.
Erstellen Sie ein MatLab-Programm, das folgende Anforderungen erfüllt:
(a) Der Integrand
1
( ) tx
f t t a
= −
+ ist als MatLab-Funktion integrand zu formulieren, entweder als inline-Funktion oder in einem m-File. Parameter sind t,x,a.
(b) Der Parameter a soll im ganzen Programm einen festen Wert haben, z.B. a = 2π2.
(c) Die Integralfunktion z x( ) ist im Intervall x∈[1,N], N = 14, als blaue Kurve zu zeichnen.
Benutzen Sie als Schrittweite für x und die MatLab-Routine quadl für die Integration. Bein Aufruf enthält quadl als ersten Parameter den Namen integrand der Funktion
0,1 h =
1
( ) tx
f t t a
= −
+ , falls diese als m-File vorliegt, muss @ vorangestellt werden. Die weiteren Parameter x,a des Integranden müssen in quadl ebenfalls angegeben werden, aber
erst ab dem sechsten Parameter. (Anmerkung: warning off schaltet evtl. Warnungen aus, die bei der Intergration mit quadl evtl. auftreten können.)
(d) Das Bild soll auch die x-Achse enthalten. Ihr Bild zeigt den tatsächlichen Funktionsverlauf, den wir im folgenden schneller durch die Rekursionsformeln für ganzzahlige
nachbilden wollen.
x= ∈n N (e) Ermitteln Sie nun durch Vorwärtsrekursion die Zahlen z nn, =2, 3,...,N und zeichnen Sie
diese als rote Kringel ins vorhandene Bild ein. Wahrscheinlich treffen nicht alle Punkte die Integralkurve genau, da sich allmählich Rundefehler fortpflanzen, obwohl der Anfangswert
mit einem relativen Fehler < eps behaftet, also praktisch exakt ist.
z1
(f) Setzen Sie zN+1=0 (kleiner Fehler!) und berechnen Sie durch Rückwärtsrekursion die Zahlen . Zeichnen Sie diese als grüne Kreuze ins vorhandene Bild ein.
Trotz des kleinen Fehlers im Anfangswert treffen die grünen Kreuze die blaue Kurve genauer als die roten Punkte. Erklärung dieses Phänomens: Die Subtraktion verursacht hier viel geringere Rundefehler; der Anfangsfehler von
, , 1,..., 2,1 z nn =N N−
1 0
zN+ = wird während der Rechnung gedämpft.
(g) Addieren Sie die maximalen Bitstellenverluste p (s.u.) bei beiden Rekursionen Schritt für Schritt auf und geben Sie die Zahlen tabellarisch in drei Spalten aus: nebeneinander n p p.
Hinweis: Bei Subtraktion x – y gehen maximal log2 x y
p x y
⎛ + ⎞
= ⎜ −⎝ ⎠⎟ Bit verloren.
Zweierlogarithmus in MatLab: log2.
(h) Daraus erkennen Sie auch die Anzahl m der verläßlichen Dezimalstellen in jedem rekursiv berechneten Folgenglied zn: 3
10
m≈ p. Geben Sie n p p m m für die beiden Rekursionen tabellarisch aus.
