Skript - Mathematische Optimierung

Nichtlineare Optimierung
Skript zur Vorlesung von
Dr. Stefan Körkel
Frühjahrssemester 2016
Universität Mannheim
c Dr. Stefan Körkel
Fakultät für Wirtschaftsinformatik und Wirtschaftsmathematik
Universität Mannheim
Version: 16. März 2016
Korrekturen und Anmerkungen bitte an [email protected]
Hinweise zur Vorlesung
Einordnung: Modul MAC 507, aber 4-stündig plus Übungen
Webseite zur Vorlesung:
http://mathopt.math.uni-mannheim.de/config/links/Vorlesung%20Nichtlineare%20Optimierung
Dozent
Dr. Stefan Körkel, Lehrstuhlvertreter Mathematische Optimierung
Sprechstunde nach den Vorlesungen oder nach Vereinbarung
Büro: Raum A5, B128
email: [email protected]
Termine
Vorlesung
Dienstags, 8:30 - 10:00, Raum B6, A101
Mittwochs, 8:30 - 10:00, Raum A5, C012
Beginn: 16.2.
Übungen
Freitags, 10:15-11:45, Raum A5, C012
Scheinkriterien
• Regelmäßige Teilnahme an den Übungen, mindestens einmal vorrechnen
• Mindestens 50% der Punkte aus den Übungsaufgaben
3
• Mindestens 3 erfolgreich bearbeitete Programmieraufgaben
• Bestehen der Klausur am Semesterende, Note ist die Klausurnote.
Mögliche Klausurtermine:
DI, 7. Juni
MI, 8. Juni
FR, 10. Juni
DI, 14. Juni
MI, 15. Juni
FR, 17. Juni
Übungen
Theoretische Aufgaben
Pro Woche ein Übungsblatt
Ausgabe der Aufgabenblätter mittwochs in der Vorlesung oder als PDF auf der Webseite
Rückgabe eine Woche später mittwochs in der Vorlesung
Praktische Programmieraufgaben
Programmiersprache freigestellt, z. B. Matlab, Octave, C, Fortran
Abgabe nach zwei Wochen als ausgedruckten Quelltext und ausgedruckte Ergebnisse.
Literaturempfehlungen
• Skript zu dieser Vorlesung: wird während des Semesters zeitnah zur Vorlesung erstellt
• Manuskript zur Vorlesung im Sommersemester 2012 an der Universität Heidelberg
• Skript von Prof. Kolb zur 2-stündigen Vorlesung 2015
• Nocedal, Wright: Numerical Optimization, Springer
• Ulbricht, Ulbricht: Nichtlineare Optimierung, Birkhäuser
• Jarre, Stoer: Optimierung, Springer
• Geiger, Kanzow: Numerische Verfahren zur Lösung unrestringierter Optimierungsverfahren, Springer
• Geiger, Kanzow: Theorie und Numerik restringierter Optimierungsaufgaben, Springer
• Bertsekas: Nonlinear Programming, Athena Scientific
Weitere Lehre
• Seminar Konvexe Optimierung
Freitags, 13:45-15:15, Raum B6, A104
Vorbesprechung und Themenvergabe: 19. Februar
• Bachelorarbeiten aus dem Bereich der nichtlinearen Optimierung
5
Inhaltsverzeichnis
1 Nichtlineare Optimierungsprobleme
7
Allgemeine Problemklasse . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
7
Lokale und globale Optima . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
9
Konvexe Optimierungsprobleme . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
9
Konvergenz von iterativen Verfahren . . . . . . . . . . . . . . . . . . . . . . . . . . . 10
Schreibweisen für Ableitungen . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 12
2 Unbeschränkte Optimierungsprobleme
14
Optimalitätsbedingungen . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 14
Abstiegsverfahren allgemein . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 15
Steilster Abstieg: Gradientenverfahren . . . . . . . . . . . . . . . . . . . . . . . . . . 16
Konvergenz des Gradientenverfahrens . . . . . . . . . . . . . . . . . . . . . . . . . . . 18
3 Konjugierte Gradienten (CG-)Methoden
23
Lineare CG-Verfahren . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 23
Nichtlineare CG-Methoden . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 32
Vergleich Gradientenverfahren – CG-Verfahren . . . . . . . . . . . . . . . . . . . . . . 39
4 Linesearch
40
Wahl der Schrittweite . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 40
Konvergenz von Linesearch-Methoden . . . . . . . . . . . . . . . . . . . . . . . . . . . 42
Linesearch-Algorithmen . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 43
6
Kapitel 1
Nichtlineare Optimierungsprobleme
Allgemeine Problemklasse
Wir betrachten in dieser Vorlesung Optimierungsprobleme der Form
min f (x)
g(x) = 0
h(x) ≥ 0
(1.1)
Bemerkung 1.1 o. B. d. A. min: durch Multiplikation der Zielfunktion mit −1 können Maximierungsprobleme in Minimierungsprobleme transformiert werden.
o. B. d. A. ≥: durch Multiplikation mit −1 können ≤-Ungleichungen in ≥-Ungleichungen transformiert werden.
Definition 1.2 f : Rn → R heißt Zielfunktion, g : Rn → Rm Gleichungsnebenbedingungen, h : Rn → Rk Ungleichungsnebenbedingungen und x ∈ Rn Optimierungsvariablen.
Wir betrachten endlich-dimensionale Probleme mit endlichen vielen Nebenbedingungen. Im
allgemeinen ist m < n.
Wir betrachten glatte Probleme: f, g, h ∈ C 2 oder C 3 , d. h. zweimal oder dreimal stetig
differenzierbar.
Wir unterscheiden
• unbeschränkte Probleme
• beschränkte Probleme
Natürlich können unbeschänkte Probleme auch mit Verfahren für beschränkte Probleme gelöst
werden. Wir werden aber auch Verfahren speziell für unbeschränkte Probleme behandeln.
Nicht behandelt werden in dieser Vorlesung:
7
8
KAPITEL 1. NICHTLINEARE OPTIMIERUNGSPROBLEME
• (Gemischt-) Ganzzahlige Probleme
• (Semi-) Infinite Probleme
• Nichtdifferenzierbare Probleme
• Optimalsteuerungsprobleme
Spezialfälle
• Lineare Optimierung: Zielfunktion und Nebenbedingungen sind linear:
min cT x
Ax = 0
x≥0
Dafür gibt es maßgeschneiderte Methoden, die in einer eigenen Vorlesung Lineare Opti”
mierung“ behandelt werden.
• Quadratische Optimierung: Die Zielfunktion ist quadratisch, die Nebenbedingungen
sind linear:
1
min xT Hx + g T x
2
Ax = b
Cx ≥ d
Treten als Teilprobleme in den Verfahren dieser Vorlesung auf.
• Eindimensionale Probleme: Eindimensionale Minimierung von
f :R⊃I→R
Treten als Hilfsprobleme für mehrdimensionale Probleme auf (Liniensuche), wichtig für
die Theorie.
• Nichtlineare Ausgleichsprobleme:
1
min kF1 (x)k22
2
F2 (x) = 0
F3 (x) ≥ 0
oder mit anderen Normen k.k1 , k.k∞ . Die spezielle Zielfunktion muss bei der Wahl des
Lösungsalgorithmus berücksichtigt werden. Siehe Kapitel ??.
KAPITEL 1. NICHTLINEARE OPTIMIERUNGSPROBLEME
9
Lokale und globale Optima
Definition 1.3 (Zulässiger Punkt) Ein Punkt x ∈ Rn heißt zulässiger Punkt des Optimierungsproblems (1.1), wenn er alle Nebenbedingungen erfüllt:
x zulässig :⇐⇒ g(x) = 0 ∧ h(x) ≥ 0.
Definition 1.4 (Lokales Minimum) Ein Punkt x∗ ∈ Rn heißt lokales Minimum des Optimierungsproblems (1.1), wenn er zulässig ist und sein Zielfunktionswert nicht größer ist als
die Zielfunktionswerte aller zulässigen Punkte in einer Umgebung von x∗ :
x∗ lokales Minimum :⇐⇒ g(x∗ ) = 0 ∧ h(x∗ ) ≥ 0 ∧
∃ Umgebung U von x∗ , so dass ∀x ∈ U
mit g(x) = 0, h(x) ≥ 0 gilt:
f (x∗ ) ≤ f (x).
Definition 1.5 (Globales Minimum) Ein Punkt x∗ ∈ Rn heißt globales Minimum des
Optimierungsproblems (1.1), wenn er zulässig ist und sein Zielfunktionswert nicht größer ist
als die Zielfunktionswerte aller zulässigen Punkte:
x∗ globales Minimum :⇐⇒ g(x∗ ) = 0 ∧ h(x∗ ) ≥ 0 ∧
f (x∗ ) ≤ f (x) ∀x mit g(x) = 0, h(x) ≥ 0.
Definition 1.6 (Striktes Minimum) Ein Minimum heißt strikt, wenn gilt:
f (x∗ ) < f (x)
f (x∗ ) < f (x)
∀x ∈ U, x zulässig, x 6= x∗
∀x, x zulässig, x 6= x∗
(lokal)
(global)
Die Bestimmung von bzw. Suche nach globalen Optima heißt globale Optimierung, Wir beschränken uns in dieser Vorlesung auf die Berechnung von lokalen Optima.
Konvexe Optimierungsprobleme
Definition 1.7 (Konvexe Menge) Eine Menge C ⊆ Rn heißt konvex, wenn die Strecke
zwischen zwei beliebigen Punkten aus C ganz in C liegt, d. h.
∀x1 , x2 ∈ C, θ ∈ [0, 1] : θx1 + (1 − θ)x2 ∈ C
TODO: Abbildung
Definition 1.8 (Konvexe Funktion) Eine Funktion f : Rn ⊃ D → R heißt konvex, wenn
D eine konvexe Menge ist und für alle x1 , x2 ∈ D und θ ∈ [0, 1] gilt:
f (θx1 + (1 − θ)x2 ) ≤ θf (x1 ) + (1 − θ)f (x2 )
10
KAPITEL 1. NICHTLINEARE OPTIMIERUNGSPROBLEME
TODO: Abbildung
Satz 1.9 Sei ein Minimierungsproblem gegeben mit einer konvexen Zielfunktion f und einer
konvexen Menge S aller zulässigen Punkte:
min f (x)
x∈S
f konvex
S konvex
Sei x∗ ein lokales Minimum. Dann ist x∗ auch ein globales Minimum.
Beweis. Sei x∗ lokales Minimum. Dann existiert R > 0, so dass f (x∗ ) ≤ f (x) für alle x ∈ S mit
kx − x∗ k < R.
Annahme: x∗ ist nicht globales Minimum, d. h. es existiert x̄ ∈ S, x̄ 6= x∗ mit f (x̄) < f (x∗ ). Es
gilt kx̄ − x∗ k ≥ R, denn sonst wäre f (x∗ ) ≤ f (x̄).
Betrachte x̂ = (1 − θ)x∗ + θx̄ mit θ =
R
2kx̄−x∗ k
Es gilt nun kx̂ − x∗ k = kθ(x̄ − x∗ )k =
R
2
∈ (0; 21 ]. Wegen der Konvexität von S ist x̂ ∈ S.
< R. Daraus folgt f (x∗ ) ≤ f (x̂).
Wegen der Konvexität von f ist aber f (x̂) ≤ (1 − θ)f (x∗ ) + θf (x̄) < f (x∗ ).
TODO: Abbildung (3D)
Mehr zur konvexen Optimierung z. B. in dem Buch
Boyd, Vandenberghe: Convex Optimization, Cambridge University Press.
Konvergenz von iterativen Verfahren
Verfahren zur Bestimmung von Lösungen, d. h. Optima nichtlinearer Optimierungsprobleme
sind im allgemeinen iterativ.
Algorithmus 1.10 (Iteratives Verfahren)
1. Wähle einen Startwert x0 . Setze j := 0.
2. Mache einen Iterationsschritt: xj+1 := X j (xj ).
3. Überprüfe ein Konvergenzkriterium und stoppe ggf.
4. Setze j := j + 1 und gehe zu 2.
Definition 1.11 (Konvergenz) Sei k.k eine Norm im Variablenraum. Sei {x0 , x1 , . . .} eine
durch das iterative Verfahren 1.10 ohne Schritt 3 berechnete Folge.
1. Die Folge {x0 , x1 , . . .} konvergiert, wenn es ein x∗ gibt mit kxj − x∗ k → 0.
2. Das iterative Verfahren konvergiert lokal gegen ein x∗ , wenn es eine Umgebung U von
x∗ gibt, so dass für alle Startwerte x0 ∈ U die Folgen {x0 , x1 , . . .} gegen x∗ konvergieren.
3. Das iterative Verfahren konvergiert global, wenn es für jeden Startwert x0 ein x∗ gibt,
so dass die Folge {x0 , x1 , . . .} gegen x∗ konvergiert.
KAPITEL 1. NICHTLINEARE OPTIMIERUNGSPROBLEME
11
Fragen zur Konvergenz
• Konvergiert ein Verfahren lokal?
• Konvergiert ein Verfahren global?
Methoden zur Erzielung globaler Konvergenz bei zunächst nur lokal konvergenten Verfahren
nennt man Globalisierungsmethoden.
Unterscheide:
lokale Optima — lokale Konvergenz
globale Optima — globale Konvergenz
Globale Optimierung — Globalisierungsmethoden
• Wie schnell ist die Konvergenz?
Definition 1.12 (Konvergenzordnung) Sei {xj } ⊂ Rn die Folge der Iterierten.
• {xj } konvergiert gegen ein x∗ ∈ Rn (mindestens) linear, falls ein c ∈ (0; 1) existiert, so
dass
kxj+1 − x∗ k ≤ ckxj − x∗ k
für alle j ∈ N hinreichend groß.
• {xj } konvergiert gegen x∗ (mindestens) superlinear, falls eine Nullfolge {cj } ⊂ R+
existiert mit
kxj+1 − x∗ k ≤ cj kxj − x∗ k
für alle j ∈ N hinreichend groß.
• Gilt xj → x∗ , so konvergiert {xj } gegen x∗ (mindestens) quadratisch, falls ein C > 0
existiert mit
kxj+1 − x∗ k ≤ Ckxj − x∗ k2
für alle j ∈ N hinreichend groß.
Beispiele:
• {xj } mit kxj − x∗ k =
j
j
∗
• {x } mit kx − x k =
• {xj } mit kxj − x∗ k =
1 j
2
j
konvergiert linear gegen x∗
1
j
konvergiert superlinear gegen x∗
j
1 2
2
konvergiert quadratisch gegen x∗
Bemerkung 1.13 (R-Konvergenz) Die Konvergenzbegriffe aus Definition 1.12 nennt man
auch Q-lineare, Q-superlineare und Q-quadratische Konvergenz. Daneben gibt es auch:
12
KAPITEL 1. NICHTLINEARE OPTIMIERUNGSPROBLEME
• R-linear :⇐⇒ ∃{aj , j ∈ N}, aj → 0 Q-linear : kxj − x∗ k ≤ aj ∀j ∈ N
• R-superlinear :⇐⇒ ∃{aj , j ∈ N}, aj → 0 Q-superlinear : kxj − x∗ k ≤ aj ∀j ∈ N
• R-quadratisch :⇐⇒ ∃{aj , j ∈ N}, aj → 0 Q-quadratisch : kxj − x∗ k ≤ aj ∀j ∈ N
(Q: quotient“, R: root“)
”
”
Dabei gilt:
• R-superlineare Konvergenz ist äquivalent zu
p
lim j kxj − x∗ k = 0
j→∞
• Aus der Q-quadratischen Konvergenz folgt stets R-quadratische Konvergenz, aber nicht
umgekehrt (linear, superlinear analog).
Schreibweisen für Ableitungen
Definition 1.14 Sei f : Rn ⊃ D → R, f ∈ C 2 .
Dann ist
• f 0 (x) =
∂f
∂x1
...
∂f
∂xn

∂f
∂x1
die Ableitung, ein Zeilenvektor,

 
• ∇f (x) = f 0 (x)T =  ...  der Gradient, ein Spaltenvektor,
∂f
∂xn


• ∇2 f (x) = ∇(f 0 )(x) = (∇f )0 (x) = 
metrische Matrix.
∂2f
∂x1 ∂x1
...
..
.
∂2f
∂xn ∂x1
∂2f
∂x1 ∂xn
..
.
...
∂2f
∂xn ∂xn


 die Hessematrix, eine sym-
Definition 1.15 Sei H eine symmetrische Matrix.
H ist positiv definit :⇐⇒ H 0 :⇐⇒ xT Hx > 0 ∀x 6= 0.
H ist positiv semidefinit :⇐⇒ H 0 :⇐⇒ xT Hx ≥ 0 ∀x.
Definition 1.16 Sei g : Rn ⊃ D → Rm , g ∈ C 1 .
Dann ist
KAPITEL 1. NICHTLINEARE OPTIMIERUNGSPROBLEME
 ∂g1
∂x1
...

• g 0 (x) =  ...
∂gm
∂x1
 ∂g1
∂x1

..  die Jacobimatrix,
. 
...
...

• ∇g(x) =  ...
∂g1
∂xn
∂g1
∂xn
∂gm
∂xn
∂gm
∂x1

..  = g 0 (x)T .
. 
...
∂gm
∂xn
13
Kapitel 2
Unbeschränkte Optimierungsprobleme
Wir betrachten Probleme der Form
min f (x)
(2.1)
mit f : Rn → R, f ∈ C 2 .
Optimalitätsbedingungen
Lemma 2.1 (Eindimensionaler Fall) Sei f : R ⊃ (a; b) → R, f ∈ C 2 . Dann gilt:
1. Wenn x∗ ∈ (a; b) ein lokales Minimum von f ist, dann ist f 0 (x∗ ) = 0 (notwendige Bedingung 1. Ordnung) und f 00 (x∗ ) ≥ 0 (notwendige Bedingung 2. Ordnung).
2. Wenn f 0 (x∗ ) = 0 und f 00 (x∗ ) > 0 für ein x∗ ∈ (a; b), dann ist x∗ ein striktes lokales
Minimum (hinreichende Bedingung 2. Ordnung).
Beweis. Übungsaufgabe.
Satz 2.2 (Notwendige Bedingung 1. und 2. Ordnung für den unbeschränkten Fall)
Sei f : Rn → R, f ∈ C 2 . Dann gilt: Wenn x∗ ein lokales Minimum von f ist, ist
1. ∇f (x∗ ) = 0 (notwendige Bedingung 1. Ordnung) und
2. ∇2 f (x∗ ) 0 (notwendige Bedingung 2. Ordnung).
Beweis. Sei p ∈ Rn beliebig. Die eindimensionale Funktion f˜(t) := f (x∗ + tp), t ∈ R hat ein
lokales Minimum in t = 0. Also gilt nach Lemma 2.1:
d
1. f˜0 (0) = dt
f (x∗ + tp)t=0 = f 0 (x∗ )p = ∇f (x∗ )T p = 0. Daraus folgt ∇f (x∗ ) = 0, da p
beliebig ist.
14
KAPITEL 2. UNBESCHRÄNKTE OPTIMIERUNGSPROBLEME
2. f˜00 (0) =
d2
f (x∗
dt2
+ tp)
15
= pT ∇2 f (x∗ )p ≥ 0. Also ist ∇2 f (x∗ ) 0.
t=0
Satz 2.3 (Hinreichende Optimalitätsbedingung für den unbeschränkten Fall) Sei f :
Rn → R, f ∈ C 2 . Sei x∗ ∈ Rn mit ∇f (x∗ ) = 0 und ∇2 f (x∗ ) 0. Dann gilt: x∗ ist ein striktes
lokales Minimum von f .
Beweis. Es existiert eine Umgebung U von x∗ mit ∇2 f (x) 0 ∀x ∈ U , denn die Eigenwerte
hängen stetig von den Matrixeinträgen ab. Entwicklung von f in eine Taylorreihe um x∗ ergibt:
1
f (x) = f (x∗ ) + ∇f (x∗ )T (x − x∗ ) + (x − x∗ )T ∇2 f (x̃)(x − x∗ )
{z
} 2|
{z
}
|
=0
>0
mit x̃ ∈ U . Daraus folgt: f (x) > f (x∗ ) für x ∈ U , x 6= x∗ .
Abstiegsverfahren allgemein
• Abstiegsverfahren sind einfache Verfahren und leicht zu implementieren.
• Sie sind eine Referenz im Vergleich mit anderen Verfahren im Hinblick auf Implementierbarkeit und Konvergenzverhalten.
Algorithmus 2.4 (Konzept der Abstiegsverfahren)
1. Wähle einen Startwert x0 ∈ Rn . j := 0.
2. Bestimme eine Richtung ∆xj ∈ Rn .
3. Bestimme eine Schrittweite αj > 0.
4. Iteriere xj+1 := xj + αj ∆xj . j := j + 1.
5. Gehe zu 2 oder stoppe bei Konvergenz“ oder Abbruch“.
”
”
Dabei wird
• ∆xj so gewählt, dass entlang xj + α∆xj , α > 0 ein Abstieg möglich ist,
• αj > 0 so gewählt, dass entlang der Abstiegsrichtung ein “qualifizierter” Abstieg erfolgt.
Die Wahl von αj kann z. B. durch Liniensuche (Linesearch) erfolgen, s. u.
Definition 2.5 (Abstiegsrichtung) Sei f : Rn → R und x ∈ Rn . Ein Vektor ∆x ∈ Rn heißt
Abstiegsrichtung von f in x, wenn es ein t̄ > 0 gibt mit
f (x + t∆x) < f (x)
für alle t ∈ (0; t̄].
16
KAPITEL 2. UNBESCHRÄNKTE OPTIMIERUNGSPROBLEME
Lemma 2.6 Sei f : Rn → R stetig differenzierbar, x ∈ Rn und ∆x ∈ Rn mit
∇f (x)T ∆x < 0.
Dann ist ∆x eine Abstiegsrichtung von f in x.
Beweis. Es gilt
f (x + h∆x) − f (x)
< 0.
h&0
h
∇f (x)T ∆x = lim
Also ist für alle hinreichend kleinen h > 0
f (x + h∆x) − f (x)
<0
h
und ∆x ist damit eine Abstiegsrichtung.
Bei Abstiegsverfahren wählt man als Richtung eine Abstiegsrichtung und eine Schrittweite, die
den Abstieg bewerkstelligt.
Definition 2.7 (Abstiegsverfahren) Seit {x0 , x1 , . . .} eine Folge von Iterierten eines Verfahrens zur Minimierung von f (x), x ∈ Rn . Das Verfahren ist ein Abstiegsverfahren, wenn
f (xj+1 ) < f (xj ) für (fast) alle j ∈ N.
Steilster Abstieg: Gradientenverfahren
Lemma 2.8 (Eigenschaften des Gradienten) Sei f : Rn → R stetig differenzierbar, x ∈
Rn und g = ∇f (x). Dann gilt:
1. g ist die Richtung des steilsten Anstiegs von f in x.
2. Sei ϕ : (a; b) → Rn eine stetig differenzierbare Kurve, die auf einer Höhenlinie von f
verläuft: f (ϕ(t)) = konst ∀t ∈ (a; b). Sei t̃ ∈ (a; b) und x̃ = ϕ(t̃). Dann ist
∇f (x̃)T ϕ0 (t̃) = 0,
d. h. der Gradient steht senkrecht auf der Höhenlinie {x : f (x) = konst}.
Beweis. Übungsaufgabe
TODO: Abbildung
KAPITEL 2. UNBESCHRÄNKTE OPTIMIERUNGSPROBLEME
17
Gradientenverfahren
• Als Richtung ∆xj nimmt man den steilsten Abstieg entlang des negativen Gradienten.
• Der Schritt soll die Länge 1 haben: k∆xj k = 1.
Taylorentwicklung (Linearisierung):
f (xj + ∆xj ) − f (xj ) ≈ ∇f (xj )T ∆xj
(2.2)
∆xj ist also die Lösung des Minimierungsproblems
min{∇f (xj )T ∆x s.t. k∆xk = 1}.
∆x
TODO: Abbildung. Siehe auch Übungsaufgabe.
√
Für die euklidische Norm kxk = xT x ist die Lösung
∆xj = −
∇f (xj )
.
k∇f (xj )k
Algorithmus 2.9 (Gradientenverfahren)
Eingabe: Startwert x0 , Abbruchgenauigkeit > 0, ᾱ zur Schrittweitensteuerung.
1. Startwert x0 , j := 0.
2. Berechne f (xj ) und ∇f (xj ).
3. Falls k∇f (xj )k ≤ : STOP.
4. ∆xj := −∇f (xj )/k∇f (xj )k.
5. Liniensuche: Bestimme αj ∈ (0; ᾱ] so, dass f (xj + αj ∆xj ) = min f (xj + α∆xj ).
α∈(0;ᾱ]
6. xj+1 := xj + αj ∆xj .
7. j := j + 1, weiter bei 2.
Bemerkung 2.10
1. Der Rechenaufwand dieses Verfahrens besteht in den Funktions- und Gradientenauswertungen.
2. Falls f (xj ) oder ∆f (xj ) schon in der Liniensuche des vorherigen Schritts berechnet wurden, entfällt ein Teil von 2.
3. Die Bestimmung der Schrittweite αj lässt sich als Methode zur Bestimmung des Gültigkeitsbereichs der Linearisierung (2.2) verstehen.
4. Statt der exakten“ Liniensuche in Schritt 5. wird in der Praxis oft eine inexakte“ Li”
”
niensuche verwendet, s. u.
18
KAPITEL 2. UNBESCHRÄNKTE OPTIMIERUNGSPROBLEME
Konvergenz des Gradientenverfahrens
Definition 2.11 (Niveaumenge) Die Menge N (x̂) := {x ∈ Rn : f (x) ≤ f (x̂)} heißt Niveaumenge von f zum Punkt x̂.
Satz 2.12 Sei x0 ∈ Rn beliebig und die Niveaumenge N (x0 ) sei kompakt. Dann existiert ein
Häufungspunkt x∗ der Folge {x0 , x1 , . . .}, die vom Gradientenverfahren (ohne Abbruch) berechnet wurde, und ∇f (x∗ ) = 0, oder es ist ∇f (xj ) = 0 für ein xj .
Beweis. Es gilt: f (xj+1 ) ≤ f (xj ) ∀j, daher ist {xj } ⊂ N (x0 ). Also existiert ein Häufungspunkt
0
x∗ ∈ N (x0 ) und eine konvergente Teilfolge xj → x∗ mit
0
0
0
0
0
0
f (x∗ ) ≤ f (xj +1 ) = f (xj − αj ∇f (xj )/k∇f (xj )k) ≤ f (xj ).
Sei ∇f (xj ) 6= 0 ∀j. Annahme: ∇f (x∗ ) 6= 0.
Dann haben wir Abstieg in x∗ : Es existiert ein α∗ > 0 und δ > 0, so dass
f (x∗ − α∗ ∇f (x∗ )/k∇f (x∗ ))k) = f (x∗ ) − α∗ ∇f (x∗ )T ∇f (x∗ )/k∇f (x∗ )k + O(kα∗ k2 ) ≤ f (x∗ ) − δ.
Weil aber f und ∇f stetig sind, gilt auch
δ
0
0
0
f (|{z}
xj −α∗ ∇f (xj ) /k ∇f (xj ) k) ≤ f (x∗ ) −
|
{z
}
{z
}
|
2
∗
→x
→∇f (x∗ )
→∇f (x∗ )
für fast alle j 0 .
Damit ist
0
0
0
f (xj − α∗ ∇f (xj )/k∇f (xj )k) ≤ f (x∗ ) −
δ
0
0
0
0
< f (x∗ ) ≤ f (xj − αj ∇f (xj )/k∇f (xj )k).
2
0
Das steht im Widerspruch zur Minimalität von αj in der Linesearch in Schritt j 0 .
TODO: Abbildung Steilster Abstieg, Zickzack-Konvergenz
Lemma 2.13 Für das Gradientenverfahren mit exakter Liniensuche gilt:
T
∆xj+1 ∆xj = 0 ∀j.
Beweis. Wegen
f (xj + αj ∆xj ) = min f (xj + α∆xj )
α
ist
d
j
j 0=
f (x + α∆x )
= ∇f (xj+1 )T ∆xj .
dα
α=αj
Im Gradientenverfahren ist
∇f (xj+1 )
j+1
∆x
=−
.
k∇f (xj+1 )k
Es folgt
∇f (xj+1 )T ∆xj
T
= 0.
∆xj+1 ∆xj = −
k∇f (xj+1 )k
KAPITEL 2. UNBESCHRÄNKTE OPTIMIERUNGSPROBLEME
19
Konvergenzgeschwindigkeit
Die Orthogonalität aufeinanderfolgender Schritte kann zu einer langsamen Zickzack-Konvergenz
in der Nähe des Minimums führen.
In der Nähe des Minimums
1
f (x) = f (x∗ ) + ∇f (x∗ )T (x − x∗ ) + (x − x∗ )T ∇2 f (x∗ )(x − x∗ ) + O(kx − x∗ k3 )
| {z }
| {z }
2
0
=0
sieht f für den Algorithmus quadratisch aus. Für die Untersuchung der Asymptotik reicht es
aus, folgendes Problem zu betrachten:
1
min xT Ax
2
mit A 0 positiv definit.
Lemma 2.14 Bei exakter Liniensuche gilt für f (x) =
∇f (xj ) = Axj :
1 T
x Ax
2
mit A positiv definit, g j :=
T
1. xj+1 = xj −
gj gj j
g
g j T Ag j


2. f (xj+1 ) = 1 − jT j
2

g g
j

 f (x )
T
T −1 j
j
j
j
g Ag
g A g
Beweis.
1.
gj
mit αj aus exakter Linesearch
j
kg k
d
d
1
0=
f (xj + α∆xj )
=
(xj + α∆xj )T A(xj + α∆xj )
dα
dα 2
α=αj
α=αj
d
1
T
T
T
T
=
= ∆xj Axj + αj ∆xj A∆xj
α∆xj Axj + α2 ∆xj A∆xj dα
2
α=αj
xj+1 = xj − αj
T
T
∆xj Axj
kg j kg j g j
kg j k3
α =−
=
=
∆xj T A∆xj
g j T Ag j
g j T Ag j
j
T
x
j+1
gj gj j
=x − T
g
g j Ag j
j
20
KAPITEL 2. UNBESCHRÄNKTE OPTIMIERUNGSPROBLEME
2.
1 2
T
T
f (xj+1 ) = f (xj + αj ∆xj ) = f (xj ) + αj ∆xj Axj + αj ∆xj A∆xj
2
jT
jT
j
g
1 2 g Ag
= f (xj ) − αj j g j + αj
kg k
2
kg j k2
2
2
2
T
jT j
jT j
g
g
g
g
gj gj
1
1
= f (xj ) −
+
= f (xj ) −
T
T
2 g j Ag j
2 g j T Ag j
g j Ag j
T
T
T
g j A−1 g j = xj AA−1 Axj = xj Axj = 2f (xj )


2
2
jT j
jT j
g g
1 g g
2f (xj )
j
·
=

f (xj+1 ) = f (xj ) − 1− 
 f (x )
T
T
T
T
2 g j Ag j
g j A−1 g j
g j Ag j g j A−1 g j
g j = Axj ,
Lemma 2.15 (Ungleichung von Kantorovich) Sei A positiv definit. Dann gilt für alle x 6=
0:
2
xT x
4λmin λmax
≥
,
T
T
−1
(x Ax) (x A x)
(λmin + λmax )2
wobei λmin und λmax der kleinste und der größte Eigenwert von A sind.
Beweis. Seien 0 < λmin = λ1 ≤ . . . ≤ λn = λmax die Eigenwerte und u1 , u2 , . . . , un paarweise
orthonormale Eigenvektoren von A. Sei x ∈ Rn , x 6= 0 beliebig, aber fest. Da die u1 , u2 , . . . , un
eine Basis des Rn bilden, existiert eine Darstellung der Form
x=
n
X
βi ui mit
i=1
n
X
βi2 > 0.
i=1
Wegen der Orthonormalität der ui gilt
2
xT x
=n
F (x) := T
P
(x Ax) (xT A−1 x)
n
P
βi2
2
i=1
λi βi2
i=1
n
P
2
λ−1
i βi
.
i=1
1
βi2
n
.
wird daraus F (x) = n
Mit γi := P
n
P
P −1
2
βj
λi γi
λi γi
j=1
Sei λ̄ :=
n
X
λi γi . Wegen γi ≥ 0 und
i=1
i=1
n
X
i=1
γi = 1 ist λ̄ eine Konvexkombination der Eigenwerte
i=1
von A. Deshalb ist λ1 ≤ λ̄ ≤ λn .
Wir betrachten im R2 die Punkte Pi :=
λi , λ−1
i
und Q :=
λ̄,
n
X
i=1
!
λ−1
i γi
. Die Punkte
P1 , . . . , Pn liegen auf dem Graphen der Funktion r : R+ → R : λ 7→ r(λ) = λ−1 . Wegen r00 (λ) >
KAPITEL 2. UNBESCHRÄNKTE OPTIMIERUNGSPROBLEME
21
0 für alle λ > 0 ist r strikt konvex. TODO: Abbildung. Also liegen die Punkte P2 , . . . , Pn−1
n
X
nicht oberhalb der Geraden durch P1 und Pn . Der Punkt Q liegt wegen Q =
γi Pi und
γi ≥ 0,
n
X
i=1
γi = 1 in der konvexen Hülle von {P1 , . . . , Pn } und somit ebenfalls nicht oberhalb
i=1
der Geraden g durch P1 und Pn :
λ1 + λn − λ
,
g(λ) =
λ1 λn
g(λ̄) ≥
n
X
λ−1
i γi
i=1
Daraus folgt:
1
F (x) =
λ̄
n
P
λ−1
i γi
≥
λ1 λn
4λ1 λn
λ1 λn
≥ min
.
=
λ1 ≤λ≤λn λ (λ1 + λn − λ)
λ̄ λ1 + λn − λ̄
(λ1 + λn )2
i=1
Satz 2.16 (Konvergenz des Gradientenverfahrens) Das Gradientenverfahren konvergiert
linear
2
λmax − λmin
j+1
f (x ) ≤
f (xj )
λmax + λmix
und es gilt
1
λmax − λmin
j+1
j+1
2
kxj kA
kx kA := kA x k2 ≤
λmax + λmin
Falls also λmin λmax , ist die Konvergenzrate nahe bei 1 und das Verfahren konvergiert nur
sehr langsam:
λmax − λmin
2
κ(A) − 1
=1−
=
λmax + λmin
κ(A) + 1
κ(A) + 1
λmax
mit κ(A) =
Kondition von A.
λmin
Beweis. Aus Lemma 2.14 und der Ungleichung von Kantorovich folgt:
(λmax − λmin )2
4λmax λmin
j
j
j+1
f (x ) =
f (x ) ≤ 1 −
2
2 f (x )
(λmax + λmin )
(λmax + λmin )
1
1
1
1
1
f (xj+1 ) = (xj+1 )T Axj+1 = (xj+1 )T A 2 A 2 xj+1 = kxj+1 k2A
2
2
2
p
p
λ
−
λ
λmax − λmin
max
min
j+1
j+1
j
kx kA = 2f (x ) ≤
2f (x ) =
kxj kA .
λmax + λmin
λmax + λmin
Bemerkung 2.17 Aus ∇f (x∗ ) = 0 folgt
1
f (x) = f (x∗ ) + (x − x∗ )T ∇2 f (x∗ )(x − x∗ ) + O(kx − x∗ k3 )
2
Die Aussage von Satz 2.16 gilt allgemein, wenn man für λmax und λmin den größten und kleinsten Eigenwert von ∇2 f (x∗ ) nimmt.
22
KAPITEL 2. UNBESCHRÄNKTE OPTIMIERUNGSPROBLEME
Backtracking-Linesearch
Linesearch: Minimiere f entlang der Abstiegsrichtung ∆xj :
f (xj + αj ∆xj ) = min f (xj + α∆xj )
α>0
Exakte Linesearch: Bestimme das (globale) Minimum.
Inexakte Linesearch: Bestimme ein αj , das bestimmte Bedingungen erfüllt, z. B.
• Sufficient Decrease (d. h. Schritt nicht zu groß):
f (xj + αj ∆xj ) ≤ f (xj ) + c1 · αj · ∇f (xj )T ∆xj
für eine Konstante c1 ∈ (0; 1). D. h. der realisierte Abstieg ist proportional zu αj und zur
Richtungsableitung ∇f (xj )T ∆xj < 0.
• Schritt nicht zu klein.
Eine einfache Methode, diese beiden Bedingungen zu erfüllen, ist die Backtracking-Linesearch.
Algorithmus 2.18 (Backtracking-Linesearch)
• Wähle α > 0, ρ, c ∈ (0; 1).
• Wiederhole bis f (xj + α∆xj ) ≤ f (xj ) + c · α · ∇f (xj )T ∆xj :
α := ρ · α.
• αj := α.
Mehr zu Linesearch-Methoden in Kapitel 4.
Kapitel 3
Konjugierte Gradienten
(CG-)Methoden
Lineare CG-Verfahren
Sei A ∈ Rn×n positiv definit.
Betrachte das quadratische Minimierungsproblem
1
min xT Ax − bT x =: φ(x)
2
(3.1)
Der Gradient der Zielfunktion ist
∇φ(x) = Ax − b =: r(x) ( Residuum“)
”
Im Minimum x∗ gilt:
Ax∗ = b
(3.2)
Umgekehrt: Für x∗ mit Ax∗ = b ist φ(x∗ ) minimal, den die hinreichende Bedingung ist erfüllt.
Das quadratische Problem (3.1) und das lineare Gleichungssystem (3.2) sind also äquivalent.
Wir wollen Problem (3.1) bzw. (3.2) nun iterativ lösen ohne vollständiges Abspeichern oder
gar Zerlegen von A, vielmehr durch Auswertung nur von Matrix-Vektor-Produkten mit A.
Damit ist das Verfahren
1. geeignet für große lineare Gleichungssysteme bzw. QPs,
2. anwendbar auf nichtlineare Optimierungsprobleme.
Ansatz für das Verfahren: Verwende konjugierte“ Richtungen.
”
Definition 3.1 (Konjugiertheit) Eine Menge von Vektoren {p0 , p1 , . . . , pl } heißt konjugiert
bzgl. der symmetrischen, positiv definiten Matrix A, falls pi 6= 0, i = 0, . . . , l und
T
pi Apj = 0 für alle i 6= j.
23
(3.3)
24
KAPITEL 3. KONJUGIERTE GRADIENTEN (CG-)METHODEN
Bemerkung 3.2 Konjugierte Vektoren sind linear unabhängig.
Beweis. Übungsaufgabe.
Algorithmus 3.3 (Konjugierte-Richtungen-Methode)
1. Wähle einen Startwert x0 ∈ Rn und eine Menge {p0 , p1 , . . . , pn−1 } von konjugierten Richtungen.
2. Für j = 0, 1, . . . , n − 1 berechne
rj := Axj − b
und
xj+1 := xj + αj pj
mit
(3.4)
T
r j pj
.
α =− T
pj Apj
j
Bemerkung 3.4 α∗ = −
(3.5)
(Ax − b)T p
ist das Minimum der eindimensionalen Funktion φ(x +
pT Ap
αp).
Beweis. Aus
d
φ(x + αp)
= ∇φ(x + α∗ p)T p = (Ax − b)T p + α∗ pT Ap = 0
dα
∗
α=α
folgt
α∗ = −
(Ax − b)T p
.
pT Ap
Satz 3.5 Für beliebiges x0 ∈ Rn konvergiert die Folge {xj }, die durch die Konjugierte-Richtungen-Methode berechnet wird, in höchstens n Schritten gegen die Lösung des linearen Systems
(3.2).
Beweis.
• Weil die Richtungen pi , i = 0, . . . , n − 1 linear unabhängig sind, spannen sie den gesamten
Rn auf. Deshalb gibt es σi , i = 0, . . . , n − 1, so dass x∗ − x0 = σ0 p0 + . . . + σn−1 pn−1 .
T
Multipliziere mit pj A und wende die Konjugiertheitsbedingung (3.3) an:
T
T
pj A(x∗ − x0 ) = σj pj Apj .
KAPITEL 3. KONJUGIERTE GRADIENTEN (CG-)METHODEN
25
• Wenn xj durch Algorithmus 3.3 berechnet wurde, ist xj = x0 + α0 p0 + . . . + αj−1 pj−1 und
T
pj A(xj − x0 ) = 0. Daher folgt
T
T
T
T
T
T
pj A(x∗ − x0 ) = pj A(x∗ − xj ) + pj A(xj − x0 ) = pj (b − Axj ) = −pj rj = αj (pj Apj ).
Daraus folgt σj = αj und x∗ = x0 + α0 p0 + . . . + αn−1 pn−1 .
Beispiel 3.6 A Diagonalmatrix
Die Höhenlinien von φ sind Ellipsen, deren Achsen parallel zu den Koordinatenachsen sind:
TODO: Abbildung achsenparallele Ellipsen
Das Minimum kann durch eindimensionale Minimierungen entlang der Koordinatenachsen bestimmt werden.
Wenn A nicht diagonal ist, sind die Achsen nicht mehr entlang der Koordinatenachsen ausgerichtet.
TODO: Abbildung Ellipsen
Das Minimum wir dann nicht in n Schritten entlang der Koordinatenachsen erreicht.
Aber wir können das Problem transformieren: x̂ = S −1 x mit S = p0 p1 · · · pn−1 ∈ Rn×n
mit konjugierten Richtungen p0 , p1 , . . . , pn−1 bzgl. A:
1
φ̂(x̂) = φ(x) = φ(S x̂) = x̂T (S T AS)x̂ − (S T b)T x̂.
2
Wegen der Konjugiertheit ist S T AS diagonal und n Schritte entlang der x̂-Koordinatenachsen
erreichen x∗ . Diese korrespondieren zu n Schritten entlang pi im x-Raum.
Für diagonales A bestimmt jeder Schritt eine Komponente von x∗ . Dies kann man verallgemeinern. Dazu benutzen wir
Bemerkung 3.7
rj+1 = Axj+1 − b = A(xj + αj pj ) − b = Axj − b + αj Apj = rj + αj Apj .
(3.6)
Satz 3.8 (Expanding Subspace Minimierung) Sei x0 ∈ Rn ein beliebiger Startwert und
sei die Folge {xj } durch Algorithmus 3.3 berechnet. Dann gilt
T
rj pi = 0 für i = 0, . . . , j − 1
(3.7)
und xj minimiert φ(x) = 12 xT Ax − bT x über
x0 + Span{p0 , . . . , pj−1 }.
(3.8)
26
KAPITEL 3. KONJUGIERTE GRADIENTEN (CG-)METHODEN
Beweis.
• Zeige: x̃ minimiert φ über (3.8) genau dann, wenn
r(x̃)T pi = 0 für alle i = 0, . . . , j − 1.
Sei
σ := σ0 . . . σj−1
T
und h(σ) := φ(x0 + σ0 p0 + . . . + σj−1 pj−1 ).
h(σ) ist strikt konvex quadratisch, daher hat es ein eindeutiges Minimum σ ∗ mit
∂h ∗
(σ ) = 0,
∂σi
i = 0, . . . , j − 1.
Nach der Kettenregel ist
0=
∂h ∗
∗
(σ ) = ∇φ(x0 + σ0∗ p0 + . . . σj−1
pj−1 )T pi ,
∂σi
i = 0, . . . , j − 1.
∗
pj−1 erfüllt damit r(x̃)T pi = 0, i = 0, . . . , j − 1.
x̃ := x0 + σ0∗ p0 + . . . σj−1
• Zeige (3.7) durch Induktion über j.
Für j = 1 gilt x1 = x0 + α0 p0 minimiert φ entlang p0 . Daraus folgt
T
∇φ(x1 )T p0 = r1 p0 = 0.
T
Gelte rj−1 pi = 0 für i = 0, . . . , j − 2. Dann ist nach (3.6) rj = rj−1 + αj−1 Apj−1 und
nach Definition von αj−1
T
T
T
pj−1 rj = pj−1 rj−1 + αj−1 pj−1 Apj−1 = 0.
Für pi , i = 0, . . . , j − 2 gilt mit (3.6)
T
pi r j =
T
pi rj−1
| {z }
+αj−1
=0 (Induktionsannahme)
T
pi Apj−1
| {z }
=0
=0 (Konjugiertheit)
Induktionsschluß.
Bis jetzt: beliebige Wahl konjugierter Richtungen.
Ansätze zur Bestimmung konjugierter Richtungen:
• Eigenvektoren von A → Berechnung viel zu teuer.
• Modifizierte Gram-Schmidt-Orthogonalisierung → benötigt Abspeichern aller Richtungen, sehr aufwändig (siehe Übungsaufgabe).
KAPITEL 3. KONJUGIERTE GRADIENTEN (CG-)METHODEN
27
Konjugierte Gradienten Methode
Generiere eine Menge von konjugierten Vektoren, wobei pj nur aus pj−1 bestimmt wird. Die
alten Vektoren p0 , . . . , pj−2 müssen nicht gespeichert werden.
pj ist eine Linearkombination aus −rj = −∇φ(xj ) (Richtung des steilsten Abstiegs) und pj−1 :
pj = −rj + β j pj−1 .
β j ergibt sich aus der Konjugiertheit:
T
0=p
j−1 T
j
Ap = −p
j−1 T
j j−1 T
j
Ar + β p
Ap
j−1
rj Apj−1
⇒β =
.
pj−1 T Apj−1
j
Für p0 wählt man den steilsten Abstieg in x0 :
p0 = −r0 .
In jedem Schritt macht man außerdem eine eindimensionale Minimierung entlang pj zur Bestimmung der Schrittweite αj , siehe (3.5).
Algorithmus 3.9 (CG-Verfahren, Version 1)
1. Wähle Startwert x0 , setze r0 := Ax0 − b, p0 := −r0 , j := 0.
2. Solange rj 6= 0:
Schrittweite:
T
r j pj
α := − T
pj Apj
j
neue Iterierte:
xj+1 := xj + αj pj
neuer Gradient:
rj+1 := Axj+1 − b
Faktor für Linearkombination:
T
β j+1 :=
rj+1 Apj
pj T Apj
neue Richtung:
pj+1 := −rj+1 + β j+1 pj
j := j + 1
Definition 3.10 (Krylov-Unterraum)
κ(r0 ; j) := Span{r0 , Ar0 , . . . , Aj r0 }
ist der Krylov-Unterraum vom Grad j für r0 zu A. (Hier und im folgenden: A hoch j“.)
”
28
KAPITEL 3. KONJUGIERTE GRADIENTEN (CG-)METHODEN
Satz 3.11 Sei die j-te Iterierte des CG-Verfahrens (Algorithmus 3.9) nicht der Lösungspunkt
x∗ . Dann gilt:
T
1. rj ri = 0 für i = 0, 1, . . . , j − 1,
2. Span{r0 , r1 , . . . , rj } = Span{r0 , Ar0 , . . . , Aj r0 },
3. Span{p0 , p1 , . . . , pj } = Span{r0 , Ar0 , . . . , Aj r0 },
T
4. pj Api = 0 für i = 0, 1, . . . , j − 1.
Beweis. Induktion über j für 2., 3., 4.
Induktionsanfang: j = 0: 2., 3. trivial. 4. gilt für j = 1 nach Konstruktion von β 1 .
Induktionsannahme: Gelte 2., 3., 4. für j.
Induktionsschritt: Zeige: 2., 3., 4. gelten für j + 1.
• rj ∈ Span{r0 , Ar0 , . . . , Aj r0 } wegen 2. für j, pj ∈ Span{r0 , Ar0 , . . . , Aj r0 } wegen 3. für j.
Apj ∈ Span{Ar0 , . . . , Aj+1 r0 }. rj+1 := rj + αj Apj wegen (3.6). Daraus folgt
rj+1 ∈ Span{r0 , Ar0 , . . . , Aj+1 r0 } und Span{r0 , . . . , rj+1 } ⊆ Span{r0 , Ar0 , . . . , Aj+1 r0 }.
• Aj+1 r0 = A(Aj r0 ) ∈ Span{Ap0 , Ap1 , . . . , Apj } wegen 3. für j.
i+1
i
Api = r αi−r für i = 0, 1, . . . , j wegen (3.6). Also ist Aj+1 r0 ∈ Span{r0 , . . . , rj+1 }.
Daraus folgt Span{r0 , . . . , Aj+1 r0 } ⊆ Span{r0 , . . . , rj+1 }. Daraus folgt 2. für j + 1.
• Span{p0 , . . . , pj , pj+1 } = Span{p0 , . . . , pj , rj+1 } wegen Berechn. der neuen Richtung,
= Span{r0 , Ar0 , . . . , Aj r0 , rj+1 } wegen 3. für j,
= Span{r0 , r1 , . . . , rj , rj+1 } wegen 2. für j,
= Span{r0 , Ar0 , . . . , Aj+1 r0 } wegen 2. für j + 1.
Daraus folgt 3. für j + 1.
T
T
T
• pj+1 = −rj+1 + β j+1 pj . Daraus folgt pj+1 Api = −rj+1 Api + β j+1 pj Api = 0 für i = j
nach Konstruktion von β j+1 .
Für i ≤ j − 1: Nach Induktionsannahme sind p0 , . . . , pj konjugiert. Nach Satz 3.8 ist
T
rj+1 pi = 0 für i = 0, . . . , j. Aus 3. folgt für i = 0, . . . , j−1: Api ∈ A Span{r0 , Ar0 , . . . , Ai r0 }
= Span{Ar0 , A2 r0 , . . . , Ai+1 r0 } ⊂ Span{p0 , p1 , . . . , pi+1 }.
T
Daraus folgt rj+1 Api = 0 für i = 0, . . . , j − 1.
T
T
T
Also ist pj+1 Api = − rj+1 Api +β j+1 pj Api = 0 für i = 0, . . . , j − 1.
| {z }
| {z }
=0
=0 (∗)
((*) gilt Induktionsannahme)
Induktionsschluss.
Zu zeigen bleibt 1.
T
rj pi = 0 für i = 0, . . . , j − 1, j = 1, . . . , n − 1 nach (3.7).
ri = −pi + β i pi−1 . Daraus folgt ri ∈ Span{pi , pi−1 } für i = 1, . . . , j − 1.
T
Daraus folgt rj ri = 0 für i = 1, . . . , j − 1.
T
T
rj r0 = −rj p0 = 0 nach Definition von p0 .
Wichtig im Beweis: p0 = −r0 .
KAPITEL 3. KONJUGIERTE GRADIENTEN (CG-)METHODEN
29
Korollar 3.12 Für A ∈ Rn×n konvergiert die Folge {xj } der durch das CG-Verfahren berechneten Iterierten gegen die Lösung x∗ in höchstens n Schritten.
Beweis. Wegen der Konjugiertheit der pj kann Satz 3.5 angewendet werden.
Bemerkung 3.13 Zum Namen: Konjugierte Gradienten“ ist eine falsche Benennung. Die
”
Gradienten sind orthogonal. Die Suchrichtungen sind konjugiert.
Umformungen:
T
αj = −
T
T
r j pj
rj (rj − β j pj−1 )
rj rj
=
=
pj T Apj
pj T Apj
pj T Apj
und nach (3.6)
T
j
j
α Ap = r
j+1
−r
j
=⇒
β
j+1
T
rj+1 rj+1
rj+1 Apj
=
=
pj T Apj
rj T rj
Algorithmus 3.14 (CG-Verfahren, Standardform)
1. Wähle Startwert x0 , setze r0 := Ax0 − b, p0 := −r0 , j := 0.
2. Solange rj 6= 0:
T
rj rj
pj T Apj
:= xj + αj pj
αj := −
xj+1
rj+1 := rj + αj Apj
T
rj+1 rj+1
β
:=
rj T rj
j+1
p
:= −rj+1 + β j+1 pj
j := j + 1
j+1
Bemerkung 3.15
1. Die Vektoren x, r und p werden nur in der aktuellen Iteration benötigt und können sofort
mit neuen Werten überschrieben werden.
2. Der Hauptrechenaufwand besteht in dem Matrix-Vektor-Produkt Api . Zusätzlich sind die
T
T
inneren Produkte pj (Apj ) und rj+1 rj+1 und drei Additionen von Vektoren auszuwerten.
30
KAPITEL 3. KONJUGIERTE GRADIENTEN (CG-)METHODEN
Konvergenzrate des CG-Verfahrens
Nach Satz 3.11 ist
xj+1 = xo + α0 p0 + α1 p1 + . . . + αj pj
= x0 + γ 0 r0 + γ 1 Ar0 + . . . γ j Aj r0
= x0 + γ 0 I + γ 1 A + . . . γ j Aj r0
für Konstanten γ j . (Hier und im folgenden: A hoch j“.)
”
Sei Pj∗ (.) das Polynom vom Grad j mit Koeffizienten γ i , i = 0, . . . , j:
Pj∗ (A) = γ 0 I + γ 1 A + . . . + γ j Aj
Dann ist
xj+1 = x0 + Pj∗ (A)r0 .
(3.9)
Betrachte die Norm kzk2A := z T Az. Sei x∗ Minimum von φ(x) = 21 xT Ax − bT x. Dann ist
1
1
kx − x∗ k2A = (x − x∗ )T A(x − x∗ ) = φ(x) − φ(x∗ ),
2
2
denn 21 (x − x∗ )T A(x − x∗ ) − φ(x) + φ(x∗ ) = (x − x∗ )T (Ax∗ − b) = 0.
xj+1 minimiert φ und damit auch kx − x∗ k2A über die Menge
x0 + Span{p0 , . . . , pj } nach Satz 3.8
= x0 + Span{r0 , Ar0 , . . . , Aj r0 } nach Satz 3.11. Das heißt Pj∗ minimiert
min kx0 + Pj (A)r0 − x∗ kA
Pj
(3.10)
(Minimum über alle Polynome vom Grad j).
Es ist
r0 = Ax0 − b = Ax0 − Ax∗ = A(x0 − x∗ ),
xj+1 − x∗ = x0 + Pj∗ (A)r0 − x∗ = (I + Pj∗ (A)A)(x0 − x∗ ).
Seien 0 < λ1 ≤ λ2 ≤ . . . ≤ λn die Eigenwerte von A und v1 , . . . , vn zugehörige orthonormale
Eigenvektoren:
n
X
A=
λi vi viT .
i=1
In der Basis von Eigenvektoren sei
0
∗
x −x =
n
X
ξi vi .
i=1
Aus Avi = λi vi folgt Pj (A)vi = Pj (λi )vi , i = 1, . . . , n.
KAPITEL 3. KONJUGIERTE GRADIENTEN (CG-)METHODEN
Also:
x
j+1
∗
− x = (I +
Pj∗ (A)A)
n
X
n
X
ξi vi =
(1 + Pj∗ (λi )λi )ξi vi .
i=1
Sei
kzk2A
T
= z Az =
n
X
31
i=1
λi (viT z)2 . Dann ist
i=1
kxj+1 − x∗ k2A =
n
X
viT
λi
i=1
=
n
X
n
X
!2
(1 + Pj∗ (λk )λk )ξk vk
k=1
λi (1 + Pj∗ (λi )λi )2 ξi2
(viT vk = δik )
i=1
= min
Pj
n
X
λi (1 + Pj (λi )λi )2 ξi2
i=1
≤ min max (1 + Pj (λi )λi )2
n
X
Pj i=1,...,n
!
λk ξk2
k=1
0
= min max (1 + Pj (λi )λi ) kx − x∗ k2A
2
Pj i=1,...,n
Also ist
min max (1 + Pj (λi )λi )2 ≥ 0
(3.11)
Pj i=1,...,n
die Konvergenzrate des CG-Verfahrens.
Satz 3.16 Wenn A nur l verschiedene Eigenwerte hat, terminiert die CG-Iteration in höchstens
l Iterationen.
Beweis. Seien τ1 < τ2 < . . . < τl die verschiedenen Werte der Eigenwerte. Definiere ein Polynom
Ql (λ) :=
(−1)l
(λ − τ1 ) . . . (λ − τl ).
τ1 τ2 . . . τl
Es ist Ql (τi ) = 0, i = 1, . . . , l und Ql (0) = 1. Also ist Ql (λ) − 1 ein Polynom l-ten Grades mit
Nullstelle λ = 0.
Polynomdivision: P̄l−1 (λ) = (Ql (λ) − 1)/λ ist ein Polynom (l − 1)-ten Grades. Sei j := l − 1 in
(3.11):
0 ≤ min max (1 + Pl−1 (τi )τi )2 ≤ max (1 + P̄l−1 (τi )τi )2 = max Ql (τi )2 = 0.
Pl−1 i=1,...,l
i=1,...,l
i=1,...,l
Daraus folgt
kxl − x∗ k2A = 0 und xl = x∗ .
Satz 3.17 Wenn A die Eigenwerte λ1 ≤ . . . ≤ λn hat, gilt
2
λn−j − λ1
j+1
∗ 2
kx − x kA ≤
kx0 − x∗ k2A .
λn−j + λ1
32
KAPITEL 3. KONJUGIERTE GRADIENTEN (CG-)METHODEN
Beweis. Siehe Buch von Luenberger, Ye: Linear and Nonlinear Programming, Springer.
Beispiel 3.18 (Geclusterte Eigenwerte) A habe m große und n−m kleine Eigenwerte ≈ 1.
Nach m + 1 CG-Schritten:
kxm+1 − x∗ k2A ≈ εkx0 − x∗ k2A
(ε klein).
TODO: Abbildung
Die Kondition von A ist κ(A) = kAk2 kA−1 k2 =
λn
.
λ1
Satz 3.19
kxj − x∗ kA ≤ 2
!j
p
κ(A) − 1
p
kx0 − x∗ kA
κ(A) + 1
Präkonditionierung
Transformiere das lineare System, um die Eigenwertverteilung zu verändern. Das ergibt eine
Beschleunigung der CG-Methode.
x̂ = Cx, C regulär.
1
φ̂(x̂) = φ(x) = φ(C −1 x̂) = x̂T (C −T AC −1 )x̂ − (C −T b)T x̂
2
bzw.
(C −T AC −1 )x̂ = C −T b.
Die Konvergenzrate hängt von den Eigenwerten von C −T AC −1 ab. Wähle C so, dass die Kondition von C −T AC −1 viel kleiner als die von A ist.
Nichtlineare CG-Methoden
Wir betrachten nun
min f (x),
wobei f eine allgemeine konvexe oder eine allgemeine nichtlineare C 2 -Funktion ist.
KAPITEL 3. KONJUGIERTE GRADIENTEN (CG-)METHODEN
33
Die Fletcher-Reeves-Methode
Zwei einfache Änderungen gegenüber Algorithmus 3.14:
• Linesearch statt Minimierung einer quadratischen Funktion zur Bestimmung der Schrittweite αj ,
• ersetze das Residuum r durch den Gradienten ∇f der nichtlinearen Funktion.
Algorithmus 3.20 (Fletcher-Reeves)
1. Startwert x0 .
Berechne f 0 := f (x0 ) und ∇f 0 := ∇f (x0 ).
Setze p0 := −∇f 0 , j := 0.
2. Solange ∇f j 6= 0:
Berechne αj durch Linesearch entlang xj + αpj .
Berechne xj+1 := xj + αj pj .
Berechne f j+1 := f (xj+1 ) und ∇f j+1 := ∇f (xj+1 ).
T
Berechne βFj+1
R :=
∇f j+1 ∇f j+1
.
∇f j T ∇f j
j
Berechne pj+1 := −∇f j+1 + βFj+1
Rp .
Setze j := j + 1.
(f j wird in der Linesearch benötigt.)
Für f (x) = 12 xT Ax − bT x, A 0, und α exaktes Minimum ist dies das CG-Verfahren aus
Algorithmus 3.14.
In jeder Iteration wird nur die Zielfunktion und ihr Gradient ausgewertet. Keine MatrixOperationen sind nötig, es müssen nur wenige Vektoren gespeichert werden.
Das Verfahren ist also geeignet für große Probleme.
Bestimmung von αj
T
T
∇f j pj = −k∇f j k2 + βFj R ∇f j pj−1 < 0: Abstieg
T
Exakte Linesearch: αj−1 minimiert f entlang xj−1 + αpj−1 . Daraus folgt: ∇f j pj−1 = 0, daher
T
ist ∇f j pj < 0.
T
Inexakte Linesearch: vermeide ∇f j pj > 0 durch Bedingungen an die Schrittweite αj :
T
• f (xj + αj pj ) ≤ f (xj ) + c1 αj ∇f j pj
T
• |∇f (xj + αj pj )T pj | ≤ −c2 ∇f j pj
(Starke Wolfe-Bedingungen) mit 0 < c1 < c2 < 12 . Siehe Kapitel 4.
34
KAPITEL 3. KONJUGIERTE GRADIENTEN (CG-)METHODEN
Varianten
Andere Wahl von β j+1 .
Polak-Ribière (PR):
T
∇f j+1 (∇f j+1 − ∇f j )
k∇f j k2
Wenn f strikt konvex quadratisch und die Linesearch exakt ist, gilt
βPj+1
R =
j+1
βPj+1
R = βF R ,
T
weil dann ∇f j+1 ∇f j = 0 (Satz 3.11).
Für PR können die starken Wolfe-Bedingungen keinen Abstieg garantieren, jedoch für
PR+:
j+1
βPj+1
R+ = max{βP R , 0}.
In der Praxis ist PR robuster und effizienter als Fletcher-Reeves (FR).
Hestens-Stiefel (HS):
T
j+1
βHS
∇f j+1 (∇f j+1 − ∇f j )
=
(∇f j+1 − ∇f j )T pj
stimmt wegen Satz 3.8 für f strikt konvex quadratisch und exakte Linesearch mit βFj+1
R überein.
Sei
j
1
Z
∇2 f (xj + τ αj pj )dτ
Ḡ :=
0
die durchschnittliche Hessematrix von f auf der Strecke zwischen xj und xj+1 . Dann gilt:
∇f j+1 − ∇f j = αj Ḡj pj
und für
j+1 j
pj+1 = −∇f j+1 + βHS
p
und
T
j+1
βHS
∇f j+1 αj Ḡj pj
=
(αj Ḡj pj )T pj
folgt
T
T
T
j+1 j
pj+1 Ḡj pj = −∇f j+1 Ḡj pj + βHS
p Ḡj pj = 0,
d. h. die Konjugiertheit von pj+1 und pj bzgl. Ḡj .
Später sehen wir: globale Konvergenz für |β j | ≤ βFj R . Daher verwendet man die
Modifiziere PR-Methode (FR–PR):
Für j ≥ 2:
 j
j
j

−βF R , falls βP R < −βF R ,
β j = βPj R ,
falls |βPj R | ≤ βFj R ,

 j
βF R ,
falls βPj R > βFj R .
D. h. verwende das robustere und effizientere PR, solange es globale Konvergenz garantiert,
ansonsten FR.
Andere Varianten:
KAPITEL 3. KONJUGIERTE GRADIENTEN (CG-)METHODEN
k∇f j+1 k2
,
(∇f j+1 − ∇f j )T pj
T
j 2
∇f j+1
j
j kŷ k
= ŷ − 2p j T j
(ŷ ) p
(ŷ j )T pj
35
• β j+1 =
• β j+1
mit ŷ j = ∇f j+1 − ∇f j .
Neustart-Strategien
Neustart alle n Iterationen durch Wahl von β j = 0, d. h. Steilster-Abstieg-Schritt in dieser
Iteration.
Bemerkung: In der Nähe eines strikten lokalen Minimums ist f so gut wie strikt konvex quadratisch. Wenn die Iterierten konvergieren, kommen sie irgendwann in diese Nähe. Dann wird
irgendwann der Algorithmus neugestartet und verhält sich ab dann fast wie das lineare CGVerfahren. Wir können zwar nicht endliche Terminierung in n Schritten erwarten, aber wesentlichen Fortschritt hin zur Lösung.
In der Praxis sind solche Neustarts im allgemeinen nicht so wichtig, da wir hoffen, dass das
Verfahren für große n in wesentlich weniger als n Schritten konvergiert“.
”
Andere Neustart-Strategie: Neustart, wenn zwei folgende Gradienten weit entfernt von
orthogonal sind:
T
∇f j ∇f j−1
≥ ν ≈ 0.1.
k∇f j k2
Verhalten der Fletcher-Reeves-Methode
Annahme 3.21 Sei f zweimal stetig differenzierbar und die Niveaumenge N (x) = {x : f (x) ≤
f (x0 )} sei beschränkt.
Dann existiert eine Schrittlänge αj , die die starken Wolfe-Bedingungen erfüllt, siehe Lemma
4.2 in Kapitel 4.
Lemma 3.22 Sei Algorithmus 3.20 mit einer Liniensuche implementiert, die die starken WolfeBedingungen erfüllt mit 0 < c2 < 21 . Dann erzeugt das Verfahren Abstiegsrichtungen pj , die die
folgenden Ungleichungen erfüllen:
T
1
∇f j pj
2c2 − 1
−
≤
≤
j
2
1 − c2
k∇f k
1 − c2
für alle j = 0, 1, . . . .
2ξ − 1
ist auf 0; 12 monoton steigend und t(0) = −1, t
1−ξ
2c1 − 1
Daher ist −1 <
< 0 wegen c2 ∈ 0; 21 .
1 − c2
Beweis. Die Funktion t(ξ) :=
(3.12)
1
2
= 0.
Wenn (3.12) gezeigt ist, folgt sofort die Abstiegseigenschaft. Zeige (3.12) durch Induktion.
36
KAPITEL 3. KONJUGIERTE GRADIENTEN (CG-)METHODEN
0
(2*x-1)/(1-x)
-1/(1-x)
-0.5
-1
-1.5
-2
0
0.1
0.2
0.3
0.4
0.5
Abbildung 3.1: Rechte und linke Seite in (3.12) für x = c2 ∈ 0; 21 .
Für j = 0 ist der mittlere Term gleich −1 wegen p0 = −∇f 0 und beide Ungleichungen sind
erfüllt.
Gelte (3.12) für j ≥ 1. Aus der zweiten Wolfe-Bedingung
T
T
T
c2 ∇f j pj ≤ ∇f j+1 pj ≤ −c2 ∇f j pj ,
j+1
der Definition von βFj+1
in Algorithmus 3.20 folgt
R und der Berechnung von p
T
T
T
T
j+1
∇f j+1 pj
pj
∇f j+1 pj+1
∇f j pj
j+1 ∇f
≤
−1
+
=
−1
+
β
·
=
−1 + c2
FR
k∇f j k2
k∇f j k2
k∇f j+1 k2
k∇f j+1 k2
T
∇f j pj
≤ −1 − c2
.
k∇f j k2
Nach Induktionsannahme ist
T
−1
∇f j pj
≤
1 − c2
k∇f j k2
und daher
T
−1
−1
∇f j+1 pj+1
−1
2c2 − 1
= −1 + c2
≤
≤ −1 − c2
=
.
j+1
2
1 − c2
1 − c2
k∇f k
1 − c2
1 − c2
Daraus folgt (3.12) für j + 1.
Sei θj der Winkel zwischen pj und der Richtung des steilsten Abstiegs −∇f j :
T
−∇f j pj
cos θ =
.
k∇f j kkpj k
j
Wenn θj ≈ 90◦ , d. h. cos θj ≈ 0, ist pj eine schlechte“ Suchrichtung (siehe Kapitel 4).
”
Aus (3.12) folgt:
1 − 2c2 k∇f j k
1
k∇f j k
j
·
≤
cos
θ
≤
·
1 − c2
kpj k
1 − c2 kpj k
für alle j = 0, 1, . . . .
(3.13)
KAPITEL 3. KONJUGIERTE GRADIENTEN (CG-)METHODEN
37
Daraus folgt: cos θj ≈ 0 genau dann, wenn k∇f j k kpj k.
Wenn pj fast orthogonal zum Gradienten ist, ist es wahrscheinlich, dass man einen Schritt mit
sehr kleiner Schrittweite machen muss, um Abstieg zu bekommen. D. h.
xj+1 ≈ xj ,
∇f j+1 ≈ ∇f j ,
βFj+1
R ≈ 1.
Daraus folgt
pj+1 ≈ pj
und k∇f j+1 k kpj+1 k,
d. h. pj+1 ist wieder eine schlechte Suchrichtung.
TODO: Abbildung
Für Polak-Ribière gilt bei einer Suchrichtung mit cos θj ≈ 0 und ∇f j+1 ≈ ∇f j : βPj+1
R ≈ 0, d. h.
j+1
j+1
p
≈ −∇f
ist eine gute Suchrichtung. PR macht bei schlechten Suchrichtungen sozusagen
einen Neustart. Dies gilt auch für PR+, HS und FR–PR.
FR sollte nicht ohne Neustart-Strategie implementiert werden.
Globale Konvergenz
Annahme 3.23
1. Die Niveaumenge N (x0 ) := {x : f (x) ≤ f (x0 )} ist beschränkt.
2. Es gibt eine offene Umgebung U von N (x0 ), in der die Zielfunktion f Lipschitz-stetig
differenzierbar ist.
Bemerkung 3.24 Unter der Annahme 3.23 existiert eine Konstante γ̄ > 0, so dass k∇f (x)k ≤
γ̄ für alle x ∈ N (x0 ).
Wir benutzen nun Zoutendijks Theorem, siehe Satz 4.3: Für eine Linesearch der Form
xj+1 = xj + αj pj , wobei pj eine Abstiegsrichtung ist und αj die starken Wolfe-Bedingungen
erfüllt, gilt:
∞
X
cos2 θj · k∇f j k2 < ∞.
(3.14)
j=0
Werde das nichtlineare CG-Verfahren in den Iterationen j1 , j2 , . . . neugestartet. In diesen
Iterationen ist cos θj = 1. Dann folgt
X
k∇f j k2 =
j=j1 ,j2 ,...
X
cos2 θj · k∇f j k2 ≤
j=j1 ,j2 ,...
∞
X
cos2 θj · k∇f j k2 < ∞.
j=0
Wenn ji+1 − ji = n̄, ist die Folge {ji , i = 0, 1, . . .} unendlich und
lim k∇f ji k = 0 bzw.
i→∞
Für nicht-neugestartete Verfahren gilt:
lim inf k∇f j k = 0.
j→∞
(3.15)
38
KAPITEL 3. KONJUGIERTE GRADIENTEN (CG-)METHODEN
Satz 3.25 (Al-Baali) Es gelten die Annahmen 3.23, Algorithmus 3.20 sei mit einer Linesearch implementiert, die die starken Wolfe-Bedingungen erfüllt mit 0 < c1 < c2 < 12 . Dann
gilt
lim inf k∇f j k = 0.
j→∞
Beweis. Durch Widerspruch. Annahme: es existiert eine Konstante γ > 0, so dass k∇f j k ≥ γ
für alle j hinreichend groß.
Setze linke Ungleichung aus (3.13) in die Zoutendijk-Bedingung (3.14) ein:
2
∞
∞
X
X
k∇f j k4 1 − 2c2
2 j
j 2
∞>
.
cos θ · k∇f k ≥
kpj k2
1 − c2
j=0
j=0
(3.16)
Aus der zweiten Wolfe-Bedingung und der linken Ungleichung in (3.12) folgt
T
T
|∇f j pj−1 | ≤ −c2 ∇f j−1 pj−1 ≤
Daraus folgt mit pj = −∇f j + βFj R pj−1 und βFj R =
c2
k∇f j−1 k2 .
1 − c2
k∇f j k2
in Algorithmus 3.20:
k∇f j−1 k2
2
T
kpj k2 ≤ k∇f j k2 + 2βFj R |∇f j pj−1 | + βFj R kpj−1 k2
2
2c2 j
≤ k∇f j k2 +
βF R k∇f j−1 k2 + βFj R kpj−1 k2
1 − c2
2
1 + c2
=
k∇f j k2 + βFj R kpj−1 k2 .
1 − c2
Sei c3 :=
1 + c2
≥ 1.
1 − c2
kpj k2 ≤ c3 k∇f j k2 + βFj R
j 4
= c3 k∇f k
j
X
2 c3 k∇f j−1 k2 + βFj−1
R
2 c3 k∇f j−2 k2 + . . . + βF1 R
k∇f k k−2 ,
k=0
wobei
βFj R
2
βFj−1
R
2
2
. . . βFj−k
=
R
k∇f j k4
k∇f j−k−1 k4
und p0 = −∇f 0 . Es gilt also
kpj k2 ≤
und
∞
X
j=1
c3 · γ̄ 4
·j
γ2
∞
X1
1
≥
γ
4
kpj k2
j
j=1
für eine Konstante γ4 > 0. Andererseits folgt aus (3.16)
∞
X
j=0
∞
1
1 X k∇f j k4
≤
<∞
kpj k2
γ 4 j=0 kpj k2
.
2
kp0 k2 . . .
KAPITEL 3. KONJUGIERTE GRADIENTEN (CG-)METHODEN
39
Dieses globale Konvergenzresultat kann auf beliebige Wahl von β j mit |β j | ≤ βFj R ausgeweitet
werden.
Im allgemeinen: Falls Konstanten c4 , c5 > 0 existieren, so dass für j = 1, 2, . . .
cos θj ≥ c4
k∇f j k
kpj k
und
k∇f j k
≥ c5 > 0,
kpj k
folgt aus der Zoutendijk-Bedingung (3.14), dass
lim k∇f j k = 0.
j→∞
Dies gilt für Polak-Ribière unter den Annahmen, dass f strikt konvex ist und die verwendete
Linesearch exakt.
Für allgemeine nichtkonvexe Funktionen f kann man kein globales Konvergenzresultat für PR
beweisen.
Satz 3.26 Betrachte die Polak-Ribière-Methode mit einer exakten Linesearch. Es existiert eine
zweimal stetig differenzierbare Zielfunktion f : R3 → R und ein Startpunkt x0 ∈ R3 , so dass
die Folge {k∇f j k} von Null weg beschränkt ist.
Beweis: Powell, nicht konstruktiv.
Vergleich Gradientenverfahren – CG-Verfahren
Gradientenverfahren
CG-Verfahren
Mit qualifizierter Linesearch
Abstiegsverfahren
Abstiegsverfahren
Suchrichtungen
steilster Abstieg
Kombination aus
steilstem Abstieg und
vorheriger Richtung
sind orthogonal
sind konjugiert
Verhalten
Zick-Zack
Krylov-Raum-Eigenschaft
Für quadratische Probleme
im allgemeinen keine
endliche Terminierung
linear
2
λmax − λmin
λmax + λmin
endliche Terminierung
Lokale Konvergenz
Rate in k.kA -Norm
linear
2
λn−j − λ1
λn−j + λ1
gut für geclusterte
Eigenwerte, kann
vorkonditioniert werden
Globale Konvergenz
global konvergent
global konvergent
Implementierung
beide geeignet für große Probleme, da Matrix-frei
und Speicherung nur weniger Vektoren nötig
Kapitel 4
Linesearch
In jeder Iteration eines Abstiegsverfahrens:
• Suchrichtung: pj
• Schrittweite: αj > 0
• Iterationsvorschrift: xj+1 = xj + αj pj
Bemerkung 4.1 Falls ∇f (xj )T pj < 0, ist pj eine Abstiegsrichtung.
Wahl der Schrittweite
f sei Zielfunktion eines unbeschränkten Minimierungsproblems.
Ziele:
• Wesentliche Reduktion von f durch den Schritt.
• Effiziente Bestimmung der Schrittweite.
Betrachte
φ(α) := f (xj + αpj ),
α > 0.
Es ist φ0 (0) = ∇f (xj )T pj < 0.
Exakte Linesearch: Bestimme globales Minimum der Funktion φ(α), α > 0. Im allgemeinen
zu teuer.
Inexakte Linesearch: Untersuche eine Folge von Kandidaten für α, stoppe, wenn bestimmte
Bedingungen erfüllt sind.
Einfachste Bedingung: Abstieg in f :
f (xj + αj pj ) < f (xj ).
40
KAPITEL 4. LINESEARCH
41
Abstieg allein reicht nicht aus, um Konvergenz zu erzielen:
TODO: Abbildung
Minimum x∗ mit f (x∗ ) = −1. Iterierte xj mit f (xj ) = 5j , j = 1, 2, . . ..
Abstieg in jeder Iteration, aber Grenzwert = 0 und nicht −1.
Wir benötigen hinreichenden Abstieg (sufficient decrease).
Die Wolfe-Bedingungen
1. αj sollte Sufficient Decrease ergeben:
f (xj + αj pj ) ≤ f (xj ) + c1 α∇f (xj )T pj
(4.1)
für ein c1 ∈ (0; 1).
Der Abstieg in f ist mindestens proportional zu α und ∇f (xj )T pj .
(4.1) heißt auch Armijo-Bedingung.
TODO: Abbildung
In der Praxis: c1 klein, z. B. c1 = 10−4 .
Sufficient Decrease bewirkt, dass die Schritte nicht zu groß sind.
2. Für kleine Werte von α ist (4.1) erfüllt. Um Fortschritt zu erzielen, muss man auch
gewährleisten, dass die Schritte nicht zu klein sind.
Zweite Bedingung:
∇f (xj + αj pj )T pj ≥ c2 ∇f (xj )T pj
(4.2)
(Curvature Condition) mit c2 ∈ (c1 ; 1).
TODO: Abbildung
Bedeutet: φ0 (αj ) ist größer als c2 mal φ0 (0). Falls φ0 (α) 0, können wir α vergrößern, d.
h. akzeptieren es nicht. Falls aber φ0 (α) . 0 oder > 0, können wir keine Abnahme von f
mehr erwarten.
Typische Werte: c2 = 0.9 für Newton- oder Quasi-Newton-Verfahren, c2 = 0.1 für nichtlineares CG-Verfahren.
(4.1) und (4.2) zusammen heißen Wolfe-Bedingungen.
Starke Wolfe-Bedingungen: (4.1) und
|∇f (xj + αj pj )T pj | ≤ c2 |∇f (xj )T pj |
(4.3)
mit 0 < c1 < c2 < 1.
Wir erlauben nicht, dass die Ableitung zu positiv wird. Daher schließen wir Punkte aus, die
weit von stationären Punkten entfernt sind.
42
KAPITEL 4. LINESEARCH
Lemma 4.2 Sei f : Rn → R steig differenzierbar. Sei pj eine Abstiegsrichtung in xj . Sei f
entlang {xj +αpj , α > 0} nach unten beschränkt. Falls 0 < c1 < c2 < 1, existieren Intervalle von
Schrittweiten α, die die Wolfe-Bedingungen (4.1) und (4.2) bzw. die starken Wolfe-Bedingungen
(4.1) und (4.3) erfüllen.
Beweis. φ(α) := f (xj + αpj ) ist nach unten beschränkt für alle α > 0. Wegen 0 < c1 < c2 < 1
ist die Halbgerade l(α) := f (xj ) + αc1 ∇f (xj )T pj , α > 0 nach unten unbeschränkt und muss
daher mindestens einmal den Graphen von φ schneiden. Sei α0 der kleinste Schnittpunkt für
α > 0, d. h.
f (xj + α0 pj ) = f (xj ) + α0 c1 ∇f (xj )T pj .
Für alle 0 < α < α0 gilt die Bedingung (4.1). Aus dem Mittelwertsatz folgt, dass ein α00 ∈ (0; α0 )
existiert mit
f (xj + α0 pj ) − f (xj ) = α0 ∇f (xj + α00 pj )T pj .
Also ist
∇f (xj + α00 pj )T pj = c1 ∇f (xj )T pj > c2 ∇f (xj )T pj
wegen c1 < c2 . Daher erfüllt α00 auch (4.2). Wegen stetiger Differenzierbarkeit von f gilt das
auch für ein Intervall. Weil 0 > ∇f (xj + α00 pj )T pj , gilt auch die starke Wolfe-Bedingung (4.3).
TODO: Abbildung
Goldstein-Bedingung: bewirkt ebenfalls, dass die Schrittweite nicht zu klein ist:
f (xj ) + (1 − c)αj ∇f (xj )T pj ≤ f (xj + αj pj ) ≤ f (xj ) + cαj ∇f (xj )T pj
(4.4)
mit 0 < c < 12 .
Konvergenz von Linesearch-Methoden
Wir betrachten den Winkel θj zwischen der Suchrichtung pj und der Richtung des steilsten
Abstiegs −∇f (xj ):
−∇f (xj )T pj
.
cos θj =
k∇f (xj )kkpj k
Satz 4.3 (Zoutendijk) Betrachte eine Iteration der Form xj+1 = xj + αj pj , wobei pj eine
Abstiegsrichtung ist und αj die Wolfe-Bedingungen (4.1) und (4.2) erfüllt. Sei f in Rn nach
unten beschränkt und stetig differenzierbar in einer offenen Menge U , die die Niveaumenge
N (x0 ) := {x : f (x) ≤ f (x0 )} enthält, x0 sei der Startpunkt der Iteration. Sei ∇f Lipschitzstetig in U , d. h. ∃L > 0, so dass
k∇f (x) − ∇f (x̄)k ≤ Lkx − x̄k
∀x, x̄ ∈ U.
Dann gilt:
∞
X
j=0
cos2 (θj ) · k∇f (xj )k2 < ∞.
(4.5)
KAPITEL 4. LINESEARCH
43
Beweis. Aus (4.2) folgt
T
∇f (xj+1 ) − ∇f (xj ) pj ≥ (c2 − 1)∇f (xj )T pj
Aus der Lipschitz-Bedingung für ∇f folgt
1
∇f (xj+1 ) − ∇f (xj ) − Lαj pj 2
2
2
T
1
1
2
= ∇f (xj+1 ) − ∇f (xj ) − Lαj ∇f (xj+1 ) − ∇f (xj ) pj + L2 αj kpj k2
2
2
T
1
1
2
2
≤ L2 αj kpj k2 − Lαj ∇f (xj+1 ) − ∇f (xj ) pj + L2 αj kpj k2 .
2
2
0≤
Daraus folgt
T
∇f (xj+1 ) − ∇f (xj ) pj ≤ Lαj kpj k2 .
Also
αj ≥
c2 − 1 ∇f (xj )T pj
·
.
L
kpj k2
Mit (4.1) erhalten wir
1 − c2 ∇f (xj )T pj
j+1
j
·
f (x ) ≤ f (x ) − c1
L
kpj k2
mit c := c1
2
= f (xj ) − c cos2 (θj )k∇f (xj )k2
1 − c2
. Daraus folgt
L
f (xj+1 ) ≤ f (x0 ) − c
j
X
cos2 (θk )k∇f (xk )k2 .
k=0
Weil f nach unten beschränkt ist, ist f (x0 ) − f (xj+1 ) kleiner als eine positive Konstante für
alle j. Deshalb gilt auch
∞
X
cos2 (θk )k∇f (xk )k2 < ∞.
k=0
Für die Goldstein-Bedingungen (4.4) können ähnliche Ergebnisse gezeigt werden.
Die Voraussetzungen von Satz 4.3 sind nicht zu restriktiv:
• Falls f nicht nach unten beschränkt ist, ist das Minimierungsproblem nicht wohldefiniert.
• Die Glattheitsvoraussetzung wird auch für viele lokale Konvergenzbeweise benötigt.
Linesearch-Algorithmen
Sei
φ(α) := f (xj + αpj ).
44
KAPITEL 4. LINESEARCH
Bestimme αj , das eine der Bedingungen von oben erfüllt.
Für quadratische Probleme f (x) = 21 xT Ax − bT x ist
αj = −
(Ax − b)T pj
∇f (xj )T pj
=
−
.
pj T Apj
pj T Apj
Sei αj die Schrittweite im j-ten Schritt. Im folgenden: α0 , α1 , . . . Testschrittweiten zur Bestimmung von αj .
Backtracking-Linesearch
Für i = 0, 1, . . .: αi = α0 · ρi ( ρ hoch i“) mit 0 < ρ < 1, solange bis die Sufficient-Decrease”
Bedingung erfüllt ist:
φ(αi ) ≤ φ(0) + c1 αi φ0 (0).
(4.6)
Durch Konstruktion des Verfahrens gilt außerdem eine Schritt-nicht-zu-klein-Bedingung.
Erweiterung davon:
Interpolationsmethode
Algorithmus 4.4 (Interpolationsmethode)
• Sei ein Startwert α0 gegeben. Wenn α0 (4.6) erfüllt, stoppe. Ansonsten enthält (0; α0 )
Schrittweiten, die (4.6) erfüllen.
• Berechne die quadratische Interpolation φq von φ(0), φ0 (0) und φ(α0 ):
φ(α0 ) − φ(0) − α0 φ0 (0)
α2 + φ0 (0)α + φ(0)
φq (α) =
α02
( α hoch 2“), diese hat das Minimum
”
φ0 (0)α02
α1 = −
2 (φ(α0 ) − φ(0) − α0 φ0 (0))
Falls α1 (4.6) erfüllt, stoppe. Sonst:
• Berechne die kubische Interpolation φc von φ(0), φ0 (0), φ(α0 ) und φ(α1 ):
φc (α) = aα3 + bα2 + αφ0 (0) + φ(0)
( α hoch 3“ und α hoch 2“) mit
”
”
2
1
a
α0 −α12
φ(α1 ) − φ(0) − φ0 (0)α1
= 2 2
.
b
φ(α0 ) − φ(0) − φ0 (0)α0
α0 α1 (α1 − α0 ) −α03 α13
Das Minimum α2 von φc liegt im Intervall [0; α1 ]:
p
−b + b2 − 3aφ0 (0)
α2 =
.
3a
Falls α2 (4.6) erfüllt, stoppe.
KAPITEL 4. LINESEARCH
45
• Ansonsten iteriere mit kubischen Interpolationen von φ(0), φ0 (0) und den letzten zwei
Werten φ(αi−2 ) und φ(αi−1 ), bis ein αi gefunden ist, das (4.6) erfüllt.
Falls dabei αi zu nahe bei αi−1 ist oder viel kleiner als αi−1 , setze αi := αi−1 /2. ( Safeguard”
Strategie“, um zu erreichen, dass der Algorithmus vorankommt und die αi nicht zu klein
werden.)

Download Report