NumerikI - UZHuser.math.uzh.ch/berner/files/numerik.pdf · 1 ZahlendarstellungundRundungsfehler...

69
Numerik I Mitschri der Vorlesung von Prof. Dr. Chipot Tobias Berner Universität Zürich 2007

Transcript of NumerikI - UZHuser.math.uzh.ch/berner/files/numerik.pdf · 1 ZahlendarstellungundRundungsfehler...

Numerik I

Mitschri der Vorlesung vonProf. Dr. Chipot

Tobias BernerUniversität Zürich2007

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