4.6 Berechnung von Eigenwerten

4.6 Berechnung von Eigenwerten

In diesem Abschnitt befassen wir uns mit dem Eigenwertproblem: zu gegebener Matrix A∈R^n×n sind die Eigenwerte (und gegebenenfalls Eigenvektoren) gesucht. Wir erinnern an Definition 4.9 und formulieren

Satz 4.70 (Eigenwert). Es sei A∈R^n×n. Dann gilt:

• Die Matrix A hat genau n Eigenwerte, ihrer Vielfachheit nach gezählt.

• Die Eigenwerte sind Nullstellen des charakteristischen Polynoms:

det(A−λI) = 0.

• Reelle Matrizen können komplexe Eigenwerte haben, komplexe Eigenwerte treten stets als konjugierte Paare λ,¯λauf.

• Eigenvektoren zu unterschiedlichen Eigenwerten sind linear unabhängig.

• Fallsnlinear unabhängige Eigenvektoren existieren, so existiert eine reguläre Matrix S ∈R^n×n, so dassS⁻¹AS =D eine Diagonalmatrix ist. Die Spaltenvektoren von S sind die Eigenvektoren und die Diagonaleinträge von D die Eigenwerte.

• Symmetrische Matrizen A =A^T haben nur reelle Eigenwerte. Es existiert eine Or- thonormalbasis von Eigenvektoren und eine Diagonalisierung Q^TAQ=D mit einer orthogonalen Matrix.

• Bei Dreiecksmatrizen stehen die Eigenwerte auf der Diagonalen.

Zu einem Eigenwert λ sind die Eigenvektoren als Lösung des homogenen linearen Glei- chungssystems bestimmt:

(A−λI)w= 0.

Umgekehrt gilt:

Definition 4.71 (Rayleigh-Quotient). Sei A ∈ R^n×n sowie w ∈ Rⁿ ein Eigenvektor.

Dann ist durch den Rayleigh-Quotienten der zugehörige Eigenwert gegeben:

λ= (Aw, w)2

kwk²₂ .

Mit Hilfe dieser Definition folgt eine einfache Schranke für die Eigenwerte:

|λ| ≤sup

w6=0

(Aw, w)₂

kwk²₂ ≤ kAk2kwk²₂

kwk²₂ =kAk2.

Wir haben in Abschnitt 4.1 bereits gesehen, dass diese Schranke für die Beträge der Ei- genwerte in jeder Matrixnorm mit verträglicher Vektornorm gilt.

4.6.1 Konditionierung der Eigenwertaufgabe

Bevor wir auf konkrete Verfahren zur Eigenwertberechnung eingehen, analysieren wir die Kondition der Aufgabe, d.h. die Abhängigkeit der Eigenwerte von der Störung der Matrix.

Hierfür benötigen wir zunächst einen Hilfsatz:

Hilfsatz 4.72. Es seien A, B ∈R^n×n beliebige Matrizen. Dann gilt für jeden Eigenwert λ von A, der nicht zugleich Eigenwert von B für jede natürliche Matrizennorm die Ab- schätzung:

k(λI−B)⁻¹(A−B)k ≥1.

Beweis: Es sei w6= 0 ein Eigenvektor zuλ. Dann gilt:

(A−B)w= (λI−B)w.

Wennλkein Eigenwert vonB ist, so ist die Matrix (λI−B) regulär, d.h. es folgt:

(λI−B)⁻¹(A−B)w=w.

Und schließlich durch Normbilden:

1 = k(λI−B)⁻¹(A−B)wk

kwk ≤sup

x

k(λI−B)⁻¹(A−B)xk

kxk =k(λI−B)⁻¹(A−B)k.

Auf dieser Basis erhalten wir ein einfaches Kriterium zur Eingrenzung der Eigenwerte einer Matrix:

Satz 4.73 (Gerschgorin-Kreise). Es sei A∈R^n×n. Alle EigenwerteλvonA liegen in der Vereinigung der Gerschgorin-Kreise:

K_i={z∈C: |z−a_ii| ≤

n

X

k=1,k6=i

|aik|}, i= 1, . . . , n.

Angenommen zu den Indexmengen I_m = {i1, . . . , i_m} und I_m^′ = {1, . . . , n} \I_m seien die Vereinigungen Um = ^S_i∈ImKi und U_m^′ = ^S_i∈I′

mKi disjunkt. Dann liegen genau m Eigenwerte (ihrer algebraischen Vielfachheit nach gezählt) in U_m und n−m Eigenwerte in U_m^′ .

Beweis:(i)Es seiD∈R^n×ndie DiagonalmatrixD= diag(a_ii). Weiter seiλein Eigenwert vonAmitλ6=a_ii. (In diesem Fall wäre die Aussage des Satzes trivial erfüllt). Hilfsatz 4.72 besagt bei Wahl der maximalen Zeilensummennorm:

1≤ k(λI−D)⁻¹(A−D)k∞= max

i







n

X

j=1,j6=i

|(λ−a_ii)⁻¹a_ij|







= max

i







|λ−a_ii|⁻¹

n

X

j=1,j6=i

|aij|





 .

D.h. es existiert zu jedem λein Indexi∈ {1, . . . , n}so dass gilt:

|λ−a_ii| ≤

n

X

j=1,j6=i

|aij|.

(ii) Es seien durchImundI_m^′ Indexmengen mit oben genannter Eigenschaft gegeben. Wir definieren die Matrix

A_s:=D+s(A−D),

mitA0=Dund A1 =A. Entsprechend definieren wir die Vereinigungen:

U_m,s:= ^[

i∈Im

K_i(A_s), U_m,s^′ := ^[

i∈I^′_m

K_i(A_s), K_i(A_s) ={z∈C, |z−a_ii| ≤s^X

j6=i

|aij|}.

Aufgrund der stetigen Abhängigkeit der Kreisradien von s gilt U_m,s∩U_m,s^′ = ∅ für s ∈ [0,1]. Im Fall s = 0 gilt A₀ = D und jeder Eigenwert von A₀ liegt im Mittelpunkt (also λi =aii) des trivialen Kreises mit Radius Null. Das Ergebnis folgt nun durch die stetige

Abhängigkeit der Eigenwerte von s.

Wir betrachten hierzu ein Beispiel:

Beispiel 4.74 (Gerschgorin-Kreise). Wir betrachten die Matrix

A=







2 0.1 −0.5 0.2 3 0.5

−0.4 0.1 5





,

mit den Eigenwerten spr(A)≈ {1.91,3.01,5.08}. Eine erste Schranke liefern verschiedene mit Vektornormen verträgliche Matrixnormen von A:

kAk_∞= 5.5, kAk₁ = 6, kAk₂≈5.1.

Die Gerschgorin-Kreise von A sind:

K₁ ={z∈C, |z−2| ≤0.6}, K₂ ={z∈C, |z−3| ≤0.7}, K₃ ={z∈C, |z−5| ≤0.5}.

Diese Abschätzung kann verfeinert werden, daAund A^T die gleichen Eigenwerte besitzen.

Die Radien der Gerschgorin-Kreise können auch als Summe der Spaltenbeträge berechnet werden. Zusammen ergibt sich:

K₁ ={z∈C, |z−2| ≤0.6}, K₂ ={z∈C, |z−3| ≤0.2}, K₃ ={z∈C, |z−5| ≤0.5}.

Alle drei Kreise sind disjunkt.

Wir können nun den allgemeinen Stabilitätssatz für das Eigenwertproblem beweisen:

Satz 4.75 (Stabilität des Eigenwertproblems). Es seiA∈R^n×n eine Matrix mitnlinear unabhängigen Eigenvektoren w₁, . . . , w_n. Durch A˜ = A +δA sei eine beliebig gestörte Matrix gegeben. Dann existiert zu jedem Eigenwert λ( ˜A) von A˜=A+δA ein Eigenwert λ(A) von A so dass mit der Matrix W := (w₁, . . . , w_n) gilt:

|λ(A)−λ( ˜A)| ≤cond2(W)kδAk.

Beweis:Es gilt füri= 1, . . . , ndie BeziehungAwi=λi(A)wi, oder in Matrixschreibweise AW =Wdiag(λ_i(A)), also

W⁻¹AW = diag(λ_i(A)).

Wir betrachten nun einen Eigenwert ˜λ=λ( ˜A). Falls ˜λauch Eigenwert vonA ist, so gilt die Behauptung. Also sei ˜λnun kein Eigenwert vonA. Dann folgt:

k(˜λI−A)⁻¹k2 =kW⁻¹[˜λI−diag(˜λ_i(A))]⁻¹Wk2≤cond₂(W)k[˜λI−diag(λ_i(A))]⁻¹k2. Für die (symmetrische) Diagonalmatrix ˜λI−diag(λi(A)) gilt

k[˜λI−diag(λ_i(A))]⁻¹k₂= max

i=1,...,n|˜λ−λ_i(A)|⁻¹. Mit Hilfsatz 4.72 folgt das gewünschte Ergebnis:

1≤ k[˜λI−A)⁻¹δAk2≤ k[˜λI−A)⁻¹k kδAk2≤cond2(W) max

i=1,...,n|˜λ−λi(A)|⁻¹₂ kδAk2. Die Konditionierung des Eigenwertproblems einer MatrixAhängt von der Konditionszahl der Matrix der Eigenvektoren ab. Für symmetrische (hermitesche) Matrizen existiert eine Orthonormalbasis von Eigenvektoren mit cond₂(W) = 1. Für solche Matrizen ist das Eigenwertproblem gut konditioniert. Für allgemeine Matrizen kann das Eigenwertproblem beliebig schlecht konditioniert sein.

4.6.2 Direkte Methode zur Eigenwertberechnung

Die Eigenwerte einer Matrix A können prinzipiell als Nullstellen des charakteristischen Polynoms χ_A(z) = det(zI−A) berechnet werden. Die Berechnung der Nullstellen kann zum Beispiel mit dem Newton-Verfahren geschehen, Startwerte können mit Hilfe der Gerschgorin-Kreise bestimmt werden.

In Kapitel 2 haben die Berechnung von Nullstellen eines Polynoms jedoch als ein sehr schlecht konditioniertes Problem kennengelernt. Wir betrachten hierzu ein Beispiel.

Beispiel 4.76 (Direkte Berechnung von Eigenwerten). Es sei A∈R^5×5 eine Matrix mit den Eigenwertenλ_i = 1, i= 1, . . . ,5. Dann gilt:

χ_A(z) =

5

Y

i=1

(z−i) =z⁵−15z⁴+ 85x³−225x²+ 274z−120.

Der Koeffizient −55 vor z⁹ sei mit einem relativen Fehler von 0.1% gestört:

˜

χ_A(z) =z⁵−0.999·15z⁴+ 85x³−225x²+ 274z−120.

Dieses gestörte Polynom hat die Nullstellen (d.h. Eigenwerte):

λ1≈0.999, λ2 ≈2.05, λ3≈2.76, λ_4/5 ≈4.59±0.430i.

Die Eigenwerte können also nur mit einem (ab λ₃) wesentlichen Fehler bestimmt werden.

Es gilt etwa:

|4.59 + 0.43i−5|

5 ≈0.1,

d.h. der Fehler in den Eigenwerten beträgt 10%, eine Fehlerverstärkung von 100.

Das Aufstellen des charakteristischen Polynoms erfordert eine Vielzahl schlecht konditio- nierter Additionen sowie Multiplikationen. Für allgemeine Matrizen verbietet sich dieses direkte Verfahren. Lediglich für spezielle Matrizen wie Tridiagonalsysteme lassen die Ei- genwerte bestimmen, ohne das zunächst die Koeffizienten des charakteristischen Polynoms explizit berechnet werden müssen.

4.6.3 Iterative Verfahren zur Eigenwertberechnung

Das vorangegangene Beispiel widerspricht scheinbar zunächst der (auf jeden bei hermite- schen Matrizen) guten Konditionierung der Eigenwertsuche. Wir leiten nun stabile nume- rische Verfahren her, die die Eigenwerte nicht mehr direkt berechnen, sondern sie iterativ approximieren.

Zur Herleitung eines ersten einfachen Iterationsverfahren machen wir die folgende Beob- achtung: angenommen, zu gegebener Matrix A∈R^n×n existiere eine Basis aus Eigenvek- toren{w1, . . . , wn}, d.h. die Matrix sei diagonalisierbar. Weiter gelte|λn|>|λn−1| ≥ · · · ≥

|λ1|>0. Dann seix∈Rⁿ in Basisdarstellung:

x⁽⁰⁾ =

n

X

j=1

α_jw_j.

Wir nehmen an, dassαn6= 0, dass also die Komponente zum Eigenvektor des betragsmäßig größten Eigenwertes nicht trivial ist. Wir definieren die Iteration

x⁽ⁱ⁺¹⁾ = Ax⁽ⁱ⁾ kAx⁽ⁱ⁾k. Dann gilt

x⁽ⁱ⁺¹⁾= A_kAx^Ax(i−1)_(i−1)k

kA_kAx^Ax(i−1)_(i−1)kk = A²x⁽ⁱ⁻¹⁾ kA²x⁽ⁱ⁻¹⁾k

kAx⁽ⁱ⁻¹⁾k

kAx⁽ⁱ⁻¹⁾k =· · ·= Aⁱ⁺¹x⁽⁰⁾

kAⁱ⁺¹x⁽⁰⁾k. (4.10)

Weiter gilt mit der Basisdarstellung von x⁽⁰⁾:

Aⁱx⁽⁰⁾=

n

X

j=1

α_jλⁱ_jw_j =α_nλⁱ_n



w_n+

n−1

X

j=1

α_j α_n

λⁱ_j λⁱ_n



. (4.11)

Es gilt|λj/λn|<1 für j < n, daher folgt:

Aⁱx⁽⁰⁾ =

n

X

j=1

α_jλⁱ_jw_j =α_nλⁱ_n(w_n+o(1)),

Hieraus folgt durch Normierung:

Aⁱx⁽⁰⁾

kAⁱx⁽⁰⁾k = α_nλⁱ_n

|αnλⁱ_n|

! w_n

kwnk +o(1)→span{wn} (i→ ∞).

Die Iteration läuft in den Raum, der durchw_n aufgespannt wird. Für einen Vektorw, der Vielfaches eines Eigenvektors istw=sw_n gilt:

Aw=sAw_n=sλ_nw_n=λ_nw.

Diese vektorwertige Gleichung gilt in jeder Komponente, kann daher nach dem Eigenwert aufgelöst werden:

λ_n= [Aw]_k w_k .

Wir fassen zusammen:

Satz 4.77 (Potenzmethode nach von Mises). Es sei A∈R^n×n eine Matrix n linear un- abhängigen Eigenvektoren {w1, . . . , wn}. Der betragsmäßig größte Eigenwert sei separiert

|λ_n|>|λ_n−1| ≥ · · · ≥ |λ₁|. Es sei x⁽⁰⁾ ∈Rⁿ ein Startwert mit nichttrivialer Komponente in Bezug auf w_n. Für einen beliebigen Indexk∈ {1, . . . , n} konvergiert die Iteration

˜

x⁽ⁱ⁾=Ax⁽ⁱ⁻¹⁾, x⁽ⁱ⁾ := x˜⁽ⁱ⁾

kx˜⁽ⁱ⁾k, λ⁽ⁱ⁾:= x˜⁽ⁱ⁾_k x⁽ⁱ⁻¹⁾_k , gegen den betragsmäßig größten Eigenwert:

|λ⁽ⁱ⁺¹⁾−λ_n|=O

λ_n−1 λn

i!

, i→ ∞.

Beweis: Wir knüpfen an der Vorbereitung des Beweises (4.10) an. Es gilt:

λ⁽ⁱ⁾ = x˜⁽ⁱ⁾_k

x⁽ⁱ⁻¹⁾_k = [Ax⁽ⁱ⁻¹⁾]_k

x⁽ⁱ⁻¹⁾_k = [Aⁱx⁽⁰⁾]_k [Aⁱ⁻¹x⁽⁰⁾]_k.

Weiter, mit (4.11) gilt:

λ⁽ⁱ⁾=

a_nλⁱ_n

[w_n]_k+^Pn−1_{j=1 α}^αj

n

λⁱ_j λⁱn[w_j]_k

anλⁱ⁻¹n

[wn]_k+^Pn−1_{j=1 α}^αj

n

λⁱ_j⁻¹ λⁱn⁻¹[wj]_k

=λ_n [w_n]_k+^Pn−1_{j=1 α}^αj

n

λⁱ_j λⁱ_n[w_j]_k [wn]_k+^Pn−1_{j=1 α}^αk

n

λⁱ_j⁻¹ λⁱn⁻¹[wj]_k

=λn

[w_n]_k+^α_αⁿ⁻¹

n +o(1) ^λ_λⁿ⁻¹

n

i[w_n]_k [w_n]_k+^α_αⁿ⁻¹

n +o(1) ^λ_λⁿ⁻¹

n

i−1[w_n]_k

=λn 1 +O

λn−1

λ_n

i−1!!

.

Die letzte Abschätzung nutzt die Reihenentwicklung:

1 +xⁿ⁺¹

1 +xⁿ = 1 +O(|x|ⁿ).

Die Potenzmethode ist eine sehr einfache und numerisch stabile Iteration. Sie kann al- lerdings nur den betragsmäßig größten Eigenwert ermitteln. Der Konvergenzbeweis kann verallgemeinert werden auf Matrizen, deren größter Eigenwert mehrfach vorkommt. Kon- vergenz gilt jedoch nur dann, wenn aus |λn|=|λi| auch λ_n =λ_i folgt. Es darf also nicht zwei verschiedene Eigenwerte geben, die beide den gleichen, größten Betrag annehmen.

Dies schließt zum Beispiel den Fall zweier komplex konjugierter Eigenwerten λund ¯λals beträgsgrößte aus.

Bemerkung 4.78 (Potenzmethode bei symmetrischen Matrizen). Die Potenzmethode kann im Fall symmetrischer Matrizen durch Verwenden des Rayleigh-Quotienten verbes- sert werden. Die Iteration

λ⁽ⁱ⁾= (Ax⁽ⁱ⁻¹⁾, x⁽ⁱ⁻¹⁾)2

(x⁽ⁱ⁻¹⁾, x⁽ⁱ⁻¹⁾)₂ = (˜x⁽ⁱ⁾, x⁽ⁱ⁻¹⁾)2, liefert die Fehlerabschätzung:

λ⁽ⁱ⁾=λ_n+O

λ_n−1 λn

2i! .

Beispiel 4.79 (Potenzmethode nach von Mises). Es sei:

A=







2 1 2

−1 2 1

1 2 4





,

mit den Eigenwerten:

λ1≈0.80131, λ2= 2.2865, λ3= 4.9122.

Wir starten die Potenzmethode mit x⁽⁰⁾ = (1,1,1)^T, und wählen zur Normierung die Maximumsnorm. Weiter wählen wir als Index k= 1:

i= 1 : x˜⁽¹⁾=Ax⁽⁰⁾ =





 5 2 7





, x⁽¹⁾= x˜⁽¹⁾ kx˜⁽¹⁾k_∞ ≈





 0.714 0.286

1





, λ⁽¹⁾ = x˜⁽¹⁾₁ x⁽⁰⁾₁ ≈5,

i= 2 : x˜⁽²⁾=Ax⁽¹⁾ =





 3.71 0.857

5.29





, x⁽²⁾= x˜⁽²⁾ kx˜⁽²⁾k∞

≈





 0.702 0.162

1





, λ⁽²⁾ = x˜⁽²⁾₁

x⁽¹⁾₁ ≈5.19,

i= 3 : x˜⁽³⁾=Ax⁽²⁾ =





 3.57 0.124

1





, x⁽³⁾= x˜⁽³⁾ kx˜⁽³⁾k∞

≈





 0.710 0.124

1





, λ⁽³⁾ = x˜⁽³⁾₁

x⁽²⁾₁ ≈5.077,

i= 4 : x˜⁽⁴⁾=Ax⁽³⁾ =





 3.54 0.538

4.96





, x⁽⁴⁾= x˜⁽⁴⁾ kx˜⁽⁴⁾k∞

≈





 0.715 0.108

1





, λ⁽⁴⁾ = x˜⁽⁴⁾₁

x⁽³⁾₁ ≈4.999,

i= 5 : x˜⁽⁵⁾=Ax⁽⁴⁾ =





 3.54 0.502

4.93





, x⁽⁵⁾= x˜⁽⁵⁾ kx˜⁽⁵⁾k∞ ≈





 0.717 0.102

1





, λ⁽⁵⁾ = x˜⁽⁵⁾₁

x⁽⁴⁾₁ ≈4.950,

i= 6 : x˜⁽⁶⁾=Ax⁽⁵⁾ =





 3.54 0.486

4.92





, x⁽⁶⁾= x˜⁽⁶⁾ kx˜⁽⁶⁾k∞

≈





 0.719 0.0988

1





, λ⁽⁶⁾ = x˜⁽⁶⁾₁

x⁽⁵⁾₁ ≈4.930.