numeriki - uzhuser.math.uzh.ch/berner/files/numerik.pdf · 1 zahlendarstellungundrundungsfehler...
Post on 31-Mar-2019
213 Views
Preview:
TRANSCRIPT
Inhaltsverzeichnis
1 Zahlendarstellung und Rundungsfehler 41.1 Zahlendarstellung . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 4
1.1.1 Basiswechsel . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 51.2 Maschinenzahlen . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 6
1.2.1 Gleitpunktdarstellug . . . . . . . . . . . . . . . . . . . . . . . . . . 71.2.2 Rundungsfehler . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 81.2.3 Fehlerfortp anzung bei Operationen . . . . . . . . . . . . . . . . 8
1.3 Kondition . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 91.3.1 Die Kondition einer Funktion . . . . . . . . . . . . . . . . . . . . . 91.3.2 Die Kondition eines Algorithmus . . . . . . . . . . . . . . . . . . . 151.3.3 Globaler Fehler einer Berechnung . . . . . . . . . . . . . . . . . . 16
2 Auswertung elementarer Funktionen 172.1 Gewöhnliche Polynome . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 17
2.1.1 Auswertung von P (ξ) . . . . . . . . . . . . . . . . . . . . . . . . . 172.2 Trigonometrische Polynome . . . . . . . . . . . . . . . . . . . . . . . . . . 21
3 Iterative Verfahren 243.1 Vorbereitung . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 243.2 Algorithmus . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 24
3.2.1 Bisektions Methode . . . . . . . . . . . . . . . . . . . . . . . . . . 243.2.2 Fixpunkt Methode . . . . . . . . . . . . . . . . . . . . . . . . . . . 253.2.3 Newton Methode . . . . . . . . . . . . . . . . . . . . . . . . . . . . 263.2.4 Sekanten Methode . . . . . . . . . . . . . . . . . . . . . . . . . . . 293.2.5 Geschwindigkeit der Algorithmen . . . . . . . . . . . . . . . . . . 32
3.3 Konvergenzbeschleunigung . . . . . . . . . . . . . . . . . . . . . . . . . . . 343.3.1 Die ∆2-Methode von Aitken . . . . . . . . . . . . . . . . . . . . . 34
3.4 Algebraische Gleichungen . . . . . . . . . . . . . . . . . . . . . . . . . . . 36
4 Approximation und Interpolation 414.1 Polynom-Interpolation . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 41
4.1.1 Lagrange & Newton Interpolation . . . . . . . . . . . . . . . . . . 414.1.2 Hermite Interpolation . . . . . . . . . . . . . . . . . . . . . . . . . 49
4.2 Spline Funktionen . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 514.2.1 Lineare Splines . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 524.2.2 kubische Splines . . . . . . . . . . . . . . . . . . . . . . . . . . . . 53
4.3 Ausgleichung im quadratischen Mittel . . . . . . . . . . . . . . . . . . . . 56
2
INHALTSVERZEICHNIS
5 Numerische Integration 595.1 Einführung . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 595.2 Interpolarische Formeln . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 595.3 Zusammengsetzte Formeln . . . . . . . . . . . . . . . . . . . . . . . . . . . 635.4 Konvergenzuntersuchungen . . . . . . . . . . . . . . . . . . . . . . . . . . 635.5 Gauss-Quadratur . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 645.6 Varia . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 67
Index 68
3
1 Zahlendarstellung und Rundungsfehler
1.1. Zahlendarstellung
De nition 1.1 Sei N ∈ N, N ≥ 2. N ist die Basis der Darstellung.Sei x ∈ R. Die Darstellung von x im N -er System ist eine Folge
±a0a1 . . . ap.ap+1ap+2 . . . ai ∈ {0,1, . . . ,N − 1}.
Notation(±a0a1 . . . ap.ap+1ap+2 . . . )N = ±a0Np + a1Np−1 + ⋅ ⋅ ⋅ + apN0 + ap+1
N+ ap+2
N2 + . . .wobei N0 = 1.
Beispiel 1.2 (+11.52)10 = 1 ⋅ 10 + 1 + 510+ 2
102
Satz 1.3 Sei x ∈ R, dann gilt
1 + x + ⋅ ⋅ ⋅ + xn = 1 − xn+1
1 − x.
Sei x ∈ R, ∣x∣ < 1, dann gilt+∞∑k=0
xk = 1
1 − x.
B(1 − x)(1 + x + ⋅ ⋅ ⋅ + xn) = 1 + x + x2 + ⋅ ⋅ ⋅ + xn − x − x2 ⋅ ⋅ ⋅ − xn+1 = 1 − xn+1
d.h. 1 + x + ⋅ ⋅ ⋅ + xn = 1−xn+1
1−x
∑+∞k=0 xk = limn→+∞ 1 + x + ⋅ ⋅ ⋅ + xn = limn→+∞1−xn+1
1−x = 11−x
Bemerkung 1.5 Die Reihe (Summe) ap+1N+ ap+2
N2 + ⋅ ⋅ ⋅ +ap+kNk konvergiert. (wobei ap+k ∈
{0, . . . ,N − 1})Denn 0 ≤ uk = ap+k
Nk ≤ N−1Nk ≤ N
Nk = 1Nk−1 . Weiter ist N eine ganze Zahl ≥ 2, und so-
mit ist 1N≤ 1
2< 1. Das heisst, die Reihe 1
Nk−1 konvergiert und somit konvergiert auch uk .
Bemerkung 1.6 Die Darstellung einer Zahl ist nicht eindeutig bestimmt.Beispiel 1.7
• 0.999 . . . = 910+ 9
102+ ⋅ ⋅ ⋅ + 9
10k= 9
10(1 + 1
10+ 1
102+ . . . ) = 9
101
1− 110
= 1,
• 0.25999 . . . = 0.26,
• (0.111 . . .)2 = 12+ 1
22+ 1
23+ ⋅ ⋅ ⋅ = 1
2(1 + 1
2+ 1
22+ . . . ) = 1
21
1− 12
= 1.
4
1.1 Zahlendarstellung
1.1.1. Basiswechsel
N Bezeichnung einer Zahl in Basis N2 Binärzahl8 Oktalzahl10 Dezimalzahl16 Hexadezimalzahl
Basiswechsel vomN er ins 10er-System
Beispiel 1.8
1. Was ist (1111)2 im Dezimalsystem:(1111)2 = 1 ⋅ 20 + 1 ⋅ 21 + 1 ⋅ 22 + 1 ⋅ 23 = 1 + 2 + 4 + 8 = 15,
2. (231)4 = 1 ⋅ 40 + 3 ⋅ 41 + 2 ⋅ 42 = 1 + 12 + 32 = (45)10.
3. Die Ziffern des Hexadezimalsystems sind 0,1, . . . ,9,A,B,C,D,E,F .(A1D)16 = 13 + 1 ⋅ 161 + 10 ⋅ 162 = 29 + 2560 = (2589)10,
4. (0.111)2 = 12+ 1
4+ 1
8= 4+2+1
8= 7
8.
Basiswechsel vom 10er insN er-System
• Finden des ganzen Teils:x = a0 + a1N + ⋅ ⋅ ⋅ + apNp = a0 +Nq1q1 = a1 + a2N + . . .
• Finden des dezimalen Teils:Sei x ∈ (0,1) x = a1
N+ a2
N2 + . . .Nx = a1 + a2
N+ . . .
a1 = [Nx] = der ganze Teil von NxNx − a1 = Nx − [Nx] = a2
N+ . . .
Wir setzen f(x) = Nx − [Nx]a2 = [Nf(x)], f(f(x)) = f2(x) = Nf(x) − [Nf(x)] = a3
N+ . . .
Beispiel 1.9
1. (24)10 = (11000)2
2. (24)10 = (33)7
3. (0.9375)10 = (0.1111)2.Denn 0.9375 = a0
2+ a1
22+ . . . a0 = [2 ⋅ 0.9375] = [1.875] = 1
a1 = [(2 ⋅ 0.9375 − [1.875]) ⋅ 2] = [2 ⋅ 0.875] = [1.75] = 1a2 = [2 ⋅ 0.75] = [1.5] = 1a3 = [2 ⋅ 0.5] = [1] = 1.
De nition 1.10 x besitzt eine periodische Entwicklung, fallsx = ±a0a1 . . . ap.ap+1ap+2 . . . b1 . . . bq b1 . . . bq . . .
5
Zahlendarstellung und Rundungsfehler
Satz 1.11 Sei x ∈ Rx besitzt eine periodische Entwicklung⇐⇒ x ∈ Q.
B
“⇒” Sei x mit einer periodischen Entwicklung.x = ±a0a1 . . . ap.ap+1ap+2 . . . ap+kb1 . . . bqb1 . . . bq . . .A = ±a0a1 . . . ap.ap+1ap+2 . . . ap+k ∈ Qx = A + b1
Nk+1 + b2Nk+2 + ⋅ ⋅ ⋅ + bq
Nk+q + b1Nk+q+1 + ⋅ ⋅ ⋅ + bq
Nk+2q + . . .
b1Nk+1 + b2
Nk+2 + ⋅ ⋅ ⋅ + bqNq + b1
Nk+q+1 + ⋅ ⋅ ⋅ + bqNk+2q + . . .
= 1Nk ( b1N + ⋅ ⋅ ⋅ +
bqNq ) + 1
NkNq ( b1N + ⋅ ⋅ ⋅ +bqNq ) + 1
NkN2q ( b1N + ⋅ ⋅ ⋅ +bqNq )
= ( b1N+ ⋅ ⋅ ⋅ + bq
Nq ) 1Nk (1 + 1
Nq + 1N2q + . . . )
= ( b1N+ ⋅ ⋅ ⋅ + bq
Nq ) 1Nq
11− 1
Nq∈ Q.
“⇐” Sei x ∈ Q, x > 0. Also x = [x] + abwobei a < b.
Es genügt zu zeigen, dass abeine periodische Entwickung besitzt.
y = ab= a1
N+ a2
N2 + . . .a1 = [Ny], Ny = [Ny] + a2
N+ . . .
f(y) = a2
N+ . . .
Nf(y) = a2 + a3
N+ . . .
a2 = [Nf(y)]f(f(y)) = Nf(y) − [Nf(y)] = a3
N+ . . .
a1a2 ⋅ ⋅ ⋅ = [Ny][Nf(y)][Nf2(y)] . . .Ny = Na
b, Na = qb + r1, 0 ≤ r1 < b
Nab= Ny = q + r1
b[Ny] = q Ny − [Ny] = r1
b0 ≤ r1 < b
f(y) = r1b
a1a2 ⋅ ⋅ ⋅ = [Ny][Nf(y)][Nf2(y)] ⋅ ⋅ ⋅ = [q][N r1b][N r2
b] . . .
f(y) = r1b, f2(y) = r2
b, . . . , fk(y) = rk
bDie ri nehmen eine endliche Anzahl von Werten an. Sei rk die erste Zahl mit fürein p, rk−p = rk. Dann ist rk−p, rk−p+1, . . . , rk−1 die Periode der Zahl y.
Beispiel 1.13 Darstellung von 127
im Dualsystem.127= 1 + 5
757= a1
2+ a2
22+ . . .
1 + 37= 10
7= a1 + a2
2+ . . . 3
7= a2
2+ a3
22+ . . .
67= a2 + a3
2+ . . .
67= a3
2+ . . .
1 + 57= 12
7= a3 + . . . also: 12
7= 1.101101101 ⋅ ⋅ ⋅ = 1.101
1.2. Maschinenzahlen
Sei RM die Menge der Maschinenzahlen, d.h. die Menge der Zahlen die im Computerbenutzt werden können.
6
1.2 Maschinenzahlen
1.2.1. Gleitpunktdarstellung ( oating point representation)
Jede Zahl ist als ±. a1a2 . . . atMantisse
± b1 . . . bsExponent
gespeichert.
Beispiel 1.14Zahl t s Maschinenzahl12.306 3 1 +.123 + 2-0.000120121 4 2 −.1201 − 3.
Die kleinste und die grösste Zahl die man speichern kann
Wir arbeiten im N er-System.t = Anzahl der Ziffern der gespeicherten Zahl.s = Anzahl der Ziffern des Exponenten.
Satz 1.15 Sei t und s gegeben.
• Die kleinste Zahl, die man speichenr kann, istN−Ns
• Die grösste Zahl, die man speichern kann, ist (1 − 1
N t)NNs−1
B
• Die kleinste Zahl: .1Nem wobei em ist das Minimum des Exponenten.em = −(N − 1)(N − 1)(N − 1)(N − 1) . . . (N − 1) (Folge von Ziffern)
= −(N − 1) + (N − 1)N + ⋅ ⋅ ⋅ + (N − 1)Ns−1
= −(N − 1)(1 +N + ⋅ ⋅ ⋅ +Ns−1)= (N − 1) 1−N
s
1−N= 1 −Ns.
Also ist die kleinste Zahl .1Nem = 1.Nem−1 = N−Ns
• Die grösste Zahl ist: (N − 1)(N − 1) . . . (N − 1) ⋅NeM wobei eM der grösste Ex-ponent ist.A = N−1
N+ N−1
N2 + ⋅ ⋅ ⋅ + N−1Nt = N−1
N(1 + 1
N+ ⋅ ⋅ ⋅ + 1
Nt−1 ) = 1 − 1Nt
eM = (N − 1)(N − 1) . . . (N − 1)= N − 1 + (N − 1)N + ⋅ ⋅ ⋅ + (N − 1)Ns−1
= (N − 1)(1 +N + ⋅ ⋅ ⋅ +Ns−1)= Ns − 1
also grösste Zahl: (1 − 1Nt )NNs−1
RM = {N−Ns
≤ .a1 . . . atNe ≤ (1 − 1Nt )NNs−1}
Beispiel 1.17 t = 6, e ∈ [−9,9]23= a1
10+ a2
102+ . . .
203= a1 + a2
10+ . . . a1 = 6
23= a2
10+ . . .
23= 0.66666666
7
Zahlendarstellung und Rundungsfehler
0.666667 ist näher an 23als 0.666666 und 2
3wird als
23= 0.666667 gespeichert.
1.2.2. Rundungsfehler
Sei x eine Zahl. In der Maschine arbeitet man mit einer Approximation rd(x), d.h. dienächstgelegene Maschinenzahl die erfüllt:
∣x − rd(x)∣ ≤ ∣x − y∣ ∀y ∈ RM .
De nition 1.18
• Absoluter Fehler: ea = ea(x) = x − rd(x),
• Relativer Fehler: er = er(x) =x − rd(x)
x, wobei x ≠ 0.
Beispiel 1.19 t = 2, e ∈ [−9,9].x = 1.04 10−9 = 0.104 ⋅ 10−8 rd(x) = 0.1 10−8.y = 9.89 108 = 0.989 ⋅ 109 rd(y) = 0.99 109.ea(x) = x − rd(x) = 0.004 ⋅ 10−8 = 0.4 ⋅ 10−10.ea(y) = y − rd(y) = −0.001 ⋅ 109 = −0.1 ⋅ 107.er(x) = 0.4⋅10−10
0.104⋅10−8 =0.4⋅10−20.104
= 0.040.104
10−1 ∣er(x)∣ ≤ 1210−1.
er(y) = 0.1⋅107.989⋅109 =
0.10.989
10−2 ≤ 1210−1.
Satz 1.20 Seim =min{∣x∣ ∣ x ∈ RM},M =max{∣x∣ ∣ x ∈ RM}.Seim ≤ x ≤M .Dann gilt
∣x − rd(x)x
∣ ≤ 1
2N1−t
(d.h. der relative Fehler ist von x unabhängig).
B Sei x wie oben.x = ±.a1 . . . atat+1 . . . Ne,rd(x) = ±.a1 . . . at ⋅Ne (ev at = at + 1 je nach at+1).∣x − rd(x)∣ ≤ 1
2N−tNe,
∣er(x)∣ = ∣x−rd(x)x∣ ≤ 1
2N−tNe
∣x∣ ≤12N−tNe
Ne/N =12N1−t
da ∣x∣ ≥ 0.1 ⋅Ne = 1N⋅Ne.
Beispiel 1.22 Für N = 2∣εx∣ = ∣er(x)∣ ≤ 1
221−t = 2−t = eps =Maschinengenauigkeit. (precision of the computer)
1.2.3. Fehlerfortp anzung bei Operationen
Sei ○ eine Operation (+,−, ⋅, /).Für x, y gilt fl(x ○ y) = rd(x ○ y) ≃ x ○ y = x ○ y(1 − εx○y), ∣εx○y ∣ ≤ 1
2N1−t.
8
1.3 Kondition
Multiplikation
rd(x) = x(1 − εx), rd(y) = y(1 − εy).rd(x) rd(y) = x(1 − εx)y(1 − εy) = xy(1 − εx − εy − εxεy) ≃ xy(1 − εx − εy)
Division
rd(x) = x(1 − εx), rd(y) = y(1 − εy).rd(x)rd(y) =
xy( 1−εx1−εy )
11−εy = 1 + εy + ε
2y + . . .
Ô⇒ rd(x)rd(y) =
xy(1 − εx)(1 + εy + ε2y + . . . ) = x
y(1 − εx + εy + . . . ) ≃ x
y(1 − (εx − εy))
Also: ε xy≅ εx − εy .
Addition
rd(x) + rd(y) = x(1 − εx) + y(1 − εy)= x + y − xεx − yεy= (x + y)(1 − x
x+y εx −y
x+y εy).Also: εx+y = x
x+y εx +y
x+y εy kann sehr gross sein für x + y nah von 0.Falls x ⋅ y > 0 dann 0 ≤ x
x+y ≤ 1, 0 ≤ yx+y ≤ 1 Ô⇒∣εx+y ∣ ≤ ∣εx∣ + ∣εy ∣
Beispiel 1.23 (Auslöschung (Cancellation error)) t = 7, N = 10,g=“garbage” (unsichere Stelle).x = 0.100011g, y = −0.100010g′x + y = 0.1g′′10−5.Es ist gut die Subtraktion zu vermeiden. Hier sind verschiedene Beispiele√x + δ −
√x = (
√x + δ −
√x) ⋅
√x+δ+
√x√
x+δ+√x= δ√
x+δ+√x
f(x + δ) − f(x) = f ′(x + θδ)δ, θ ∈ (0,1)√x + δ −
√(x) = δ
2√x+θδ
cos(a + b) = cosa cos b − sina sin bcos(a − b) = cosa cos b + sina sin bcos(a + b) − cos(a − b) = −2 sina sin bcosp − cos q = −2 sin(p+q
2) sin(p−q
2), wobei p = a + b, q = a − b.
cos(x + δ) − cos(x) = −2 sin(x + δ2) sin δ
2
1.3. Kondition
1.3.1. Die Kondition einer Funktion
x ∈ Rm f→ y = f(x) ∈ Rn
rd(x) → f(rd(x)) ?≃ y = f(x)
• m = n = 1
9
Zahlendarstellung und Rundungsfehler
f(rd(x)) = f(x(1 − εx))= f(x − xεx)= f(x) + f ′(x − θxεx)(−xεx)≃ f(x) + f ′(x)(−xεx)≅ f(x)(1 − xf ′(x)
f(x) εx)
εf(x) = xf ′(x)f(x) εx.
De nition 1.24Kond(f)(x) = ∣x ⋅ f
′(x)f(x)
∣
Konditionszahl der Funktion f an der Stelle x. (x ≠ 0, f(x) ≠ 0).
Diese Zahl misst, wie stabil die Funktion f bezüglich kleinen Perturbationen ist.
Beispiel 1.25 f(x) = xn, f ′(x) = nxn−1,Ô⇒Kond(f)(x) = ∣x⋅n⋅x
n−1
xn ∣ = n.
Für die Fälle x = 0 oder f(x) = 0 haben wir das Folgende.
• x = 0, f(x) ≠ 0.f(0 + e0) ≅ f(0) + f ′(0)e0 = f(0)(1 + f ′(0)
f(0) e0),
Kond(f)(0) = ∣ f′(0)f(0) ∣.
• x ≠ 0, f(x) = 0.f(rd(x)) = f(x(1 − εx)) = f(x − xεx) ≅ f(x) − f ′(x)xεx = −f ′(x)xεx,Kond(f)(x) = ∣xf ′(x)∣.
• x = 0, f(x) = 0.f(0 + e0) ≅ f(0) + f ′(0)e0 = f ′(0)e0,Kond(f)(0) = ∣f ′(0)∣.
Wir betrachten jetzt den Fall einer Funktion von Rn nach Rn.Sei f ∶ Rm → Rn, x = (x1, . . . , xm)↦ (f1(x1, . . . , xm), . . . , fn(x1, . . . , xm))Sei x ∈ Rm, xi → fj(x1, . . . , xm), ist eine Funktion von R nach R mit Konditionszahl
γij(x) =xi
∂fj∂xi(x)
fj(x).
Wir haben dann n ⋅m Konditionszahlen und wir können de nieren:
Kond(f)(x) = ∥(γij(x))∥,
für ∥ ∥ eine Norm der Matrizen.
10
1.3 Kondition
De nition 1.26 Sei V ein Vektorraum über R oder C.Eine Norm ist eine AbbildungV → R+ = [0,+∞)x↦ ∥x∥mit den Eigenschaen
1. ∥x∥ = 0 ⇐⇒ x = 0,
2. ∥λx∥ = ∣λ∣∥x∥ ∀λ ∈ R bzw C, ∀x ∈ V ,
3. ∥x + y∥ ≤ ∥x∥ + ∥y∥ ∀x, y ∈ V .
Beispiel 1.27
1. Rn euklidische Norm:
∥x∥2 = (n
∑i=1
x2i)
12
,
n = 2 ∣x∣2 = (x21 + x2
2)12 .
2. 1 ≤ p < +∞
∣x∣p = (n
∑i=1∣xi∣p)
1p
.
Für p = +∞ ∣x∣∞ =maxi=1...m ∣xi∣ = limp→∞ ∣x∣p.
3. Für C das selbe, wobei ∣z∣2 = zz für z ∈ C.
A = (aij)i=1...n,j=1...m A ist eine n ×m Matrix.Die Matrizen bilden einen Vektorraum über R.
Beispiel 1.28
∥A∥p =⎛⎝∑i,j∣aij ∣p
⎞⎠
1p
.
De nition 1.29 Sei ∥ ∥ eine Norm auf Rm, ∣ ∣ eine Norm auf Rn.Dann
∥∣A∣∥ = supx≠0
∣Ax∣∥x∥
ist eine Norm.
Satz 1.30 Es gilt
∥∣A∣∥∞ = supx≠0
∣Ax∣∞∥x∥∞
= maxi=1...n
m
∑j=1∣aij ∣.
B
11
Zahlendarstellung und Rundungsfehler
1.
∣Ax∣∞ = maxi=1...n
∣(Ax)i∣
= maxi=1...n
∣m
∑j=1
aijxj ∣
≤ maxi=1...n
m
∑j=1∣aij ∣∣xj ∣
≤ maxi=1...n
m
∑j=1∣aij ∣∥x∥∞
≤ ∥x∥∞ maxi=1...n
m
∑j=1∣aij ∣,
∣Ax∣∞∥x∥∞
≤ maxi=1...n
m
∑j=1∣aij ∣.
2. ∥∣A∥∣∞ ≥maxi=1...n∑mj=1 ∣aij ∣
Wir suchen x0 ∈ Rm mit ∣Ax0∣∞∥x0∥∞ = ∥∣A∥∣∞
(es wird implizieren, supx≠0∣Ax0∣∞∥x0∥∞ ≥
∣Ax0∣∞∥x0∥∞ ≥ ∥∣A∥∣∞).
∃i0 ∈ {1, . . . , n}mit maxi=1,...,n∑mj=1 ∣aij ∣ = ∑
mj=1 ∣ai0j ∣
ai0j . . . ai0m,
x0,j ∶= sign(ai0j) = {+1 ai0j ≥ 0−1 ai0j < 0
(Ax0)i0 = ∑mj=1 ai0jx0,j = ∑m
j=1 ∣ai0j ∣∥x0∥∞ =maxj=1...m ∣x0,j ∣ = 1
∣Ax0∣∞∥x0∥∞
= ∣Ax0∣∞
= maxi=1...m
∣(Ax0)i∣
≥ ∣(Ax0)i0 ∣
=m
∑j=1∣ai0j ∣
= maxi=1...n
m
∑j=1∣aij .∣
De nition 1.32 Eine Norm ist eine Matrix Norm falls gilt
∥AB∥ ≤ ∥A∥∥B∥.
Sei x ein Vektor.Der relative Fehler über x ist ∥∆x∥
∥x∥ wobei ∆x = (∆x1, . . . ,∆xm).
12
1.3 Kondition
Sei f ∶ Rm → Rn, x↦ f(x) = y.x +∆x↦ f(x +∆x) = y +∆y, ∥∆x∥
∥x∥ ↝∥∆y∥∥y∥ .
Für alle j gilt:fj(x +∆x) = fj(x) +∑m
i=1∂fj∂xi(x)(∆x)i + ⋅ ⋅ ⋅ ≅ fj(x) +∑m
i=1∂fj∂xi(x)(∆x)i,
Das heisst, man hat
f(x +∆x) =⎛⎜⎝
f1(x +∆x)⋮
fn(x +∆x)
⎞⎟⎠
= (fj(x) +m
∑i=1
∂fj
∂xi(x)(∆x)i)
= f(x) + (m
∑i=1
∂fj
∂xi(x)(∆x)i)
f(x +∆x) = f(x) + (∇f)(∆x)wobei∇f ist die Jacobi Matrix gegeben mit∇f(x) = (∂fj
∂xi)ji= (aji)
∥∆y∥∞ = ∥f(x +∆x) − f(x)∥∞ = ∥(∇f)∆x∥∞ ≤ ∥(∇f)∥∞∥∆x∥∞.⇐⇒ ∥∆y∥
∥y∥ ≤∥x∥∞∥(∇f)∥∞∥∆x∥∞
∥f∥∞∥x∥∞ .Dann können wir setzen
Kond(f)(x) = ∥x∥∞∥(∇f)(x)∥∞∥f(x)∥∞
so dass
∥∆y∥∞∥y∥∞ ≤ Kond(f)(x) ∥∆x∥∞
∥x∥∞ .
Beispiel 1.33 n =m = 2, f ∶ R2 → R2.
f(x) = f ( x1
x2) = (
1x1+ 1
x21x1− 1
x2
) = ( f1(x)f2(x)
)
Notation: ∂∂xi= ∂xi
γ11 =x1∂x1
f1f1
= −x2
x1+x2γ12 =
x1∂x1f2
f2= x2
x1−x2
γ21 =x2∂x2
f1f1
= −x1
x1+x2γ22 =
x2∂x2f2
f2= x1
x2−x1
(schlecht für x1 ≃ −x2).
13
Zahlendarstellung und Rundungsfehler
(∇f) = ( ∂x1f1 ∂x2f1∂x1f2 ∂x2f2
) =⎛⎝− 1
x21− 1
x22
− 1x21
1x22
⎞⎠
∥∣∇f∥∣∞ = max{ 1x21+ 1
x22, ∣ 1
x22− 1
x21∣}
= max{x21+x
22
x21x
22,∣x2
2−x21∣
x21x
22}
= 1x21x
22max{x2
1 + x22, ∣x2
2 − x21∣}
= x21+x
22
x21x
22.
∥x∥∞ = max{∣x1∣, ∣x2∣}.∥f(x)∥∞ = max{∣ 1
x1+ 1
x2∣, ∣ 1
x2− 1
x1∣}
= max{ ∣x1+x2∣∣x1x2∣ ,
∣x2−x1∣∣x1x2∣ }
= 1∣x1∣∣x2∣ max{∣x1 + x2∣, ∣x1 − x2∣}
Kond(f)(x) = ∥x∥∞∥∇f∥∞∥f∥∞ = max{∣x1∣,∣x2∣}∣x1∣∣x2∣
x21+x
22
max{∣x1+x2∣,∣x1−x2∣} .
x1 ≃ x2: Kond(f)(x) ≅ ∣x1∣2x21
∣x1∣22∣x1∣ ≃ 1,x1 ≃ −x2 Kond(f)(x) ≃ 1.In beiden Fällen ist diese Globale Kondition gut. Es war nicht der Fall mit der γij .
Was passiert falls y(x) = f(x) ist die Lösung einer algebraischen Gleichung ist?
Beispiel 1.34 (Algebraische Gleichungen)P (y) = a0 + a1y + a2y2 + ⋅ ⋅ ⋅ + an−1yn−1 + yn = 0.Wir suchen y(a0, . . . , an−1) Nullstellen dieser Gleichung.
Wir nehmen an, dass y bezüglich a differenzierbar ist.∂
∂ai(a0 + a1y + a2y2 + ⋅ ⋅ ⋅ + an−1yn−1 + yn) = 0,
a1∂y∂ai(a)+ a22y(a) ∂y
∂ai+ ⋅ ⋅ ⋅ + aiiy(a)i−1 ∂y
∂ai(a)+ yi(a)+ ⋅ ⋅ ⋅ + nyn−1(a) ∂y
∂ai(a) = 0,
∂y∂ai(a1 + 2a2y + 3a3y2 + ⋅ ⋅ ⋅ + nyn−1) + yi(a) = 0,
∂y∂ai(a)P ′(y(a)) + yi(a) = 0,
∂y∂ai(a) = − −y
i(a)P ′(y(a)) , (P ′(y(a)) ≠ 0).
(Kond y)(a) = ∥ −aiyi(a)
y(a)P ′(y(a))∥ =1
∣y(a)∣∣P ′(y(a))∣∥aiyi(a)∥ = ∑
n−1i=0 ∣ai∣∣y(a)∣i
∣y(a)∣∣P ′(y(a))∣ .
Zum Beispiel nehmen wirP (y) = ∏n
i=1(y − i) = a0 + a1y + ⋅ ⋅ ⋅ + an−1yn−1 + ynP (−1) = ∏n
i=1(−1 − i) = (−2)(−3) . . . (−(n + 1))= a0 + a1(−1) + . . . an−1(−1)n−1 + (−1)n
∣P (−1)∣ = (n + 1)! ≤ ∣a0∣ + ∣a1∣ + ⋅ ⋅ ⋅ + ∣an−1∣ + 1∑n−1
i=0 ∣ai∣ ≥ (n + 1)! − 1.P ′(1) = [(y − 1)(y − 2) . . . (y − n)]′(1)
= (y − 1)[(y − 2) . . . (y − n)]′(1) + (y − 2) . . . (y − n)(1)= (−1)n−1(n − 1)!
Kond(1) = ∑n−1i=0 ∣ai∣P ′(1) ≥
(n+1)!−1(n−1)! ≃ (n + 1)n ∼ n
2.
Beispiel 1.35 (Lineare Systeme) Sei A eine n ×m Matrix.
S ∶ x = Ay x ∈ Rn
14
1.3 Kondition
Gesucht ist y = y(x) ∈ Rn.Eine einzige Lösung existiert ⇐⇒ A ist invertierbar.y(x) = A−1x(Kond y)(x) = ∥x∥∥∇y(x)∥∥y(x)∥ = ∥x∥∥∇y(x)∥∥A−1x∥A−1 = (γij)
y(x) =⎛⎜⎝
y1(x)⋮
yn(x)
⎞⎟⎠=⎛⎜⎝
∑nj=1 γ1jxj
⋮∑n
j=1 γnjxj
⎞⎟⎠
∇y(x) = ( ∂yi
∂xj) = (γij) = A−1, (∂yi(x)
∂xj= γij).
(Kond y)(x) = ∥x∥∥A−1∥
∥A−1x∥ =∥Ay∥∥A−1∥∥y∥ ≤ ∥A∥∥A−1∥.
Kond(S) ∶= ∥A∥∥A−1∥, wobei ∥A∥ eine Operatornorm ist.
in Dimension 1: ∥A∥∥A−1∥ = ∣A∣ ∣ 1A∣ = 1,
A−1 rd(x)−A−1xA−1x
= rdx−xx
.
in Dimension 4: A =⎛⎜⎜⎜⎝
10 7 8 77 5 6 58 6 10 97 5 9 10
⎞⎟⎟⎟⎠.
A ist invertierbar, denn A =⎛⎜⎜⎜⎝
0 1 0 11 1 0 10 0 0 11 1 1 0
⎞⎟⎟⎟⎠.
Aij = Aij mod 2, detA = detA = 1 ≠ 0.y ∶= (1,1,1,1)⊺ x = Ay = (32,23,33,31)⊺δx ∶= (0.1,−0.1,0.1,0.1)⊺, Ay = x + δxy = (9.2,−12.6,4.5,−1.1)⊺y − y = (8.2,−13.6,3.5,−2.1)⊺
max ∣ yi−yi
yi∣ ≥ 2.1, ∣ (δx)i
xi∣ ≤ ∣0.1∣
23= 1
230,
230 ⋅ 2.1 ist die Fehlermultiplikation.
Beispiel 1.36 (Hermitesche Matrix)
Hn =⎛⎜⎜⎜⎝
1 1/2 1/3 ⋯ 1/n1/2 1/3 1⋮
1/n 1 . . . 1/(n − 1)
⎞⎟⎟⎟⎠
n = 10, Kond(Hn) = 1.6 ⋅ 1013.Das System ist schlecht konditioniert.
1.3.2. Die Kondition eines Algorithmus
f ∶ Rm → Rn
Wir möchten y = f(x) berechnen.Ein Algorithmus zur Berechnung von f(x) ist eine Funktion fA ∶ Rm
M → RnM .
Annahme: ∀x ∈ RmM ∃xA ∈ Rm s.d. fA(x) = f(xA) (*).
xA ist nicht zwingend eindeutig bestimmt.
15
Zahlendarstellung und Rundungsfehler
KondA(x) =minxA die (*) erfüllen
∥x−xA∥∥x∥eps ≤
∥x−xA∥∥x∥eps ∀xA, s.d. (*) gilt.
Beispiel 1.37 f(x) = x1x2 . . . xn, x = (x1, . . . , xn).Algorithmus A: p1 = x1, und pk = xk(1 − εk)pk−1, k = 2, . . . ,m.
fA(x) = x1(x2(1 − ε2)) . . . xn(1 − εn) = f(xA)mit xA = (x1, x2(1 − ε2), . . . , xn(1 − εn)).Es folgt x − xA = (0,−ε2x2, . . . ,−εnxn).
Wir wählen eine p-Norm.∥x∥p = (∑n
i=1 ∣xi∣p)1/p ,∥x − xA∥p = ∥(0,−ε2x2, . . . , εnxn)∥p
≤ (∑ ∣εi∣p∣xi∣p)1/p≤ eps∥x∥p.
Und (kondA)(x) ≤ eps∥x∥peps∥x∥p = 1.
Der Algorithmus ist gut konditioniert.
1.3.3. Globaler Fehler einer Berechnung
Wir werten f(x) aus, f eine Funktion von Rm nach Rn.Wir benutzen einen Algorithmus fA.
Der globale Fehler ist ∥fA(rd(x))−f(x)∥∥f(x)∥ .Wir setzen: x∗ = rd(x), y∗ = rd(y), y∗A = fA(x∗), y = f(x).
Wir nehmen an, dass∃x∗A mit fA(x∗) = f(x∗A).
Dann (KondA)(x) = ∥x∗A−x
∗∥∥x∗∥ eps−1.
Der globale Fehler ist: ∥y∗A−y∥∥y∥ ≤ ∥y
∗A−y
∗∥∥y∥ + ∥y
∗−y∥∥y∥ .
Sei ε = ∥x∗−x∥∥x∥ .
∥y∗−y∥∥y∥ = ∥f(x∗)−f(x)∥
∥f(x)∥
≤ (Kond f)(x) ⋅ ∥x∗−x∥∥x∥
= (Kond f)(x)ε∥y∗A−y
∗∥∥y∥ ≃ ∥y∗A−y
∗∥∥y∗∥
= ∥fA(x∗)−f(x∗)∥∥f(x∗)∥
= ∥f(x∗A)−f(x∗)∥
∥f(x∗)∥
≤ (Kond f)(x) ∥x∗A−x
∗∥∥x∗∥
≃ (Kond f)(x)(Kond)(A)eps,und∥y∗A−y∥∥y∥ ≤ (Kond f)(x)(ε +Kond(A)eps)
16
2 Auswertung elementarer Funktionen
2.1. Gewöhnliche Polynome
De nition 2.1 Ein Polynom ist eine Funktion des Types
P (x) = a0 + a1x + a2x2 + ⋅ ⋅ ⋅ + anxn,
wobei ai ∈ R, bzw C.Falls an ≠ 0, ist P vom GradmΠn = {P (x) = a0 + a1x + ⋅ ⋅ ⋅ + anxn∣ ai ∈ K = R oder C}.
Bemerkung 2.2 Πn ist ein Vektorraum der Dimension n + 1Eine Basis ist (1, x, . . . , xn).
2.1.1. Auswertung von P (ξ)
Um P (ξ) auszuwerten brauchen wir:anξ
n = anξξ . . . ξ n Multiplikationen.an−1ξ
n−1 = an−1ξξ . . . ξ n − 1 Multiplikationen.⋮a1ξ 1 Multilplikation
+n Additionen.Alles zusammen brauchenwir alson+(n−1)+⋅ ⋅ ⋅+1+n = (n+1)n
2+n ≃ n2 Operationen
für n gross.
Horner-Algorithmus
P (ξ) = a0 + a1ξ + a2ξ2 + ⋅ ⋅ ⋅ + anξn
= a0 + ξ(a1 + ξ(a2 + ⋅ ⋅ ⋅ + ξ(an − 1 + ξan)) . . . ).
Mit dieser Methode brauchen wir nMultiplikationen und n Additionen. Also insgesamt2n Operationen. 2n≪ n2
2und die Anzahl der Operationen ist viel kleiner für n gross.
Horner-Schema
bn−1 = anbi = ai+1 + ξbi+1 i = n − 2, . . . ,0
P (ξ) = a0 + ξb0
17
Auswertung elementarer Funktionen
De nition 2.3 Wir setzen
Pn(x) = P (x) = a0 + a1x + ⋅ ⋅ ⋅ + anxn =n
∑i=0
aixi
Pn−1(x) = b0 + b1x + ⋅ ⋅ ⋅ + bn−1xn−1 =n−1∑i=0
bixi
De nition 2.4 ξ ist eine Nullstelle von P falls P (ξ) = 0(ξ ∈ R bzw ξ ∈ C)
Beispiel 2.5 1+x2 besitzt keine Nullstelle inR aber es besitzt zwei Nullstellen inC, denn(1 + x2) = (x + i)(x − i).
Satz 2.6 Sei P ein Polynom
P (ξ) = 0 ⇐⇒ P (x) = (x − ξ)Q(x)
mitQ ein Polynom.
B
“⇐” klar.
“⇒” P (x) = P (ξ) + (x−ξ)1!
P ′(ξ) + ⋅ ⋅ ⋅ + (x−ξ)nP (n)(ξ)n!
ist Polynom von Grad n.P (ξ) = 0 Ô⇒ P (x) = (x−ξ)
1!P ′(ξ) + ⋅ ⋅ ⋅ + (x−ξ)
nP (n)(ξ)n!
= (x − ξ)Q(x).
De nition 2.8 Sei ξ eine Nullstelle von P .Die Vielfachheit von ξ ist k ⇐⇒ P (x) = (x − ξ)kQ(x)mit Q(ξ) ≠ 0.
Bemerkung 2.9 Leibniz-Formel
(f ⋅ g)(q) =q
∑i=0
f (i)g(q−1)(qi).
Satz 2.10 Sei ξ eine Nullstelle von P (ein Polynom).Die Vielfachheit von ξ ist k ⇐⇒ P (ξ) = P ′(ξ) = ⋅ ⋅ ⋅ = P (k−1)(ξ) = 0 undP (k)(ξ) ≠ 0.
B
“⇒” Die Vielfachheit von ξ ist k Ô⇒ P (x) = (x − ξ)kQ(x)mit Q(ξ) ≠ 0P (q)(x) = ((x − ξ)kQ(x))(q) = ∑q
i=0 (qi)((x − ξ)k)(i)Q(q−i)(x).
(x − ξ)k′ = k(x − ξ)k−1 und((x− ξ)k)i = k . . . (k− i+1)(x− ξ)k−i i ≤ k und ((x− ξ)k)(i) = 0wenn i > k.q = 1, . . . , k − 1 < kP (q)(x) = ∑q
i=0 (qi)k(k − 1) . . . (k − i + 1)(x − ξ)k−iQ(q−i)
18
2.1 Gewöhnliche Polynome
i ≤ q < k Ô⇒ k − i > 0P (q)(ξ) = ∑q
i=0 (qi)k . . . (k − i + 1)(ξ − ξ)k−iQ(q−i)(ξ) = 0.
P (k)(x) = ∑ki=0 (
ki)k(k − 1) . . . (k − i + 1)(x − ξ)k−iQ(k−i)(x)
= (kk)k . . .1Q(x) +∑k−1
i=0 (ki)k . . . (k − i + 1)(x − ξ)k−iQ(k−i)(x)
P (k)(ξ) = k!Q(ξ) ≠ 0.
“⇐” P (ξ) = P ′(ξ) = ⋅ ⋅ ⋅ = P (k−1)(ξ) = 0, P (k)(ξ) ≠ 0TaylorÔ⇒ P (a) = P (ξ) + (x − ξ)P ′(ξ) + . . .
+ (x−ξ)k−1P (k−1)(ξ)(k−1)! + (x−ξ)
kP (k)(ξ)k!
+ ⋅ ⋅ ⋅ + (x−ξ)nP (n)(ξ)n!
,
P (x) = (x−ξ)kP (k)(ξ)k!
+ ⋅ ⋅ ⋅ + (x−ξ)nP (n)(ξ)n!
,P (x) = (x − ξ)k(P
(k)
k!+ (x − ξ)R(x)) = (x − ξ)kQ(x), Q(ξ) = P (k)(ξ)
k!≠ 0.
Satz 2.12 Sei P ein Polynom vom Grad ≤ n. Falls P n+ 1 Nullstellen besitzt, so ist P ≡ 0.d.h. a0 = a1 = ⋅ ⋅ ⋅ = an = 0.
Bemerkung 2.13 Eine k-fache Nullstelle ist oben als k Nullstellen gez”ahlt.B P besitzt n + 1 Nullstellen, d.h.P (x) = (x − ξ1) . . . (x − ξn+1)Q(x) = (xn+1 + . . . )Q(x) Ô⇒ Q ≡ 0sonst kriegen wir in P einen Term akx
k mit k > nP (ξ) ≡ 0 ∀ξ Ô⇒ ak ≡ 0.Bemerkung 2.15 P (x) = a0 + a1x + ⋅ ⋅ ⋅ + akxk . . .
P = 0 ⇐⇒ a0 = a1 = ⋅ ⋅ ⋅ = an = 0P = 0 ⇐⇒ P (ξ) = 0 ∀ξ ∈ R bzw C
Sind ”aquivalent, dennP (ξ) = 0 ∀ξ Ô⇒ P ′(ξ) = P ′′(ξ) = ⋅ ⋅ ⋅ = P (k)(ξ) = 0 ∀sÔ⇒ P (k)(0) = 0 = akk! Ô⇒ ak = 0.
Satz 2.16 Sei Pn ∈ Πn vom Grad n. Sei Pn−1 das Polynom des HS-Verfahrens.Dann gilt
Pn(x) = Pn−1(x)(x − ξ) + Pn(ξ).
BPn−1(x)(x − ξ) + Pn(ξ) = ∑n−1
i=0 bixi(x − ξ) + Pn(ξ)
= ∑n−1i=1 xi(bi−1 − ξbi) − ξb0 + bn−1xn + Pn(ξ)
= ∑ni=0 aix
i (a0 = Pn(ξ) − ξb0)= Pn(x).
Satz 2.18 Sei Pn wie oben. Es gilt
P (k)n (ξ) = k!Pn−k(ξ)
(wobei Pn−k ist das Polynom des HS-Verfahrens für Pn−k+1.)
19
Auswertung elementarer Funktionen
B Mittels Induktion.Aus Satz 2.16:Pn(x) = Pn−1(x)(x − ξ) + Pn(ξ)P ′n(x) = P ′n−1(x)(x − ξ) + Pn−1(x)P ′n(ξ) = Pn−1(ξ)P ′′n (x) = P ′′n−1(x)(x − ξ) + 2P ′n−1(x)P ′′n (ξ) = 2P ′n−1(ξ) = 2Pn−2(ξ)P(k)n (x) = P
(k)n−1(x − ξ) + kP
(k−1)n−1 (x)
Nach Induktion ist die Formel gültig fürP(k−1)n (x) = P (k−1)n−1 (x)(x − ξ) + (k − 1)P
(k−2)n−1 (x)
P(k−1)n (x)′ = P
(k−1)n−1 (x)′(x − ξ) + P
(k−1)n−1 (x) + (k − 1)P
(k−2)n−1 (x)′
P(k)n = P
(k)n−1(x)(x − ξ) + kP
(k−1)n−1 (x)
P(k)n (ξ) = kP
(k−1)n−1 (ξ)
= k(k − 1)P (k−2)n−2 (ξ)= k(k − 1)(k − 2) ⋅ . . . ⋅ 1Pn−k(ξ)= k!Pn−k(ξ).
Bemerkung 2.20 Die Berechnung von Pn(ξ) benötigt 2n OperationenP ′(ξ)1!= Pn−1(ξ) 2(n − 1)
⋮P (k)
k!= Pn−k(ξ) 2(n − k)
Um P (ξ), P ′(ξ), . . . , k!P (k!)(ξ) zu berechnen brauchen wir2(n + n − 1 + ⋅ ⋅ ⋅ + n − k)2((k + 1)n − (1 + ⋅ ⋅ ⋅ + k)) = 2((k + 1)n − (k + 1)k
2) = (k + 1)(2n − k
2) Operationen.
Beispiel 2.21 P (x) = 2 + x + x2 + x3, ξ = 2bn−1 = anbi = ai+1 + ξbi+1 i = n − 2, . . . ,0P (ξ) = a0 + ξb01 1 1 2 Pn
1 3 7 16=Pn(2) Pn−11 5 17=Pn−1(2) = P ′(2)
1 7=Pn−2(2) = P ′′(2)2!
P ′′(2) = 14.
P (x) = a0 + ⋅ ⋅ ⋅ + anxn
Algorithmus um P (k)(ξ)/k! in ak zu speichern:
a[j] := a_jfor k := 0 to n-1 do
for j := n-1 to k doa[j] := a[j] + ξ * a[j+1]
endend
20
2.2 Trigonometrische Polynome
an an−1 a1 a0an = bn−1 bn−2 b0 Pn(ξ) k = 0
P ′n(ξ)1!
Pn(ξ) k = 1⋮
P (n)(ξ)n!
2.2. Trigonometrische Polynome (FFT - Fast Fourier Transform)
De nition 2.22 Ein Trigonometrisches Polynom ist eine Funktion des Types
tN(x) =a02+
N
∑j=1
aj cos(jx) + bj sin(jx)
N ist der “Grad” des Polynomes.
Falls ai = 0 ∀i ist tN ein reines Sinus-Polynom.Falls bi = 0 ∀i ist tN ein reines Cosinus-Polynom.
Bemerkung 2.23 tN ist periodisch mit Periode 2πtN(x + 2π) = tN(x) ∀x ∈ RBemerkung 2.24 Sei f periodisch mit Periode p, falls f auf [a, a + p] de niert ist, so istf auf ganz R eindeutig bestimmt.
De nition 2.25 Sei xk = 2kπN
k = 0, . . . ,N − 1tN ∣{xk} heisst diskretes Fourier Polynom.
Bemerkung 2.26 Um die tN(xk) zu berechnen brauchen wir 2N ⋅N = 2N2 Multipli-kationen.Berechnung von tN(xk)mit complexen ExponentialenSei i mit i2 = −1.exp(jxi) = cos(jx) + i sin(jx) = ∑∞k=0
(jxi)kk!
. Wir setzenγj = aj − ibj .
γj ⋅ exp(jxki) = (aj − ibj)(cos(jxk) + i sin(jxk))= aj cos(jxk) + bj sin(jxk) + i(. . . ).
aj cos(jxk) + bj sin(jxk) = Re(γj exp(jxki))= 1
2(γj exp(jxki) + γj exp(−jxki))
= 12(γj exp(jxki) + γj exp((N − j)xki)), da
exp((N − j)xki) = exp(−jxki +Nxki)= exp(−jxki + 2πki)= cos(−jxk + 2πk) + i sin(−jxk + 2πk)= cos(jxk) + i sin(−jxk)= exp(−jxki).
tN(xk) = a0
2+∑N−1
i=1 aj cos(jxk) + bj sin(jxk) + aN cos(Nxk) + bN sin(Nxk)= a0
2+∑N−1
j=112(γj exp(jxki) + γj exp((N − j)xki)) + aN
= a0
2+ aN +∑N−1
j=1 cj exp(jxki)= ∑N−1
j=0 cj exp(jxki),
21
Auswertung elementarer Funktionen
wobei cj =γj+γN−j
2für j ≠ 0, c0 = a0
2+ aN
Der FFT-AlgorithmusWir setzen wN = exp( 2πiN
).exp(jxki) = exp( j2kπiN
) = wjkN und
tN(xk) = ∑N−1j=0 cjw
jkN . Sei FN die Matrix de niert mit
FN = (wjkN )j,k=0...N−1
FN =
⎛⎜⎜⎜⎜⎝
1 1 1 . . . 11 wN w2
N wN−1N
⋮ ⋮1 wN−1
N w2(N−1)N . . . w
(N−1)2N
⎞⎟⎟⎟⎟⎠
.
Wir bemerken, dassws
N = wrN ⇐⇒ r = s + lN l ∈ N.
tN(xk) = ∑N−1j=0 cjw
jkN ∀k = 0, . . . ,N − 1, cj =
γj+γN−j2
, γj = aj − ibj .
Es gilt⎛⎜⎝
tN(x0)⋮
tN(xN−1)
⎞⎟⎠= FN
⎛⎜⎝
c0⋮
cN−1
⎞⎟⎠,
da die i-te Komponente von FNC ist∑N−1k=0 wik
N ck = ∑N−1k=0 ckw
kiN = tN(xi).
Die Anzahl der Multiplikationen ist N2. Wir betrachten den Fall N = 2Mw2
N = (exp 2πiN)2 = exp 4πi
2M= exp 2πi
M= wM und
tN(x2l) = ∑N−1j=0 cjw
2ljN
= ∑N−1j=0 cjw
ljM
= ∑M−1j=0 (cj + cj+M)w
ljM
da wl(j+M)M = wlj
MwlMM = wlj
M .
tN(x2l+1) = ∑N−1j=0 cjw
(2l+1)jN = ∑N−1
j=0 cjwljMwj
N = ∑M−1j=0 (cj − cj+M)w
ljMwj
N
da wl(j+M)M w
(j+M)N = wlj
MwjNwM
N = −wljMW j
N .
Wir haben jetzttN(x2l) = ∑M−1
j=0 (cj + cj+M)wljM , l = 0, . . . ,M − 1,
tN(x2l+1) = ∑M−1j=0 (cj − cj+M)w
jNwlj
M , l = 0, . . . ,M − 1.
Sei N = 2nSei ∆N die Anzahl der Multiplikationen der FFT-Methode.∆N =M + 2∆M
22
2.2 Trigonometrische Polynome
∆2n = 2n−1 + 2∆2n−1
= 2n−1 + 2(2n−2 + 2∆2n−2)= n2n−1
= n2n
2
= nN2
= N lnN2 ln2
≃ N lnN.
Es ist viel besser als N2.
23
3 Iterative Verfahren
Sei f ∶ R→ R, (oder auch Rn → Rm). Man versucht die Gleichung f(x) = 0 zu lösen.
3.1. Vorbereitung
Satz 3.1 Sei f ∶ [a, b]→ [a, b] eine stetige Funktion.Falls f(a) ⋅ f(b) < 0 so exisitert x ∈ (a, b)mit f(x) = 0.Falls f streng monoton ist, ist die Lösung dieser Gleichung eindeutig bestimmt.
B f(a) ⋅ f(b) < 0, d.h. f(a), f(b) haben nicht das selbe Vorzeichen.Aus dem Zwischenwertsatz folgt ∃x ∈ (a, b)mit f(x) = 0.Beispiel 3.3 f(x) = x + 0.5 cosx − 1f(0) = −0.5 < 0, f(π) = π − 1.5 > 0∃x ∈ (0, π)mit f(x) = 0.f ′(x) = 1 − 0.5 sinx > 0 also ist f monoton.Ô⇒ die Nullstelle ist eindeutig bestimmt.
3.2. Algorithmus
3.2.1. Bisektions Methode
f(a)f(b) ≤ 0 Ô⇒ ∃x ∈ [a, b]mit f(x) = 0.Wir berechnen f(a)f(a+b
2) ≤ 0 dann ist in [a, a+b
2] eine Lösung.
Sonst f(a)f(a+b2) > 0 und
f(b)f(a+b2) = f(b)f(a+b
2) ⋅ f(a)f(
a+b2 )
f(a)f( a+b2 )= f( a+b2 )
2f(a)f(b)f(a)f( a+b2 )
≤ 0.Programm
a_0 = a;b_0 = b;für n = 0,1,2, ... tunfalls f(a_n)f((a_n+b_n)/2) <= 0dann a_(n+1) = a_n, b_(n+1) = (a_n+b_n)/2sonst a_(n+1) = (a_n+b_n)/2, b_(n+1) = b_n
24
3.2 Algorithmus
Satz 3.4 Die obige Folge konvergiert gegen eine Lösung der Gleichung f(x) = 0.
B an+1 ≥ an ∀n und bn+1 ≤ bn ∀n. Zusätzlich gilta ≤ an, bn ≤ b und beide Folgen sind beschränkt, beide konvergieren.∃α,β ∈ [a, b]mit an → α, bn → βbn+1 − an+1 = 1
2(bn − an) = ⋅ ⋅ ⋅ = 1
2n+1(b0 − a0) = 1
2n+1(b − a)→ 0
Ô⇒ bn − an → 0 Ô⇒ β − α = 0.Bemerkung 3.6 Wir können das Verfahren beenden, wenn 1
2n(b − a) ≤ ε, wobei ε die
Gewünschte Genauigkeit ist.
3.2.2. Fixpunkt Methode
De nition 3.7 Ein Fixpunkt ist ein Punkt x mit f(x) = x.
Bemerkung 3.8 Die Gleichung f(x) = 0 kann als eine Fixpunkt Gleichung geschriebenwerden.z.B. f(x) = 0 ⇐⇒ x + f(x) = x, g(x) = f(x) + x.Programm
x_0fuer n=1,2,... x_(n+1) = g(x_n)
Satz 3.9 Falls die obige Folge konvergiert und g stetig ist, konvergiert diese Folge gegen einenFixpunkt von g.
B xn → x∞ (die Folge konvergiert gegen x∞).xn+1 = g(xn)→ g(x∞)Ô⇒ x∞ = g(x∞) und der Limes ist ein Fixpunkt.Bemerkung 3.11 Der Satz gilt auch für g ∶ Rn → Rn.
Satz 3.12 Sei g ∶ [a, b]→ [a, b] differenzierbar mit ∣g′(x)∣ ≤ k < 1.Dann besitzt g einen einzigen Fixpunkt in [a, b] und die Folgex0 fest, xn+1 = g(xn) , n ≥ 0konvergiert gegen diesen Fixpunkt.
B g ∶ [a, b]→ [a, b].Betrachte x − g(x)Ô⇒ a − g(a) ≤ 0 und b − g(b) ≥ 0.x − g ist stetig Ô⇒ ∃x ∈ [a, b]mit g(x) = x.(x − g(x))′ = 1 − g′(x) > 0 Ô⇒ die Nullstelle von x − g(x) ist eindeutig bestimmt.Sei x∞ der einzige Fixpunkt von g auf [a, b]x0 gegeben in [a, b], xn+1 = g(xn). Dann gilt∣xn − x∞∣ = g(xn−1) − g(x∞) = g′(ξ)(xn−1 − x∞)mit ξ ∈ (xn−1, x∞) ⊂ [a, b].∣xn−x∞∣ = ∣g′(ξ)∣∣xn−1−x∞∣ ≤ k∣xn−1−x∞∣ ≤ ⋅ ⋅ ⋅ ≤ kn∣x0−x∞∣→ 0, wenn n→ +∞.Beispiel 3.14 g ∶ [0, b]→ [0, b], g(x) = x2. g′(x) = 2x.Also wenn b < 1
2, dann ∣g′(x)∣ ≤ 2b < 1 auf [0, b].
Wenn x0 > 1, rekursiv xn+1 = g(xn) = x2n > xn > 1 und (xn) divergiert.
25
Iterative Verfahren
De nition 3.15 Sei D ⊂ Rn abgeschlossen.g ∶D →D so, dass es ein k < 1 gibt mit
∣g(x) − g(y)∣ ≤ k∣x − y∣
∀x, y ∈D.Ein solches g heisst Kontraktion.
Satz 3.16 SeiD ⊂ Rn abgeschlossen und g ∶D →D eine Kontraktion.Dann hat g einen eindeutigen Fixpunkt x∞ ∈D und die Folgex0 ∈D, xn+1 = g(xn) ∀n ≥ 0 konvergiert gegen x∞
BEindeutigkeit:Annahme: Seien x∞, x
′∞ Fixpunkte
Ô⇒ g(x∞) = x∞, g(x′∞) = x′∞Ô⇒ ∥g(x∞) − g(x′∞)∥ = ∥x∞ − x′∞∥ ≤ k∥x∞ − x′∞∥Ô⇒ ∥x∞ − x′∞∥ = 0 da k < 1Ô⇒ x∞ = x′∞.
Existenz:Ziel: Wir wollen zeigen, dass (xn)n∈N eine Cauchy Folge ist.∥x2 − x1∥ = ∥g(x1) − g(x0)∥ ≤ k∥x1 − x0∥,∥x3 − x2∥ = ∥g(x2) − g(x1)∥ ≤ k∥x2 − x1∥ ≤ k2∥x1 − x0∥Rekursiv folgt ∥xn+1 − xn∥ ≤ kn∥x1 − x0∥.
Sei ε > 0. Für n ≥ 0, p ∈ N, p ≠ 0∥xn+p − xn∥ = ∥xn+p − xn+p−1 + xn+p−1 − xn+p−2 + ⋅ ⋅ ⋅ + xn+1 − xn∥
≤ ∥xn+p − xn+p−1∥ + ⋅ ⋅ ⋅ + ∥xn+1 − xn∥≤ (kn+p−1 + ⋅ ⋅ ⋅ + kn)∥x1 − x0∥= kn 1−kp
1−kDa k < 1 ∃n0 so dass kn0
1−k ∥x1 − x0∥ < ε.Ô⇒ ∀n ≥ n0, ∀p ∈ N, p ≠ 0Ô⇒ ∥xn+p − xn∥ < kn
1−k ∥x1 − x0∥ ≤ kn0
1−k ∥x1 − x0∥ < ε.Das heisst (xn) ist eine Cauchy FolgeÔ⇒ ∃x∞ ∈ Rn so dass xn → x∞Ô⇒ x∞ ∈D (da D abgeschlossen ist).Bemerkung 3.18
1. Satz 3.16 gilt auch, wenn wir Rn durch einen Banachraum ersetzen.
2. Fehlerabschätzung: ∥xn − x∞∥ ≤ kn
1−k ∥x1 − x0∥
3.2.3. Newton Methode
Sei x mit f(x) = 0. Aus der Taylorschen Entwicklung folgt:
f(x + h) = f(x) + f ′(x)h + o(h).
26
3.2 Algorithmus
Für x = xn, h = xn+1 − xn gilt:f(xn+1) = f(xn + xn+1 − xn) = f(xn) + f ′(xn)(xn+1 − xn) + o(xn+1 − xn)
Berechne xn+1 als Lösung von f(xn) + f ′(xn)(xn+1 − xn) = 0Ô⇒ xn+1 = xn − f(xn)
f ′(xn) .
Gleichung einer Tangente an f in (xn, f(xn)):y = y(x) = f(xn) + f ′(xn)(x − xn)Program
x_0 gegebenfor n = 0,1,..x_(n+1) = x_n - f(x_n)/f'(x_n)
Satz 3.19 Sei f ∶ I → R, f ∈ C2(I), mit I ⊂ R Intervall.Sei x∞ ∈ I so, dass f(x∞) = 0 und f ′(x∞) ≠ 0Dann ist die Folge
x0 ∈ (x∞ − ε, x∞ + ε), xn+1 = xn −f(xn)f ′(xn)
∀n ≥ 0
für genug kleine ε wohlde niert und konvergiert gegen x∞.
B Sei g(x) = x − f(x)f ′(x)
∣g′(x)∣ = ∣1 − (f′(x))2−f(x)f ′′(x)(f ′(x))2 ∣ = ∣ f(x)f
′′(x)(f ′(x))2) ∣ ≤
12für x ∈ [x∞ − ε, x∞ + ε] und ε klein
genug.
Da f , f ′, f ′′ stetig und f ′(x∞) ≠ 0 ∃r, so dass ∣f ′(x)∣ ≥ ∣f′(x0)∣2
∀x ∈ [x∞ − r, x∞ + r] und ∃M so dass ∣f ′′(x)∣ ≤M ∀x ∈ [x0 − r, x0 + r].
Da f stetig auf x∞ und f(x∞) = 0 Ô⇒ ∃ε > 0, ε ≤ r, so dass∣f(x)∣ < ∣f
′(x∞)∣28M
auf [x∞ − ε, x∞ + ε].
Also gilt auf [x∞ − ε, x∞ + ε] ⊂ [x∞ − r, x∞ + r]∣g′(x)∣ = ∣ f(x)f
′′(x)f ′(x)2 ∣ = ∣f(x)∣
∣f ′′(x)∣f ′(x)2 < ∣f(x)∣
2Mf ′(x∞)2 <
f ′(x∞)28M
M
( f ′(x∞)2 )
2 = 12.
∀x ∈ [x∞ − ε, x∞ + ε] gilt∣g(x) − x∞∣ = ∣g(x) − g(x∞)∣
= ∣g′(x)(x − x∞)∣≤ ∣g′(c)∣∣x − x∞∣≤ 1
2∣x − x0∣
< 12ε
< εc ist zwischen x und x∞, also c ∈ [x∞ − ε, x∞ + ε].
Daher ist g ∶ [x∞ − ε, x∞ + ε] → [x∞ − ε, x∞ + ε] eine C1-Funktion und ∣g′(x)∣ ≤ 12
auf (x∞ − ε, x∞ + ε)Ô⇒ (Satz 3.16) xn → x∞.
27
Iterative Verfahren
Satz 3.21 Sei f ∶ Rn → Rn so, dass f ∈ C2(Rn) und sei x∞ ∈ Rn so, dass f(x∞) = 0und∇f(x∞) invertierbar. Dann existiert ε > 0 so dass die Folgex0 ∈ Dε ∶= {x ∈ Rn∣ ∥x − x∞∥ ≤ ε},xn+1 = xn − [∇f(xn)]−1f(xn)
wohlde niert ist und gegen x∞ konvergiert.
B g(x) = x − [∇f(x)]−1f(x).De niere L ∶ V (x∞)→Mn×n, L(x) = [∇f(x)]−1, wobei V (x∞) eine Umgebung vonx∞ ist, auf der∇f invertierbar ist.
g(x) = x −L(x)f(x) ⇐⇒ gi(x) = xi −∑nk=1Lik(x)fk(x) ∀i ∈ {1,2, . . . , n}.
∂gi∂xj= δij −∑n
k=1 (∂Lik
∂xj(x)fk(x) +Lik(x)∂fk∂xj
(x)) ∀i, j ∈ {1,2, . . . , n}.
∇g(x) = In −L(x)∇f(x) −M(x), wobei M(x) = (∑nk=1
∂Lik
∂xj(x)fk(x))
i,j∈{1,...,n}.
Daher∇g(x) = −M(x) und∇g(x∞) = −M(x∞) = 0 undM ∶ V (x∞)→Mn×n ist stetig, da L auf V (x∞) C1.
Sei ∥ ⋅ ∥ die Operatornorm auf Mn×n assoziiert zur Euklidischen Norm in Rn.Wir haben ∥x∥↦ ∥M(x)∥ ist stetig und ∥M(x∞)∥ = 0.Ô⇒ ∃ε > 0 so, dass ∥M(x)∥ = ∥∇g(x)∥ ≤ 1
2∀x ∈Dε.
Behauptung: g ist eine Kontraktion auf Dε.Beweis: ∥g(x) − g(y)∥ = ∥ ∫
10
ddtg(x + t(y − x))dt∥
= ∥ ∫10 ∇g(x + t(y − x))(y − x)dt∥
≤ ∫10 ∥∇g(x + t(x − y))(y − x)∥dt
≤ ∫10 ∥∇g(x + t(x − y))∥dt∥y − x∥
∀x, y ∈Dε
∥g(y) − g(x)∥ ≤ ∫10 ∥∇g(x + t(y − x)∥dt∥y − x∥ ≤ ∫
10
12∥y − x∥ = 1
2∥y − x∥,
da ∥∇g(x + t(y − x))∥ ≤ 12∀t ∈ [0,1].
∀x ∈Dε ∥g(x) − x∞∥ = ∥g(x) − g(x∞)∥ ≤ 12∥x − x∞∥ ≤ 1
2ε ≤ εÔ⇒g(x) ∈Dε.
Somit ist g ∶ Dε → Dε eine Kontraktion, Dε ist kompakt Ô⇒ (xn) konvergiert gegenden eindeutigen Fixpunkt von g in Dε, also gegen x∞.Bemerkung 3.23 if x, y ∈Dε Ô⇒ x + t(y − x) ∈Dε ∀t ∈ [0,1]since∥x + t(y − x) − x∞∥ = ∥(1 − t)(x − x∞) + t(y − x∞)∥
≤ (1 − t)∥x − x∞∥ + t∥y − x∞∥≤ (1 − t)ε + tε= ε
Beispiel 3.24 f ∶ C→ C, f(z) = ez − zz = x + iy ez = exeiy = ex(cos y + i sin y)f(z) − z = (ex cos y − x) + i(ex sin y − y)
f(z) = 0⇐⇒ { ex cos y − x = 0ex sin y − y = 0
f ∶ R2 → R2
f(x, y) = (ex cos y − x, ex sin y − y)
28
3.2 Algorithmus
∇f(x, y) = (ex cos y − 1 −ex sin yex sin y ex cos y − 1)
det(∇f(x, y)) = (ex cos y − 1)2 + e2x sin2 y= e2x cos2 y − 2ex cos y + 1 + e2x sin2 y= e2x − 2ex cos y + 1 = (ex − cos y)2 + sin2 y > 0
3.2.4. Sekanten Methode
Die Newton Methode ist gegeben mit:xn+1 = xn − f(xn)
f ′(xn) .Ist f nicht differenzierbar, so kann man f ′(xn)mit f(xn−1)−f(xn)
xn−1−xnersetzen.
Die Methode lautetx−1, x0, x−1 ≠ x0 gegebenxn+1 = xn − xn−1−xn
f(xn−1)−f(xn)f(xn) n ≥ 0
Die Gerade durch die Punkte (xn−1, f(xn−1)), (xn, f(xn)) ist:
t↦ (xn, f(xn)) + t(xn−1 − xn, f(xn−1) − f(xn))
f(xn) + t(f(xn−1) − f(xn)) = 0 Ô⇒ t = − f(xn)f(xn−1)−f(xn)
Ô⇒ der Schnittpunkt mit y = 0 ist (xn − f(xn) xn−1−xn
f(xn−1)−f(xn) ,0).
Program
x_(-1), x_(0) gegeben, x_(-1) != x_(0)if x_(n-1) = x_n thenstop
elsex_(n+1) = x_n - (x_(n-1)-x_(n))/(f(x_(n-1))-f(x_(n)))*f(x_(n))
Satz 3.25 Sei f ∶ Ix∞ → R, mit Ix∞ offenes Intervall inR, das x∞ enthält. f ∈ C2(Ix∞).Sei weiter f(x∞) = 0, f ′(x∞) ≠ 0Dann ist für ε klein genug die Folgex−1, x0 ∈ [x∞ − ε, x∞ + ε], x−1 ≠ x0,xn+1 = xn − xn−1−xn
f(xn−1)−f(xn)f(xn) wenn xn+1 ≠ xn und xn+p = xn ∀p ∈ N wennxn−1 = xn
wohlde niert und konvergiert gegen x∞.
Bemerkung 3.26 Die Folge ist de niert für f(xn−1) ≠ f(xn)Falls xn−1, xn nah von x∞ sind, giltf(xn−1) − f(xn) = f ′(ηn)(xn−1 − xn), ηn ∈ (xn−1, xn)f ′(ηn) ≠ 0 für xn−1, xn nah von x∞. (f ′(x∞) ≠ 0)Für xn−1, xn nah von x∞, f(xn−1) − f(xn) ≠ 0 ⇐⇒ xn−1 ≠ xn.Sonst xn−1 = xn Ô⇒ xn = xn−1 − f(xn−1) (xn−1−xn−2)
f(xn−1)−f(xn−2)(n ist die “erste” Zahl mit xn = xn−1)Ô⇒ f(xn−1) = 0 Ô⇒ f(xn) = 0 Ô⇒ wir haben die Nullstelle erreicht.B Wir nehmen an, dass xn ≠ xn−1 ∀nfalls xn, xn−1 nah von x∞ sind, ist die Folge de niert.
29
Iterative Verfahren
f ′(x∞) ≠ 0∃δ > 0 mit ∣f ′(x)∣ ≥m > 0, ∣f ′′(x)∣ ≤M ∀x ∈ [x∞ − δ, x∞ + δ]und x∞ ist die einzige Nullstelle von f in [x∞ − δ, x∞ + δ].
Sei ε = δ ∧ mM=min(δ, m
M)
Behauptung:Für x−1, x0 ∈ [x∞ − ε, x∞ + ε], gilt xn ∈ [x∞ − ε, x∞ + ε] ∀n und xn → x∞Beweis:Nach Induktionsannahme x−1, x0, . . . , xn ∈ [x∞ − ε, x∞ + ε]xn+1
?∈ [x∞ − ε, x∞ + ε]
• xn = xn−1 Ô⇒ f(xn) = 0 xn = x∞Wir haben den Limes erreicht.
• xn ≠ xn−1 Ô⇒ f(xn) ≠ f(xn−1)xn+1 = xn − f(xn) (xn−1−xn)
f(xn−1)−f(xn)
= xn(f(xn−1)−f(xn))−f(xn)(xn−1−xn)f(xn−1)−f(xn)
= xnf(xn−1)−f(xn)xn−1f(xn−1)−f(xn) .
xn+1 − x∞ = xnf(xn−1)−xn−1f(xn)−x∞(f(xn−1)−f(xn))f(xn−1)−f(xn)
= (xn−x∞)f(xn−1)−(xn−1−x∞)f(xn)f(xn−1)−f(xn)
= (xn−x∞)(f(xn−1)−f(x∞))−(xn−1−x∞)(f(xn)−f(x∞))f(xn−1)−f(xn)
= (xn−x∞)(xn−1−x∞)f(xn−1)−f(xn) (
f(xn−1)−f(x∞)xn−1−x∞ − f(xn)−f(x∞)
xn−x∞ ).
Wir setzen g(x) = f(x)−f(x∞)x−x∞
xn+1 − x∞ = (xn−x∞)(xn−1−x∞)f(xn−1)−f(xn) (g(xn−1) − g(xn))
= (xn − x∞)(xn−1 − x∞)g′(ξn) (xn−1−xn)f(xn−1)−f(xn)
für ξn ∈ (xn−1, xn) ⊂ [x∞ − ε, x∞ + ε].
xn+1 − x∞ = (xn − x∞)(xn−1 − x∞)g′(ζn) xn−1−xn
f ′(ξn)(xn−1−xn)
= (xn − x∞)(xn−1 − x∞) g′(ξn)
f ′(ζn)für ξn, ζn ∈ (xn−1, xn) ⊂ [x∞ − ε, x∞ + ε].
Es gilt g′(x) = f ′(x)(x−x∞)−(f(x)−f(x∞))(x−x∞)2 .
Da f(x∞) = f(x) + f ′(x)(x∞ − x) + f ′′(γ)2!(x∞ − x)2
Ô⇒f(x∞) − f(x) + f ′(x)(x − x∞) = f ′′(γ)2!(x∞ − x)2,
gilt g′(x) = f ′′(γ)2
, γ ∈ (x∞, x)g′(ξn) = f ′′(γn)
2, γn ∈ (x∞, ξn) ⊂ [x∞ − ε, x∞ + ε] ⊂ [x∞ − δ, x∞ + δ]
∣g′(ξn)∣ ≤ ∣f′′(γn)∣2≤ M
2.
∣xn+1 − x∞∣ = ∣xn − x∞∣∣xn−1 − x∞∣ ∣g′(ξn)∣∣f ′(ζn)∣
≤ ∣xn − x∞∣∣xn−1 − x∞∣ M2∣f ′(ζn)∣
≤ ∣xn − x∞∣∣xn−1 − x∞∣ M2m .
30
3.2 Algorithmus
∣xn+1 − x∞∣ ≤ 12∣xn − x∞∣∣xn−1 − x∞∣Mm , xn, xn−1 ∈ [x∞ − ε, x∞ + ε]
∣xn+1 − x∞∣ ≤ 12∣xn − x∞∣εM
m∣xn+1 − x∞∣ ≤ 1
2∣xn − x∞∣ ≤ 1
2ε < ε.
Zusätzlich gilt:∣xn+1 − x∞∣ ≤ 1
2∣xn − x∞∣ ≤ ( 12)
2 ∣xn−1 − x∞∣≤ (1
2)n ∣x1 − x∞∣→ 0 für n→ +∞.
31
Iterative Verfahren
3.2.5. Geschwindigkeit der Algorithmen
De nition 3.28 Sei xn eine Folge. DieOrdnung der Konvergenz dieser Folge gegen x∞ist die Zahl p mit
lim∣xn+1 − x∞∣∣xn − x∞∣p
= C ≠ 0.
Die Ordnung der Konvergenz dieser Folge ist mindestens pfalls ∣xn+1 − x∞∣ ≤ C ∣xn − x∞∣p , C ≠ 0, gilt für n gross.Beispiel 3.29 p die Ordnung der Konvergenz ist mindestens 10.∣xn − x∞∣ ≤ 1
10
∣xn+1 − x∞∣ ≤ C ∣xn − x∞∣10 ≤ C ( 110)10
∣xn+2 − x∞∣ ≤ C ∣xn+1 − x∞∣10 ≤ C11 ( 110)100.
Zum Beispiel:
1. Bisektionsmethode∣bn − an∣ = 1
2∣bn−1 − an−1∣
die Konvergenz der Länge der Intervalle ist linear.
2. Fixpunkt-Algorithmusx = g(x)x0 fest, xn+1 = g(xn)∣xn+1 − xn∣ = ∣g(xn) − g(x∞)∣ = ∣g′(η)∣∣xn − x∞∣da g(x) = g(x∞) + g′(η)(x − x∞).Die Konvergenz ist mindestens linear.Für g′(x∞) ≠ 0 ist die Folge linear konvergent, da∣xn+1−x∞∣∣xn−x∞∣ = g
′(ηn), ηn ∈ (x∞, xn), g′(ηn)→ g′(x∞) für n→ +∞
g′(x∞) = 0Ô⇒C = 0.∣xn+1 − x∞∣ = ∣g(xn) − g(x∞)∣ = ∣g
′′(ηn)∣2!∣xn − x∞∣2 da
g(xn) = g(x∞) + g′(x∞)(x − x∞) + g′′(ηn)2!(xn − x∞)2,
g(xn) − g(x∞) = g′′(ηn)2!(xn − x∞)2.
Die Konvergenz ist mindestens quadratisch p = 2 falls g′′(x∞) ≠ 0.
Zum Beispiel für die Newton Methode:
xn+1 = xn −f(xn)f ′(xn)
= g(xn)
mit g′(x∞) = f ′′ff ′2(x∞) = 0. Die Methode ist mindestens quadratisch.
3. Sekanten MehthodeDie Ordnung dieses Verfahrens ist mindestensp = 1+
√5
2> 1+
√4
2= 3
2= 1.5.
32
3.2 Algorithmus
p ist die positive Nullstelle vonp2 − p − 1, dap2 − p − 1 = (p − 1
2)2 − 1
4− 1 = (p − 1
2)2 − 5
4= (p − 1
2−√52)(p − 1
2+√52).
Es folgt von oben:xn+1−x∞ = (xn−x∞)(xn−1−x∞) 12
f ′′(ξn)f ′(ζn) , ξn, ζn ∈ [x∞−ε, x∞+ε], ξn, ζn →
x∞.
∣xn+1 − x∞∣ = ∣(xn − x∞)∣∣(xn−1 − x∞)∣Cn, Cn → C .
Sei yn = ∣xn−x∞∣∣xn−1−x∞∣p mit p die positive Nullstelle von p2 − p − 1 = 0
p2 = p + 1 oder p = 1 + 1p.
yn+1 = ∣xn+1−x∞∣∣xn−x∞∣p =
∣xn−x∞∣∣xn−1−x∞∣∣xn−x∞∣p Cn = ∣xn−1−x∞∣
∣xn−x∞∣p−1Cn = Cny− 1
pn , da
y− 1
pn = ( ∣xn−x∞∣
∣xn−1−x∞∣p )− 1
p = ∣xn−1−x∞∣
∣xn−x∞∣1p= ∣xn−1−x∞∣∣xn−x∞∣p−1 .
Wir betrachten die Folge y0 gegeben; yn+1 = Cny− 1
pn .
yn konvergiert (siehe unten). Sei l der Limes. Dann giltl = Cl−
1p , l1+
1p = C , lp = C , l = C
1p .
an = ln ynln(yn+1) = ln(Cny
− 1p
n ) = lnCn − 1pln yn
an+1 + 1pan = lnCn = Ln, a0 fest.
Behauptung: Die Folge an konvergiert (i.e. die Folge yn = ean konvergiert, unddann die Ordnung der Konvergenz von xn ist mindestens p = 1+
√5
2).
Beweis:an+1 = Ln − 1
pan
a1 = L0 − 1pa0
a2 = L1 − 1p(L0 − 1
pa0) = L1 − 1
pL0 + 1
p2 a0a3 = L2 − 1
p(L1 − 1
pL0 + 1
p2 a0) = L2 − 1pL1 + 1
p2L0 − 1p3 a0
an = ∑n−1k=0 (− 1
p)kLn−k−1 + (− 1
p)na0.
Wir setzen An−1 = ∑n−1k=0 (− 1
p)kLn−k−1.
Ln = lnCn konvergiert Ô⇒ ∃L mit ∣Ln∣ ≤ L ∀n, Ln → l
An ist eine Cauchy Folge, d.h.∣An+m −An∣ ≤ ε′ für n gross genug, ∀m.
33
Iterative Verfahren
Beweis:∣An+m −An∣ = ∑n+m
k=0 (− 1p)kLn+m−k −∑n
k=0 (− 1p)kLn−k
= ∑⌊n/2⌋k=0 (−1p)k(Ln+m−k −Ln−k) +∑k>⌊n/2⌋ (− 1
p)kLn+m−k
−∑k>⌊n/2⌋ (− 1p)kLn−k
Da Ln konvergiert gilt
- Ln ist beschränkt Ô⇒ ∣Ln∣ ≤M ∀n- Ln ist eine Cauchy Folge.∣Lk+m −Lk ∣ ≤ ε für k > k0(ε)
Wir wählen n2≥ k0 dann
∣An+m −An∣ ≤ ∑⌊n/2⌋k=01pk ∣Ln+m−k −Ln−k ∣ + 2M ∑k>⌊n/2⌋
1pk
≤ 11−1/pε + 2Mε für n gross genug.
Bemerkung 3.30 Um die Verfahren zu stoppen ∣xn − xn−1∣ ≤ ε, ∣f(xn)∣ ≤ ε
3.3. Konvergenzbeschleunigung
De nition 3.31 Sei xn eine Folge, die gegen x∞ konvergiert. Die Folge xn konvergiertschneller, falls
limn→∞
∣xn − x∞∣∣xn − x∞∣
= 0.
3.3.1. Die∆2-Methode von Aitken
Sei xn eine Folge mit xn+1 − x∞ = k(xn − x∞) k ∈ (0,1)Wir habenxn+1 − x∞ = k(xn − x∞),xn − x∞ = k(xn−1 − x∞),xn+1 − xn = k(xn − xn−1),xn+1 − kxn = (1 − k)x∞,
x∞ = 11−k (xn+1 − kxn) = xn+1 + ( 1
1−k − 1)xn+1 − kxn
1−k = xn+1 + k1−k (xn+1 − xn).
k = xn+1−xn
xn−xn−1und
1 − k = 1 − xn+1−xn
xn−xn−1= (xn−xn−1)−(xn+1−xn)
xn−xn−1.
Dann folgtx∞ = xn+1 + (xn+1−xn)2
xn−xn−1/ (xn−xn−1)−(xn+1−xn)
(xn−xn−1) = xn+1 − (xn+1−xn)2(xn+1−xn)−(xn−xn−1) .
∆xn = xn+1 − xn
∆(∆(xn−1)) =∆2xn−1 =∆(xn − xn−1) = (xn+1 − xn) − (xn − xn−1).
x∞ = xn+1 − (∆xn)2∆2xn−1
(x∞ ist bestimmt durch xn−1, xn, xn+1).
34
3.3 Konvergenzbeschleunigung
Satz 3.32 Sei xn ≠ x∞ eine Folge mit
limn→∞
xn+1 − x∞x∞ − x∞
= k ∈ [−1,1).
Dann ist die Folge
xn = xn+1 −(∆xn)2
∆2xn−1
so, dasslimn→∞
xn − x∞xn − x∞
= 0.
B Sei kn = xn+1−x∞xn−x∞ , kn → k.
xn+1 − x∞ = kn(xn − x∞)
∆2xn−1 = (xn+1 − xn) − (xn − xn−1)= (xn+1 − x∞) + (x∞ − xn) − (xn − x∞) − (x∞ − xn−1)= (xn+1 − x∞) − 2(xn − x∞) + (xn−1 − x∞)= (knkn−1 − 2kn−1 + 1)(xn−1 − x∞).
∆xn = xn+1 − xn
= (xn+1 − x∞) + (x∞ − xn)= kn(xn − x∞) − (xn − x∞)= (kn − 1)(xn − x∞).
xn − x∞ = (xn+1 − x∞) − (kn−1)2(xn−x∞)2(knkn−1−2kn−1+1)(xn−1−x∞)
= (xn − x∞) (kn − (kn−1−1)2(knkn−1−2kn−1+1)
xn−x∞(xn−1−x∞)) .
xn−x∞xn−x∞ = (kn − (kn−1)2kn−1
(knkn−1−2kn−1+1)→ (k −(k−1)2k(k2−2k+1)) = k − k = 0.
Anwendung Algorithmus von Steffensenxn+1 = g(xn)x0 gegeben
Program
Für n = 0,1,.., tuny_n = g(x_n)z_n = g(y_n)x_n+1 = z_n-((z_n-y_n)^2)/(z_n-2y_n+x_n)
x0 festXn+1 = g(Xn),
xn+1 = g(g(xn)) − (g(g(xn))−g(xn))2g(g(xn))−2g(xn)+xn
= G(xn),
xn+1 =Xn+2 − (Xn+2−Xn+1)2Xn+2−2Xn+1+Xn
=Xn+2 − (∆Xn+1)2∆2Xn
.
35
Iterative Verfahren
3.4. Algebraische Gleichungen
Eine algebraische Gleichung ist eine Gleichung P (x) = 0, mitP (x) = a0 + a1x + a2x2 + ⋅ ⋅ ⋅ + anxn.
Satz 3.34 Sei z ∈ C eine Nullstelle von P . Dann gilt
∣z∣ ≤ 1 + maxk=0...n−1
∣ak ∣∣an∣
Beispiel 3.35 P (x) = 1 + 2x + 7x4
∣z∣ ≤ 1 + 27.
Hilfssatz 3.36 (Gerschgoring) SeiA = (Aij) eine n × nMatrix. SeiDi = {x ∈ C∣ ∣x −Aii∣ ≤ ∑j≠i ∣Aij ∣}.Alle Eigenwerte vonA liegen in ∪ni=1Di.
B Sei λ ein Eigenwert von A.Sei v mit v ≠ 0, Av = λv, v = (v1, . . . , vn)⊺.Sei i s.d. ∣vi∣ =maxk=1...n ∣vk ∣.∑n
j=1Aijvj = λviÔ⇒∑j≠iAijvj +Aiivi = λvi, und∣(λ −Aii)∣∣vi∣ = ∣∑j≠iAijvj ∣
≤ ∑j≠i ∣Aij ∣∣vj ∣.Es folgt:∣λ −Aii∣ ≤ ∑j≠i ∣Aij ∣ ∣vj ∣
∣vi∣≤ ∑j≠i ∣Aij ∣.
λ Eigenwert von A ⇐⇒ λ ist ein Eigenwert von A⊺,denn:λ Eigenwert ⇐⇒ det(A − λI) = 0
⇐⇒ det((A − λI)⊺) = 0⇐⇒ det(A⊺ − λI) = 0.
Seien A =⎛⎜⎝
. . .Ai1 Ai2 . . . Aii . . . Ain
. . .
⎞⎟⎠,
Di = {x ∈ C∣ ∣x −Aii∣ ≤ ∑j≠i ∣Aij ∣},
D′i = {x ∈ C∣ ∣x −Aii∣ ≤ ∑j≠i ∣Aji∣}.
Die Eigenwerte von A sind in ∪iDi, ∪iD′i enthalten.
36
3.4 Algebraische Gleichungen
Beispiel 3.38
A =⎛⎜⎜⎜⎝
1 2 3 45 6 7 88 7 6 54 3 2 1
⎞⎟⎟⎟⎠
D1 = {x∣ ∣x − 1∣ ≤ 9} = D4,D2 = {x∣ ∣x − 6∣ ≤ 20} = D3.
Behauptung: D1 ⊂D2
Beweis: Sei x ∈D1: ∣x − 6∣ = ∣x − 1 − 5∣ ≤ ∣x − 1∣ + 5 ≤ 9 + 5 < 20.
D′1 = {x∣ ∣x − 1∣ ≤ 17},D′2 = {x∣ ∣x − 6∣ ≤ 12}.
Behauptung: D′2 ⊂D′1.Beweis: ∣x − 1∣ = ∣x − 6 + 5∣ ≤ ∣x − 6∣ + 5 ≤ 12 + 5 = 17.
Ô⇒ λ ∈D′1 ∩D2.B (B S .) Es gilt
∆n = det
⎛⎜⎜⎜⎜⎜⎝
−an−1an− x −an−2
an. . . . . . − a0
an
1 −x 0 . . . 00 1 0⋮0 . . . 1 −x
⎞⎟⎟⎟⎟⎟⎠
= (−1)n
anP (x)
Nach Entwicklung bezüglich der letzten Spalte gilt∆n = (−1)n+1 −a0
an⋅ 1 + (−1)2n(−x)∆n−1
= (−1)n a0
an− x∆n−1
= (−1)n a0
a1+ (−1)n a1
anx + x2∆n−2, da
∆n−1 = (−1)n−1 a1
an− x∆n−2 und die Formel folgt.
z Nullstelle von P ⇐⇒ z ist Eigenwerte von A =⎛⎜⎜⎜⎝
−an−1an
−an−2an
. . . − a0
an
1 0⋮0 0 1
⎞⎟⎟⎟⎠.
Aus dem Satz von Gerschgorin folgtz ∈ ∪iDi = ∪n−2k=1{z∣ ∣z∣ ≤ 1 +
∣ak ∣∣an∣} ∪ {z∣ ∣z +
an−1an∣ ≤ 1} ∪ {z∣ ∣z∣ ≤ ∣ a0
an∣}.
z ist s.d.∣z∣ ≤ 1 + ∣ak ∣
∣an∣ für ein k ∈ {1, . . . , n − 2}∣z∣ ≤ ∣a0∣
∣an∣
∣z + an−1an∣ ≤ 1 Ô⇒ ∣z∣ = ∣z + an−1
an− an−1
an∣ ≤ ∣z + an−1
an∣ + ∣an−1
an∣ ≤ 1 + ∣an−1∣
∣an∣
Ô⇒∣z∣ ≤maxk=0,...,n−1 1 + ∣ak ∣∣an∣ .
Bemerkung 3.40 Es gilt auch ∣z∣ ≤max{(1 +maxk=1...,n−1∣ak ∣∣an∣),
∣a0∣∣an∣}
Newton Methodex0 in einem Kreis (indem die Lösung liegt), xn+1 = xn − P (xn)
P ′(xn) .Wir hoffen, dass xn gegen eine Lösung von P (x) = 0 konvergiert.Programm
37
Iterative Verfahren
für p=0,1,..* Auswertung von P(x_p), P'(x_p) *b_n-1 = a_nc_n-1 = a_nfür i=n-2 ... 1b_i = a_i + x_p * b_i+1c_i = b_i + x_p * c_i+1
x_p+1 = x_p - (a_0 + b_0 * x_p) / c_1
Bemerkung 3.41 Probleme sind zu erwarten im Fall P (x) = (x − x∞)Q(x), Q(x∞) =0.Bairstow MethodeSei P ∈ R[x] (Polynom mit reellen Koeffizienten).
Sei λ ∈ C eine Nullstelle von P , dann ist λ auch eine Nullstelle.Denn: λ Nullstelle: a0 + a1λ + a2λ2 + ⋅ ⋅ ⋅ + anλn = 0Ô⇒ a0 + a1λ + ⋅ ⋅ ⋅ + anλn = 0Ô⇒ a0 + a1λ + ⋅ ⋅ ⋅ + anλn
Ô⇒ a0 + a1λ + ⋅ ⋅ ⋅ + anλn= 0.
Dann giltP (x) = (x − λ)(x − λ)Q(x)
= (x2 − (λ + λ) + λλ)Q(x)= (x2 − 2Re(λ) + ∣λ∣2)Q(x)= (x2 − rx − s)Q(x).
Sei P ein PolynomP (x) = (x2 − rx − s)Q(x) +A(r, s)(x − r) +B(r, s).Wir suchen r, s mit A(r, s) = 0, B(r, s) = 0oder falls F die Abbildung
(rs) F→ (A(r, s)
B(r, s))
ist, F (r, s) = 0.
Die Lösung dieser Gleichung mit der Newton Methode lautet
(rn+1sn+1) = (rn
sn) − (∇F (rn, sn))−1F (rn, sn)
r0, s0 gegeben.
∇F = (Ar As
Br Bs)mit zum Beispiel
Ar(r, s) = Ar = ∂A∂r(r, s),
As(r, s) = As = ∂A∂s(r, s).
Es folgt
(rn+1sn+1) = (rn
sn) − (Ar(rn, sn) As(rn, sn)
Br(rn, sn) Bs(rn, sn))−1
(A(rn, sn)B(rn, sn)
).
38
3.4 Algebraische Gleichungen
Wir versuchen A, B, Ar , As, Br , Bs zu berechnen.Sei Q = b2 + b3x + ⋅ ⋅ ⋅ + bnxn−2, A = b1, B = b0, es gilt(a0 + a1x + ⋅ ⋅ ⋅ + anxn) = (x2 − rx − s)(b2 + b3x + ⋅ ⋅ ⋅ + bnxn−2) + b1(x − r) + b0.
Wir identi zieren die verschiedenen Potenzen von x:xn ∶ an = bnxn−1 ∶ an−1 = bn−1 − rbnxn−2 ∶ an−2 = bn−2 − rbn−1 − sbnxk ∶ ak = bk − rbk+1 − sbk+2x2 ∶ a2 = b2 − rb3 − sb4x ∶ a1 = b1 − rb2 − sb3
∶ a0 = b0 − rb1 − sb2.
Es folgtbn = anbn−1 = an−1 + rbnbk = ak + rbk+1 + sbk+2b1 = a1 + rb2 + sb3b0 = a0 + rb1 + sb2
A(r, s) = b1 = ein Polynom in r und s. B(r, s) = b0 = ein Polynom in r und s.
Wir setzen ck = (bk)r = ∂bk∂r
, dk = (bk)s = ∂bk∂s
Also: Ar = c1, As = d1, Br = c0, Bs = d0
Nach Ableitung der obigen Gleichung kommt:(bn−1)r = (rbn)r = bn + r(bn)r = bn(bk)r = bk+1 + r(bk+1)r + s(bk+2)r.
cn = 0cn−1 = bn = ancn−2 = bn−1 + rcn−1ck = bk+1 + rck+1 + sck+2c1 = b2 + rc2 + sc3c0 = b1 + rc1 + sc2.
dn = 0dn−1 = 0dn−2 = bn = andk = bk+2 + rdk+1 + sdk+2.
Es folgt, dass dk = ck+1 gilt.
Also: Ar = c1, As = d1 = c2, Br = c0, Bs = d0 = c1
(Ar As
Br Bs) = (c1 c2
c0 c1)
(c1 c2c0 c1
)−1
= 1c21−c0c2
( c1 −c2−c0 c1
)
39
Iterative Verfahren
(rp+1sp+1) = (rp
sp) − 1
c21−c0c2( c1 −c2−c0 c1
)(b1b0)
Bairstow-Algorithmus r0, s0 gegeben. ai gegeben.
Für p = 0,1, ...b[n] = a[n]b[n-1] = a[n-1]+r[p]*b[n]c[n] = 0c[n-1] = b[n]Für k = n-2,..0 tunb[k] = a[k] + r[p]*b[k+1] + s[p]*b[k+2]c[k] = b[k+1] + r[p]*c[k+1] + s[p]*c[k+2]
r[p+1] = r[p] - (c[1]*b[1] - c[2]*b[0])/(c[1]*c[1] - c[0]*c[2])s[p+1] = s[p] - (c[1]*b[0] - c[0]*b[1])/(c[1]*c[1] - c[0]*c[2])
bis stop
Am Ende kriegen wir r, s. Mit die Nullstelle von x2 − rx − s = 0 sind zwei Nullstellenvon P .
40
4 Approximation und Interpolation
De nition 4.1 Seien (xi, yi), i = 0, . . . , n Punkte aus R2.Eine Funktion f
interpoliert diese Daten, genau dann wennf(xi) = yi ∀i = 0, . . . , n,
approximiert diese Daten, genau dann wenn∣f(xi) − yi∣ ≤ ε ∀i = 0, . . . , n.
(xi) sind die Knoten oder Gitterpunkte.h =maxi=0,...,n−1 ∣xi+1 − xi∣ heisst die Maschenweite dieses Gitters.
4.1. Polynom-Interpolation
4.1.1. Lagrange & Newton Interpolation
De nition 4.2 Πn = {a0 + a1x + ⋅ ⋅ ⋅ + anxn} = die Menge der Polynome von Grad≤ n.
Satz 4.3 Sei (xi, yi) ∈ R2 , xi ≠ xj ∀i ≠ j.∃!Pn ∈ Πn mit Pn(xi) = yi ∀i = 0, . . . , n.Pn heisst das Lagrangesche Interpolations Polynom, (L.I.P.).
B
1. Existenzli(x) =∏n
j=0j≠i
x−xj
xi−xj∈ Πn
li(xj) = δij wobei δij das Kronecker-Symbol ist.
li(x), i = 0, . . . , n heisst die Lagrangesche Basis.
Wir setzen Pn(x) = ∑ni=0 yili(x) ∈ Πn.
Pn(xj) = ∑ni=0 yili(xj) = yj lj(xj) = yj ∀j = 0, . . . , n.
Bemerkung: Die n + 1 Polynome li(x) bilden eine Basis von Πn ⇐⇒ dieli linear unabhängig sind.∑n
i=0 αili(x) = 0 Ô⇒ αi = 0
41
Approximation und Interpolation
2. EindeutigkeitSeien P , Q zwei interpolierende Polynome. P (xi) = Q(xi) = yi ∀i impliziert(P −Q)(xi) = 0 ∀i = 0, . . . , n. Das heisstP −Q ∈ Πn und P −Q besitzt n + 1 Nullstellen Ô⇒ P = Q
Beispiel 4.5 n = 1, (x0, y0), (x1, y1).l0(x) = x−x1
x0−x1, l1(x) = x−x0
x1−x0.
P1(x) = y0l0(x) + y1l1(x) = y0 x−x1
x0−x1+ y1 x−x0
x1−x0
= y0(x−x0+x0−x1)
x0−x1+ y1 x−x0
x1−x0
= y0 + (x − x0) ( y1−y0
x1−x0) Newton Form des L.I.P.
Im Allgemeinen giltPn(x) = a0 + a1(x − x0) + ⋅ ⋅ ⋅ + an(x − x0) ⋅ . . . ⋅ (x − xn−1)für ai ∈ R. Diese Form heisst die Newtonsche Form des Polynoms.Diese Behauptung folgt aus
Satz 4.6 1, (x − x0), . . . , (x − x0) ⋅ (x − x1) ⋅ . . . ⋅ (x − xn−1) ist eine Basis vonΠn.
B Die obigen Polynome bilden eine Basis ⇐⇒ sie sind linear unabhängig.Sei αi mit α0 +α1(x− x0)+ ⋅ ⋅ ⋅ +αn(x− x0)(x− x1) ⋅ . . . ⋅ (x− xn−1) = 0 (∀x ∈ R).x = x0 Ô⇒ α0 = 0.Ô⇒ α1(x − x0) + α2(x − x0)(x − x1) + ⋅ ⋅ ⋅ + an(x − x0) . . . (x − xn−1) = 0 ∀x.Ô⇒ (x − x0)(α1 + α2(x − x1) + ⋅ ⋅ ⋅ + an(x − x1) . . . (x − xn−1)) = 0 ∀x.Ô⇒ α1 + α2(x − x1) + ⋅ ⋅ ⋅ + an(x − x1) . . . (x − xn−1) = 0 ∀x.(erst ∀x ≠ x0, dann falls wir den Limes nehmen für alle x).Ô⇒ . . .Ô⇒ αi = 0 ∀i.
De nition 4.8 Seien (xi, xi), (xi+1, yi+1), . . . , (xi+k, yi+k) k + 1 Punkte.∃!Pk ∈ Πk mit Pk(xi+j) = yi+j ∀j = 0, . . . , k,
Pk = akxk + . . . ein Polynom von Grad < k undak heissen die k − te dividierte Differenzen der Daten(xi, yi), . . . , (xi+k, yi+k).
Bezeichnung: ak = f[xi, xi+1, . . . , xi+k].
Beispiel 4.9 k = 1, (x0, y0), (x1, y1),P1(x) = y0 + (x − x0) y1−y0
x1−x0,
f[x0, x1] = y1−y0
x1−x0.
Satz 4.10 (Newtonsche Form des L.I.P.) Seien (x0, y0), . . . , (xn, yn) Punkte aus R2.Sei Pn das Lagrange Interpolationspolynom der Daten xi ≠ xj ∀i ≠ j.Dann gilt:
Pn(x) =n
∑j=0
f[x0, . . . , xj](x − x0) ⋅ . . . ⋅ (x − xj−1)
für j = 0 f[x0] = y0.
42
4.1 Polynom-Interpolation
B Zu zeigen ist:Pn(x) = f[x0] + f[x0, x1](x − x0) + f[x0, x1, x2](x0 − x1)(x − x1)
+ ⋅ ⋅ ⋅ + f[x0, . . . , xn](x − x0)(x − x1) . . . (x − xn−1).Pn(x) = f[x0, x1, . . . , xn−1]xn +Q(x)mit gradQ(x) ≤ n − 1,
= f[x0, x1, . . . , xn−1](x − x0) . . . (x − xn−1) +R(x)mit gradR(x) ≤ n − 1.R(x0) = Pn(x0) = y0R(x1) = Pn(x1) = y1
⋮R(xn−1) = Pn(xn−1) = yn−1Ô⇒R(x) = Pn−1(x)das PolynomausΠn−1, welches dieDaten (x0, y0), . . . , (xn−1, yn−1)interpoliert.Nach Induktion folgt:Pn−1(x) = R(x)
= f[x0, . . . , xn−1](x − x0) ⋅ . . . + ⋅ ⋅ ⋅ + f[x0, x1](x − x0) + f[x0].
Satz 4.12 (Rekursionsformel) Es gilt∀i f[xi] = yi, ∀k,∀x0, . . . , xk
f[x0, . . . , xk] =f[x1, . . . , xk] − f[x0, . . . , xk−1]
xk − x0.
f[x0, . . . , xk] ist unabhängig von der Ordnung der Punkte.
B Sei Pk ∈ Πk das Interpolationspolynom der Daten (x0, y0), . . . , (xk, yk)Pk ist unabhängig von der Ordnung der Daten.Pk(x) = f[x0, . . . , xk]xk + . . .und f[x0, . . . , xk] ist unabhängig von Ordnung der Daten.
k = 0.(xi, yi), P0(x) = yi = yi1 = y1x0 der Koeffizient der höchsten Ordnung f[xi] ist yi,d.h. f[xi] = yi ∀i.
k = 1.(x0, y0)„ (x1, y1) P1(x) = y0 + (x − x0) y1−y0
x1−x0
Ô⇒ f[x0, x1] = y1−y0
x1−x0= f[x1]−f[x0]
x1−x0, d.h. das ist die Formel für k = 1.
Induktion über kSeien (x0, y0), . . . , (xk, yk)Sei Pk−1 das Lagrange Interpolationspolynom der Daten (x0, y0), . . . , (xk−1, yk−1).Sei Qk−1 das Lagrange Interpolationspolynom der Daten (x1, y1), . . . , (xk, yk).
Sei R(x) = x−xk
x0−xkPk−1 + x−x0
xk−x0Qk−1 ∈ Πk.
R(x0) = Pk−1(x0) = y0.R(xk) = Qk−1(xk) = yk .Für j ≠ 0, kR(xj) = xj−xk
x0−xkPk−1(xj) + xj−x0
xk−x0Qk−1(xj) = yj(xj−xk
x0−xk+ xj−x0
xk−x0) = yj xk−x0
xk−x0= yj .
Also: R(x) = Pk(x) das Lagrange Polynom der Daten, d.h.
43
Approximation und Interpolation
Pk(x) = x−xk
x0−xkPk−1 + x−x0
xk−x0Qk−1(x).
In Pk(x), f[x0, . . . , xk] ist der Koeffizient der höchsten Ordnung. So,Pk−1(x) = f[x0, . . . , xk−1]xk−1 + Polynom von Grad < k − 1,
Qk−1(x) = f[x1, . . . , xk]xk−1 + Polynom von Grad < k − 1.
Der Koeffizient höchster Ordnung in R(x) ist dann1
x0−xkf[x0, . . . , xk−1] + 1
xk−x0f[x1, . . . , xk] = f[x1,...,xk]−f[x0,...,xk−1]
xk−x0.
Dann haben wir: f[x0, . . . , xk] = f[x1,...,xk]−f[x0,...,xk−1]xk−x0
.
Bemerkung 4.14 Das Schemader dividiertenDifferenzen ist
f[x0] = y0f[x0, x1] = f[x1]−f[x0]
x1−x0= y1−y0
x1−x0
f[x1] = y1 f[x0, x1, x2] = f[x1,x2]−f[x0,x1]x2−x0
f[x1, x2] = y2−y1
x2−x1
f[x2] = y2⋮ ⋮ ⋮
Beispiel 4.15 Gesucht: Lagrange InterpolationspolynomderDaten (0,1), (1,2), (2,1), (3,0)1
12 −1−1 1/3
1 0−1
0
P3(x) = f[x0]+f[x0, x1](x − x0)+f[x0, x1, x2](x − x0)(x − x1)+f[x0, x1, x2, x3](x − x0)(x − x1)(x − x2)
= 1 + x − x(x − 1) + 13x(x − 1)(x − 2).
Seien (xi, yi), i = 0, . . . , n die Interpolationsdaten. Die progressive Form des L.I.P. istPn(x) = f[x0] + (x − x0)f[x0, x1] ⋅ ⋅ ⋅ + (x − x0) . . . (x − xn−1)f[x0, . . . , xn].
Falls wir xn, . . . , x0 betrachten, giltPn(x) = f[xn] + (x − xn)f[xn, xn−1] + ⋅ ⋅ ⋅ + (x − xn) . . . (x − x1)f[xn, . . . , x0]
Pn(x) = ∑nj=0 f[x0, . . . , xj](x − x0) . . . (x − xj−1)
= ∑0j=n f[xn, . . . , xj](x − xn) . . . (x − xj+1)
Algorithmus (Berechnung von Pn(x)).yi ist in di gespeichert.*wir speichern f[xn, . . . , xi]*
for k=1 to n dofor i=0 to n-k dod[i] = (d[i+1]-d[i])/(x[i+1]-x[i]);
end;end;
*Horner Algorithmus*
44
4.1 Polynom-Interpolation
P = d[0]for i=1 to n doP = d[i]+(x-x[i])*PP = P[n](x)
end;
Wir verfolgen den Algorithmus:d0 = y0 . . . , dn = ynk = 1 für i = 0 bis n − 1 tund0 = d1−d0
x1−x0, d1 = d2−d1
x2−x1, …, dn−1 = dn−dn−1
xn−xn−1
k = 2 für i = 0, . . . , n − 2 tund0 = d1−d0
x2−x0, …, dn−2 = dn−1−dn−2
xn−xn−2s
d0 = y0f[x0, x1]
y1
y2 f[x1, x2]⋮
dn−2 = yn−2f[xn−1, xn−2]
dn−1 = yn−1 f[xn, xn−1, xn−2]f[xn, xn−1]
dn = yn
Am Ende sind die f[xn, . . . , xk] in dk gespeichert.
Satz 4.16 (Fehlerabschätzung) Sei f ∈ Cn+1([a, b]). Seien x0, . . . , xn ∈]a, b[
∣f(x) − Pn(x)∣ ≤Mn+1
(n + 1)!∣ωn(x)∣ ∀x ∈ [a, b]
Wobei ωn(x) =∏ni=0(x − xi),Mn+1 =maxx∈[a,b] ∣f (n+1)(x)∣.
B Pn(x) ist das L.I.P. der Daten (xi, f(xi))i=0,...,n. Sei Q das L.I.P. der Daten(x0, . . . , xn, x).
Q(x) = Pn(x) + f[x0, . . . , xn, x](x − x0) . . . (x − xn),
Q(x) − Pn(x) = f[x0, . . . , xn, x]ωn(x),
f(x) − Pn(x) = f[x0, . . . , xn, x]ωn(x).f − Pn ist gleich 0 an n + 1 Punkten x0, . . . , xn.
f(xi) = f(xi+1) Ô⇒ ∃ξi ∈]xi, xi+1[mit f ′(ξi) = 0, so
45
Approximation und Interpolation
f ′ − P ′n ist gleich 0 an n Punkten aus ]a, b[,f ′′ − P ′′n ist gleich 0 an n − 1 Punkten aus ]a, b[,Ô⇒ fn − P (n)n ist gleich 0 an 1 Punkt aus ]a, b[.
d.h. ∃ξ ∈]a, b[mit P (n)n (ξ) = f (n)(ξ). AberPn(x) = f[x0, . . . , xn]xn + . . . und es folgt, dassP(n)n = n!f[x0 − xn].
∃ξ ∈]a, b[mit f[x0, . . . , xn] = f(n)(ξ)n!
∣f(x) − Pn(x)∣ = ∣f(n+1)(ξ)∣(n+1)! ∣ωn(x)∣ ≤maxξ∈[a,b]
∣f(n+1)(ξ)∣(n+1)! ∣ωn(x)∣.
Behauptung: Das L.I.P. ist im Allgemeinen keine gute Approximation einer Funktion.
f(x)p(x)
Die blaue Funktion kann beliebig weit von p sein.
Der Fall xi = x0 + ih (xi+1 = xi + h)xs = x0 + sh, s ∈ Rfs = f(xs).Wir de nieren die niten Differenzen mit ∆ifs = fs für i = 0∆ifs =∆i−1fs+1 −∆i−1fs, i ≠ 0
Zum Beispiel ∆fs = fs+1 − fs = f(xs+1) − f(xs) (eine nite Differenz).
∆2fs =∆fs+1 −∆fs = (fs+2 − fs+1) − (fs+1 − fs) = fs+2 − 2fs+1 + fs.Wir habenxs+2 = x0 + (s + 2)h = x0 + sh + 2h = xs + 2h,xs+1 = xs + h,xs = xs, und∆2fsh2 = f(xs+2h)−2f(xs+h)+f(xs)
h2 ≃ f ′′(xs).
46
4.1 Polynom-Interpolation
Es folgt aus:f(xs + 2h) = f(xs) + 2hf ′(xs) + (2h)2 f ′′(xs+2θh)
2!
f(xs + h) = f(xs) + hf ′(xs) + h2 f ′′(xs+θ′h)2!
∆2fsh2 = (f(xs) + 2hf ′(xs) + (2h)2 f ′′(xs+2θh)
2!− 2f(xs) − 2hf ′(xs)
−2h2 f ′′(xs+θ′h)2!
+ f(xs))/h2
= 2h2f ′′(xs+2θh)−h2f ′′(xs+θ′h)h2 .
So, limh→0∆2fsh2 = f ′′(xs).
Satz 4.18 Es gilt
f[xk, . . . , xk+i] =∆ifki!hi
.
B i = 0f[xk] = fk = ∆0fk
0!h0 .Induktionsannahme: f[xk, . . . , xk+i−1] = ∆i−1fk
(i−1)!hi−1 ∀k.f[xk, . . . , xk+i] = f[xk+1,...,xk+i]−f[xk,...,xk+i−1]
xk+i−xkmit
f[xk+1, . . . , xk+i] = f[xk+1, . . . , xk+1+(i−1)] = ∆i−1fk+1(i−1)!hi−1 und
f[xk, . . . , xk+i−1] = ∆i−1fk(i−1)!hi−1 .
Ô⇒ f[xk, . . . , xk+i] = ∆i−1fk+1−∆i−1fk(i−1)!hi−1(xk+i−xk) =
∆ifk(i−1)!hi−1ih
= ∆ifki!hi .
Satz 4.20 (Progressive Newton Formel) Es gilt
Pn(x) = Pn(x0 + s(x)h) =n
∑i=0(s(x)
i)∆if0
(x = x0 + s(x)h, d.h. s(x) = x−x0
h)
Für ganze Zahlen (si) = s!
i!(s−i)! =s(s−1)...(s−i+1)(s−i)...1)
i!(s−i)(s−i−1)...1 = s(s−1)...(s−i+1)i!
,und wir benützen die selbe Formel für alle s(x), d.h.(s(x)
i) = s(x)(s(x)−1)...(s(x)−i+1)
i!B Pn(x) = ∑n
i=0 f[x0, . . . , xi](x − x0) . . . (x − xi−1)= ∑n
i=0∆if0i!hi (x − x0) . . . (x − xi−1).
f[x0, . . . , xi] = ∆if0i!hi und x = x0 + s(x)h = x0 + sh,
x − x0 = x0 + sh − x0 = sh,x − xk = x0 + sh − (x0 + kh) = (s − k)h.
Es folgtPn(x) = ∑n
i=0∆if0i!hi sh(s − 1)h . . . (s − i + 1)h
= ∑ni=0
∆if0i!
s(s − 1) . . . (s − i + 1)= ∑n
i=0 (si)∆if0.
∣Pn(x) − f(x)∣ ≤max[a,b]∣f(n+1)(ξ)∣(n+1)! ∏
ni=0(x − xi)
Für welche Knoten ist ∣ω(x)∣ am kleinsten auf [a, b]? Wir möchten zeigen, dass es für dieTschebyscheff-Knoten ist.
47
Approximation und Interpolation
De nition 4.22 (Tschebyscheff-Knoten)
xj =1
2(a + b + (b − a) cos (2j + 1)
2(n + 1)π)
Hilfssatz 4.23 Wir können a = −1, b = 1 wählen.B Wir betrachten die folgende Abbildungφ ∶ [a, b] → [c, d], x ↦ c + (x − a) ( d−c
b−a). Es ist eine bijektive Abbildung von [a, b] auf[c, d].yi = φ(xi) = c + (xi − a) ( d−cb−a),y = φ(x) = c + (x − a) ( d−c
b−a),y − yi = (x − xi) d−cb−a .ω(y) =∏n
i=0(y − yi) =∏ni=0(x − xi) ( d−cb−a)
n+1
Ô⇒ ω(x) = ( b−ad−c)
n+1ω(y).
So max[a,b] ω(x) = ( b−ad−c)n+1
max[c,d] ω(y).Wir wählen c = −1, d = +1,max[a,b] ω(x) = ( b−a2 )
n+1max [−1,1]ω(y),
und das Problem ist das Maximum der ω auf [−1,1] zu nden.d.h. betrachte das Problem auf [−1,1]. Welche yi minimieren das max[−1,1] ω(y)?yi = cos( 2i+1
2(n+1)π).
Wirwerden sehen, dass Für y ∈ [−1,1] setzenwirTn(y) = cos(narccos y) n = 0,1, . . .(Tschebischeff Polynome).
Satz 4.25 Es gilt
1. Tn+1(yi) = 0 ∀i = 0, . . . , n,
2. Tn+1(y) = 2yTn(y) − Tn−1(y), T0(y) = 1, T1(y) = y,
3. Tn ∈ Πn d.h. ist ein Polynom von Grad ≤ n,
4. Tn+1(y) = 2n∏ni=0(y − yi),
5. max[−1,1] ∣Tn+1(y)∣ = 1, Tn+1(cos kπn+1) = (−1)
k , k = 0, . . . , n + 1.
B
1. Tn+1(yi) = cos((n + 1)arccos(yi))= cos((n + 1) 2i+1
2(n+1)π)= cos((2i + 1)π
2)
= 0.
2. Tn+1(y) = cos((n + 1)arccos y)= cos(narccos y + arccos y)= cos(narccos y) cos(arccos y) − sin(narccos y) sin(arccos y)= yTn(y) − sin(narccos y) sin(arccos y).
cos(a − b) − cos(a + b) = 2 sina sin b und es folgt
48
4.1 Polynom-Interpolation
Tn+1(y) = yTn(y) − 12(Tn−1(y) − Tn+1(y)),
12Tn+1(y) = yTn(y) − 1
2Tn−1(y)
Tn+1(y) = 2yTn(y) − Tn−1(y), T0 = 1, T1 = y.
3. Tn ∈ Πn
T0 ∈ Π0, T1 ∈ Π1, . . . ,Πn−1 ∈ Πn−1 (nach Induktionsannahme).Tn = 2yTn−1 − Tn−2 ∈ Πn.Der Koeffizient des Terms der grössten Ordnung in Tn ist 2n−1
4. Es folgt aus 1), dass Tn+1(y) = c∏ni=0(y − yi) = 2n∏
ni=0(y − yi).
5. maxy∈[−1,1] ∣Tn+1(y)∣ =max[−1,1] ∣ cos((n + 1)arccos y)∣und das Maximum ist erreicht für (n + 1)arccos y = kπ⇐⇒ y = cos kπ
n+1Es gilt Tn+1(cos kπ
n+1) = coskπ = (−1)k
Satz 4.27 max[a,b] ∣ω(x)∣ =max[a,b] ∣∏ni=0(x − xi)∣ ist minimal für
xi =1
2(a + b + (b − a) cos( 2j + 1
2(n + 1)π) .
B ωT =∏ni=0(x − xi), wobei xi die Tschebischeff-Knoten sind.
Sei ω =∏ni=0(x − zi), für andere Knoten zi.
Falls der Satz falsch ist, gilt max[a,b] ∣ω∣ <max[a,b] ∣ωT ∣ (für Knoten zi)
Sei P = ω − ωT ∈ Πn.ωT = cn+1Tn+1, wobei cn+1 eine Konstante ist.Die Extremwerte von ωT sind cn+1 ⋅ (Extremwerte von Tn+1) = cn+1(−1)kSeien ej , ej+1 zwei aufeinanderfolgende Extremwerte, j = 0, . . . , n + 1.ωT (ej) +wT (ej+1) = 0,∣ωT (ej)∣ = ∣ωT (ej + 1)∣ =max ∣ωT ∣Ô⇒ P (ej)P (ej+1) = (ω − ωT (ej))(ω − ωT (ej+1)) < 0und P besitzt n + 1 Nullstellen , P ∈ Πn Ô⇒ P = 0 Ô⇒ ω = ωT .
4.1.2. Hermite Interpolation
Beispiel 4.29 Das Problem der Strassenlaterne.Gegeben: α, (x0, h0), (x1, h1).Wir suchen p ∈ Π2 , p(x) = a + bx + cx2 mitp(x1) = h1, p(x0) = h0
p′(x0) = −tangα = h′0.Wir suchen p von der Form p(x) = A +B(x − x0) +C(x − x0)2.A = h0, p′(x) = B + 2C(x − x0) Ô⇒ B = h′0.p(x) = h0 + h′0(x − x0) +C(x − x0)2Ô⇒ p(x1) = h1 = h0 + h′0(x − x0) +C(x1 − x0)2Ô⇒ h1−h0
x1−x0= h′0 +C(x1 − x0)
Ô⇒ C = (h1−h0
x1−x0− h′0)/(x1 − x0).
49
Approximation und Interpolation
VerallgemeinerungGesucht ist ein Polynom P mit
(B)
⎧⎪⎪⎪⎪⎪⎨⎪⎪⎪⎪⎪⎩
P (x0) = f(x0), P ′(x0) = f ′(x0), . . . , P (i0)(x0) = f (i0)(x0),P (x1) = f(x1), P ′(x1) = f ′(x1), . . . , P (i1)(x1) = f (i1)(x1),⋮P (xn) = f(xn), P ′(xn) = f ′(xn), . . . , P (in)(xn) = f (in)(xn).
Zum Beispiel: Sei f eine Funktion, genug o differenzierbar.Gesucht ist P ein Polynom mitP (xj) = f(xj), und P (i)(xj) = f (i)(xj) für j = 0, . . . , n, i = 1, . . . , ij .
De nition 4.30 P heisst das Hermitesche Interpolationspolynom von f .
Beispiel 4.31 Seien a, b gegeben, m = a+b2
(a, f0, f ′0), (m,f1), (b, f2, f ′2).Gesucht ist P mitP (a) = f0, P ′(a) = f ′0, P (m) = f1, P (b) = f2, P ′(b) = f ′2.
(ξ0, ξ1, ξ2, ξ3, ξ4) = (a, a,m, b, b),P (x) = f[ξ0] + f[ξ0, ξ1](x − ξ0) + ⋅ ⋅ ⋅ + f[ξ0, ξ1, . . . , ξ4](x − ξ0) . . . (x − ξ3).
f[x0, x1] = f(x1)−f(x0)x1−x0
→ f ′(x0), für x1 → x0.
De nition 4.32
f[ k + 1-malxj , . . . , xj] =
f (k)(xj)k!
.
Wir bezeichnen mit (ξ0, . . . , ξN) = (i0 + 1-mal
x0, . . . , x0,i1 + 1-mal
x1, . . . , x1, . . . ,in + 1-mal
xn, . . . , xn) das erwei-tertes Gitter, N = i0 + i1 + ⋅ ⋅ ⋅ + in + n.
Satz 4.33 ∃!P ∈ ΠN mit P (i)(xj) = f (i)(xj) ∀i = 0, . . . , ij , ∀j = 0, . . . , n, N =∑n
j=0 ij + nund gilt
P (x) =N
∑j=0
f[ξ0, . . . , ξj]j−1∏k=0(x − ξk)
(die f[ξ0, . . . , ξj] sind durch Induktion de niert).
B
1. EindeutigkeitSeienP,Q s.d. (B) gilt. (P−Q)(x0) = 0, (P−Q)′(x0) = 0,…, (P−Q)(i0)(x0) =0 Ô⇒ x0 Nullstelle mit Vielfachheit i0 + 1Ô⇒ (P −Q) = (x − x0)i0+1Qx1 ist eine Nullstelle von P −Q mit Vielfachheit i1 + 1 …
P −Q = (x − x0)i0+1 ⋅ . . . ⋅ (x − xn)in+1CÔ⇒P −Q ∈ ΠN Ô⇒ C = 0 Ô⇒ P −Q = 0.
50
4.2 Spline Funktionen
2. ExistenzWir suchen N + 1 Unbekannte (die Koeffizienten des Polynomes P ).Die Gleichung (B) ist ein lineares System der Koeffizienten. Das heisst, gesuchtsind die ai mit P (x) = ∑N
i=0 aixi,
∑Ni=0 aix
i0 = f(x0),
⋮∑N
i=0 aiαi = f (in)(x0).
(B) besitzt eine einzige Lösung ⇐⇒ das homogene System (B) besitzt nur dieNull-Lösung. Es folgt aus der Eindeutigkeit.
P (x) = ∑Nj=0 f[ξ0, . . . , ξj]∏
j−1k=0(x − ξk).
folgt nach Induktion (siehe Beweis für die Newtonsche Forme des L.I.P.).
Satz 4.35 Sei f ∈ CN+1([a, b]) es gilt
f(x) − P (x) = fN+1(ζ)(N + 1)!
N
∏k=0(x − ξk)
wobei ζ ∈ [a, b].
4.2. Spline Funktionen
∣f(x) − Pn(x)∣ ≤ Mn+1n+1 ∏
nk=0 ∣x − xk ∣, x ∈ [a, b],
≤ Mn+1n+1 hn+1,
Mn+1 =max[a,b] ∣fn∣, b − a = h.
De nition 4.36 Sei t0 < t1 < ⋅ ⋅ ⋅ < tn = b eine Unterteilung von [a, b]. Eine Spline-funktion ist eine stetige Funktion S ∶ [a, b]→ R, so dassS∣[tj ,tj+1] ∈ Πk−1 ∀j = 0, . . . , n − 1S∣[tj ,tj+1] ist die Einschränkung von S auf [tj , tj+1],k heisst die Ordnung der Splinefunktion, k ≥ 1.k = 1: “konstante Funktionen”k = 2: lineare Splinefunktion.
De nition 4.37 Wir bezeichnen mit Sk,κ die Menge der Splinefunktionen, mit S ∈Cκ([a, b]), S∣[tj ,tj+1] ∈ Πk−1 ∀j = 0, . . . , n − 1 (κ ≤ k − 1)
Satz 4.38 Sk,κ ist ein Vektorraum, und
dimSk,κ = nk − (n − 1)(κ + 1).
51
Approximation und Interpolation
B Wir suchen nk Koeffizienten. Eine Spline Funktion ist eine Folge (P1, . . . , Pn), mit Pi ∈ Πk−1 ∀i undP1(t1) = P2(t1), . . . , P (κ)1 (t1) = P (κ)2 (t1),⋮Pn−1(tn−1) = Pn(tn−1), . . . , P (κ)n−1(tn−1) = P
(κ)n (tn−1).
RRRRRRRRRRRRRRR
(n−1)(κ+1) Bedingungen
Um S ∈ Sk,κ zu nden, haben wir nk Unbekannte und (n − 1)(κ + 1) BedingungenÔ⇒ dimSk,κ = nk − (n − 1)(κ + 1) falls die linearen Bedingungen unabhängig sind.
Sei lij die lineare Form de niert mitlij(P ) = P (j)i+1(ti) − P
(j)i (ti)
für i = 1, . . . , n, j = 0, . . . , κ, P = (P1, . . . , Pn) ∈ (Πk−1)n = Vlij ∈ V ′ = L(V,R) (Dualraum).
Behauptung: lij , i = 1, . . . , n, j = 0, . . . , κ sind linear unabhängig.Beweis:∑i,j γij lij = 0
?Ô⇒ γij = 0.∑i,j γij lij = 0 Ô⇒ ∑i,j γij lij(P ) = 0 ∀P ∈ (Πk−1)n.
Nach der Wahl von P mit lij(P ) = δij (siehe Hermitesche Interpolation) kriegen wirγij = 0.dimSk,κ = nk − (n − 1)(κ + 1)
• κ = k − 1: dimSk,k−1 = nk − (n − 1)k = kIn diesem Fall sind alle Pi gleich P einem Polynom.
• Die grösst mögliche Glattheit einer Splinefunktion ist für κ = k − 2Sk,k−2S4,2 gewöhnlicher Spline.
4.2.1. Lineare Splines
k = 2, κ = k − 2 = 0Sk,k−2 = Sk
Wir interpolieren Daten (ti, fi), fi = f(ti).Sei Sj(t) die Splinefunktion zwischen tj , tj+1.
∆fj = fj+1 − fj , ∆tj = tj+1 − tjSj(t) = fj + ∆fj
∆tj(t − tj) = fj + fj+1−fj
tj+1−tj (t − tj), j = 0, . . . , n − 1.
Satz 4.40 (Feherabschätzung) Sei f ∈ C2([a, b]) . Es gilt
∣S(x) − f(x)∣ ≤ h2
8maxx∈[a,b]
∣f (2)(x)∣.
B Sj(x) ist das Lagrange Interpolationspolynom von Grad < 2 auf (tj , tj+1) undes gilt∣Sj(x) − f(x)∣ ≤maxx∈[tj ,tj+1]
∣f(2)(x)∣2!∣(x − tj)(x − tj+1)∣.
52
4.2 Spline Funktionen
Sei x↦ x2 − (tj + tj+1)x + tjtj+1 = (x − tj)(tj+1 − x) =∶ F (x)Diese Funktion ist maximal für F ′(x) = 0 = −2x + (tj + tj+1), d.h. x = tj+tj+1
2.
F ( tj+tj+12) = ( tj+1−tj
2)2 und
∣Sj(x) − f(x)∣ ≤ max[tj ,tj+1]∣f(2)∣2!
(tj+1−tj)2
4
≤ max[a,b] ∣f (2)∣h2
8∀x ∈ [tj , tj+1].
Was passiert für f nur stetig?
Sj(x) = f(tj) + f(tj+1)−f(tj)tj+1−tj (x − tj), x ∈ [tj , tj+1].
Sj(x) − f(x) = f(tj) − f(x) + f(tj+1−f(tj)tj+1−tj (x − tj)
= (f(tj)−f(x))(tj+1−tj)+(f(tj+1)−f(tj))(x−tj)tj+1−tj
∣Sj(x) − f(x)∣ = ∣(f(tj)−f(x))∣(tj+1−x)+∣(f(tj+1)−f(x))∣(x−tj)tj+1−tj
≤ ω(f,h)(tj+1−x)+ω(f,h)(x−tj)tj+1−tj
= ω(f, h)→ 0 h→ 0.
De nition 4.42 Sei f eine Funktion
ω(f, h) = sup{∣f(x) − f(y)∣ ∣x, y ∈ [a, b], ∣x − y∣ ≤ h}
heisst der Stetigkeit-Modul von f .(ω(f, h)→ 0 für h→ 0 (falls f stetig).)
4.2.2. kubische Splines
k = 4, κ = k − 2 = 2, S4,2 = S4
Seien (ti, f(ti)) gegeben. Im die Spline zu nden haben wir 4n Unbekannte, die 4 Ko-effizienten der n Polynome.
S(tj) = f(tj) ∀j = 0, . . . , n n + 1 Bedingungen.Für C2-Verbindungen haben wir 3(n − 1) Bedingungen.Zusammen: n + 1 + 3(n − 1) = 4n − 2 < 4n Bedingungen. Es fehlen zwei.
De nition 4.43 Die natürlichen Splinefunktionen sind die Funktionen mit:
S ∈ S4
S(tj) = f(tj) ∀j = 0, . . . , n
S′′(t0) = S′′(tn) = 0.
Beispiel 4.44 Finden S ∈ S4 mit S(−1) = 1, S(0) = 2, S(1) = −1, S′′(−1) = S′′(1) =0.
S(x) = S1(x) = a + bx + cx2 + dx3 auf (−1,0),S(x) = S2(x) = a′ + b′x + c′x2 + d′x3 auf (0,1).
53
Approximation und Interpolation
S′′1 (x) = 2c + 6dx,S′′2 (x) = 2c′ + 6d′x,S′′1 (0) = S′′2 (0) Ô⇒ c = c′.
S′1(x) = b + 2cx + 3dx2,S′2(x) = b′ + 2c′x + 3d′x2,S′1(0) = S′2(0) Ô⇒ b = b′.
S1(0) = S2(0) Ô⇒ a = a′.
S′′1 (−1) = S′′2 (1) = 0,2c − 6d = 0,2c′ + 6d′ = 0,Ô⇒ d′ = −d, c = 3d.
Wir suchen a, b, c, dS(0) = S1(0) = 2 Ô⇒ a′ = a = 2,S1(−1) = 1Ô⇒ a − b + c − d = 1, c = 3d, a = 2,S2(1) = −1Ô⇒ a′ + b′ + c′ + d′ = −1, a′ = a, b′ = b, c′ = c, d′ = −d,2 − b + 2d = 1, 2 + b + 2d = −1 Ô⇒ b = b′ = −1s, d = −1, c = −3, a = 2.
also S1(x) = 2 − x − 3x2 − x3 und S2(x) = 2 − x − 3x2 + x3.
Finden der natürlichen kubischen Splines (Existenz und Eindeutigkeit)Die natürliche Splinefunktion ist bestimmt falls wir zi = S′′(ti), i = 0, . . . , n kennen.
S′′(t0) = S′′(tn) = 0. Seizj = S′′(tj), j = 1, . . . , n − 1.
S′′(t) = a + bt auf tj , tj+1. Es folgtS′′(t) = zj+1(t−tj)+zj(tj+1−t)
tj+1−tj = zj+1 (t−tj)tj+1−tj + zj(tj+1−t)tj+1−tj und nach Integration
S(t) = zj+1 (t−tj)3
6(tj+1−tj) + zj(tj+1−t)3
6(tj+1−tj) +Cj(t − tj) +Dj(tj+1 − t).
Wir habenS(tj) = yj (= f(tj)),S(tj+1) = yj+1.
Sei hj ∶= tj+1 − tj
S(tj) = yj ⇐⇒ zj(hj)3
6hj+Djhj = yj ⇐⇒ Dj = yj
hj− zjhj
6,
S(tj+1) = yj+1 ⇐⇒ zj+1h3j
6hj+ cjhj = yj+1 ⇐⇒ Cj = yj+1
hj− zj+1hj
6.
So, S(t) = zj+1 (t−tj)3
6hj+ zj (tj+1−t)
3
6hj+ ( yj+1
hj+1− zj+1hj
6)(t − tj) + ( yj
hj− zjhj
6)(tj+1 − t),
t ∈ (tj , tj+1).
Wir brauchen jetzt die BedingungenS′j(tj) = S′j−1(tj), j = 1, . . . , n − 1, mit Sj die Einschränkung von S auf [tj , tj+1].
54
4.2 Spline Funktionen
S′j(t) = zj+1(t−tj)2
2hj− zj (tj+1−t)
2
2hj+ (yj+1
hj− zj+1hj
6) − ( yj
hj− zjhj
6),
S′j−1(t) = zj+1(t−tj−1)2
2hj−1− zj−1 (tj−t)
2
2hj−1+ ( yj
hj−1− zjhj−1
6) − ( yj−1
hj−1− zj−1hj−1
6),
und es folgt−zj hj
2+ (yj+1
hj− zj+1hj
6) − ( yj
hj− zjhj
6) = zj
hj−12+ ( yj
hj−1− zjhj−1
6)
−( yj−1hj−1− zj−1hj−1
6) ∀j = 1, . . . , n − 1
Sei bj ∶= yj+1−yj
hj. Dann gilt
−zj hj
2− zj+1 hj
6+ zj hj
6+ bj = zj hj−1
2− zj hj−1
6+ zj−1 hj−1
6+ bj−1,
bj − bj−1 = zj+1 hj
6− zj
6(hj + hj−1) + zj
2(hj + hj−1) + zj−1 hj−1
6,
zj+1hj + 2zj(hj + hj−1) + zj−1hj−1 = 6(bj − bj−1) j = 1, . . . , n − 1,z0 = zn = 0.Hilfssatz 4.45 SeiA = (aij) eine n × nMatrix. Falls ∣aii∣ > ∑j≠i ∣aij ∣ istA invertierbar.(Diagonal dominant).
Korollar 4.46 (S) besitzt eine einzige Lösung.B Die Matrix des Systems ist Diagonal dominant.∣ajj ∣ = 2(hj + hj−1)∑k≠j ∣ajk ∣ = hj + hj−1Ô⇒∣ajj∣ > ∑k≠j ∣ajk ∣.B (B H .) z.z. Ax = 0 Ô⇒ x = 0.Ax = 0 ⇐⇒ ∑n
j=1 aijxi = 0 ∀i = 1, . . . , n. Fallsx ≠ 0, sei m mit ∣xm∣ =maxj ∣xj ∣ ≠ 0.
∑nj=1 amjxj = 0 ⇐⇒ ∣amm∣∣xm∣ = ∣ −∑j≠m amjxj ∣ ≤ ∑j≠m ∣amj ∣∣xj ∣Ô⇒ ∣amm∣ ≤ ∑j≠m ∣amj ∣ ∣xj ∣
∣xm∣ ≤ ∑j≠m ∣amj ∣ – unmöglich.
Bemerkung 4.49 ∣aii∣ ≥ ∑j≠i ∣aij ∣ /Ô⇒ A invertierbar.
Beispiel: A = (1 11 1
).
Satz 4.50 (Eigenschaen der Splines)Der natrürliche kubische Spline ist die einzigeC2-Funktion mit S(ti) = fi, die das Integral
∫b
av′′(x)2dx
minimiert.D.h. S ∈M = {v ∈ C2([a, b])mit v(ti) = fi ∀i = 0, . . . , n, v′′(a) = v′′(b) = 0} und
∫b
aS′′
2dx ≤ ∫
b
av′′
2dx ∀v ∈M.
B Sei v ∈M , u = v − S ⇐⇒ v = u + S∫
ba v′′
2dx = ∫
ba (u + S)
′′2dx = ∫ba u′′
2dx + ∫
ba S′′
2dx + 2 ∫
ba u′′S′′dx
v(ti) = S(ti) = fi.
55
Approximation und Interpolation
Wir haben∫
ba S′′u′′dx = ∑n−1
i=0 ∫ti+1ti
S′′u′′dx
= ∑n−1i=0 ∫
ti+1ti(S′′u′)′ − S′′′u′dx
= ∑n−1i=0 ∫
ti+1ti(S′′u′)′dx −∑n−1
i=0 ∫ti+1ti
S′′′u′dx
∑n−1i=0 ∫
ti+1ti(S′′u′)′dx = ∑n−1
i=0 S′′u′(ti+1) − S′′u′(ti)= S′′u′(ti) − S′′u′(t0) + S′′u′(t2) − S′′u′(t1) + S′′u′(t3) − S′′u′(t2) + . . .= S′′u′(tn) − S′′u′(t0)= S′′u′(b) − S′′u′(a)= 0, (S, u sind C2).
∫ba S′′u′′ = ∑n
i=0 − ∫ti+1ti
S′′′u′dx
= −∑ni=0 ∫
ti+1ti(S′′′u)′ − S(4)udx (auf (ti, ti+1), S ∈ Π3Ô⇒S(4) = 0)
= ∑n−1i=1 S′′′u(ti+1) − S′′′u(ti) = 0 (u(ti) = 0 ∀i).
Es gilt∫
ba v′′
2dx = ∫
ba u′′
2 + ∫ba S′′
2dx > ∫
ba S′′
2dx für u′′ ≠ 0.
(u′′ = 0 Ô⇒ u = Ai + tBi auf (ti, ti++)u = v − S, u(ti) = v(ti) − S(ti) = fi − fi = 0, u(ti+1) = 0Ai + tiBi = 0, Ai + ti+1Bi = 0 Ô⇒ Bi = 0 = Ai.)
∫ba v′′
2dx = ∫
ba S′′
2dx ⇐⇒ u′′ = 0 auf (ti, ti+1) ⇐⇒ u = 0 ⇐⇒ v = S
Für v ≠ S gilt:∫
ba v′′
2dx > ∫
ba S′′
2dx
S ist die einzige Funktion, die ∫ba v′′
2dx auf M minimiert.
4.3. Ausgleichung im quadratischen Mittel
Beispiel 4.52 Daten (xi, yi), i = 1, . . . ,5. Gesucht C1,C2, oder die Gerade x ↦ C1 +C2x, so dass gilt:
5
∑i=1((C1 +C2xi) − yi)2 ≤
5
∑i=1(λ1 + λ2xi − yi)2 ∀λ1, λ2 ∈ R.
Wir setzen
A =
⎛⎜⎜⎜⎜⎜⎝
1 x1
1 x2
1 x3
1 x4
1 x5
⎞⎟⎟⎟⎟⎟⎠
, Ac = A(c1c2) =⎛⎜⎝
c1 + c2x1
⋮c1 + c2x5
⎞⎟⎠.
∥x∥ = (∑5i=1 x
2i )1/2 (euklidische Norm in R5), c = (c1
c2), y =
⎛⎜⎝
y1⋮y5
⎞⎟⎠.
Es gilt∑5i=1(c1 + c2xi − yi)2 = ∥Ac − y∥2, und das Problem lautet:
nden c ∈ R2 mit ∥Ac − y∥2 ≤ ∥Aλ − y∥2 ∀λ ∈ R2.
Wir suchen die beste lineare Kombination von 1 und x um die Daten zu approximierenim quadratischen Mittel.
56
4.3 Ausgleichung im quadratischen Mittel
Seien Daten (xi, yi), i = 1, . . . , p. Seien Θ1,Θ2, . . . ,Θn n gegebene Funktionen. (z.B.1, x, n = 2).Wir suchen C = (c1, . . . , cn)⊺ ∈ Rn, s.d. die Funktion
n
∑i=1
ciΘi
die “beste” Approximation im Sinn der quadratischen Mittel ist. Sei A die p × n Matrix
A =⎛⎜⎝
Θ1(x1) Θ2(x1) . . . Θn(x1)⋮ ⋮
Θ1(xp) Θ2(xp) . . . Θn(xp)
⎞⎟⎠.
Sei ∥x∥ = (∑pi=1 x
2i )1/2 die euklidische Norm in Rp.
Wir suchen c mitp
∑i=1(c1Θ(xi)+⋅ ⋅ ⋅+cnΘ(xi)−yi)2 ≤
p
∑i=1(λ1Θ1(xi)+⋅ ⋅ ⋅+λnΘn(xi)−yi)2 ∀λ ∈ Rn,
d.h. ∥Ac − y∥2 ≤ ∥Aλ − y∥2 ∀λ ∈ Rp.
Satz 4.53 Sei cmit∥Ac − y∥2 ≤ ∥Aλ − y∥2 ∀λ ∈ Rn,
dann giltA⊺Ac = A⊺y
(die beiden Aussagen sind äquivalent)
Falls die Spalten von A unabhängig sind, ist c eindeutig bestimmt und gegeben durch c =(A⊺A)−1A⊺y.
B Sei Θi = (Θ1(x1), . . . ,Θi(xp))⊺ die i-te Spalte von A.V = [Θ1, . . . ,Θn] der erzuegte Vektorraum.
∥Ac− y∥2 ≤ ∥Aλ− y∥2 ∀λ ∈ Rn ⇐⇒ Ac ist die orthogonale Projektion von y auf V ,d.h. der einzige Vektor mit Ac − y�V
Ac ist der einzige Vektor mit ⟨Θi,Ac − y⟩ = 0 ∀i = 1, . . . , n⇐⇒ A⊺(Ac − y) = 0 ⇐⇒ A⊺Ac = A⊺y.
Falls die Spalten von A linear unabhängig sind, ist (A⊺A) invertierbar.Wir müssen zeigenA⊺Ax = 0 Ô⇒ x = 0A⊺Ax = 0 Ô⇒ ⟨A⊺Ax,x⟩ = 0 Ô⇒ ⟨Ax,Ax⟩ = 0Ô⇒ Ax = 0 ⇐⇒ x1Θ1 + . . . xnΘn = 0Ô⇒ x = 0 (die Spalten sind linear unabhängig).Beispiel 4.55 Finden Sie den Ausgleich im quadratischen Mittel der Daten.x 0 1 2 4y 0 1 3 5
57
Approximation und Interpolation
A =⎛⎜⎜⎜⎝
1 01 11 21 4
⎞⎟⎟⎟⎠
y =⎛⎜⎜⎜⎝
0135
⎞⎟⎟⎟⎠
y = c1 + c2x, c ∶ A⊺Ac = A⊺y
A⊺A = (4 77 21
), A⊺Ac = (4 77 21
)(c1c2) = A⊺y = ( 9
27)
4c1 + 7c2 = 97c1 + 21c2 = 27
c1 = 0, c2 = 97
y = c1 + c2x + c3x2
Θ1 = 1, Θ2 = x, Θ3 = x2
A =⎛⎜⎜⎜⎝
1 0 01 1 11 2 41 4 16
⎞⎟⎟⎟⎠, A⊺A =
⎛⎜⎝
4 7 217 21 7321 73 237
⎞⎟⎠, A⊺y =
⎛⎜⎝
92793
⎞⎟⎠
…
58
5 Numerische Integration
5.1. Einführung
Wir möchten das Integral I(f) = ∫ba f(x)dx berechnen.
Sei F eine Stammfunktion von f , dann gilt ∫ba f(x)dx = F (b) − F (a).
∃f wobei die Stammfunktion unberechbar ist.z.B. f(x) = e−x
2
, f(x) = sin(x2), f(x) = e−√
x2+ln(1+x)…
Sei f(x) ein Polynom, dann ist ∫ba f(x)dx = F (b) − F (a) berechbar.
Bemerkung 5.1 I(f) = ∫ba f(x)dx Fläche unter der Kurve.
5.2. Interpolarische Formeln
Sei a ≤ x0 < ⋅ ⋅ ⋅ ≤ xn ≤ b eine Unterteilung von (a, b).Sei f ∶ [a, b]→ R.Sei li(x) =∏j≠i
x−xj
xi−xj, i = 0, . . . , n die Lagransche Basis. Das L.I.P. ist gegeben mit
P (x) = ∑ni=0 f(xi)li(x).
Dann können wir de nieren:In(f) = I(p) = ∫
ba ∑
ni=0 f(xi)li(x) = ∑n
i=0 f(xi) ∫ba li(x)dx = ∑n
i=0 f(xi)ai. (0)(0) ist eine “Integrationsformel”.Beispiel 5.2
1. n = 0, x0 = a+b2
, l0 = 1,I0(f) = f(x0) ∫
ba l0(x)dx = f(a+b2 )(b − a).
Diese Integrationsformel heisst die Mittelpunktsformel.
2. n = 1, x0 = a, x1 = b, l0(x) = x−x1
x0−x1= x−b
a−b , l1(x) =x−x0
x1−x0= x−a
b−a ,I1(f) = f(x0) ∫
ba l0dx + f(x1) ∫
ba lidx
= f(a) ∫ba
x−ba−bdx + f(b) ∫
ba
x−ab−a dx
= f(a) 1a−b
(x−b)22∣ba + f(b) 1
b−a(x−a)2
2∣ba
= f(a) b−a2+ f(b) b−a
2
= f(a)+f(b)2
(b − a).Diese Integrationsformel heisst die Trapezformel.
3. n = 2, x0 = a, x1 = a+b2
, x2 = b,l0 = (x−x1)(x−x2)
(x0−x1)(x0−x2) = 2(x− a+b
2 )(x−b)(b−a)2 = 2(x−x1)(x−x2)
(b−a)2 .
59
Numerische Integration
∫ba l0(x)dx = ∫
x2
x0
2(x−x1)(x−x2)(b−a)2
= 2(b−a)2 ∫
x1
x0((x − x2) + (x2 − x1))(x − x2)
= 2(b−a)2 (
13(x − x2)3∣x2
x0+ 1
2(x2 − x1)(x − x2)2∣x2
x0)
= 2(b−a)2 (−
13(a − b)3 − 1
2(b − a)2 (b−a)
2)
= (b − a)2( 13− 1
4)
= (b−a)6
,
∫ba l1(x)dx = 4
6(b − a),
∫ba l2(x)dx = b−a
6.
I2(f) = f(a) b−a6 + f(a+b2)4 b−a
6+ f(b) b−a
6= b−a
6(f(a) + 4f(a+b
2) + f(b)).
Diese Integrationsformel heisst die Simpsonsformel/Keplersche Fassregel.
4. n = 3, x0 = a, x1 = a + b−a3
, x2 = a + 2(b−a)3
, x3 = b,I3(f) = b−a
8(f(a) + 3f(x1) + 3f(x2) + f(b)).
Diese Integrationsformel heisst die Drei-Achtel-Regel/Pulcherima.
De nition 5.3 Sei a ≤ x0 < x1 < ⋅ ⋅ ⋅ < xn ≤ b
In(f) =n
∑i=0
f(xi)ai
mit ai = ∫ba li(x)dx, li(x) =∏j≠i
x−xj
xi−xjheisst eine Newton-Cotes-Formel.
Für x0 = a, xn = b ist die Formel “abgeschlossen”.Für x0 > a, xn < b ist die Formel “offen”.
De nition 5.4 (Ordnung einer Integrationsformel) Sei In eine Integrationsformel.Die Ordnung von In ist k⇐⇒ In(p) = I(p) ∀p ∈ Πk−1
(I(p) = ∫ba p(x)dx)
In ist “interpolarisch” falls die Ordnung n + 1 ist.
Satz 5.5 Die Formel (0) ist interpolarisch.(In(f) = ∑n
i=0 f(xi) ∫ba li(x)dx = I(Pn(x))).
B Wir müssen zeigen, dass In(p) = I(p) ∀p ∈ Πn gilt.Sei p ∈ Πn. Das Lagrange Interpolationspolynom von p ist p.Es folgt In(p) = I(p).
Beispiel 5.7 Behauptung: Die Mittelpunktsformel besitzt die Ordnung 2.I0(f) = f(a+b2 )(b − a).
I0(c) = c(b − a) = ∫ba cdx = I(c) für alle Konstanten c.
I0(x − a+b2) = 0(b − a) = 0 und
60
5.2 Interpolarische Formeln
∫ba (x −
a+b2) = (x−
a+b2 )
2
2∣ba = 0.
I0(p) = I(p) ∀p ∈ Π1.
Beispiel 5.8 Die Ordnung der Simpsonsregel ist 4.I2(f) = (b−a)6
(f(a) + 4f(a+b2) + f(b)).
Für k > 0 giltI2((x − a+b
2)k) = (b−a)
6((a−b
2)k + ( b−a
2)k) = b−a
6( b−a
2)k((−1)k + 1k).
k = 0:I2(1) = (b − a) = I(1) = ∫
ba dx = (b − a).
k = 1:I(x − a+b
2) = ∫
ba (x −
a+b2) = (x−
a+b2 )
2
2∣ba = 0 = I2(x − a+b
2).
k = 2:I2((x − a+b
2)2) = (b−a)
3
12,
I((x − a+b2)2) = ∫
ba (x −
a+b2)2 = 1
3(x − a+b
2)3)∣ba =
(b−a)312= I2((x − a+b
2)2).
k = 3:I2((x − a+b
2)3) = 0 = I((x − a+b
2)3).
k = 4:I2((x − a+b
2)4) = b−a
6( (b−a)
4
162 = (b−a)
5
6⋅8
I((x − a+b2)4) = 1
5(x − a+b
2)5∣ba = 2
5( b−a
2)5 = (b−a)
5
5⋅16 ≠ I2((x −a+b2)4).
I2((x − a+b2)k) = I((x − a+b
2)k) k = 0,1,2,3.
Satz 5.9 Sei In(f) = ∑ni=0 f(xi)ai eine interpolarische Integrationsformel.
Sei f ∈ Cn+1([a, b]). Dann gilt
∣In(f) − I(f)∣ ≤maxx∈[a,b] ∣f (n+1)(x)∣
(n + 1)! ∫b
a
n
∏i=0∣x − xi∣dx.
B Sei pn das L.I.P. der Daten (xi, f(xi))i=0,...,n.In(f) = ∑n
i=0 f(xi)ai = ∑ni=0 pn(xi)ai = In(pn) = I(pn).
∣In(f) − I(f)∣ = ∣I(pn) − I(f)∣ = ∣I(pn − f)∣ = ∣ ∫ba (pn − f)(x)dx∣ und
∣f − pn(x)∣ ≤max[a,b]∣f(n+1)(x)∣(n+1)! ∏
ni=0 ∣(x − xi)∣.
Also ∣In(f) − I(f)∣ ≤ ∫ba ∣f − pn∣dx ≤max[a,b]
f(n+1)
(n+1)! ∫ba ∏
ni=0 ∣x − xi∣.
Hilfssatz 5.11 Seienx0, . . . , x2m, 2m+1Knotenmit symmetrischer Lage gegeben, d.h.xj−a = b − x2m−j ∀j = 0, . . . ,m.Dann gilt ∫
ba ∏
2mj=0(x − xj)dx = 0
61
Numerische Integration
B ω(x) =∏nj=0(x − xj), y = a + b − x.
ω(a + b − y) = ∏2mj=0(a + b − y − xj)
= ∏2mj=0(x2m−j − y)
= −∏2mj=0(y − x2m−j)
= −ω(y).Es folgt:∫
ba ω(x)dx = ∫
ab ω(a + b − y)(−)dy = ∫
ba ω(a + b − y)dy = − ∫
ba ω(y)dy,
Ô⇒ 2 ∫ba ω = 0 Ô⇒ ∫
ba ω = 0.
Korollar 5.13 Seien x0, . . . , x2m, 2m + 1 Knoten mit symmetrischer Lage gegeben. Seix2m+1 ein anderer Knoten. Sei pk das L.I.P. der Daten (x0, f(x0)), . . . , (xk, f(xk)).Dann gilt I(p2m+1) = I(p2m).B P2m+1 = P2m + f[x0, . . . , x2m+1]∏2m
i=0(x − xi). Aus dem Hilfssatz folgt:I(p2m+1) = I(p2m) + f[x0, . . . , x2m+1]I(∏2m
i=0(x − xi)) = I(p2m+1).Korollar 5.15 Seien x0, . . . , x2m, 2m + 1 Knoten mit symmetrischer Lage gegeben. Seif ∈ C2m+2([a, b]). Sei In eine interpolarische Integrationsformel.Dann gilt ∣In(f) − I(f)∣ ≤max[a,b]
∣f(2m+2)(x)∣(2m+2)! ∫
ba ∏
2m+1i=0 ∣(x − xj)∣dx
(∀x2m+1 ∈ [a, b]).B Sei p2m+1 das L.I.P. der Knoten x0, . . . , x2m+1, n = 2m + 1In(f) = In(p2m) = I(p2m) = I(p2m+1).Es folgtIn(f) − I(f) = I(p2m+1) − I(f) = I(f − p2m+1),und∣In(f) − I(f)∣ = ∣ ∫
ba (f − p2m+1)∣
≤ ∫ba ∣f − p2m+1∣
≤ max[a,b] ∣f(2m+2)(x)∣(2m+2)! ∫
ba ∏
2m+1j=0 ∣x − xj ∣.
Anwendungen
1. Mittelpunktsregelm = n = 0.∣I0(f) − I(f)∣ ≤ max[a,b]
∣f(2)∣2! ∫
ba ∣x − x0∣∣x − x1∣
≤ max ∣f(2)∣2! ∫
ba (x −
a+b2)2
= max[a,b]∣f ′′∣2
13(x − a+b
2)3∣ba
= max[a,b] ∣f ′′∣ 13(b−a2)3 = 1
24max ∣f ′′∣(b − a)3,
für f ∈ C2([a, b]).
2. Simpsonsregelm = 1.∣I2(f) − I(f)∣ ≤ max[a,b]
∣f(4)∣4! ∫
ba ∏
3i=0(x − xi)
= max f(4)
24 ∫ba (x −
a+b2)2∣(x − a)(x − b)∣.
∫ba (x −
a+b2)2∣(x − a)(b − x)∣ = 1
3(x − a+b
2)3(x − a)(b − x)∣ba
+ 13 ∫
ba (x −
a+b2)3(2x − (a + b))
= 23 ∫
ba (x −
a+b2)4
= 2315(x − a+b
2)5∣ba = 4
15( b−a
2)5 = (b−a)
5
120.
Also: ∣I2(f) − I(f)∣ ≤max 12880
max[a,b] ∣f (4)∣(b − a)5.
62
5.3 Zusammengsetzte Formeln
3. Trapezregel∣I1(f) − I(f)∣ ≤ max[a,b]
∣f ′′∣2! ∫
ba ∣x − a∣∣x − b∣
= max ∣f ′′∣2 ∫
ba (x − a)(b − x)
= max ∣f ′′∣2( (x−a)
2
2(b − x)∣ba + ∫
ba(x−a)2
2)
= max ∣f ′′∣2
(x−a)36∣ba
= 112
max[a,b] ∣f ′′∣(b − a)3.
5.3. Zusammengsetzte Formeln
Sei J = [a, b]. Wir teilen J in L Intervalle, Ji. Dann giltI(f) = ∫
ba f(x)dx = ∑L
i=1 ∫Jif(x)dx
und für ∫Jif(x)dx wenden wir eine Integrationsformel an.
Bemerkung 5.17 Wir werden Jl = [a + (l − 1)h, a + lh]mit h = b−aL
, l = 1, . . . , Lwählen und auf Jl die selben Integrationsformeln betrachten.
1. Zusammengesetzte MittelpunktsformelIM(f) = h∑L
l=1 f(xl), xl = a + (l − 12)h.
∣IM(f) − I(f)∣ ≤ ∑Ll=1 ∣(hf(xl) − ∫
a+lha+(l−1)h f(x)∣
≤ ∑Ll=1maxJl
∣f ′′∣h3
24
= h2
24max[a,b] ∣f ′′∣Lh
= h2
24max[a,b] ∣f ′′∣(b − a).
2. Zusammengesetzte TrapezformelIT (f) = ∑L
l=112(f(a + (l − 1)h) + f(a + lh))h
= h2(f(a) + 2∑L−1
l=1 f(a + lh) + f(b)).
∣IT (f) − I(f)∣ ≤ ∑Ll=1 ∣
hf(a+(l−1)h)+f(a+lh)2
− ∫a+lha+(l−1)h f(x)∣
≤ ∑Ll=1
maxJl∣f ′′∣
12h3
≤ max[a,b]∣f ′′∣12(b − a)h2.
3. Die Zusammengesette SimpsonsregelIS(f) = ∑L
l=116(f(a + (l − 1)h) + 4f(a + (l − 1
2)h) + f(a + lh))
= 16(f(a) + 2∑L−1
l=1 f(a + lh) + 4∑Ll=1 f(a + (l − 1
2)h) + f(b)).
∣IS(f) − I(f)∣ ≤ ∑Ll=1 ∣ISloc − ∫
a+lha+(l−1)h f(x)∣
≤ ∑Ll=1
12880
maxJl∣f (4)∣h5
≤ max[a,b]∣f(4)∣(b−a)
2880h4.
(ISloc ist die lokale Simpsonsformel).
5.4. Konvergenzuntersuchungen
Sei In(f) = ∑nj=0 a
nj f(xj) eine Integrationsformel.
Annahme: ∀p Polynom, In(p)→ I(p) = ∫ba p(x)dx für n→ +∞.
63
Numerische Integration
Die Annahme gilt für die obigen zusammengesetzten Integrationsformeln (IM , IT , IS).
Satz 5.18 Falls ∃C eine Konstante mit ∑nj=0 ∣anj ∣ ≤ C dann konvergiert die Integrations-
formel, d.h. In(f)→ I(f) ∀f ∈ C([a, b]).
B Sei f ∈ C([a, b]), sei ε > 0, ∃p ein Polynom mit ∣f(x) − p(x)∣ ≤ ε ∀x ∈ [a, b](Weierstrass).∣In(f) − I(f)∣ ≤ ∣In(f) − In(p)∣ + ∣In(p) − I(p)∣ + ∣I(p) − I(f)∣∣In(f) − In(p)∣ = ∣∑n
j=0 anj (f(xj) − p(xj))∣
≤ ∑nj=0 ∣anj ∣∣f(xj) − p(xj)∣
≤ ε∑nj=0 ∣anj ∣
≤ εC.
∣I(p) − I(f)∣ = ∣ ∫ba (f(x) − p(x))∣
≤ ∫ba ∣f − p∣ ≤ (b − a)ε.
∣In(f) − I(f)∣ ≤ Cε + (b − a)ε + ∣In(p) − I(p)∣≤ (c + (b − a) + 1)ε für n gross, siehe Annahme.
Korollar 5.20 Fallsanj ≥ 0∀j, ndannkonvergiert die Integrationsformel∀f ∈ C([a, b]).B ∑n
j=1 ∣anj ∣ = ∑nj=1 a
nj = In(1) → I(1) und dann In(1) is beschränkt und die
Integrationsformel konvergiert.Korollar 5.22 Die IT , IM , IS Integrationsformel konvergieren.B anj ≥ 0In(p)→ I(p) (Siehe Fehlerabschätzung für Cα-Funktionen).
5.5. Gauss-Quadratur
Wir versuchen die Ordnung einer Integrationsformel zu maximieren.Wir betrachten die IntegrationsformelnIn(f) = ∑n
j=0 ajf(xj) für aj = ∫ba lj(x)dx. (∗)
Satz 5.24 Die Ordnung von (∗) ist höchstens 2n + 2.
B Sei ω2 =∏ni=0(x − xi)2, ω2 ∈ Π2n+2
In(ω2) = ∑nj=0 ajω
2(xj) = 0 < I(ω2) = ∫ba ω2dx.
Satz 5.26 Die Integrationsformel (∗) besitzt die Ordnung 2n + 2 genau dann wenn∫
ba ωpdx = 0 ∀p ∈ Πn.(ω(x) =∏n
i=0(x − xi)).
B
1. Die Ordnung von (∗) ist 2n + 2 d.h. In(q) = I(q) ∀q ∈ Π2n+1.Es gilt ωp ∈ Π2n+1 ∀p ∈ Πn
Ô⇒ 0 = In(ωp) = I(ωp) = ∫ba ωpdx.
2. ∫ba ωpdx = 0 ∀p ∈ Πn.
Sei q ∈ Π2n+1∃Q,R: q = ωQ +R, Q,R ∈ Πn
In(q) = In(ωQ) + In(R) = In(R).Die Formel (∗) ist interpolarisch und es folgtIn(q) = I(R) = I(ωQ) + I(R) = I(q).
64
5.5 Gauss-Quadratur
Satz 5.28 Sei q ∈ Πn+1 mit ∫ba qpdx = 0 ∀p ∈ Πn, dann
q besitzt n + 1 verschiedene Nullstellen in (a, b).
B q = q0∏Ni=1(x − ti)ri ,
q0 ist ein Polynom, mit q0 ≠ 0 auf (a, b)ti sind die Nullstellen von q in (a, b)mit Vielfachheit ri.
Seien (ξ1, . . . , ξM) die ti mit ri ungerade.Sei p =∏M
i=1(x − ξi).
Behauptung: M = n + 1Falls nicht M ≤ n, p ∈ Πn
0 = ∫ba qp = ∫
ba q0∏(x−tj)rj ∏M
i=1(x−ξi) = ∫ba q0∏(x−tj)r
′j , r′j gerade. Dieses letzte
Integral ist positiv, und wir kriegen einen Widerspruch.
Es folgt, dass M = n + 1 und q = C∏n+1j=1 (x − tj), ti ≠ tj ∀i ≠ j.
Satz 5.30 Sei ω0(x) = 1. Für n = 0,1,2, . . . ∃!ωn+1 ∈ Πn+1 mit∫
ba ωn+1p = 0 ∀p ∈ Πn
und ωn+1(b) = 1.
B Man muss n + 2 Koeffizienten (von ωn+1) nden.Sei {p0, . . . , pn} eine Basis von Πn.Wir suchen n + 2 Unbekannte, die(S) ∫
ba ωn+1pj = 0 ∀j = 0, . . . , n, ωn+1(b) = 1
lösen.
Bezüglich der Koeffizienten ist (S) ein lineares System von n + 2 Gleichungen.
(S) ist eind. lösbar⇐⇒ ∫ba ωn+1pj = 0 ∀j = 0, . . . , nωn+1(b) = 0
}Ô⇒ ωn+1 = 0
ωn+1 ∈ Πn+1
∫ba ωn+1p = 0 ∀p ∈ Πn Ô⇒ ωn+1 = C∏n+1
i=1 (x − xi) siehe oben.ωn+1(b) = 0 = C∏n+1
i=1 (b − xi) Ô⇒ C = 0.
De nition 5.32 Die ωn heissen Legendre-Polynome (eindeutig bestimmt). Die Null-stellen von ωn heissen die Gauss-Knoten. Und für diese Knoten ist die obige Integrati-onsformel von Ordnung 2n + 2.
Bemerkung 5.33 Die Integrationsformel In(f) = ∑nj=0 ajf(xj),
aj = ∫ba lj , x0, . . . , xn Gauss-Knoten, konvergiert.
Denn: In(p) = I(p) für n > gradp, d.h. die obige Annahme gilt, undaj = In(l2j ) = I(l2j ) ≥ 0.Beispiel 5.34
n = 0.ω1(x) = (x − a+b
2) 2b−a , da
ω1(b) = 1, ∫ba ω1c = 0.
I0(f) = a0f(a+b2 ) = (b − a)f(a+b2)Mittelpunktsformel.
65
Numerische Integration
n = 1.ω2(x) = c(x − x0)(x − x1).Wir suchen x0, x1 mit einer symmetrischen Lage,d.h. xM − x0 = x1 − xM , xM = a+b
2.
Das heissta + b − x0 − x1 = 0.Setzex0 = b − λ(b − a) = λa + (1 − λ)b.x1 = a + b − x0 = (1 − λ)a + λb.ω2 = c(x − x0)(x − x1) = c(x − b + λ(b − a))(x − a − λ(b − a)).Wir müssen haben∫
ba ω2 = 0 = ∫
ba xω2
d.h.0 = ∫
ba (x − b + λ(b − a))(x − a − λ(b − a))dx
= ∫ba (x − b)(x − a) + λ(b − a)
2 − λ2(b − a)2dx= − 1
6(b − a)3 + λ(b − a)3 − λ2(b − a)3
= −(b − a)3(λ2 − λ + 16)
da∫
ba (x − a)(x − b) =
12(x − a)2(x − b)∣ba − ∫
ba
12(x − a)2
= − 16(x − a)3∣ba
= − 16(b − a)3.
λ2 − λ + 16= 0,
λ = 12± 1
2√3.
λ = 12+ 1
2√3,
1 − λ = 12− 1
2√3,
x0 = ( 12 +1
2√3)a − ( 1
2− 1
2√3)b,
x1 = ( 12 −1
2√3)a + ( 1
2+ 1
2√3)b.
Für diese symmetrische Lage∫
ba (x −
a+b2)(x − x0)(x − x1) = 0,
d.h. gilt∫
ba x(x − x0)(x − x1) = ∫ a+b
2(x − x0)(x − x1) = 0, und
∫ ω(x)p = 0 ∀p ∈ Π1.
Die Konstante c für ω2 ist so dassc(b − x0)(b − x1) = 1 Ô⇒ c = 1
(b−x0)(b−x1) .
a0 = ∫ba l0 = ∫
ba
x−x1
x0−x1= b−a
2= a1 und die Integrationsformel ist
I1(f) = (b−a)2(f(x0) + f(x1)) ist von Ordnung 4.
Satz 5.35 Sei f ∈ C2n+2([a, b]) dann gilt∣In(f) − I(f)∣ ≤max[a,b]
∣f(2n+2)∣(2n+2)! ∫
ba ω2(x)dx,
ω(x) =∏ni=0(x − xi).
66
5.6 Varia
B Sei p das Hermitsche Interpolationspolynom der Datenf(x0), f ′(x0), . . . , f(xn), f ′(xn). Es giltf − p = f(2n+2)(ξ)
(2n+2)! ω2
p ∈ Π2n+1
∣In(f) − I(f)∣ = ∣In(p) − I(f)∣= ∣I(p) − I(f)∣= ∣I(p − f)∣≤ ∫
ba ∣f − p∣
≤ max[a,b]∣f(2n+2)∣(2n+2)! ∫
ba ω2dx.
5.6. Varia
• Interpolation (Beschleunigung der Konvergenz)Wir nehmen an, dassIh(f) = I(f) + c1hk1 + c2hk2 + . . . , 0 < k1 < k2 < . . . . Dann giltIh
2(f) = I(f) + c1
2k1hk1 + c2
2k2hk2 + . . .
12k1
Ih(f) − Ih2(f) = ( 1
2k1− 1)I(f) + c′hk2
1
2k1Ih−Ih
2(f)
1
2k1−1 = I(f) + c′′hk2 .
Diese letzte Formel verbessert die Konvergenz, da k2 > k1.
• Integration singulärer FunktionenEs gilt ∫
10
dx√x= ∫
10 x−
12 dx = 2
√x∣10 = 2.
Falls f singulär ist, könnenwir die “Singularität” subtrahieren und f als (f −s)+sschreiben. Dann ist f − s regulär und für s berechbar sind wir fertig.
Zum Beispiel: f singulär. ∫10 f = ∫
10 (f − s) + s = ∫
10 (f − s) + ∫
10 s
∫10
cosx√xdx = ∫
10
cosx−1√x
dx + ∫10
1√xdx, cosx−1√
x∼ − ∣x∣
√x
2ist stetig.
67
Index
A
Absoluter Fehler . . . . . . . . . . . . . . . . . . . . . . 8Approximation . . . . . . . . . . . . . . . . . . . . . . 41Auslöschung . . . . . . . . . . . . . . . . . . . . . . . . . 9
B
Basis der Darstellung . . . . . . . . . . . . . . . . . . 4
D
diskretes Fourier Polynom . . . . . . . . . . . . 21dividierte Differenzen . . . . . . . . . . . . . . . . 42Drei-Achtel-Regel . . . . . . . . . . . . . . . . . . . 60
E
erweitertes Gitter . . . . . . . . . . . . . . . . . . . . 50Exponent . . . . . . . . . . . . . . . . . . . . . . . . . . . . 7
F
FehlerAbsoluter . . . . . . . . . . . . . . . . . . . . . . . . 8Relativer . . . . . . . . . . . . . . . . . . . . . . . . 8
Fixpunkt . . . . . . . . . . . . . . . . . . . . . . . . . . . . 25
G
Gauss-Knoten . . . . . . . . . . . . . . . . . . . . . . . 65Gitterpunkte . . . . . . . . . . . . . . . . . . . . . . . . 41Grad . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 17
H
Hermitesche Interpolationspolynom . . . 50Horner Algorithmus . . . . . . . . . . . . . . . . . 17Horner Schema . . . . . . . . . . . . . . . . . . . . . . 17
I
interpolarisch . . . . . . . . . . . . . . . . . . . . . . . 60Interpolation . . . . . . . . . . . . . . . . . . . . . . . . 41
K
Keplersche Fassregel . . . . . . . . . . . . . . . . . 60Knoten . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 41Konditionszahl . . . . . . . . . . . . . . . . . . . . . . 10Kontraktion . . . . . . . . . . . . . . . . . . . . . . . . . 26
L
Lagrangesche Basis . . . . . . . . . . . . . . . . . . 41Lagrangesche Interpolations Polynom. .41Legendre-Polynome. . . . . . . . . . . . . . . . . . 65Leibniz-Formel . . . . . . . . . . . . . . . . . . . . . . 18
M
Mantisse . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 7Maschenweite . . . . . . . . . . . . . . . . . . . . . . . 41Maschinengenauigkeit . . . . . . . . . . . . . . . . . 8Matrix Norm . . . . . . . . . . . . . . . . . . . . . . . . 12Mittelpunktsformel . . . . . . . . . . . . . . . . . . 59
N
natürlichen Splinefunktionen . . . . . . . . . 53Newton-Cotes-Formel . . . . . . . . . . . . . . . 60Newtonsche Form . . . . . . . . . . . . . . . . . . . 42Norm . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 11Nullstelle . . . . . . . . . . . . . . . . . . . . . . . . . . . 18
Vielfachheit . . . . . . . . . . . . . . . . . . . . . 18
O
Ordnung . . . . . . . . . . . . . . . . . . . . . . . . 51, 60Ordnung der Konvergenz . . . . . . . . . . . . . 32
P
Periode . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 21Polynom . . . . . . . . . . . . . . . . . . . . . . . . . . . . 17
Grad . . . . . . . . . . . . . . . . . . . . . . . . . . . 17trigonometrisches . . . . . . . . . . . . . . . 21
Pulcherima . . . . . . . . . . . . . . . . . . . . . . . . . 60
68
INDEX
R
Relativer Fehler . . . . . . . . . . . . . . . . . . . . . . . 8
S
Simpsonsformel . . . . . . . . . . . . . . . . . . . . . 60Splinefunktion . . . . . . . . . . . . . . . . . . . . . . 51Stetigkeit-Modul . . . . . . . . . . . . . . . . . . . . . 53
T
Trapezformel . . . . . . . . . . . . . . . . . . . . . . . . 59Trigonometrisches Polynom . . . . . . . . . . 21
Periode . . . . . . . . . . . . . . . . . . . . . . . . 21
V
Vielfachheit . . . . . . . . . . . . . . . . . . . . . . . . . 18
69
top related