Nichtlineare Optimierung

Skript zur Vorlesung im Frühjahrsemester 2011

Helmut Harbrecht Stand: 25. Mai 2011

3

Vorwort Diese Mitschrift kann und soll nicht ganz den Wortlaut der Vorlesung wiedergeben. Sie soll als Lernhilfe dienen und das Nacharbeiten des Inhalts der Vorlesung erleichtern. Neben den unten genannten Büchern, diente mir auch das Vorlesungsskript Numerik nichtlinearer Optimierung von Gerhard Starke (Uni Hannover) als fruchtbare Quelle. Literatur zur Vorlesung: — C. Geiger und C. Kanzow: Numerische Verfahren zur Lösung unrestringierter Optimierungsaufgaben, Springer-Verlag — C. Geiger und C. Kanzow: Theorie und Numerik restringierter Optimierungsaufgaben, Springer-Verlag — F. Jarre und J. Stoer: Optimierung, Springer-Verlag

4

Inhaltsverzeichnis

Inhaltsverzeichnis 1 Grundlagen 1.1 Einführung . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 1.2 Optimalitätskriterien . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 1.3 Konvexität . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .

5 5 5 6

2 Gradientenverfahren 10 2.1 Idealisierte Variante . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 10 2.2 Praktische Variante . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 12 3 Newton-Verfahren 3.1 Lokales Newton-Verfahren . . . 3.2 Globalisiertes Newton-Verfahren 3.3 Inexaktes Newton-Verfahren . . 3.4 Trust-Region-Verfahren . . . . . 4 Quasi-Newton-Verfahren

. . . .

. . . .

. . . .

. . . .

. . . .

. . . .

. . . .

. . . .

. . . .

. . . .

. . . .

. . . .

. . . .

. . . .

. . . .

. . . .

. . . .

. . . .

. . . .

. . . .

. . . .

. . . .

. . . .

. . . .

15 15 20 22 24 31

5 Verfahren der konjugierten Gradienten 37 5.1 CG-Verfahren für lineare Gleichungssysteme . . . . . . . . . . . . . . . . . 37 5.2 Nichtlineares CG-Verfahren . . . . . . . . . . . . . . . . . . . . . . . . . . 40 5.3 Modifiziertes Verfahren von Polak und Ribière . . . . . . . . . . . . . . . . 43 6 Nichtlineare Ausgleichsprobleme 48 6.1 Gauß-Newton-Verfahren . . . . . . . . . . . . . . . . . . . . . . . . . . . . 48 6.2 Levenberg-Marquardt-Verfahren . . . . . . . . . . . . . . . . . . . . . . . . 51 7 Optimierungsprobleme mit Nebenbedingungen 56 7.1 Optimalitätsbedingungen erster Ordnung . . . . . . . . . . . . . . . . . . . 56 7.2 Optimalitätsbedingungen zweiter Ordnung . . . . . . . . . . . . . . . . . . 62 8 Projiziertes Gradientenverfahren 65 8.1 Konvergenzeigenschaften . . . . . . . . . . . . . . . . . . . . . . . . . . . . 65 8.2 Affine Nebenbedingungen . . . . . . . . . . . . . . . . . . . . . . . . . . . 73 9 SQP-Verfahren 76 9.1 Quadratische Minimierungsprobleme mit affinen Nebenbedingungen . . . . 76 9.2 Bestimmung aktiver Nebenbedingungen . . . . . . . . . . . . . . . . . . . . 78

5

1. Grundlagen 1.1 Einführung Optimierungsaufgaben treten in zahlreichen Anwendungsproblemen in den Natur- und Ingenieurwissenschaften, der Wirtschaft oder der Industrie auf. Beispielsweise versuchen Transportunternehmen, die Fahrt- oder Flugkosten zu minimieren und dabei sicherzustellen, dass alle Aufträge ausgeführt werden. Ebenso führt die numerische Simulation vieler physikalischer Vorgänge in den Naturwissenschaften auf Optimierungsprobleme, da das zugrundeliegende mathematische Modell oftmals auf dem Prinzip der Energieminimierung beruht. Unter einem endlichdimensionalen Minimierungsproblem wird die folgende Aufgabe verstanden: Gegeben sei eine Zielfunktion f : Rn → R. Gesucht ist ein Punkt x⋆ ∈ R, so dass f (x⋆ ) ≤ f (x) für alle x ∈ Rn . Dabei ist es ausreichend, sich nur mit Minimierungsproblemen zu beschäftigen, da ein Maximierungsproblem für f immer einem Minimierungsproblem für −f entspricht. Definition 1.1 Es sei f : Rn → R. Ein Punkt x⋆ ∈ Rn heißt globales Minimum, falls gilt f (x⋆ ) ≤ f (x) für alle x ∈ Rn . Das Minimum ist ein lokales Minimum, wenn es eine Umgebung U ⊂ Rn von x⋆ gibt, so dass f (x⋆ ) ≤ f (x) für alle x ∈ U .

Das Minimim heißt strikt, wenn im Fall x 6= x⋆ jeweils die strenge Ungleichung f (x⋆ ) < f (x) gilt.

In der Regel ist es mit vertretbarem Aufwand nur möglich, ein lokales Minimum von f in einer Umgebung eines Startwertes x0 zu bestimmen.

1.2 Optimalitätskriterien Um ein lokales Minimum numerisch zu finden, versucht man iterativ die Gleichung ∇f (x) = 0 zu lösen.

6

Kapitel 1. Grundlagen

Definition 1.2 Seien U ⊂ Rn eine offene Menge und f : U → R eine stetig differenzierbare Funktion. Ein Punkt x⋆ ∈ U heißt stationärer Punkt, falls gilt ∇f (x⋆ ) = 0. Wir wiederholen einige bekannte Eigenschaften lokaler Minima aus der Analysis: Satz 1.3 (notwendige Bedingung 1. Ordnung) Ist x⋆ ein lokales Minimum von f und ist f stetig differenzierbar in einer Umgebung von x⋆ , dann gilt ∇f (x⋆ ) = 0. Der Punkt x⋆ ist also ein stationärer Punkt.

Satz 1.4 (notwendige Bedingung 2. Ordnung) Ist x⋆ ein lokales Minimum von f und ist die Hesse-Matrix ∇2 f stetig in einer Umgebung von x⋆ , dann gilt ∇f (x⋆ ) = 0 und ∇2 f (x⋆ ) ist eine positiv semidefinite Matrix.

Satz 1.5 (hinreichende Bedingung 2. Ordnung) Die Hesse-Matrix ∇2 f sei stetig in einer Umgebung von x⋆ mit ∇f (x⋆ ) = 0. Ist ∇2 f (x⋆ ) eine positiv definite Matrix, dann ist x⋆ ein striktes lokales Minimum.

1.3 Konvexität Wir wenden uns einem wichtigen und in der Praxis oft auftretenden Spezialfall zu, bei dem wir mit einem lokalen zugleich ein globales Minimum gefunden haben. Dazu sei angemerkt, dass eine Menge D ⊂ Rn konvex ist, falls aus x, y ∈ D auch λx + (1 − λ)y ∈ D folgt für alle λ ∈ (0, 1). Definition 1.6 Es sei D ⊂ Rn eine konvexe Menge. Die Funktion f : D → R heißt konvex auf D, wenn für alle λ ∈ (0, 1) und alle x, y ∈ D gilt
f λx + (1 − λ)y ≤ λf (x) + (1 − λ)f (y).

Gilt für x 6= y sogar stets die strikte Ungleichung, dann heißt die Funktion strikt konvex. Gibt es ein µ > 0, so dass
f λx + (1 − λ)y + µλ(1 − λ)kx − yk22 ≤ λf (x) + (1 − λ)f (y) für alle λ ∈ (0, 1) und alle x, y ∈ D, dann heißt die Funktion f gleichmäßig konvex.

Beispiele 1.7 1. Die Gerade f (x) := x ist konvex auf R, aber nicht strikt konvex.

1.3. Konvexität

7

2. Die Exponentialfunktion f (x) := exp(x) ist strikt konvex auf R, dort aber nicht gleichmäßig konvex. 3. Die Parabel f (x) := x2 ist gleichmäßig konvex auf R. Hingegen ist die sehr ähnlich aussehende Funktion f (x) := x4 zwar strikt konvex auf R, aber nicht gleichmäßig konvex. △ 6

f f (y) λf (x) + (1 − λ)f (y) f (x)

f λx + (1 − λ)y

•




•

•

• -

x

λx + (1 − λ)y

y

Bei einer eindimensionalen konvexen Funktion liegt die Verbindungslinie zweier Punkte oberhalb des Graphen. Bemerkung Sei f : Rn → R eine quadratische Funktion, das heißt 1 f (x) = xT Ax + bT x + c 2 mit einer symmetrischen Matrix A ∈ Rn×n , b ∈ Rn und c ∈ R. Die Funktion f ist genau dann konvex, wenn A positiv semidefinit ist. Ist die Matrix A sogar positiv definit, so ist f sogar gleichmäßig konvex. △ Satz 1.8 Seien D ⊂ Rn eine offene und konvexe Menge und f : D → R stetig differenzierbar. Die Funktion f ist genau dann konvex auf D, wenn für alle x, y ∈ D gilt f (x) − f (y) ≥ ∇f (y)T (x − y).

(1.1)

Ist diese Ungleichung strikt für alle x 6= y, dann ist f sogar strikt konvex. Die Funktion f ist genau dann gleichmäßig konvex, wenn ein µ > 0 existiert, so dass f (x) − f (y) ≥ ∇f (y)T (x − y) + µkx − yk22

(1.2)

für alle x, y ∈ D. Beweis. Es gelte zuächst (1.2). Für x, y ∈ D und beliebiges λ ∈ (0, 1) ergibt sich dann mit z := λx + (1 − λ)y f (x) − f (z) ≥ ∇f (z)T (x − z) + µkx − zk22 ,

f (y) − f (z) ≥ ∇f (z)T (y − z) + µky − zk22 .

8

Kapitel 1. Grundlagen

Multipliziert man diese Gleichungen mit λ beziehungsweise 1 − λ und addiert sie anschließend, dann folgt
λf (x) + (1 − λ)f (y) − f λx + (1 − λ)y ≥ 2µλ(1 − λ)kx − yk22 ,

das heißt, f ist gleichmäßig konvex. Sei f nun als gleichmäßig konvex auf D vorausgesetzt. Für alle x, y ∈ D und λ ∈ (0, 1) gilt dann mit einem µ > 0

f y + λ(x − y) = f λx + (1 − λ)y ≤ λf (x) + (1 − λ)f (y) − µλ(1 − λ)kx − yk22 und daher
f y + λ(x − y) − f (y) ≤ f (x) − f (y) − µ(1 − λ)kx − yk22 . λ Aufgrund der stetigen Differenzierbarkeit von f folgt somit für λ → 0+ ∇f (y)T (x − y) ≤ f (x) − f (y) − µkx − yk22 ,

dies bedeutet, es gilt (1.2). Da der soeben geführte Beweis auch im Fall µ = 0 seine Gültigkeit behält, folgt die Äquivalenz von (1.1) zur Konvexität von f . Es verbleibt zu zeigen, dass die strikte Konvexität von f die strikte Ungleichung f (x) − f (y) > ∇f (y)T (x − y) für alle x, y ∈ D mit x 6= y impliziert. Als strikt konvexe Funktion ist f insbesondere konvex, das heißt, es gilt (1.1). Für   1 1 1 y z := (x + y) = x + 1 − 2 2 2 ergibt sich daher  ∇f (y)T (x − y) = 2∇f (y)T (z − y) ≤ 2 f (z) − f (y) .

(1.3)

Ist x 6= y, dann folgt wegen der strikten Konvexität

1 1 f (z) < f (x) + f (y). 2 2 Dies eingesetzt in (1.3) liefert die Behauptung ∇f (y)T (x − y) < f (x) − f (y).

Satz 1.9 Die Funktion f : D ⊂ Rn → R sei konvex. Dann ist jedes lokale Minimum x⋆ auch ein globales Minimum von f . Ist f zusätzlich differenzierbar, so ist jeder stationäre Punkt x⋆ ein globales Minimum.

1.3. Konvexität

9

Beweis. Angenommen, der Punkt x⋆ ist ein lokales, aber kein globales Minimum. Dann gibt es einen Punkt y⋆ ∈ Rn mit f (y⋆ ) < f (x⋆ ). Für alle x = λx⋆ + (1 − λ)y⋆ ,

λ ∈ (0, 1)

(1.4)

gilt aufgrund der Konvexität f (x) ≤ λf (x⋆ ) + (1 − λ)f (y⋆ ) < f (x⋆ ). Da in jeder Umgebung von x⋆ Punkte der Form (1.4) liegen, steht dies im Widerspruch zur Annahme, dass x⋆ ein lokales Minimum ist. Folglich ist jedes lokale Minimum auch ein globales Minimum. Wir zeigen nun die zweite Aussage. Dazu sei f differenzierbar vorausgesetzt und x⋆ ein stationärer Punkt. Wir führen den Beweis wieder per Widerspruch und nehmen an, dass x⋆ kein lokales Minimum ist. Dann können wir ein y⋆ wie oben wählen und erhalten aufgrund der Konvexität
d ∇f (x⋆ )T (y⋆ − x⋆ ) = f x⋆ + t(y⋆ − x⋆ ) ds t=0
⋆ ⋆ ⋆ f x + t(y − x ) − f (x⋆ ) = lim t→0+ t ⋆ (1 − t)f x ) + tf (y⋆ ) − f (x⋆ ) ≤ lim t→0+ t ⋆ ⋆ = f (y ) − f (x ) < 0. Deshalb ist ∇f (x⋆ ) 6= 0 und folglich ist x⋆ kein stationärer Punkt.

10

Kapitel 2. Gradientenverfahren

2. Gradientenverfahren 2.1 Idealisierte Variante Im folgenden setzen wir stets voraus, dass f : D ⊂ Rn → R stetig differenzierbar ist. Zunächst wollen wir das Gradientenverfahren betrachten, das auch Verfahren des steilsten Abstiegs genannt wird. Die Idee dabei ist, die Iterierte xk in Richtung des Antigradienten −∇f (xk ) aufzudatieren xk+1 := xk − αk ∇f (xk ), so dass f (xk+1 ) < f (xk ) ist. Dass dies im Fall ∇f (xk ) 6= 0 immer möglich ist, zeigt uns das nächste Lemma.

Lemma 2.1 Vorausgesetzt es ist ∇f (xk ) 6= 0, dann gibt es ein δ > 0, so dass die Funktion
ϕ(α) = f xk − α∇f (xk ) für alle 0 ≤ α ≤ δ streng monoton fällt. Insbesondere gilt
f xk − δ∇f (xk ) < ϕ(0) = f (xk ).

Beweis. Die Funktion ϕ ist stetig differenzierbar und es gilt
d ′ f xk − α∇f (xk ) = −k∇f (xk )k22 < 0. ϕ (0) = dα α=0

Aus Stetigkeitsgründen folgt die Existenz eines δ > 0 mit ϕ′ (α) < 0 für alle 0 ≤ α ≤ δ und damit die Behauptung. Algorithmus 2.2 (idealisiertes Gradientenverfahren) input: Funktion f : Rn → R und Startnäherung x0 ∈ Rn output: Folge von Iterierten {xk }k∈N ➀ Initialisierung: setze k := 0 ➁ berechne den Antigradienten dk = −∇f (xk ) ➂ löse αk = argmin f (xk + αdk ) α∈R

➃ setze xk+1 := xk + αk dk

(2.1)

2.1. Idealisierte Variante

11

➄ erhöhe k := k + 1 und gehe nach ➁

Satz 2.3 Die Funktion f : D ⊂ Rn → R sei gleichmäßig konvex. Weiter sei f differenzierbar mit Lipschitz-stetigem Gradienten: k∇f (x) − ∇f (y)k2 ≤ Lkx − yk2

für alle x, y ∈ D.

Dann konvergiert das idealisierte Gradientenverfahren für beliebige Startnäherungen x0 ∈ D gegen das eindeutige globale Minimum x⋆ und es gilt   µ2  ⋆ f (xk+1 ) − f (x ) ≤ 1 − 2 f (xk ) − f (x⋆ ) , k = 1, 2, . . . . L

Beweis. Aufgrund von (2.1) gilt für ein beliebiges β > 0, dass
f (xk+1 ) ≤ f xk − β∇f (xk ) Z 1
T = f (xk ) − β ∇f xk − tβ∇f (xk ) ∇f (xk ) dt 0 Z 1 
T 2 = f (xk ) − βk∇f (xk )k2 − β ∇f xk − tβ∇f (xk ) − ∇f (xk ) ∇f (xk ) dt 0 Z 1




∇f xk − tβ∇f (xk ) − ∇f (xk ) k∇f (xk )k2 dt ≤ f (xk ) − βk∇f (xk )k22 + β {z }2 0 |   2L k∇f (xk )k22 . ≤ f (xk ) − β − β 2

≤tβLk∇f (xk )k2

Für β = 1/L schließen wir also f (xk+1 ) − f (x⋆ ) ≤ f (xk ) − f (x⋆ ) −

1 k∇f (xk ) − ∇f (x⋆ ) k22 . | {z } 2L

(2.2)

=0

Die gleichmäßige Konvexität impliziert

 T µkxk − x⋆ k22 ≤ ∇f (xk ) − ∇f (x⋆ ) (xk − x⋆ ) ≤ k∇f (xk ) − ∇f (x⋆ )k2 kxk − x⋆ k2 . Dies eingesetzt in (2.2) liefert f (xk+1 ) − f (x⋆ ) ≤ f (xk ) − f (x⋆ ) −

µ2 kxk − x⋆ k22 . 2L

12

Kapitel 2. Gradientenverfahren

Weiterhin erhalten wir aus der Lipschitz-Bedingung Z 1
T ⋆ f (xk ) − f (x ) = ∇f x⋆ + t(xk − x⋆ ) (xk − x⋆ ) dt 0 Z 1 
T = ∇f x⋆ + t(xk − x⋆ ) − ∇f (x⋆ ) (xk − x⋆ ) dt Z0 1




∇f x⋆ + t(xk − x⋆ ) − ∇f (x⋆ ) kxk − x⋆ k2 dt ≤ {z }2 0 | ≤tLkxk −x⋆ k2

L kxk − x⋆ k22 , 2 womit sich die gewünschte Abschätzung ergibt. =

Bemerkung Im obigen Beweis haben wir L kxk − x⋆ k22 2 gezeigt. Umgekehrt folgt aus der gleichmäßigen Konvexität aber auch Z 1
T ⋆ f (xk ) − f (x ) ≥ ∇f x⋆ + t(xk − x⋆ ) (xk − x⋆ ) dt 0 Z 1 ≥ tµkxk − x⋆ k22 dt f (xk ) − f (x⋆ ) ≤

0

µ = kxk − x⋆ k22 . 2 Daher folgt aus Satz 2.3, dass der Fehler des idealisierten Gradientenverfahrens 2 f (xk ) − f (x⋆ ) kxk − x⋆ k22 ≤ µ  k 2 µ2  ≤ 1− 2 f (x0 ) − f (x⋆ ) µ L  k L µ2 ≤ 1 − 2 kx0 − x⋆ k22 . µ L Wir erhalten demnach eine lineare Konvergenzordnung mit ρ =

p

kxk − x⋆ k2 ≤ cρk kx0 − x⋆ k2 1 − µ2 /L2 < 1.

△

2.2 Praktische Variante Um einen praktikablen Algorithmus für das Gradientenverfahren zu erhalten, müssen wir noch das eindimensionale Minimierungsproblem (2.1) lösen. Die Bedingung d f (xk + αdk ) = ∇f (xk + αdk )T dk = 0 dα stellt eine nichtlineare Gleichung für α ∈ R dar, deren exaktes Lösen viel zu teuer ist. Daher halbiert man die Schrittweite ausgehend von αk = 1 solange, bis ein Abstieg erzielt wurde.

2.2. Praktische Variante

13

Bemerkung Ist f (x) = 21 xT Ax + bT x + c eine quadratische Funktion, dann kann (2.1) exakt gelöst werden: dT dk , wobei dk := b − Axk . αk = k dk Adk In diesem Fall ist das idealisierte Gradientenverfahren leicht durchführbar und konvergiert gemäß Satz 2.3 linear. Der Kontraktionsfaktor ρ hängt dabei stark von der Kondition der Matrix A ab, denn es gilt cond2 (A) − 1 ρ= . cond2 (A) + 1 △ Algorithmus 2.4 (Gradientenverfahren) input: Funktion f : Rn → R und Startnäherung x0 ∈ Rn output: Folge von Iterierten {xk }k∈N ➀ Initialisierung: wähle σ ∈ (0, 1) und setze k := 0 ➁ berechne den Antigradienten dk := −∇f (xk ) und setze αk := 1 ➂ solange f (xk + αk dk ) > f (xk ) − σαk kdk k22

(2.3)

setze αk := αk /2 ➃ setze xk+1 := xk + αk dk ➄ erhöhe k := k + 1 und gehe nach ➁

Bemerkung Die Liniensuche ➂ wird als Armijo-Schrittweitenregel bezeichnet. Sie garantiert nicht nur Abbruchbedingung f (xk+1) < f (xk ), sondern die Armijo-GoldsteinBedingung
f xk − αk ∇f (xk ) ≤ f (xk ) − σαk k∇f (xk )k22 . (2.4)
Dabei bricht die Liniensuche mit einem αk > 0 ab, denn für ϕ(α) = f xk − α∇f (xk ) gilt nämlich
f xk − α∇f (xk ) = ϕ(α) = ϕ(0) + αϕ′ (0) + O(α) = f (xk ) − αk∇f (xk )k22 + O(α).

△

Der nachfolgende Satz liefert ein globales Konvergenzresultat für das Gradientenverfahren unter der Verwendung der Armijo-Schrittweitenregel. Er benötigt nur die Lipschitzstetigkeit von f , aber keine Konvexität. Satz 2.5 Es sei D ⊂ Rn eine offene Menge, in der f stetig differenzierbar, nach unten beschränkt und ∇f zudem Lipschitz-stetig ist. Ferner sei neben x0 auch die gesamte Niveaumenge {x ∈ Rn : f (x) ≤ f (x0 )} in D enthalten. Dann gilt für die Iterierten {xk }k≥0 aus Algorithmus 2.4 k→∞ ∇f (xk ) −→ 0.

14

Kapitel 2. Gradientenverfahren

Beweis. Da D die gesamte Niveaumenge enthält, ist sichergestellt, dass die Iterierten {xk }k≥0 alle in D enthalten sind. Nach Konstruktion ist dann die Folge {f (xk )}k≥0 monoton fallend und nach unten beschränkt. Daher folgt aus der Armijo-Goldstein-Bedingung (2.4), dass f (x0 ) ≥ f (x1 ) + σα0 k∇f (x0 )k22 ≥ · · · ≥ f (xk+1 ) + σ

k X ℓ=0

αℓ k∇f (xℓ )k22 ≥ 0.

Da die Reihe auf der rechten Seite notwendigerweise für k → ∞ konvergent ist, folgt k→∞

αk k∇f (xk )k22 −→ 0.

(2.5)

Es verbleibt zu zeigen, dass αk > ε für ein ε > 0. Für festes k ist aufgrund von Algorithmus 2.4 αk = 1 oder die Armijo-Goldstein-Bedingung ist für 2αk verletzt:
2σαk k∇f (xk )k22 > f (xk ) − f xk − 2αk ∇f (xk ) = 2αk k∇f (xk )k22 − R2 (xk , 2αk ) mit dem Taylor-Restglied  
T R2 (xk , 2αk ) = 2αk k∇f (xk )k22 − ∇f xk − ξ∇f (xk ) ∇f (xk ) ,

ξ ∈ (0, 2αk ).

Da ∇f Lipschitz-stetig ist, ist das Taylor-Restglied R2 (xk , 2αk ) durch ναk2 k∇f (xk )k22 für ein geeignetes ν > 0 beschränkt. Daher ist σαk2 k∇f (xk )k22 > 2αk k∇f (xk )k22 − 2ναk k∇f (xk )k22 = 2αk (1 − ν)k∇f (xk )k22 , dies bedeutet

2(1 − ν) =: ε > 0. σ Somit bleibt αk für alle k ≥ 0 größer als min{1, ε} und daher folgt die Behauptung aus (2.5). αk >

Beachte: Satz 2.5 besagt nicht, dass die Folge {xk }k≥0 selber konvergiert. Selbst wenn die Folge {xk }k≥0 konvergiert, braucht der Grenzwert darüberhinaus kein Minimum von f zu sein. Beispiel 2.6 Gegeben sei die Funktion f (ξ, η) = ξ 2 + (η 2 − 1)2 + ξ 2 (η 2 − 1)2 → min . Das Minimum von f ist 0 und wird offensichtlich für ξ = 0 und η = ±1 angenommen. Hat eine Iterierte xk von Algorithmus 2.4 die Form xk = [ξk , 0]T , dann gilt     4ξk 2ξ + 2ξ(η 2 − 1)2 . = ∇f (xk ) = 0 2η(η 2 − 1)(1 + ξ 2 ) (ξ,η)=(ξk ,0)

Daher hat die nächste Iterierte zwangsläufig wieder die Form xk+1 = [ξk+1 , 0]T und nach Satz 2.5 konvergiert ∇f (xk ) = [4ξk , 0]T gegen Null. Deshalb streben auch ξk und xk gegen Null für k → ∞. Dennoch ist [0, 0]T lediglich ein Sattelpunkt von f , da f (0, η) für η = 0 ein lokales Maximum aufweist. △

15

3. Newton-Verfahren 3.1 Lokales Newton-Verfahren Sei f : Rn → R zweimal stetig differenzierbar. Beim Newton-Verfahren wird das Minimierungsproblem f (x) → min gelöst, indem sukzessive die quadratischen Näherungen 1 qk (x) := f (xk ) + ∇f (xk )T (x − xk ) + (x − xk )T ∇2 f (xk )(x − xk ) 2 minimiert. Ist die Hesse-Matrix ∇2 f (xk ) positiv definit, so ist die neue Iterierte xk+1 genau die Lösung der Gleichung !

∇qk (x) = ∇f (xk ) + ∇2 f (xk )(x − xk ) = 0. Hieraus ergibt sich die Iterationsvorschrift
−1 xk+1 = xk − ∇2 f (xk ) ∇f (xk ),

k = 0, 1, 2, . . . .

(3.1)

Algorithmus 3.1 (Newton-Verfahren) input: Funktion f : Rn → R und Startnäherung x0 ∈ Rn output: Folge von Iterierten {xk }k∈N ➀ Initialisierung: setze k := 0 ➁ berechne die Newton-Richtung durch Lösen des linearen Gleichungssystems ∇2 f (xk )dk = −∇f (xk )

(3.2)

➂ setze xk+1 := xk + dk ➃ erhöhe k := k + 1 und gehe nach ➁ Der mächtigste Satz zum Newton-Verfahren ist der Satz von Newton-Kantorovich. Er liefert nicht nur die Konvergenz des Newton-Verfahrens, sondern auch die Existenz eines stationären Punkts. Satz 3.2 (Newton-Kantorovich) Sei D ⊂ Rn offen und konvex und f : D ⊂ Rn → R zweimal stetig differenzierbar mit k∇2 f (x) − ∇2 f (y)k2 ≤ Lkx − yk2

16

Kapitel 3. Newton-Verfahren f’(x)

x0 x2

x1

Abbildung 3.1: Geometrische Interpretation des Newton-Verfahrens.

für alle x, y ∈ D. Für ein gegebenes D nehmen wir an, dass die Hesse-Matrix

2 x0 ∈
−1 2

∇ f (x0 ) invertierbar ist mit α := ∇ f (x0 2 . Ferner gelte und für


−1 β := ∇2 f (x0 ) ∇f (x0 ) 2 ≤

1 2αL

 p 1  1 ± 1 − 2αβL ≥ 0 αL n − sei S := {x ∈ R : kx − x0 k2 ≤ γ } ⊂ D. Dann sind die Iterierten des Newton-Verfahrens (3.1) wohldefiniert, liegen alle in S und konvergieren gegen ein x⋆ ∈ S mit ∇f (x⋆ ) = 0, welches eindeutig ist in D ∩ {x ∈ Rn : kx − x0 k2 ≤ γ + }. Insbesondere ist die Konvergenz quadratisch, falls 2αβL < 1 ist. Dies bedeutet, es existiert ein c > 0 so dass gilt γ ± :=

kx⋆ − xk+1 k2 ≤ ckx⋆ − xk k22 . Bevor wir den Beweis dieses Satzes erbringen, stellen wir drei hilfreiche Lemmata zur Verfügung. Dabei gelten stets die Voraussetzungen des Satzes von Newton-Kantorovich. Lemma 3.3 Sei {yk }k≥0 eine Folge im Rn und {tk }k≥0 eine Folge nichtnegativer, monoton wachsender, reeller Zahlen mit tk → t⋆ < ∞. Gilt kyk+1 − yk k2 ≤ tk+1 − tk ,

k = 0, 1, . . . ,

dann existiert ein eindeutig bestimmtes y⋆ ∈ Rn mit yk → y⋆ und ky⋆ − yk k2 ≤ t⋆ − tk ,

k = 0, 1, . . . .

Dies bedeutet, die Folge {tk }k≥0 majorisiert {yk }k≥0 .

3.1. Lokales Newton-Verfahren

17

Beweis. Es gilt

k+p−1

X

(yℓ+1 − yℓ ) kyk+p − yk k2 =

2

ℓ=k k+p−1

≤

X ℓ=k

kyℓ+1 − yℓ k2

k+p−1

≤

X

|

ℓ=k

tℓ+1 − tℓ {z

=tk+p −tk

≤ t⋆ − tk .

}

Folglich ist {yk }k≥0 eine Cauchy-Folge im Rn , woraus die Behauptung folgt. Lemma 3.4 Für alle x ∈ T := {x ∈ Rn : kx − x0 k2 < 1/(αL)} ist die Hesse-Matrix ∇2 f (x) invertierbar mit

2


∇ f (x) −1 ≤ 2

Ist außerdem N(x) := x − ∇2 f (x)


−1

α . 1 − αLkx − x0 k2

(3.3)

∇f (x) ∈ T , dann gilt




N N(x) − N(x) ≤ 1 αLkx − N(x)k2 . 2 2 1 − αLkx0 − N(x)k2 Beweis. Für x ∈ T gilt

2






∇ f (x0 ) −1 ∇2 f (x0 ) − ∇2 f (x) ≤ ∇2 f (x0 ) −1 k∇2 f (x0 ) − ∇2 f (x) < 1. 2 | {z }2 | {z }2 ≤α

≤Lkx0 −xk2

Die Satz von der Neumann-Reihe liefert daher die Existenz von

∞  2
−1 2
k 
−1 2
−1 X ∇ f (x0 ) ∇ f (x0 ) − ∇2 f (x) . I − ∇2 f (x0 ) ∇ f (x0 ) − ∇2 f (x) = k=0

Speziell folgt die Anschätzung ∞ X

  k



I − ∇2 f (x0 ) −1 ∇2 f (x0 )−∇2 f (x) −1 ≤ αLkx −xk = 0 2 2 k=0

1 . 1 − αLkx0 − xk2

Die Abschätzung (3.3) ergibt sich nun aus der Identität ∇2 f (x)


−1


−1 2
−1 2 = I − ∇2 f (x0 ) ∇ f (x0 ) − ∇2 f (x) ∇ f (x0 ).

18

Kapitel 3. Newton-Verfahren

Die zweite Aussage des Lemmas sieht man wie folgt. Mit dem soeben gezeigten gilt




−1



2

N N(x) − N(x) = ∇f N(x) − N(x)

N(x) − ∇ f N(x) 2 2



−1

2

≤ ∇ f N(x)

∇f N(x) 2 2


α

∇f N(x) . ≤ 2 1 − αLkN(x) − x0 k2


Wir müssen also nur noch ∇f N(x) 2 abschätzen. Aufgrund der Identität
∇f (x) − ∇2 f (x) N(x) − x = 0 schließen wir






∇f N(x) = ∇f N(x) − ∇f (x) − ∇2 f (x) N(x) − x . 2 2 Mit

∇f (y) − ∇f (x) = und y := N(x) folgt

Z

1

0


∇2 f x + t(y − x) (y − x) dt

Z 1

 2

2

∇f N(x) = ∇ f (x + t(y − x) − ∇ f (x) (y − x) dt

2 0 2 Z 1


≤ k∇2 f (x + t(y − x) − ∇2 f (x) 2 ky − xk2 dt {z } 0 | ≤Ltky−xk2

L ≤ kN(x) − xk22 . 2

Damit ergibt sich schließlich die Behauptung. Lemma 3.5 Die Folge der Newton-Iterierten {xk }k≥0 ist wohldefiniert und durch die Folge {tk }k≥0 mit t2 β − αL 2 k t0 := 0, tk+1 := , k = 1, 2, . . . (3.4) 1 − αLtk majorisiert. Insbesondere konvergiert tk monoton von unten gegen t⋆ := γ − .

Beweis. Wir zeigen zunächst induktiv, dass 0 = t0 < t1 < · · · < tk < t⋆ gilt. Da für k = 0 nichts zu zeigen ist, können wir annehmen, dass die Aussage für ein k ∈ N0 bewiesen ist. Mit quadratischer Ergänzung folgt die Gleichung tk+1 − tk =

β−

αL 2 t 2 k

− tk + αLt2k 1 − αLtk

=0

z }| { αL αL 2 2 (t − t ) + t − (1 − αLt )t β − k k−1 k−1 k 2 k−1 = 2 1 − αLtk 1 αL(tk − tk−1 )2 . = 2 1 − αLtk

3.1. Lokales Newton-Verfahren

19

Wegen 1 − αLtk > 0 und tk−1 < tk schließen wir tk < tk+1 . Um tk+1 < t⋆ zu zeigen, betrachten wir das Polynom p(t) := αL t2 − t + β, welches die zwei Nullstellen γ − = t⋆ und 2 γ + besitzt. Wegen t⋆ − tk+1 =

t⋆ − αLtk t⋆ − β + 1 − αLtk

αL 2 t 2 k

=0

z }| { αL ⋆ 2 p(t⋆ ) (t − t ) − k 2 = 1 − αLtk 1 αL(t⋆ − tk )2 >0 = 2 1 − αLtk | {z } >0

ist in der Tat tk+1 < t⋆ , womit ist der Induktionsschritt k 7→ k + 1 erbracht ist. Insbesondere folgt zwingend die Konvergenz tk → t⋆ , da Fixpunkte der Iteration (3.4) der Bedingung 0 = p(t) genügen müssen und t⋆ = γ − < γ + gilt. Wir zeigen nun, dass alle xk existieren, in S liegen und jeweils kxk − xk−1 k2 ≤ tk − tk−1 erfüllen. Für k = 0 ist dies sicherlich richtig. Es möge also die Folge x1 , . . . , xk das Behauptete erfüllen. Wegen S ⊂ T folgt der Induktionsschritt k 7→ k + 1 dann aus Lemma 3.4 gemäß


kxk+1 − xk k2 = N N(xk−1 ) − N(xk−1 ) 2 1 αLkxk − xk−1 k2 ≤ 2 1 − αLkxk − x0 k2 1 αL(tk − tk−1 ) ≤ 2 1 − αLtk = tk+1 − tk und kxk+1 − x0 k2 ≤

k X ℓ=0

kxℓ+1 − xℓ k2 ≤

k X ℓ=0

tℓ+1 − tℓ = tk+1 − t0 ≤ t⋆ = γ − .

Wir sind nun in der Lage, den Satz von Newton-Kantorovich beweisen zu können, wobei wir nur die Existenz eines stationären Punktes x⋆ nachweisen und nicht dessen Eindeutigkeit. Beweis des Satzes von Newton-Kantorovich. Nach Lemma 3.5 werden die Newton-Iterierten {xk }k≥0 von {tk }k≥0 majorisiert und konvergieren gemäß Lemma 3.3 gegen ein eindeutig bestimmtes x⋆ ∈ S. Wegen k∇f (xk )k2 = k∇2 f (xk )(xk+1 − xk )k2



≤ ∇2 f (xk ) − ∇2 f (x0 ) (xk+1 − xk ) 2 + k∇2 f (x0 )(xk+1 − xk )k2  ≤ Lt⋆ + k∇2 f (x0 )k2 kxk+1 − xk k2 {z } {z }| | 0 und setze k := 0 ➁ berechne, falls möglich, die Newton-Richtung durch Lösen des linearen Gleichungssystems ∇2 f (xk )dk = −∇f (xk ) ➂ ist das Gleichungssystem nicht lösbar oder ist die Bedingung ∇f (xk )T dk ≤ −ρk∇f (xk )k22 nicht erfüllt, dann setze dk := −∇f (xk ) ➃ setze αk := 1 ➄ solange f (xk + αk dk ) > f (xk ) + σαk ∇f (xk )T dk

(3.5)

(3.6)

setze αk := αk /2 ➅ setze xk+1 := xk + αk dk ➆ erhöhe k := k + 1 und gehe nach ➁

Die Bedingungen (3.5) und (3.6) garantieren zusammen die Armijo-Goldstein-Bedingung (2.4) und somit einen hinreichenden Abstieg. Die Konvergenz ∇f (xk ) → 0 folgt daher analog zum entsprechenden Satz 3.6 über das Gradientenverfahren. Aus dem nachfol-

3.2. Globalisiertes Newton-Verfahren

21

genden Satz ergibt sich, dass das globalisierte Newton-Verfahren in der Umgebung eines Minimums zum lokalen Newton-Verfahren wird. Satz 3.7 Seien f : Rn → R zweimal stetig differenzierbar und x⋆ ∈ Rn ein stationärer Punkt, an dem die Hesse-Matrix ∇2 f (x⋆ ) positiv definit ist. Konvergiert die Folge {xk }k≥0 der Newton-Iterierten gegen x⋆ , dann exitistiert ein k0 ∈ N mit f (xk + dk ) > f (xk ) + σ∇f (xk )T ∇2 f (xk )∇f (xk )

(3.7)

für alle k ≥ k0 und jedes feste σ ∈ (0, 1/2). Beweis. Aus den Voraussetzung folgt aus Stetigkeitsgründen zunächst k→∞

(3.8)

k∇f (xk )k2 −→ 0.

Ferner ist ∇2 f (x⋆ ) invertierbar und aufgrund der Stetigkeit schließen wir auf die Existenz einer Umgebung U ⊂ Rn von x⋆ mit

2


∇ f (x) −1 ≤ c für alle x ∈ U. 2

Daher existiert ein k0 ∈ N, so dass

2


∇ f (xk ) −1 ≤ c für alle k ≥ k0 . 2

Aus der Newton-Gleichung dk := − ∇2 f (xk )


−1

∇f (xk ) folgt deshalb


−1 k→∞ kdk k2 ≤ ∇2 f (xk ) 2 k∇f (xk )k2 ≤ ck∇f (xk )k2 −→ 0.

Da mit ∇2 f (x⋆ ) auch ∇2 f (x⋆ )


−1

∇f (xk )T ∇2 f (xk )


−1

positiv definit ist, können wir k0 vergrößern, so dass

∇f (xk ) ≥ ck∇f (xk )k22

für alle k ≥ k0 .

(3.9)

Nach dem Taylorschen Satz exitiert ein ξ k auf der Verbindungsstrecke von xk und xk +dk , so dass gilt 1 f (xk + dk ) = f (xk ) + ∇f (xk )T dk + dTk ∇2 f (ξ k )dk . 2 Wegen xk → x⋆ und dk → 0 gilt auch ξk → x⋆ und folglich k→∞

k∇2 f (ξk ) − ∇2 f (xk )k2 −→ 0. Wegen σ ∈ (0, 1/2) ist daher

  1 1 c + c2 k∇2 f (ξk ) − ∇2 f (xk )k2 ≤ 0 σ− 2 2 | {z } 0 und bestimme eine Suchrichtung dk mit k∇f (xk ) + ∇2 f (xk )dk k2 ≤ ηk k∇f (xk )k2

(3.11)

➂ setze xk+1 := xk + dk ➃ erhöhe k := k + 1 und gehe nach ➁ Die Konvergenzrate des inexakten Newton-Verfahrens ergibt sich aus dem nachfolgenden Satz. Speziell gibt er uns eine Bedingung an die Wahl der Toleranzen {ηk } an, so dass das Verfahren dennoch quadratisch konvergiert.

3.3. Inexaktes Newton-Verfahren

23

Satz 3.9 Seien f : Rn → R zweimal stetig differenzierbar und x⋆ ∈ Rn ein stationärer Punkt, an dem die Hesse-Matrix ∇2 f (x⋆ ) positiv definit ist. Dann existiert eine Umgebung U ⊂ Rn von x⋆ , so dass für alle x0 ∈ U gelten: (i.) Ist ηk ≤ η für ein hinreichend kleines η, dann ist Algorithmus 3.8 wohldefiniert und die durch ihn erzeugte Folge {xk }k≥0 konvergiert mindestens linear gegen x⋆ . (ii.) Die Konvergenzrate ist superlinear, falls ηk → 0 gilt. (iii.) Die Konvergenzrate ist quadratisch, falls ηk = O(k∇f (xk )k2 ) gilt und ∇2 f lokal Lipschitz-stetig ist. Beweis. Da f zweimal stetig differenzierbar ist, ist ∇f lokal Lipschitz-stetig. Daher existiert eine Umgebung U ⊂ Rn von x⋆ und ein L > 0 derart, dass k∇f (x)k2 = k∇f (x) − ∇f (x⋆ )k2 ≤ Lkx − x⋆ k2

für alle x ∈ U.

Da ∇2 f (x⋆ ) invertierbar ist, folgt, indem wir U eventuell verkleinern, aufgrund der Stetigkeit

2


∇ f (x) −1 ≤ c für alle x ∈ U. 2 Aus der Definition der Ableitung einer Funktion schließen wir, indem wir wieder U eventuell verkleinern, dass 1 k∇f (x⋆ ) − ∇f (x) − ∇2 f (x)(x − x⋆ )k2 ≤ kx − x⋆ k2 . {z } 4c | =O(kx−x⋆ k2 )

Setze nun η := 1/(4cL) und wähle x0 ∈ U. Dann ist ∇2 f (x0 ) regulär und es lässt sich ein d0 mit (3.11) berechnen. Somit ist x1 wohldefiniert und wegen


−1
−1  2 kx1 − x⋆ k2 = x0 − x⋆ − ∇2 f (x0 ) ∇f (x0 ) + ∇2 f (x0 ) ∇ f (x0 )d0 + ∇f (x0 ) 2


−1 n ≤ ∇2 f (x0 ) 2 ∇2 f (x0 )(x0 − x⋆ ) − ∇f (x0 ) + ∇f (x⋆ ) 2

2

o

+ ∇ f (x0 )d0 + ∇f (x0 ) 2 {z } | ≤ηk∇2 f (x0 )k2



 1 ≤c + ηL kx0 − x⋆ k2 4c 1 = kx0 − x⋆ k2 2

wieder in U enthalten. Mit vollständiger Induktion schließen wir kxk − x⋆ k2 ≤

1 kx0 − x⋆ k2 , 2k

was die lineare Konvergenz xk → x⋆ nachweist. Zum Nachweis von (ii.) bemerken wir, dass analog zur obigen Ungleichungskette auch n o

kxk+1 − x⋆ k2 ≤ c ∇2 f (xk )(xk − x⋆ ) − ∇f (xk ) + ∇f (x⋆ ) 2 + ηk k∇2 f (xk )k2 n o

≤ c ∇2 f (xk )(xk − x⋆ ) − ∇f (xk ) + ∇f (x⋆ ) 2 + ηk Lkxk − x⋆ k2

24

Kapitel 3. Newton-Verfahren

bewiesen werden kann. Für ηk → 0 ist daher

kxk+1 − x⋆ k2 = O(kxk − x⋆ k2 ),

das heißt, {xk }k≥0 konvergiert superlinear gegen x⋆ . Die lokal quadratische Konvergenz wird ganz ähnlich nachgewiesen: Mit Z 1
⋆ ∇f (xk ) − ∇f (x ) = ∇2 f x⋆ + t(xk − x⋆ ) (xk − x⋆ ) dt 0

folgt mittels der lokalen Lipschitzstetigkeit

2

∇ f (xk )(xk − x⋆ ) − ∇f (xk ) + ∇f (x⋆ ) 2

Z 1

 2
⋆ ⋆ ⋆

= ∇ f (xk ) − f x + t(xk − x ) (xk − x ) dt

0 2 Z 1

2


∇ f (xk ) − ∇2 f x⋆ + t(xk − x⋆ ) kxk − x⋆ k2 dt ≤ {z }2 0 | ≤tLkxk −x⋆ k2

L kxk − x⋆ k22 . 2 Setzt man dies zusammen mit ηk = O(k∇f (xk )k2 ) = O(kxk − x⋆ k2 ) in die obige Abschätzung ein, dann ergibt sich ≤

kxk+1 − x⋆ k2 = O(kxk − x⋆ k22 ), das ist Aussage (iii.).

3.4 Trust-Region-Verfahren Die quadratische Approximation 1 f (xk + d) ≈ f (xk ) + ∇f (xk )T d + dT ∇2 f (xk )d = qk (d) 2 ist nur für kleines kdk2 genau. Man schränkt daher den Bereich ein, in dem man der quadratischen Approximation qk (d) vertraut. Dies erklärt den Begriff trust region. Im Zusammenhang mit dem Newton-Verfahren wird das Update dk beispielsweise aus min qk (d) unter der Nebenbedingung kdk2 ≤ ∆k

d∈Rn

(3.12)

bestimmt, wobei der Radius ∆k noch geeignet zu bestimmen ist. Wir haben also in jedem Schritt ein quadratisches Minimierungsproblem mit Nebenbedingungen zu lösen. Wir wenden uns der Wahl des Radius ∆k zu. Dieser wird so gesteuert, dass eine möglichst gute Übereinstimmung der quadratischen Näherung qk mit der Funktion f (xk + ·) sichergestellt ist. Liegt das Qualitätsmaß ρk =

f (xk ) − f (xk + dk ) qk (0) − qk (dk )

nahe bei 1, stimmt die tatsächliche Reduktion (Zähler) mit der vorhergesagten Reduktion (Nenner) gut überein und der Radius ∆k kann sogar vergrößert werden. Ist das Qualitätsmaß klein, dann sollte der Radius ∆k verkleinert und der Schritt wiederholt werden. Ansonsten behält man ∆k bei.

3.4. Trust-Region-Verfahren

25

Algorithmus 3.10 (Trust-Region-Verfahren) input: Funktion f : Rn → R und Startnäherung x0 ∈ Rn output: Folge von Iterierten {xk }k∈N ➀ Initialisierung: wähle 0 < ρ− < ρ+ < 1 und setze ∆0 := 1 und k := 0 ➁ bestimme dk durch näherungsweise Lösung von (3.12) und berechne ρk :=

f (xk ) − f (xk + dk ) qk (0) − qk (dk )

➂ falls ρk > ρ− setze xk+1 := xk + dk , sonst setze ∆k := ∆k /2 und gehe nach ➁ ➃ falls ρk > ρ+ setze ∆k+1 := 2∆k , sonst setze ∆k+1 := ∆k ➄ erhöhe k := k + 1 und gehe nach ➁ Ausgangspunkt zur näherungsweisen Lösung von (3.12) ist der Cauchy-Punkt, dessen Konstruktion in zwei Schritten erfolgt. 1. Bestimme dk aus min f (xk ) + ∇f (xk )T d unter der Nebenbedingung kdk2 ≤ ∆k .

d∈Rn

Dieser Richtungsvektor lässt sich explizit angeben: dk = −

∆k ∇f (xk ). k∇f (xk )k2

2. Bestimme τk aus min qk (τ dk ) unter der Nebenbedingung kτ dk k2 ≤ ∆k τ >0

und setze d⋆k := τk dk . Für die Berechnung von τk sind zwei Fälle zu unterscheiden. Falls ∇f (xk )T ∇2 f (xk )∇f (xk ) ≤ 0 gilt, so folgt τk = 1. Ansonsten ist   k∇f (xk )k32 1 . τk = min 1, ∆k ∇f (xk )T ∇2 f (xk )∇f (xk ) Unser Ziel ist es nun, eine Kombination zu finden, die für kleine Werte ∆k mit dem CauchyPunkt übereinstimmt und für größere Werte ∆k in das Newton-Verfahren übergeht. Wir stellen dazu die sogenannte Dogleg-Strategie vor und verwenden zur Abkürzung statt (3.12) 1 q(d) := f + gT d + dT Hd unter der Nebenbedingung kdk2 ≤ ∆k . 2 Für kleine Werte ∆k fällt der quadratische Term nicht ins Gewicht und wir gehen in kgk22 g, für große Werte ∆k wollen wir hingegen Richtung des Cauchy-Punktes dC := − gT Hg −1 die Newton-Richtung dN := −H g verwenden: ( τ dC , falls 0 ≤ τ ≤ 1, d(τ ) = dC + (τ − 1)(dN − dC ), falls 1 ≤ τ ≤ 2.

26

Kapitel 3. Newton-Verfahren

Lemma 3.11 Es sei H postitiv definit und g 6= 0. Dann ist (i.) die Funktion kd(τ )k2 streng monoton wachsend in τ , falls dN 6= dC , und
(ii.) die Funktion q d(τ ) monoton fallend in τ . Beweis. Beide Aussagen im Fall 0 ≤ τ ≤ 1 offensichtlich. Für 1 ≤ τ := σ + 1 ≤ 2 betrachten wir die Funktion 1 ϕ(σ) = kd(1 + σ)k22 2 1 = kdC + σ(dN − dC )k22 2 1 1 = kdC k22 + σdTC (dN − dC ) + σ 2 kdN − dC k22 2 2 und zeigen ϕ′ (σ) > 0 für alle σ ∈ (0, 1): ϕ′ (σ) = dTC (dN − dC ) + σ kdN − dC k22 {z } | >0 da dN 6=dC

dTC (dN T

− dC )   g g T gT g −1 = T g H g− T g g Hg g Hg   (gT g)2 gT g T −1 g H g 1− T = T g Hg (g Hg)(gT H−1 g) ≥ 0, >

wobei sich die letzte Ungleichung aus (gT g)2 = (gT H1/2 H−1/2 g)2 ≤ kH1/2 gk22 kH−1/2 gk22 = (gT Hg)(gT H−1 g) ergibt. Damit ist Aussage (i.) gezeigt.
Um Aussage (ii.) nachzuweisen, bilden wir ψ(σ) = q d(1 + σ) und erhalten ψ ′ (σ) = (g + HdC )T (dN − dC ) + σ (dN − dC )T H(dN − dC ) | {z } ≥0

 T ≤ g + HdC + H(dN − dC ) (dN − dC ) = (g + HdN )T (dN − dC ) = 0.

Bemerkung Aus diesem Lemma folgt, dass es genau einen Punkt gibt, an
dem der Dogleg-Pfad den Kreis mit Radius ∆k schneidet, falls kdN k2 ≥ ∆k . Da q d(τ ) monoton fällt, nimmt man diesen Schittpunkt als Update. Ist hingegen kdN k2 ≤ ∆k , dann verläuft der Pfad ganz im Innern des Trust-Region-Bereichs. Folglich wählt man als Update   falls kdN k2 ≤ ∆k , dN , dk = dS , falls kdC k2 < ∆k < kdN k2 ,   ∆k d , falls kdC k2 ≥ ∆k , kdC k2 C

3.4. Trust-Region-Verfahren

27

wobei  dS := dC + (τ ⋆ − 1)(dN − dC ) mit τ ⋆ = argmin (kdC + (τ − 1)(dN − dC )k2 − ∆k )2 . τ >0

Bei Implementierung muss man jedoch zuallererst überprüfen, ob überhaupt qk (dN ) < qk (0)

gilt. Denn ist ∇2 f (xk ) nicht positiv definit, muss dem nicht so sein. In diesem Fall ist die Newton-Richtung von vornherein zu verwerfen und der Cauchy-Punkt zu nehmen. △ Unser Ziel ist der Nachweis globaler Konvergenz des Newton-Verfahrens mit der vorgestellten Trust-Region-Technik. Ausgangspunkt dafür ist die Reduktion des Funktionals durch die Verwendung des Cauchy-Punktes. Lemma 3.12 Für den Cauchy-Punkt d⋆k gilt qk (d⋆k )

  k∇f (xk )k2 1 . ≤ qk (0) − k∇f (xk )k2 min ∆k , 2 k∇2 f (xk )k2

Beweis. Wir müssen drei Fälle unterscheiden. (i.) Falls ∇f (xk )T ∇2 f (xk )∇f (xk ) ≤ 0 folgt ∆k 1 ∆2k ∇f (xk )T ∇2 f (xk )∇f (xk ) k∇f (xk )k22 + k∇f (xk )k2 2 k∇f (xk )k22 ≤ −∆k k∇f (xk )k2   k∇f (xk )k2 1 . ≤ − k∇f (xk )k2 min ∆k , 2 k∇2 f (xk )k2

qk (d⋆k ) − qk (0) = −

(ii.) Im Fall 0 < ∇f (xk )T ∇2 f (xk )∇f (xk ) ≤ k∇f (xk )k32 /∆k ergibt sich ∆k 1 ∆2k ∇f (xk )T ∇2 f (xk )∇f (xk ) k∇f (xk )k22 + k∇f (xk )k2 2 k∇f (xk )k22 ∆2k k∇f (xk )k32 1 ≤ −∆k k∇f (xk )k2 + 2 k∇f (xk )k22 ∆k 1 = − ∆k k∇f (xk )k2 2   1 k∇f (xk )k2 ≤ − k∇f (xk )k2 min ∆k , . 2 k∇2 f (xk )k2

qk (d⋆k ) − qk (0) = −

(iii.) Ist schließlich ∇f (xk )T ∇2 f (xk )∇f (xk ) > k∇f (xk )k32 /∆k , so liegt der Cauchy-Punkt im Innern des Trust-Region-Bereichs und ist durch d⋆k = −

k∇f (xk )k22 ∇f (xk ) ∇f (xk )T ∇2 f (xk )∇f (xk )

28

Kapitel 3. Newton-Verfahren

gegeben. Somit gilt k∇f (xk )k42 ∇f (xk )T ∇2 f (xk )∇f (xk ) 1 k∇f (xk )k42 + ∇f (xk )T ∇2 f (xk )∇f (xk ) 2 (∇f (xk )T ∇2 f (xk )∇f (xk ))2 k∇f (xk )k42 1 =− 2 ∇f (xk )T ∇2 f (xk )∇f (xk ) 1 k∇f (xk )k42 ≤− 2 k∇f (xk )k22 k∇2 f (xk )k2   k∇f (xk )k2 1 . ≤ − k∇f (xk )k2 min ∆k , 2 k∇2 f (xk )k2

qk (d⋆k ) − qk (0) = −

Bemerkung Mit der Dogleg-Strategie erreicht man mindestens eine Reduktion des Funktionals qk (d) durch den Cauchy-Punkt, das heißt, es gilt   1 k∇f (xk )k2 . qk (dk ) ≤ qk (0) − k∇f (xk )k2 min ∆k , 2 k∇2 f (xk )k2 △ Satz 3.13 Auf der Niveaumenge N := {x ∈ Rn : f (x) ≤ f (x0 )} sei f zweimal stetig differenzierbar und nach unten beschränkt. Ist ∇2 f (x) beschränkt auf der gesamten Niveaumenge N, dann erfüllen die vom Trust-Region-Verfahren 3.10 mit der Dogleg-Strategie berechneten Iterierten {xk }k≥0 k→∞ ∇f (xk ) −→ 0. Beweis. Wir vollziehen den Beweis des Satzes in vier Schritten. (i.) Für das Maß ρk zur Anpassung des Trust-Region-Radius gilt f (xk ) − f (xk + dk ) − {qk (0) − qk (dk )} f (xk + dk ) − qk (dk ) |ρk − 1| = = qk (0) − qk (dk ) . qk (0) − qk (dk )

Die Taylor-Entwicklung von f um xk liefert

1 f (xk + dk ) = f (xk ) + ∇f (xk )T dk + dTk ∇2 f (xk + ξdk )dk 2 für ein ξ ∈ (0, 1). Damit ergibt sich 1 T 2 1 T 2 |qk (dk ) − f (xk + dk )| = dk ∇ f (xk )dk − dk ∇ f (xk + ξdk )dk 2 2   1 2 2 k∇ f (xk )k2 + max k∇ f (xk + tdk )k2 kdk k22 ≤ 0≤t≤1 2 2 ≤ ckdk k2 ,

3.4. Trust-Region-Verfahren

29

wobei c > 0 eine obere Schranke für ∇2 f (x) auf der Niveau-Menge N sei. Da nach Lemma 3.12 gilt   k∇f (xk )k2 1 , (3.13) qk (dk ) ≤ qk (0) − k∇f (xk )k2 min ∆k , 2 c erhalten wir folglich die Abschätzung |ρk − 1| ≤

ckdk k22

k∇f (xk )k2 min{∆k , k∇f (xc k )k2 }

≤

2c∆2k k )k2 } k∇f (xk )k2 {∆k , k∇f (x c

.

(3.14)

 (ii.) Bei erfolgreichem Iterationsschritt ist f (xk ) − f (xk + dk ) ≥ ρ− qk (0) − qk (dk ) , das heißt, gemäß (3.13) gilt   ρ− k∇f (xk )k2 f (xk ) − f (xk + dk ) ≥ . (3.15) k∇f (xk )k2 min ∆k , 2 c Da die Folge {f (xk )}k≥0 nach Konstruktion streng monoton fällt und nach unten beschränkt ist, muss gelten   k∇f (xk )k2 k→∞ min ∆k , −→ 0. (3.16) c (iii.) Wir zeigen nun, dass {k∇f (xk )k2 }k≥0 für eine unendliche Teilfolge {kℓ }ℓ≥0 gegen Null konvergiert. Angenommen, die Behauptung gilt nicht, dann folgt k∇f (xk )k2 ≥ ε > 0 für alle k ≥ K(ε). Aus (3.16) ergibt sich unmittelbar k→∞

∆k −→ 0,

weshalb (3.14) impliziert |ρk − 1| → 0 für k → ∞, beziehungsweise k→∞

ρk −→ 1. Demnach existiert ein M(ε) ≥ K(ε), so dass ρk ≥ ρ+ für alle k ≥ M(ε). Ab dem M(ε)-ten Schritt wird folglich ∆k in jedem Schritt von Algorithmus 3.10 verdoppelt, was jedoch im Widerspruch zu (3.16) steht. (iv.) Wir beweisen nun die Aussage des Satzes. Dazu nehmen wir an, dass eine Teilfolge von {k∇f (xk )k2 }k≥0 nicht gegen Null konvergiert. Nach Aussage (iii.) existiert dann ein ε > 0 und zwei Indizes ℓ < m, so dass k∇f (xℓ )k2 ≥ 2ε,

k∇f (xm )k2 ≤ ε,

k∇f (xk )k2 > ε,

k = ℓ + 1, . . . , m − 1.

Da {f (xk )}k≥0 eine Cauchy-Folge ist, kann ℓ dabei so groß gewählt werden, dass f (xℓ ) − f (xm )
0, νk ≥ 0 und berechne Hk+1 := Φ(Hk , pk , qk , γk , νk ) gemäß (4.2) ➅ erhöhe k := k + 1 und gehe nach ➁ Das Verfahren ist eindeutig durch die Wahl der Parameter γk , νk und die Minimierung in Schritt ➂ fixiert. Die Minimierung xk 7→ xk+1 und ihre Qualität kann man mit Hilfe eines Parameters σk beschreiben, der durch ∇f (xk+1)T dk = σk ∇f (xk )T dk = −σk ∇f (xk )T Hk ∇f (xk ) definiert ist. Falls dk eine Abstiegsrichtung ist, das heißt ∇f (xk )T dk < 0, dann ist σk eindeutig bestimmt. Bei exakter Liniensuche ist σk = 0 wegen ∇f (xk+1)T dk = ϕ′k (αk ) = 0,

wobei ϕk (α) := f (xk + αdk ).

Wir setzen für das folgende σk < 1

(4.3)

voraus. Falls ∇f (xk ) 6= 0 und Hk positiv definit ist, folgt aus (4.3) αk > 0 und deshalb qTk pk = αk {∇f (xk+1) − ∇f (xk )}T dk = αk (σk − 1)∇f (xk )T dk

= −αk (σk − 1)∇f (xk )T Hk ∇f (xk ) > 0, also auch qk 6= 0 und qTk Hk qk > 0. Die Matrix Hk+1 ist damit durch (4.2) wohldefiniert.

33 Die Forderung (4.3) kann nur dann nicht erfüllt werden, wenn ϕ′k (α) = ∇f (xk + αdk )T dk ≤ ∇f (xk )T dk = ϕ′k (0) < 0 für alle α ≥ 0 gilt. Dann ist aber Z f (xk + αdk ) − f (xk ) =

0

α

ϕ′k (t) dt ≤ α∇f (xk )T dk < 0 für alle α ≥ 0,

so dass f (xk + αdk ) für α → ∞ nicht nach unten beschränkt ist. Die Forderung (4.3) bedeutet also keine wesentliche Einschränkung. Damit ist bereits der erste Teil des folgenden Satzes gezeigt, der besagt, dass das Quasi-Newton-Verfahren 4.1 unsere oben aufgestellten Forderungen erfüllt. Satz 4.2 Falls im Quasi-Newton-Verfahren 4.1 die Matrix Hk für ein k ≥ 0 positiv definit ist, ∇f (xk ) 6= 0 und σk < 1 ist, dann ist für alle γk > 0, νk ≥ 0 die Matrix Hk+1 := Φ(Hk , pk , qk , γk , νk ) wohldefiniert und wieder positiv definit. Insbesondere erfüllt sie die Quasi-Newton-Gleichung Hk+1qk = pk . Beweis. Die Wohldefiniertheit von Hk+1 haben wir bereits gezeigt, so dass wir nur noch die positive Definitheit nachweisen müssen. Sei y ∈ Rn \ {0} ein beliebiger Vektor und Hk = LLT die Cholesky-Zerlegung von Hk . Mit Hilfe der Vektoren u := LT y,

v := LT qk

lässt sich yT Hk+1 y wegen (4.2) so schreiben:   vT v (pTk y)2 1 − νk 2γk νk T T y Hk+1 y = γk u u + 1 + γk νk T − γk T (uT v)2 − T (pTk y)(uT v) T pk qk pk qk v v pk qk   √ 2 T 2 T 2 T T (u v) (pk y) u v p y T k = γk u u − + T + γk νk vT v T − √ vT v pk qk pk qk vT v   (pTk y)2 (uT v)2 + . ≥ γk uT u − vT v pTk qk Die Cauchy-Schwarzsche Ungleichung ergibt uT u −

(uT v)2 ≥ 0, vT v

mit Gleichheit genau dann, wenn u = λv für ein λ 6= 0 (wegen y 6= 0). Für u 6= λv ist also yHk+1 y > 0. Für u = λv folgt aus der Nichtsingularität von Hk und L auch 0 6= y = λqk , so dass (pT y)2 = λ2 pTk qk > 0. yHk+1 y ≥ Tk pk qk Da y ∈ Rn \ {0} beliebig war, muss Hk+1 positiv definit sein. Die Quasi-Newton-Gleichung Hk+1qk = pk verifiziert man schließlich sofort mittels (4.2).

34

Kapitel 4. Quasi-Newton-Verfahren

Ein wesentliches Resultat ist, dass das Quasi-Newton-Verfahren im Fall einer quadratischen Funktion f : Rn → R das Minimum nach höchstens n Schritten liefert, sofern die Minimierung in ➂ exakt sind. Da sich jede genügend oft differenzierbare Funktion f in der Nähe ihres Minimums beliebig genau durch eine quadratische Funktion approximieren lässt, lässt diese Eigenschaft vermuten, dass das Verfahren auch bei der Anwendung auf nichtquadratische Funktionen rasch konvergiert. Satz 4.3 Sei

1 f (x) = xT Ax + bT x + c 2 eine quadratische Funktion mit einer positiv definiten Matrix A ∈ Rn×n . Wendet man das Quasi-Newton-Verfahren 4.1 zur Minimierung von f mit den Startwerten x0 und H0 an, wobei man die Minimierungen in ➂ exakt durchführt, so liefert das Verfahren Folgen {xk }k≥0 , {Hk }k≥0, {∇f (xk )}k≥0 , {pk }k≥0 und {qk }k≥0 mit den Eigenschaften: (i.) Es gibt ein kleinstes m ≤ n mit xm = x⋆ = −A−1 b, das heißt, xm ist das eindeutige Minimum von f , insbesondere gilt also ∇f (xm ) = 0. (ii.) Es ist pTk qk > 0 und pTk qℓ = pTk Apℓ = 0 für alle 0 ≤ k 6= ℓ < m. Die Vektoren pk sind demnach A-konjugiert. (iii.) Es gilt pTk ∇f (xℓ ) = 0 für alle 0 ≤ k < ℓ ≤ m. (iv.) Es ist Hℓ qk = λk,ℓ pk für alle 0 ≤ k < ℓ ≤ m mit ( γk γk+1 · · · γℓ−1 , für k < ℓ − 1, γk,ℓ := 1, für k = ℓ − 1. (v.) Falls m = n, so gilt zusätzlich Hm = Hn = PDP−1 A−1, wobei D = diag(γ0,n , γ1,n , . . . , γn−1,n ),

P = [p0 , p1 , . . . , pn−1 ].

Für γk ≡ 1 folgt Hn = A−1 . Beweis. Wir zeigen zunächst induktiv, dass die Bedigungen (ii.)–(iv.) für ein beliebiges m ≥ 0 gelten, falls für alle j < m Hj positiv definit und ∇f (xj ) 6= 0 ist. Da die Aussagen für m = 0 trivialerweise erfüllt sind, können wir annehmen, dass sie für ein beliebiges m ≥ 0 gelten. Der Induktionsschritt m 7→ m + 1 ergibt sich nun wie folgt. Da Hm positiv definit ist, folgt aus ∇f (xm ) 6= 0 sofort dm = −Hm ∇f (xm ) 6= 0 und ∇f (xm )T Hm ∇f (xm ) > 0. Weil exakt minimiert wird, ist αm die Nullstelle von  T 0 = ∇f (xm+1 )T dm = ∇f (xm ) + αm Adm dm ,

αm =

∇f (xm )Hm ∇f (xm ) , dT Adm

also pm = αm dm und

∇f (xm+1 )T pm = αm ∇f (xm+1 )T dm = 0.

(4.4)

35 Deshalb gilt  pTm qm = αm dTm ∇f (xm+1 ) − ∇f (xm ) = −αm dTm ∇f (xm )

= αm ∇f (xm )T Hm ∇f (xm ) >0 und folglich ist Hm+1 nach Satz 4.2 positiv definit. Weiter ist für k < m wegen Apk = qk (iv.)

(iii.)

pTk qm = pTk Apm = qTk pm = −αm qTk Hm ∇f (xm ) = −αm γk,mpTk ∇f (xm ) = 0.

(4.5)

Das ist der Induktionsschritt für Aussage (ii.). Weiter gilt für k < m   m X T T qj = 0 pk ∇f (xm+1 ) = pk ∇f (xk+1 ) + j=k+1

nach dem eben bewiesenen und Aussage (iii.). Zusammen mit (4.4) ergibt dies Aussage (iii.) für m + 1. Den Induktionsschritt für Aussage (iv.) sieht man wie folgt ein. Anhand von (4.2) verifiziert man sofort Hm+1 qm = pm . Wegen Aussage (ii.) für m+ 1 und der Induktionsvoraussetzung hat man ferner für k < m (ii.)

pTm qk = 0,

(iv.)

(ii.)

qTm Hm qk = γk,mqTm pm = 0,

so dass für k < m aus (4.2) folgt (iv.)

Hm+1 qk = γm Hm qk = γm γk,mpk = γk,m+1 pk . Der restliche Beweis ist nun einfach. Die Aussagen (ii.)–(iv.) können nur für m ≤ n Pm−1 richtig sein, da die Vektoren p0 , p1 , . . . , pm−1 linear unabhängig sind. Aus 0 = ℓ=0 λℓ pℓ folgt nämlich durch Multiplikation mit pTk A, k = 0, 1, . . . , m − 1, wegen Aussage (ii.) λk pTk Apk = 0, das heißt, λk = 0. Da wir bewiesen haben, dass die Aussagen (ii.)–(iv.) für beliebiges m gelten, solange ∇f (xm ) 6= 0 ist, muss es also einen ersten Index m ≤ n geben mit ∇f (xm ) = 0,

xm = −A−1 b,

dies bedeutet, es gilt Aussage (i.). Für den Fall m = n gilt wegen Aussage (iv.) zusätzlich Hn Q = PD für die Matrizen D = diag(γ0,n , γ1,n , . . . , γn−1,n ),

P = [p0 , p1 , . . . , pn−1 ],

Q = [q0 , q1 , . . . , qn−1 ].

Wegen AP = Q ergibt sich schließlich wegen der Nichtsingularität der Matrix P die Beziehung Hn = PDP−1 A−1 , Damit ist der Satz vollständig bewiesen.

36

Kapitel 4. Quasi-Newton-Verfahren

Es stellt sich nun die Frage, wie man die Parameter γk und νk wählen soll, um ein möglichst gutes Verfahren zu erhalten. Aussage (v.) aus Satz 4.3 legt die Wahl γk ≡ 1 nahe,
−1 weil dies D = I und folglich limm Hm = ∇2 f (x⋆ ) vermuten lässt, weshalb das Verfahren voraussichtlich ähnlich schnell wie ein Newton-Verfahren konvergiert. Im allgemeinen ist diese Vermutung für nichtquadratische Funktionen aber nur unter zusätzlichen Vorraussetzungen richtig. Nach praktischen Erfahrungen ist die Wahl γk ≡ 1,

νk ≡ 1 (BFGS-Verfahren)

am besten. Bemerkungen 1. Sowohl das DFP-Verfahren als auch das BFGS-Verfahren konvergieren superlinear in der Umgebung eines lokalen Minimus x⋆ , falls f : Rn → R zweimal stetig differenzierbar ist und die Hesse-Matrix in der Umgebung von x⋆ Lipschitz-stetig ist. 2. Eine andere Startmatrix H0 6= I ist denkbar, solange sie symmetrisch und positiv definit ist. 3. In der Praxis macht man gelegentlich Restarts, setzt also Hk := H0 , falls k ∈ mZ mit festem m ∈ N, beispielsweise m = 100. 4. Gerade bei großen Optimierungsproblemen stellt man die Matrix Hk nicht direkt auf, sondern berechnet sie rekursiv aus den Vektoren {(γk , νk , pk , qk )}k≥0. Damit auch bei vielen Schritten der Speicherplatz nicht überhand nimmt, speichert man nur die höchstens letzten m Vektoren. Man erlaubt also ein “Gedächtnis” von m Updates und ersetzt die unbekannte Matrix Hk−m durch H0 . Man spricht von einem Limited-Memory-Quasi-Newton-Verfahren. △

37

5. Verfahren der konjugierten Gradienten 5.1 CG-Verfahren für lineare Gleichungssysteme Es sei A ∈ Rn×n eine symmetrische und positive definite Matrix und b ∈ Rn . Das Verfahren der konjugierten Gradienten oder CG-Verfahren zur Lösung des linearen Gleichungssystems Ax = b geht davon aus, dass die Lösung x⋆ eindeutiges Minimum φ(x⋆ ) = 0 des Funktionals 1 1 1 φ(x) = (b − Ax)T A−1 (b − Ax) = xT Ax − xT b + bT A−1 b ≥ 0 2 2 2 ist. Ausgehend von einer Startnäherung x wollen wir φ in die Richtung d minimieren φ(x + αd) = φ(x) + Aus

α2 T d Ad − αdT (b − Ax) → min . α∈R 2

∂φ(x + αd) ! = αdT Ad − dT (b − Ax) = 0 ∂α

folgt daher α=

dT (b − Ax) . dT Ad

(5.1)

Lemma 5.1 Die Vektoren d0 , d1 , . . . , dk seien A-konjugiert, das heißt, es gelte dTi Adj = 0 für alle i 6= j. Ist xk = arg min φ(x) x∈x0 +span{d0 ,d1 ,...,dk−1 }

und setzt man xk+1 = xk + αk dk ,

αk =

dTk rk , dTk Adk

rk = b − Axk ,

(5.2)

so folgt xk+1 = arg

min

x∈x0 +span{d0 ,d1 ,...,dk }

φ(x).

Beweis. Die A-Konjugiertheit der Vektoren {dℓ } impliziert dTk A(xk − x0 ) = 0. Daher

38

Kapitel 5. Verfahren der konjugierten Gradienten

folgt α2 T φ(xk + αdk ) = φ(xk ) + dk Adk − αdTk (b − Axk ) 2 α2 T = φ(xk ) + dk Adk − αdTk (b − Ax0 ) 2 =: φ(xk ) + ϕ(α), das heißt, das Minimierungsproblem entkoppelt. Da nach Voraussetzung xk das Funktional φ über x0 + span{d0 , d1 . . . , dk−1} minimiert, wird das eindeutige Minimum angenommen, wenn ϕ(α) minimal ist. Dies ist aber nach (5.1) genau dann der Fall, wenn αk =

dTk (b − Ax0 ) dTk (b − Axk ) = . dTk Adk dTk Adk

(5.3)

Die Idee des CG-Verfahrens ist es nun, ausgehend von einer Startnäherung x0 , sukzessive über die konjugierten Richtungen dk zu minimieren. Die Folge der Residuen (5.2)

r0 = b − Ax0 ,

rk+1 = b − Axk+1 = rk − αk Adk ,

k ≥ 0,

(5.4)

erfüllt dann für alle ℓ < k dTℓ rk

=

dTℓ (b−Axk )

=

dTℓ

  k−1 X (5.3) b−Ax0 − αi Adi = dTℓ (b−Ax0 )−αℓ dTℓ Adℓ = 0. (5.5) i=0

Da die Richtungen dk paarweise A-konjugiert und folglich linear unabhängig sind, ergibt sich rn = 0, das heißt, das CG-Verfahren liefert die Lösung A−1 b nach höchstens n Schritten. Zu beantworten bleibt daher nur die Frage, wie die Suchrichtungen dk geschickt gewählt werden können. Lemma 5.2 Für beliebiges d0 = r0 erzeugt die Rekursion rk+1 = rk − αk Adk ,

dk+1 = rk+1 − βk dk ,

βk =

dTk Ark+1 dTk Adk

(5.6)

solange eine Folge nichtverschwindender A-konjugierter Vektoren d0 , d1 , . . . , dk+1 bis rk+1 = 0 ist. Beweis. Sei Kk (A, r0) := span{r0 , Ar0 , . . . , Ak−1 r0 }.

Wir zeigen zunächst induktiv, dass stets gilt

Kk (A, r0) = span{r0 , r1 , . . . , rk−1} = span{d0 , d1 , . . . , dk−1}. Da für k = 1 die Aussage klar ist, nehmen wir an, sie gilt für ein k ≥ 1. Dann folgt (5.6)

rk =

rk−1 −αk |{z}

∈Kk (A,r0 )

Adk |{z}

∈Kk+1 (A,r0 )

∈ Kk+1 (A, r0).

5.1. CG-Verfahren für lineare Gleichungssysteme

39

Gemäß (5.5) ist rk ⊥ span{d0 , d1 , . . . , dk−1 } = Kk (A, r0), dies bedeutet Kk (A, r0 ) ( span{r0 , r1 , . . . , rk } ⊆ Kk+1(A, r0 ). Da die Dimension von Kk+1 (A, r0) höchstens um 1 höher ist als die von Kk (A, r0 ), muss gelten Kk+1 (A, r0) = span{r0 , r1 , . . . , rk }. (5.6)

Aus rk = dk + βk−1 dk−1 folgt span{d0 , d1 , . . . , dk } = span{d0 , d1 , . . . , dk−1, rk } = span{r0 , r1 , . . . , rk−1, rk } = Kk+1 (A, r0). Insbesondere muss aus Dimensionsgründen dk 6= 0 sein. Es verbleibt, die A-Konjugiertheit zu zeigen: Angenommen, d0 , d1 , . . . , dk sind A-konjugiert. Der Induktionsschritt folgt dann aus (5.6)

(5.6)

dTk Adk+1 = dTk A(rk+1 − βk dk ) = dTk Ark+1 −

dTk Ark+1 T d Adk = 0 dTk Adk k

und für alle ℓ < k (5.6)

dTℓ Adk+1 = dTℓ Ark+1 −

dTk Ark+1 T d Adk = (Adℓ )T rk+1 = 0 dTk Adk | ℓ {z } =0

wegen Adℓ ∈ Kk+1 (A, r0) ⊥ rk+1 .

Bemerkung Der Raum Kk (A, r0) heißt Krylov-Raum. Die Iterierte xk des CG-Verfahrens minimiert demnach das Funktional φ(x) unter allen x aus dem verschobenen Krylov-Raum x0 + Kk (A, r0). △ Um den CG-Algorithmus endgültig zu formulieren, bemerken wir zunächst, dass gilt (5.3)

αk =

dTk rk (5.6) (rk − βk−1 dk−1 )T rk (5.5) krk k22 . = = T dTk Adk dTk Adk dk Adk

(5.7)

Wegen rk ∈ Kk+1 (A, r0) ⊥ rk+1 folgt ferner (5.6)

dTk Ark+1 = (Adk )T rk+1 = und damit (5.6)

βk =

1 dT Adk (5.7) (rk − rk+1 )T rk+1 = − k 2 krk+1k22 αk krk k2

krk+1 k22 dTk Ark+1 . = − krk k22 dTk Adk

Die Kombination von (5.2) und (5.6)–(5.8) liefert schließlich: Algorithmus 5.3 (CG-Verfahren) input: Matrix A ∈ Rn×n , rechte Seite b ∈ Rn und Startnäherung x0 ∈ Rn output: Folge von Iterierten {xk }k∈N ➀ Initialisierung: setze d0 = r0 := b − Ax0 und k := 0

(5.8)

40

Kapitel 5. Verfahren der konjugierten Gradienten

➁ berechne krk k22 dTk Adk xk+1 := xk + αk dk rk+1 := rk − αk Adk krk+1 k22 βk := krk k22 dk+1 := rk+1 + βk dk αk :=

➂ falls krk+1 k2 > ε erhöhe k := k + 1 und gehe nach ➁ Das CG-Verfahren wird generell als Iterationsverfahren verwendet, das heißt, man bricht die Iteration ab, falls die Residuennorm krk k2 kleiner als eine Fehlertoleranz ε ist. Pro Iterationsschritt wird nur eine Matrix-Vektor-Multiplikation benötigt. Allerdings hängt die Konvergenz des Verfahrens stark von der Kondition der Matrix ab. Man kann zeigen, dass die Iterierten {xk } des CG-Verfahrens bezüglich der Energienorm √ p kxkA := hx, xiA = xT Ax

der Fehlerabschätzung

k √ cond2 A − 1 kx0 − x⋆ kA kxk − x kA ≤ 2 √ cond2 A + 1 ⋆

genügen.

5.2 Nichtlineares CG-Verfahren In Anlehnung an das CG-Verfahren 5.3 ist das nichtlineare CG-Verfahren zur Lösung von nichtlinearen Optimierungsproblemen f (x) → min definiert. Algorithmus 5.4 (Nichtlineares CG-Verfahren) input: Funktion f : Rn → R und Startnäherung x0 ∈ Rn output: Folge von Iterierten {xk }k∈N ➀ Initialisierung: setze d0 = −∇f (x0 ) und k := 0 ➁ löse αk ≈ argmin f (xk + αdk ) α∈R

➂ berechne xk+1 := xk + αk dk βk := |

k∇f (xk+1)k22 k∇f (xk )k22 {z

}

oder

Verfahren von Fletcher und Reeves

βk := |

∇f (xk+1)T {∇f (xk+1 ) − ∇f (xk )} k∇f (xk )k22 {z } Verfahren von Polak und Ribière

dk+1 := −∇f (xk+1 ) + βk dk

5.2. Nichtlineares CG-Verfahren

41

➃ erhöhe k := k + 1 und gehe nach ➁ Bemerkung Ist f (x) = 21 xT Ax − bT x + c eine quadratische Funktion, dann fallen bei exakter Minimierung in ➁ sowohl das Verfahren von Fletcher und Reeves als auch das Verfahren von Polak und Ribière mit dem CG-Verfahren zusammen. Ersteres folgt aus ∇f (xk ) = Axk − b = −rk , zweiteres aus ∇f (xk )T ∇f (xk+1 ) = rTk rk+1 = 0. △ Lemma 5.5 Die Funktion f : D ⊂ Rn → R sei gleichmäßig konvex. Weiter sei f differenzierbar mit Lipschitz-stetigem Gradienten: k∇f (x) − ∇f (y)k2 ≤ Lkx − yk2

für alle x, y ∈ D.

Dann gilt für das Verfahren von Polak und Ribière bei exakter Liniensuche in ➁ −dTk ∇f (xk ) ≥

µ k∇f (xk )k2 kdk k2 . µ+L

Beweis. Bei exakter Liniensuche gilt im k-ten Schritt des Verfahrens von Polak und Ribière 0 = ∇f (xℓ + αℓ dℓ )T dℓ = ∇f (xℓ+1 )T dℓ

für alle ℓ ≤ k.

Daher können wir den Nenner in βk =

∇f (xk+1 )T {∇f (xk+1) − ∇f (xk )} k∇f (xk )k22

folgendemaßen umformen: k∇f (xk )k22 = (βk−1 dk−1 − dk )T ∇f (xk ) (5.9)

= −dTk ∇f (xk ) (5.9) T  = dk ∇f (xk+1 ) − ∇f (xk )  1 (xk+1 − xk )T ∇f (xk+1 ) − ∇f (xk ) . = αk

Aufgrund der gleichmäßigen Konvexität und der Lipschitz-Bedingung erhalten wir |βk | ≤

L k∇f (xk+1)k2 k∇f (xk+1)k2 kxk+1 − xk k2 L = . |αk | 2 µ kxk+1 − xk k2 µ kdk k2

Dies führt auf kdk+1 k2 ≤ k∇f (xk+1)k2 + |βk |kdk k2 ≤



 L 1+ k∇f (xk+1 )k2 , µ

(5.9)

42

Kapitel 5. Verfahren der konjugierten Gradienten

woraus dann die Behauptung folgt −

dTk+1 ∇f (xk+1 ) {∇f (xk+1) − βk dk }T ∇f (xk+1) = kdk+1 k2 k∇f (xk+1 )k2 kdk+1 k2 k∇f (xk+1)k2 k∇f (xk+1)k22 (5.9) = kdk+1 k2 k∇f (xk+1)k2 µ ≥ . µ+L

Bemerkung Die geometrische Interpretation von Lemma 5.5 ist, dass beim Verfahren von Polak und Ribière die Suchrichtung dk und die Richtung des steilsten Abstiegs −∇f (xk ) stets den Winkel θ mit cos θ > µ/(µ + L) einschließen. △ Satz 5.6 Die Funktion f : D ⊂ Rn → R sei gleichmäßig konvex. Weiter sei f differenzierbar mit Lipschitz-stetigem Gradienten: k∇f (x) − ∇f (y)k2 ≤ Lkx − yk2

für alle x, y ∈ D.

Dann konvergiert das Verfahren von Polak und Ribière mit exakter Liniensuche in ➁ für beliebige Startnäherungen x0 ∈ D gegen das eindeutige globale Minimum x⋆ und es gilt    µ4 ⋆ ⋆ f (x ) − f (x ) , k = 1, 2, . . . . f (xk+1) − f (x ) ≤ 1 − 2 k L (µ + L)2

Beweis. Der Beweis ist ähnlich zu dem von Satz 2.3. Aufgrund der Minimierungsbedingung gilt wieder für ein γ > 0 f (xk+1 ) ≤ f (xk + γdk ) Z 1 = f (xk ) + γ ∇f (xk + γtdk )T dk dt 0 Z 1  T T = f (xk ) + γ∇f (xk ) dk + γ ∇f (xk + γtdk ) − ∇f (xk ) dk dt 0 Z 1 ≤ f (xk ) + γ∇f (xk )T dk + γ k∇f (xk + γtdk ) + ∇f (xk )k2 kdk k2 dt {z } 0 | ≤tγLkdk k2

L ≤ f (xk ) + γ∇f (xk )T dk + γ 2 kdk k22 . 2

Für die Wahl γ := −

∇f (xk )T dk Lkdk k22

5.3. Modifiziertes Verfahren von Polak und Ribière

43

folgern wir mit Lemma 5.5 (∇f (xk )T dk )2 2Lkdk k22 µ2 k∇f (xk ) − ∇f (x⋆ ) k22 . ≤ f (xk ) − f (x⋆ ) − | {z } 2L(µ + L)2

f (xk+1) − f (x⋆ ) ≤ f (xk ) − f (x⋆ ) −

=0

Aus der gleichmäßigen Konvexität ergibt sich  T µkxk − x⋆ k22 ≤ ∇f (xk ) − ∇f (x⋆ ) (xk − x⋆ ) ≤ k∇f (xk ) − ∇f (x⋆ )k2 kxk − x⋆ k2 , während die Lipschitz-Stetigkeit impliziert Z 1
T L ⋆ f (xk ) − f (x ) = ∇f txk + (1 − t)x⋆ (xk − x⋆ ) dt ≤ kxk − x⋆ k22 . 2 0

Setzen wir diese beiden Abschätzungen in die obige ein, so erhalten wir das Behauptete.

Bemerkungen 1. Die Konvergenzabschätzung aus Satz 5.6 ist schlechter als die das Gradientenverfahren betreffende aus Satz 2.3. Dies spiegelt jedoch nicht die Tatsache wieder, dass das CG-Verfahren im Fall quadratischer Funktionen viel schneller konvergiert als das Gradientenverfahren. Die Ursache liegt in der Beweistechnik begründet, die das Verfahren von Polak und Ribière als Störung des Gradientenverfahrens auffasst. 2. Das Verfahren von Polak und Ribière konvergiert im allgemeinen schneller als das Verfahren von Fletcher und Reeves. 3. In der Praxis verwendet man Restarts: Wird der Winkel zwischen dem Antigradienten und der Suchrichtung zu groß, etwa −

∇f (xk )T dk 0 sofort ∇f (xk ) = 0 folgen. Es ist also nur zu zeigen, dass die Liniensuche ➁–➃ in jedem Iterationschritt k ∈ N0 erfolgreich ist. Zu diesem Zweck nehmen wir an, dass k ∈ N0 ein fester Iterationsindex mit ∇f (xk )T dk < 0 ist. Als erstes stellen wir fest, dass die Bedingung (5.10) wegen f (xk + αdk ) = f (xk ) + α∇f (xk )T dk + O(α) nach endlich vielen erfolglosen Schritten der Liniensuche immer erfüllt ist. Als nächstes zeigen wir, dass die Bedingung (5.11) ebenfalls nach endlich vielen erfolglosen Schritten stets erfüllt ist. Denn angenommen, dem ist nicht so. Dann gibt es eine Teilfolge {kℓ }ℓ>0 , so dass für jedes yℓ = xk + 2

−kℓ |∇f (xk )

T

kdk k22

dk |

dk ,

ℓ∈N

zumindest eine der beiden Bedingungen   ∇f (yℓ )T {∇f (yℓ ) − ∇f (xk )} T − ∇f (yℓ ) + dk > −γk∇f (yℓ )k22 , ∇f (yℓ ) 2 k∇f (xk )k2   T ∇f (yℓ ) {∇f (yℓ ) − ∇f (xk )} T − ∇f (yℓ ) + dk < −γk∇f (yℓ )k22 ∇f (yℓ ) k∇f (xk )k22

5.3. Modifiziertes Verfahren von Polak und Ribière

45

nicht erfüllt ist. Der Grenzübergang ℓ → ∞ liefert yℓ → xk und folglich gilt −k∇f (xk )k22 ≥ −γk∇f (xk )k22

oder

− k∇f (xk )k22 ≤ −γk∇f (xk )k22 .

Aus 0 < γ < 1 < γ folgt dann aber k∇f (xk )k2 = 0 im Widerspruch zu unserer Voraussetzung ∇f (xk )T dk < 0. Damit ist gezeigt, dass Algorithmus 5.7 wohldefiniert ist, sofern die Abstiegsbedingung ∇f (xk )T dk < 0 für alle k ∈ N0 erfüllt ist. Für k = 0 gilt sie aber nach Definition von d0 und für k > 0 folgt sie dann aus Bedingung (5.11). Lemma 5.9 Es sei D ⊂ Rn eine offene, beschränkte und konvexe Menge, in der f stetig differenzierbar, nach unten beschränkt und ∇f zudem Lipschitz-stetig ist. Ferner sei neben x0 auch die gesamte Niveaumenge N := {x ∈ Rn : f (x) ≤ f (x0 )} in D enthalten. Dann gelten die folgenden Aussagen: (i.) Alle Iterierten xk liegen in der Niveaumenge N. (ii.) Die Folge {f (xk )}k∈N ist konvergent. (iii.) Es gilt limk→∞ αk kdk k2 = 0. (iv.) Es ist αk kdk k22 ≤ γc2 , wobei c < ∞ eine obere Schranke von k∇f (x)k2 auf der Niveaumenge N sei. (v.) Es existiert eine Konstante θ > 0 mit |∇f (xk )T dk | αk ≥ θ kdk k22 für alle k ∈ N0 . Beweis. (i.) Diese Aussage ergibt sich unmittelbar aus der Bedingung (5.10). (ii.) Die Folge {f (xk )}k∈N ist streng monoton fallend und aufgrund der Voraussetzung nach unten beschränkt. Hieraus ergibt sich das Behauptete. (iii.) Aus (5.10) folgt σαk2 kdk k22 ≤ f (xk ) − f (xk+1 ) für alle k ∈ N0 . Der Grenzübergang k → ∞ liefert daher unter Berücksichtigung der schon bewiesenen Aussage (ii.) die Behauptung. (iv.) Aus den Schritten ➁–➃ von Algorithmus 5.7 folgt αk kdk k22 ≤

|∇f (xk )T dk | kdk k22 ≤ γk∇f (xk )k22 ≤ γc2 kdk k22

für alle k ∈ N0 . (v.) Zum Nachweis dieser Aussage führen wir eine Fallunterscheidung durch. Fall 1: αk = |∇f (xk )T dk |/kdk k22 . Dann ist offensichtlich |∇f (xk )T dk | αk ≥ . kdk k22 Fall 2: αk < |∇f (xk )T dk |/kdk k22 .

(5.12)

46

Kapitel 5. Verfahren der konjugierten Gradienten Dann verletzt die Schrittweite 2αk zumindest eine der Bedingungen in ➃. Der Punkt zk := xk + 2αk dk genügt also (5.10) oder (5.11) nicht. Nach Aussage (iii.) existiert ein K ∈ N, so dass zk ∈ D für alle k ≥ K. Im folgenden zeigen wir Aussage (v.) zunächst nur für solche k. Fall 2A: Der Punkt zk ∈ D verletzt (5.10). Dann gilt f (zk ) > f (xk ) − σ(2αk )2 kdk k22 . (5.13)

Aufgrund des Mittelwertsatzes existiert ein ξ k auf der Verbindungsstrecke von xk und zk , so dass f (zk ) = f (xk ) + ∇f (ξ k )T (zk − xk ) = f (xk ) + 2αk ∇f (ξk )T dk .

(5.14)

Aus (5.13) und (5.14) folgt daher  f (xk ) + 2αk ∇f (xk )T dk + 2αk ∇f (ξk )T dk − ∇f (xk )T dk > f (xk ) − σ(2αk )2 kdk k22 . Aus der Lipschitz-Stetigkeit von ∇f in D ergibt sich

2αk ∇f (xk )T dk + 2αk L kξ k − xk k2 kdk k2 > −σ(2αk )2 kdk k22 . | {z } =2αk kdk k2

Dies liefert unmittelbar

αk ≥

1 |∇f (xk )T dk | . 2(L + σ) kdk k22

(5.15)

Fall 2B: Der Punkt zk ∈ D verletzt die linke Ungleichung in (5.11). Wegen   ∇f (zk )T {∇f (zk ) − ∇f (xk )} T − ∇f (zk ) + ∇f (zk ) dk < −γk∇f (zk )k22 2 k∇f (xk )k2 ergibt sich mit Hilfe der Cauchy-Schwarzschen Ungleichung −k∇f (zk )k22 − k∇f (zk )k22

k∇f (zk ) − ∇f (xk )k2 kdk k2 < −γk∇f (zk )k22 , k∇f (xk )k22

das heißt, es ist 1+

k∇f (zk ) − ∇f (xk )k2 kdk k2 > γ. k∇f (xk )k22

Aus der Lipschitz-Stetigkeit von ∇f und der Tatsache, dass die Schrittweite αk der Bedingung (5.11) genügt, folgt daher nach kurzer Rechnung αk ≥

(γ − 1)k∇f (xk )k22 γ − 1 |∇f (xk )T dk | ≥ . 2Lkdk k22 2γL kdk k22

(5.16)

Fall 2C: Der Punkt zk ∈ D verletzt die linke Ungleichung in (5.11). Es ist   ∇f (zk )T {∇f (zk ) − ∇f (xk )} T dk > −γk∇f (zk )k22 . − ∇f (zk ) + ∇f (zk ) k∇f (xk )k22

5.3. Modifiziertes Verfahren von Polak und Ribière

47

Analog zum Fall 2B erhält man hieraus αk ≥

1 − γ |∇f (xk )T dk | . 2γL kdk k22

(5.17)

Wegen (5.12), (5.15), (5.16), (5.17) folgt Aussage (v.) mit   γ−1 1−γ 1 , , , θ := min 1, 2(L + σ) 2γL 2γL

und zwar zunächst für alle k ≥ K. Da nur endlich viele k übrigbleiben, folgt Aussage (v.) nach eventueller Verkleinerung von θ aber auch für alle k ∈ N0 .

Satz 5.10 Unter den Voraussetzungen von Lemma 5.9 gilt für die Iterierten {xk }k≥0 des modifizierten Verfahren von Polak und Ribière 5.7 k→∞

∇f (xk ) −→ 0. Beweis. Angenommen, die Aussage des Satzes ist falsch. Dann existieren ein ε > 0 und eine Teilfolge {xkℓ }, so dass k∇f (xkℓ −1 )k2 > ε für alle ℓ ∈ N. Aus der Aufdatierungsvorschrift für dkℓ und Lemma 5.9 (iv.) folgt dann kdkℓ k2 ≤ k∇f (xkℓ )k2 +

k∇f (xkℓ )k2 k∇f (xkℓ ) − ∇f (xkℓ −1 )k2 Lc3 kd k ≤ c + γ k −1 2 ℓ k∇f (xkℓ −1 )k22 ε2

für alle ℓ ∈ N. Zusammen mit Lemma 5.9 (iii.) ergibt sich hieraus lim αkℓ kdkℓ k22 = 0

ℓ→∞

und weiter aus Lemma 5.9 (v.) lim |∇f (xkℓ )T dkℓ | = 0.

ℓ→∞

Die rechte Ungleichung in (5.11) liefert daher lim k∇f (xkℓ )k2 = 0.

ℓ→∞

Weil nach Lemma 5.9 (iii.) gilt lim kxkℓ − xkℓ −1 k2 = lim αkℓ −1 kdkℓ −1 k2 = 0,

ℓ→∞

ℓ→∞

schließen wir k∇f (xkℓ −1 )k2 ≤ k∇f (xkℓ ) − ∇f (xkℓ −1 )k2 + k∇f (xkℓ )k2 ≤ Lkxkℓ − xkℓ −1 k2 + k∇f (xkℓ )k2 ℓ→∞

−→ 0.

Dies steht aber im Widerspruch zur Annahme, dass {k∇f (xkℓ −1 )k2 }ℓ>0 nicht nach Null konvergiert.

48

Kapitel 6. Nichtlineare Ausgleichsprobleme

6. Nichtlineare Ausgleichsprobleme 6.1 Gauß-Newton-Verfahren Ein nichtlineares Ausgleichsproblem liegt vor, falls zu gegebenen Daten und Funktionen     f1 (x1 , x2 , . . . , xn ) y1  f2 (x1 , x2 , . . . , xn )   y2      f(x) =  y =  ..  ∈ Rm ,  : Rn → Rm ..    .  . fm (x1 , x2 , . . . , xn ) ym dasjenige x⋆ = [x⋆1 , x⋆2 , . . . , x⋆n ]T gesucht wird, das die Minimierungsaufgabe m

X 1 φ(x) := ky − f(x)k22 = |yi − fi (x1 , x2 , . . . , xn )|2 → minn x∈R 2 i=1

(6.1)

löst. Nichtlineare Ausgleichsprobleme können im allgemeinen nur iterativ gelöst werden. Da hierzu Gradienteninformationen benötigt wird, setzen wir f als stetig differenzierbar voraus. Die Ableitung f ′ sei zusätzlich sogar Lipschitz-stetig. Natürlich kann man das nichtlineare Ausgleichsproblem (6.1) mit dem Gradientenverfahren oder dem Quasi-Newton-Verfahren zu lösen. Besser ist es jedoch, es mit dem NewtonVerfahren zu lösen, was auf die Iteration xk+1 = xk − φ′′ (xk )∇φ(xk ),

k = 0, 1, 2, . . .

führt. Hierzu wird jedoch neben dem Gradienten

∇φ(x) = − f ′ (x)


T


y − f(x) ,

auch die Hesse-Matrix benötigt, das ist φ′′ (x) = f ′ (x)


T

 ∂f1

(x)

∂x1  ∂f2 (x)  ∂x1

f ′ (x) =  

.. . ∂fm (x) ∂x1

f ′ (x) − y − f(x)


T

∂f1 (x) ∂x2 ∂f2 (x) ∂x2

... ...

.. . ∂fm (x) . . . ∂x2

 ∂f1 (x) ∂xn ∂f2 (x)   ∂xn

..  .  ∂fm (x) ∂xn

f ′′ (x).

In der Praxis will man die Berechnung des Tensors f ′′ (x) ∈ Rm×(n×n) jedoch vermei
T den. Daher vernachlässigt man den Term y − f(x) f ′′ (x) und erhält das Gauß-NewtonVerfahren.

6.1. Gauß-Newton-Verfahren

49

Zu dessen Herleitung linearisieren wir die Funktion f f(x + d) = f(x) + f ′ (x)d + O(kdk2 ). Ist nun xk eine Näherungslösung des Ausgleichsproblems (6.1), so erwartet man, dass die Optimallösung xk+1 = xk + dk des linearisierten Problems min ky − f(xk ) − f ′ (xk )dk22 = krk − f ′ (xk )dk k22 ,

d∈Rn

rk := y − f(xk )

eine bessere Lösung des Ausgleichsproblems ist. Gemäß Definition muss das Update dk die Normalengleichungen
T
T f ′ (xk ) f ′ (xk )dk = f ′ (xk ) rk

lösen. Dies führt auf folgenden Algorithmus:

Algorithmus 6.1 (Gauß-Newton–Verfahren) input: Funktion f : Rn → Rm , Datenvektor y ∈ Rm und Startnäherung x0 ∈ Rn output: Folge von Iterierten {xk }k∈N ➀ Initialisierung: setze k = 0 ➁ löse die Normalengleichungen f ′ (xk )


T

f ′ (xk )dk = f ′ (xk )

➂ setze xk+1 := xk + dk ➃ erhöhe k := k + 1 und gehe nach ➁


T

y − f(xk )




(6.2)

Satz 6.2 Sei D ⊂ Rn offen und f : D → Rm eine stetig differenzierbare
Abbildung. Das Minimierungsproblem (6.1) habe eine Lösung x⋆ ∈ D mit rang f ′ (x⋆ ) = n ≤ m. Sei
T λ > 0 der kleinste Eigenwert von f ′ (x⋆ ) f ′ (x⋆ ). Ferner gelte die Lipschitz-Bedingung kf ′ (x) − f ′ (z)k2 ≤ αkx − zk2

und

′



f (x) − f ′ (x⋆ ) T y − f(x⋆ ) ≤ βkx − x⋆ k2 2

(6.3) (6.4)

mit β < λ für alle x, z aus einer Umgebung von x⋆ . Dann existiert ein ε > 0, so dass für jeden Startvektor x0 ∈ Bε (x⋆ ) := {x ∈ Rn : kx − x⋆ k2 < ε} die Folge der Iterierten {xk }k≥0 mindestens linear gegen x⋆ konvergiert. Beweis. Nach Voraussetzung gibt es ein ε1 > 0, so dass (6.3) und (6.4) für alle x ∈ Bε1 (x⋆ ) gelten. Aus Stetigkeitsgründen folgt außerdem die Existenz von γ > 0, so dass kf ′ (x)k2 ≤ γ

für alle x ∈ Bε1 (x⋆ ).

T Wegen rang f ′ (x⋆ ) = n, ist die Matrix f ′ (x⋆ ) f ′ (x⋆ ) regulär mit



′ ⋆
T ′ ⋆ −1

= 1.

f (x ) f (x )

λ 2

(6.5)

50

Kapitel 6. Nichtlineare Ausgleichsprobleme


T Daher gibt es zu beliebigem δ > 1 ein ε2 > 0 derart, dass f ′ (x) f ′ (x) regulär ist und



′
T ′ −1

f (x) f (x)

≤ δ für alle x ∈ Bε2 (x⋆ ). (6.6)

λ 2 Wir wählen δ > 1 so, dass zusätzlich gilt

δ
0, so dass für jedes xk ∈ Bε3 (x⋆ ) folgt

kf(x⋆ ) − f(xk ) − f ′ (xk )(x⋆ − xk )k2

Z 1


′ ⋆ ⋆ ′ ⋆

= f x + t(x − x ) (x − x ) dt − f (x )(x − x ) k k k k k

0 2

Z 1




= f ′ xk + t(x⋆ − xk ) − f ′ (xk ) (x⋆ − xk ) dt

0 2 Z 1




′ ≤

f xk + t(x⋆ − xk ) − f ′ (xk ) kxk − x⋆ k2 dt 0 | {z }2 ≤αtkxk −x⋆ k2

≤

Wir setzen nun

α kxk − x⋆ k22 . 2



λ − βδ ε := min ε1 , ε2 , ε3 , αγδ ⋆ Aus (6.2) folgt dann für xk ∈ Bε (x ) dass



> 0.

xk+1 − x⋆ = xk + dk − x⋆  −1 h i
T
T

T = f ′ (xk ) f ′ (xk ) f ′ (xk ) y − f(xk ) − f ′ (xk ) f ′ (xk )(x⋆ − xk )  −1 h
T ′
T
′ = f (xk ) f (xk ) f ′ (xk ) y − f(x⋆ )
T
i + f ′ (xk ) f(x⋆ ) − f(xk ) − f ′ (xk )(x⋆ − xk ) .

Hieraus folgt dann mit (6.5) und (6.6)

 −1

h ′


T ⋆ ′ ′

f (xk ) T y − f(x⋆ ) kxk+1 − x k2 ≤ f (xk ) f (xk )

2 2

i + kf (xk )k2 kf(x⋆ ) − f(xk ) − f ′ (xk )(x⋆ − xk )k2  
T
αγ δ ′ ⋆ ⋆ 2

f (xk ) y − f(x ) + (6.8) kxk − x k2 . ≤ 2 λ 2
T
Wegen 0 = ∇φ(x⋆ ) = − f ′ (x⋆ ) y − f(x⋆ ) ergibt sich mit (6.4)

′





f (xk ) T y − f(x⋆ ) = f ′ (xk ) − f ′ (x⋆ ) T y − f(x⋆ ) ≤ βkxk − x⋆ k2 , 2 2 ′

dies bedeutet,

    δ αγ βδ λ − βδ ⋆ ⋆ kxk+1 − x k2 ≤ β+ kxk − x k2 kxk − x k2 ≤ + kxk − x⋆ k2 . | {z } λ 2 |λ {z 2λ } ≤ε≤(λ−βδ)/(αγδ) ⋆

=(λ+βδ)/(2λ)

6.2. Levenberg-Marquardt-Verfahren

51

Wegen (6.7) ist die Konstante (λ + βδ)/(2λ) < 1, das heißt, alle Iterierten {xk }k≥0 liegen in Bε (x⋆ ) und konvergieren mindestens linear gegen x⋆ . Bemerkung Die Voraussetzung (6.4) besagt im Prinzip, dass das Residuum r(x⋆ ) = y − f(x⋆ ) klein genug sein soll. Dies sieht man insbesondere, wenn man sie durch die stärkere Bedingung λ ky − f(x⋆ )k2 ≤ α ersetzt. Dann folgt nämlich (6.4):

′



f (x) − f ′ (x⋆ ) T y − f(x⋆ ) ≤ kf ′ (x) − f ′ (x⋆ )k2 ky − f(x⋆ )k2 ≤ λkx − x⋆ k2 . 2 {z }| {z } | ≤αkx−x⋆ k2

≤λ/α

△

Korollar 6.3 Zusätzlich zu den Voraussetzungen aus Satz 6.2 gelte f(x⋆ ) = y. Dann existiert ein ε > 0, so dass für jeden Startvektor x0 ∈ Bε (x⋆ ) die Folge der Iterierten {xk }k≥0 quadratisch gegen x⋆ konvergiert. Beweis. Aus Satz 6.2 folgt die lineare Konvergenz der Iterierten {xk }k≥0 gegen x⋆ . Zum Nachweis der quadratischen Konvergenz bemerken wir, dass aufgrund der Voraussetzung (6.4) mit β = 0 gilt. Daher folgt aus (6.8) kxk+1 − x⋆ k2 ≤

αγδ kxk − x⋆ k22 . 2λ

6.2 Levenberg-Marquardt-Verfahren Das Levenberg-Marquardt-Verfahren ist ein Trust-Region-Verfahren, also ein Verfahren, bei dem der Linearisierung nur im Bereich kdk k2 ≤ ∆k vertraut wird. Demnach wollen wir das restringierte Optimierungsproblem min

d∈Rn :kdk2 ≤∆k

ψ(d) =

min

d∈Rn :kdk2 ≤∆k

1 krk − f ′ (xk )dk22 , 2

(6.9)

lösen, um die Iterierte dann mit der Lösung dk aufzudatieren: xk+1 = xk + dk . Da die Menge B∆k (0) kompakt ist, existiert ein Minimum dk . Es treten dabei zwei Fälle auf: 1. Das Minimum dk liegt im Inneren der Kugel B∆k (0) und erfüllt
T
∇ψ(dk ) = f ′ (xk ) f ′ (xk )dk − rk = 0. (6.10) 2. Das Minimum liegt auf dem Rand, erfüllt also kdk k2 = ∆k . In diesem Fall muss die Höhenlinie von ψ in dk genau den Kreis ∂B∆k (0) tangieren, das heißt, der Gradient ∇ψ(dk ) zeigt in Richtung des Nullpunkts:
T
∇ψ(dk ) = f ′ (xk ) f ′ (xk )dk − rk = −λk dk für ein λk ≥ 0. (6.11)

52

Kapitel 6. Nichtlineare Ausgleichsprobleme

Der Parameter λk wird Lagrange-Parameter genannt. Da die Gleichung (6.10) als Spezialfall λk = 0 von (6.11) angesehen werden kann, erhalten wir: Lemma 6.4 Die Lösung dk des restringierten Problems (6.9) genügt der Gleichung  
T ′
T ′ f (xk ) f (xk ) + λk I dk = f ′ (xk ) rk (6.12)

für ein λk ≥ 0. Dabei ist der Wert λk genau dann positiv, wenn ∇ψ(dk ) 6= 0 und daher die Nebenbedingung kdk k2 = ∆k > 0 gilt. Der Vorteil des regularisierten Systems (6.12) liegt darin, dass es stets eindeutig lösbar ist, sofern λk > 0 ist. Die zugehörige Lösung dk =



f ′ (xk )


T

f ′ (xk ) + λk I

−1

f ′ (xk )


T

 −1
T rk = − f ′ (xk ) f ′ (xk ) + λk I ∇φ(xk )

erfüllt offenbar die Abstiegsbedingung dTk ∇φ(xk ) < 0, es sei denn, der Punkt xk ist stationär. Bemerkung Ist λk = 0, dann ergibt sich ein Gauß-Newton-Schritt, während dk für λk → ∞ der Richtung des steilsten Abstiegs entspricht. △ Wir müssen uns noch ein Kriterium überlegen, wie wir ∆k wählen. Eine neue Näherung xk+1 = xk + dk können wir anhand der Armijo-Goldstein-Bedingung (vergleiche (2.4)) bewerten: ρ=

φ(xk + dk ) − φ(xk ) 1 ky − f(xk )k22 − ky − f(xk + dk )k22 = .
T
dTk ∇φ(xk ) 2 dT f ′ (xk ) y − f(xk ) k

Wir wählen zwei Toleranzgrenzen 0 < ρ− < ρ+ < 1 und akzeptieren den Iterationsschritt, wenn ρ > ρ− gilt. Dann war der Trust-Region-Radius ∆k geeignet gewählt. Ist ρ ≥ ρ+ , so können wir ∆k sogar vergrößern. Ist hingegen ρ ≤ ρ− , dann verwerfen wir den Iterationsschritt und verkleinern ∆k . Dies führt auf den folgenden Algorithmus: Algorithmus 6.5 (Levenberg-Marquardt-Verfahren) input: Funktion f : Rn → Rm , Datenvektor y ∈ Rm und Startnäherung x0 ∈ Rn output: Folge von Iterierten {xk }k∈N ➀ Initialisierung: wähle 0 < ρ− < ρ+ < 1 und setze ∆0 := 1 und k := 0 ➁ bestimme die Lösung dk des restringierten Optimierungsproblems (6.9) ➂ berechne 1 ky − f(xk )k22 − ky − f(xk + dk )k22 ρk = (6.13)
T
2 dT f ′ (xk ) y − f(xk ) k

➃ falls ρk > ρ− setze xk+1 := xk + dk , sonst setze ∆k := ∆k /2 und gehe nach ➁ ➄ falls ρk > ρ+ setze ∆k+1 := 2∆k , sonst setze ∆k+1 := ∆k ➅ erhöhe k := k + 1 und gehe nach ➁

6.2. Levenberg-Marquardt-Verfahren

53

Satz 6.6 Es sei D ⊂ Rn eine kompakte Menge, in der f stetig differenzierbar und f ′ zudem Lipschitz-stetig ist. Ferner sei neben x0 auch die gesamte Niveaumenge {x ∈ Rn : φ(x) ≤ φ(x0 )} in D enthalten. Dann gilt für die Iterierten {xk }k≥0 aus Algorithmus 6.5 k→∞

∇φ(xk ) −→ 0. Beweis. (i) Zunächst beweisen wir eine obere Schranke für den Lagrange-Parameter λk aus (6.12). Dazu nehmen wir ohne Beschränkung der Allgemeinheit an, dass λk > 0 und daher kdk k2 = ∆k ist. Aus (6.12) folgt  


T
T
T dTk f ′ (xk ) f ′ (xk ) + λk I dk = dTk f ′ (xk ) rk ≤ kdk k2 f ′ (xk ) rk 2 .


T Da f ′ (xk ) f ′ (xk ) positiv semidefinit ist, kann die linke Seite nach unten durch λk abgeschätzt werden, so dass folgt

′


f (xk ) T rk k∇φ(xk )k2 2 = . (6.14) λk ≤ ∆k ∆k (ii) Als nächstes leiten wir eine untere Schranke für den Nenner νk in (6.13) her. Im Fall
T λk > 0 ist f ′ (xk ) f ′ (xk ) + λk I positiv definit, weshalb eine Cholesky-Zerlegung LLT existiert. Dabei gilt (6.14)


T k∇φ(xk )k2 . kLT Lk2 = kLLT k2 = f ′ (xk ) f ′ (xk ) + λk I 2 = kf ′ (xk )k22 +λk ≤ c + | {z } ∆k ≤c für alle x∈D

Setzen wir w = L−1 ∇φ(xk ), so folgt unter Beachtung von (6.12) hieraus
T

k∇φ(xk )k22 kwk22 k∇φ(xk )k22 =
T wT LT Lw ∇φ(xk ) ∇φ(xk ) kwk22 k∇φ(xk )k22 k∇φ(xk )k2
≥ ≥ min{∆k , k∇φ(xk )k2 }. (6.15) 2 1+c kwk2 c + k∇φ(xk )k2 /∆k

νk = ∇φ(xk )

Im Fall λk = 0 folgt

(LLT )−1 ∇φ(xk ) = wT w

νk = rTk f ′ (xk ) f ′ (xk )


+

rk = kPk rk k22 ,


+
wobei Pk = f ′ (xk ) f ′ (xk ) den Orthogonalprojektor auf img f ′ (xk ) bezeichnet. Wegen 

T ⊥ ′ ′ img f (xk ) = kern f (xk ) folgt daher

2


T 2
T
T 2 k∇φ(xk )k22 = f ′ (xk ) rk 2 = f ′ (xk ) Pk rk 2 ≤ f ′ (xk ) 2 kPk rk k22 ≤ cνk ,

das heißt, (6.15) ist auch im Fall λk = 0 gültig. (iii) Wir beweisen nun, dass die rechte Seite von (6.15) gegen Null konvergiert. Bei erfolgreichem Iterationsschritt ist ρk > ρ− und aus (6.13) und (6.15) folgt φ(xk ) − φ(xk+1 ) > ρ− νk ≥

ρ− k∇φ(xk )k2 min{∆k , k∇φ(xk )k2 }. 1+c

(6.16)

54

Kapitel 6. Nichtlineare Ausgleichsprobleme

Da nach Konstruktion {φ(xk )}k≥0 eine monoton fallende, nach unten beschränkte Folge ist, muss gelten k→∞ min{∆k , k∇φ(xk )k2 } −→ 0. (6.17) (iv) Als nächstes zeigen wir, dass k∇φ(xk )k2 für eine Teilfolge {kn }n∈N mit kn → ∞ gegen Null konvergiert. Angenommen, die Behauptung gilt nicht, dann folgt k∇φ(xk )k2 ≥ ε > 0 für alle k ≥ K(ε). Aus (6.17) ergibt sich damit unmittelbar k→∞

(6.18)

∆k −→ 0. Taylor-Entwicklung von ρk liefert jedoch φ(xk+1 ) − φ(xk ) dTk ∇φ(xk ) + O(kdk k22 ) = dTk ∇φ(xk ) dT ∇φ(xk )   k  2 ∆k ∆k (6.15) = 1 + O(∆k ) für k → ∞. = 1+O =1+O νk k∇φ(xk )k2 | {z }

ρk =

≥ε>0

Demnach existiert ein M(ε) ≥ K(ε), so dass ρk > ρ+ für alle k ≥ M(ε). Ab dem M(ε)-ten Schritt wird folglich ∆k in jedem Schritt von Algorithmus 6.5 verdoppelt, was jedoch im Widerspruch zu (6.18) steht. (v) Wir beweisen nun die Aussage des Satzes. Dazu nehmen wir an, dass eine Teilfolge von {k∇φ(xk )k2 }k≥0 nicht gegen Null konvergiert. Nach Aussage (iv) existiert dann ein ε > 0 und zwei Indizes ℓ < m, so dass k∇φ(xℓ )k2 ≥ 2ε,

k∇φ(xm )k2 ≤ ε,

k∇φ(xk )k2 > ε,

k = ℓ + 1, . . . , m − 1.

Da {φ(xk )}k≥0 eine Cauchy-Folge ist, kann ℓ dabei so groß gewählt werden, dass φ(xℓ ) − φ(xm )
1 eine Lipschitz-Konstante von ∇φ in D bezeichne. Wegen kxk+1 − xk k2 ≤ ∆k folgt aus (6.16) dass φ(xk ) − φ(xk+1) ≥

ερ− min{kxk+1 − xk k2 , ε}, 1+c

k = ℓ, ℓ + 1, . . . , m − 1.

Summation ergibt m−1 (6.19) ερ− X ε2 ρ− min{kxk+1 − xk k2 , ε} ≤ φ(xℓ ) − φ(xm ) < , 1 + c k=ℓ (1 + c)L

was wegen L > 1 nur erfüllt sein kann, wenn

min{kxk+1 − xk k2 , ε} = kxk+1 − xk k2 ,

k = ℓ, ℓ + 1, . . . , m − 1,

und insgesamt m−1 X k=ℓ

kxk+1 − xk k2
∆k , so wissen wir, dass das Minimum von (6.9) auf dem Rand liegt. Wir suchen dann dasjenige Tupel (λk , dk ), das (6.12) und kdk k2 = ∆k löst.
T Es bezeichne z1 ≥ · · · ≥ zn ≥ 0 die Eigenwerte von f ′ (xk ) f ′ (xk ) und {vi }ni=1 die zugehörigen orthonormalen Eigenvektoren. Entwickeln wir die rechte Seite von (6.12) in diese Eigenbasis n X
T ′ f (xk ) rk = ξi vi , i=1

dann folgt

dk (λk ) =



′

f (xk )


T

′

f (xk ) + λk I

n −1 X

ξi vi =

i=1

n X i=1

ξi vi . zi + λk

Die Forderung kdk (λk )k2 = ∆k führt auf die nichtlineare Gleichung r(λk ) :=

n X i=1

|ξi |2 ! = ∆2k . |zi + λk |2

Diese kann mit dem Hebden-Verfahren gelöst werden, einem Newton-Verfahren für die Gleichung 1 1 ! p − = 0. r(λ) ∆k Ausgehend vom Startwert λ(0) = 0 konvergiert die zugehörige Iteration   1 r 3/2 (λ(i) ) −1/2 (i) (i+1) (i) r (λ ) − , i = 0, 1, 2, . . . . λ = λ + 2 ′ (i) r (λ ) ∆k

sehr schnell gegen die Lösung λk . Die explizite Spektralzerlegung kann vermieden werden, indem man r(λ) = kdk (λ)k22 und  −3
T
T r ′ (λ) = −2rTk f ′ (xk ) f ′ (xk ) f ′ (xk ) + λI f ′ (xk ) rk = −2dk (λ)T g(λ)  
T ′ ′ mit f (xk ) f (xk ) + λI g(λ) = dk (λ)

benutzt. Für jede Hebden-Iterierte sind zwei Gleichungssysteme mit derselben Systemmatrix zu lösen. Diese entsprechen genau den Normalengleichungen zu den Ausgleichsproblemen

 ′

 ′    2    2

f (xk )

f (xk ) 0 √ rk

→ min .

√ √ → min, dk (λ) − g(λ) −

0 2 λI λI dk (λ)/ λ 2

Verwendet man die QR-Zerlegung QR = f ′ (xk ), so können letztere durch Anwendung von jeweils n(n + 1)/2 Givens-Rotationen effizient gelöst werden.

56

Kapitel 7. Optimierungsprobleme mit Nebenbedingungen

7. Optimierungsprobleme mit Nebenbedingungen 7.1 Optimalitätsbedingungen erster Ordnung Wir betrachten im folgenden Optimierungsprobleme mit Nebenbedingungen, das heißt, Probleme der Form min f (x) unter den Nebenbedingungen

x∈Rn

ci (x) = 0 für alle i ∈ JG ,

(7.1)

ci (x) ≥ 0 für alle i ∈ JU .

Dabei seien JG und JU Indexmengen für Gleichungs- und Ungleichungsnebenbedingungen. Definition 7.1 Es bezeichne G := {x ∈ Rn : ci (x) = 0 für i ∈ JG und ci (x) ≥ 0 für i ∈ JU } ⊂ Rn die Menge der zulässigen Punkte des Minimierungsproblems unter Nebenbedingungen (7.1). Dann verstehen wir unter einem lokalen Minimum des Problems (7.1) einen Punkt x⋆ ∈ G, für den für alle x ∈ U ∩ G gilt f (x⋆ ) ≤ f (x) mit einer Umgebung U ⊂ Rn von x⋆ . Als generelle Voraussetzung für das folgende seien sowohl f als auch ci für alle i ∈ JG ∪JU stetig differenzierbar im betrachteten Bereich. Definition 7.2 Es sei x⋆ ∈ G ein lokales Minimum des Minimierungsproblems (7.1). Die Folge {xk }k∈N ⊂ Rn heißt zulässige Folge für das Problem (7.1), falls folgende Eigenschaften erfüllt sind: (i.) Es gilt limk→∞ xk = x⋆ ∈ G. (ii.) Es ist xk 6= x⋆ für alle k ∈ N. (iii.) Es existiert ein K ∈ N, so dass xk ∈ G für alle k ≥ K. Die Richtung d ∈ Rn heißt Grenzrichtung einer zulässigen Folge, falls xkℓ − x⋆ =d ℓ→∞ kxkℓ − x⋆ k2 lim

(7.2)

7.1. Optimalitätsbedingungen erster Ordnung

57

für eine Teilfolge {xkℓ }ℓ∈N gilt.

Satz 7.3 Ist x⋆ ∈ G ein lokales Minimum des Minimierungsproblems (7.1), so gilt ∇f (x⋆ )T d ≥ 0 für alle Grenzrichtungen d ∈ Rn zu zulässigen Folgen mit Grenzwert x⋆ . Beweis. Angenommen, es gibt eine zulässige Folge {xk }k∈N mit ∇f (x⋆ )T d < 0 für eine Grenzrichtung d. Es sei {xkℓ }ℓ∈N eine zur Grenzrichtung gehörende Teilfolge mit (7.2). Dann gilt f (xkℓ ) = f (x⋆ ) + ∇f (x⋆ )T (xkℓ − x⋆ ) + O(kxkℓ − x⋆ k2 )

= f (x⋆ ) + kxkℓ − x⋆ k2 ∇f (x⋆ )T d + O(kxkℓ − x⋆ k2 ).

Wegen ∇f (x⋆ )d < 0 gibt es ein L ∈ N, so dass 1 f (xkℓ ) < f (x⋆ ) + kxkℓ − x⋆ k2 ∇f (x⋆ )T d 2 für alle ℓ ≥ L gilt. In jeder Umgebung von x⋆ finden wir daher ein xkℓ mit f (xkℓ ) < f (x⋆ ) im Widerspruch zur Voraussetzung. Mit Hilfe dieses Satzes können wir den Begriff des stationären Punktes auf das Minimierungsproblem mit Nebenbedingungen (7.1) verallgemeinern. Definition 7.4 Der Punkt x⋆ ∈ G wird als stationärer Punkt des Minimierungsproblems mit Nebenbedingungen (7.1) bezeichnet, falls ∇f (x⋆ )T d ≥ 0

(7.3)

für alle Grenzrichtungen d ∈ Rn zulässiger Folgen mit Grenzwert x⋆ gilt. Man beachte, dass die Bedingung (7.3) im Fall eines unrestringierten Minimierungsproblems mit der üblichen Bedingung ∇f (x⋆ ) = 0 zusammenfällt. Definition 7.5 Zu gegeben Punkt x⋆ ∈ G sei JA (x⋆ ) := JG ∪ {i ∈ JU : ci (x⋆ ) = 0}. die Menge der aktiven Nebenbedingungen. Sind die zu dieser Menge gehörigen Gradienten {∇ci (x⋆ ) : i ∈ JA (x⋆ )} linear unabhängig, so sagen wir, der Punkt x⋆ erfüllt die LICQ-Bedingung.

58

Kapitel 7. Optimierungsprobleme mit Nebenbedingungen

Bemerkung LICQ steht abkürzend für linear independence constraint qualification. △ Lemma 7.6 Sei x⋆ ein lokales Minimum des Minimierungsproblems (7.1). Ist d ∈ Rn eine Grenzrichtung einer zulässigen Folge, so gilt ∇ci (x⋆ )T d = 0 für alle i ∈ JG ,

∇ci (x⋆ )T d ≥ 0 für alle i ∈ JU ∩ JA .

Gelten umgekehrt für eine Richtung d ∈ Rn mit kdk2 = 1 diese beiden Bedingungen und ist zusätzlich die LICQ-Bedingung für x⋆ erfüllt, dann ist d eine Grenzrichtung zu x⋆ . Beweis. (i.) Es sei {xk }k∈N eine zulässige Folge, für die d eine Grenzrichtung ist. Wir erhalten, gegebenfalls durch Übergang zu einer Teilfolge, xk − x⋆ = d. k→∞ kxk − x⋆ k2 lim

Daraus folgt xk = x⋆ + kxk − x⋆ k2 d + O(kxkℓ − x⋆ k2 ). Ist i ∈ JG , so haben wir für alle k ≥ K ci (xk ) kxk − x⋆ k2 ci (x⋆ ) + kxk − x⋆ k2 ∇ci (x⋆ )T d + O(kxk − x⋆ k2 ) = kxk − x⋆ k2 ⋆ O (kxk − x k2 ) . = ∇ci (x⋆ )T d + kxk − x⋆ k2

0=

Durch den Grenzübergang k → ∞ schließen wir in diesem Fall ∇ci (x⋆ )T d = 0. Ist i ∈ JU ∪ JA (x⋆ ), dann erhalten wir für alle k ≥ K ⋆ O (kxk − x k2 ) ci (xk ) = 0≤ kxk − x⋆ k2 kxk − x⋆ k2

und analog zum obigen Vorgehen ∇ci (x⋆ )T d ≥ 0. Damit ist die erste Aussage des Satzes gezeigt. (ii.) Wir setzen A := [∇ci (x⋆ )T ]i∈JA (x⋆ ) und beachten, dass diese Matrix aufgrund der LICQ-Bedingung vollen Rang m := |JA (x⋆ )| ≤ n besitzt. Ohne Beschränkung der Allgemeinheit können wir annehmen, dass alle Nebenbedingungen aktiv sind, also JG ∪ JU = JA (x⋆ ) ist. Zur Matrix A ∈ Rm×n können wir eine Matrix B ∈ Rn×(n−m) mit vollem Rand n − m finden, so dass AB = 0 gilt. Dabei bilden die Spalten von B eine Basis des Kerns von A. Für t > 0 und d ∈ Rn mit kdk2 = 1 definieren wir die Abbildung φ : Rn × R → Rn durch 

c(x) − tAd φ(x, t) = BT (x − x⋆ − td)



mit c(x) = [ci (x)]i∈JA (x⋆ )

7.1. Optimalitätsbedingungen erster Ordnung

59

und betrachten die Lösung des Gleichungssystems φ(x, t) = 0. Für t = 0 ist x = x⋆ die Lösung dieses Gleichungssystems und die Jakobi-Matrix     A ∇c(x⋆ )T ⋆ ∈ Rn×n = ∇x φ(x , 0) = BT BT

ist nichtsingulär. Nach dem Satz über implizite Funktionen besitzt dieses Gleichungssystem also für hinreichend kleines t > 0 eine eindeutige Lösung x = x(t). Sei {tk }k∈N ⊂ (0, 1] eine Folge mit tk → 0 und für k ≥ K sei xk die eindeutige L¨soung von φ(x, tk ) = 0. Wir zeigen, dass {xk }k∈N eine zulässige Folge mit Grenzrichtung d ist. Zunächst folgern wir xk ∈ G für alle k ∈ N aus i ∈ JG ⇒ ci (xk ) = tk ∇ci (x⋆ )T d = 0,

i ∈ JU ∩ JA ⇒ ci (xk ) = tk ∇ci (x⋆ )T d ≥ 0.

Aufgrund des Satzes über implizite Funktionen hängt die Lösung xk des Gleichungssystems φ(x, t) = 0 stetig von tk ab. Dies impliziert lim xk = x⋆ .

k→∞

Als nächstes zeigen wir, dass xk 6= x⋆ gilt. Angenommen, es ist xk = x⋆ für ein k ∈ N, dann würde gelten   c(x⋆ ) − tk Ad ⋆ = 0. φ(x , tk ) = BT (x⋆ − x⋆ − tk d) Da alle Nebenbedingungen als aktiv angenommen wurden, gilt c(x⋆ ) = 0 und es folgt   A d = 0. BT

Daraus ergibt sich d = 0 im Widerspruch zu kdk2 = 1. Folglich ist {xk }k∈N eine zulässige Folge. Wir müssen noch nachweisen, dass d eine Grenzrichtung ist. 0 = φ(xk , tk )   c(xk ) − tk Ad = BT (xk − x⋆ − tk d)   ⋆ c(x ) + A(xk − x⋆ ) + O(kxk − x⋆ k2 ) − tk Ad = BT (xk − x⋆ − tk d)   A (xk − x⋆ − tk d) + O(kxk − x⋆ k2 ). = BT

Durch die spezielle Wahl von

dk := erhalten wir hieraus

 lim dk −

xk − x⋆ kxk − x⋆ k2

tk d k→∞ kxk − x⋆ k2 Wegen kdk k2 = 1 für alle k ∈ N und kdk2 = 1 folgt



tk =1 k→∞ kxk − x⋆ k2

= 0.

lim

und damit limk→∞ dk = d. Dies beweist die zweite Aussage des Satzes.

60

Kapitel 7. Optimierungsprobleme mit Nebenbedingungen

Definition 7.7 Zu einem Punkt x⋆ ∈ G und aktiven Nebenbedingungen JA (x⋆ ) bezeichnen wir K(x⋆ ) := {d ∈ Rn : ∇ci (x⋆ )T d = 0 für i ∈ JG und

∇ci (x⋆ )T d ≥ 0 für i ∈ JU ∩ JA (x⋆ )}

als Linearisierungskegel an x⋆ .

Lemma 7.8 Genau dann gibt es keine Richtung d ∈ K(x⋆ ) mit ∇f (x⋆ )T d < 0, wenn ein Vektor λ = [λi ]i∈JA (x⋆ ) existiert mit ∇f (x⋆ ) =

X

i∈JA (x⋆ )

λi ∇ci (x⋆ ) = AT λ

und λi ≥ 0 für i ∈ JU ∩ JA (x⋆ ). Beweis. (i.) Wir betrachten den Kegel   X n ⋆ ⋆ b := y ∈ R : y = K λi ∇ci (x ), λi ≥ 0 für i ∈ JU ∩ JA (x ) . i∈JA (x⋆ )

b offensichtlich abgeschlossen ist, lautet die Aussage des Lemmas also: Da K

b ∇f (x⋆ )T d ≥ 0 für alle Richtungen d ∈ K(x⋆ ) ⇐⇒ ∇f (x⋆ ) ∈ K.

b und d ∈ K(x⋆ ), dann gilt (ii.) Seien ∇f (x⋆ ) ∈ K X ∇f (x⋆ )T d = λi ∇ci (x⋆ )T d i∈JA (x⋆ )

=

X

i∈JG

≥0

λi ∇ci (x⋆ )T d +

X

i∈JU ∩JA (x⋆ )

λi ∇ci (x⋆ )T d

wegen λi ≥ 0 und ∇ci (x⋆ )T d ≥ 0 für alle i ∈ JU ∩ JA (x⋆ ). b dann sei y b der Punkt mit minimalen b ∈ K (iii.) Angenommen, es gilt ∇f (x⋆ ) 6∈ K, Abstand zu ∇f (x⋆ ), das heißt, die Lösung des Minimierungsproblems b min ky − ∇f (x⋆ )k2 unter der Nebenbedingung y ∈ K.

y∈Rn

b ein Kegel ist, liegt mit y b Insbesondere muss y b auch tb b also Da K y für alle t ≥ 0 in K. erfüllen: d ⋆ 2 y − ∇f (x )k2 0 = ktb dt t=1 2 T ⋆ = 2tkb yk2 − 2b y ∇f (x ) t=1
b − ∇f (x⋆ ) . = 2b yT y

7.1. Optimalitätsbedingungen erster Ordnung

61

b folgt aus der Konvexität von K, b dass Für beliebiges y ∈ K Deshalb ist

b) − ∇f (x⋆ )k22 ≥ kb kb y + θ(y − y y − ∇f (x⋆ )k22

was für θ → 0 auf

für alle θ ∈ (0, 1).


b )T y b − ∇f (x⋆ ) + θ2 ky − y bk22 ≥ 0, 2θ(y − y



b )T y b − ∇f (x⋆ ) = yT y b − ∇f (x⋆ ) 0 ≤ (y − y

(7.4)

führt. b −∇f (x⋆ ) 6= 0 die Bedingungen d ∈ K(x⋆ ) und ∇f (x⋆ )T d < 0 Wir zeigen nun, dass d := y erfüllt. Zweiteres folgt sofort aus
T b − ∇f (x⋆ ) y b − dT d = −kdk22 < 0. dT ∇f (x⋆ ) = dT (b y − d) = y

Ersteres sieht man hingegen wie folgt ein. Durch geeignete Wahl von λi , i ∈ JA (x⋆ ), erreicht man b und − ∇ci (x⋆ ) ∈ K, b i ∈ JG =⇒ ∇ci (x⋆ ) ∈ K b i ∈ JU ∩ JA (x⋆ ) =⇒ ∇ci (x⋆ ) ∈ K.

Setzen wir in (7.4) die spezielle Wahl y = ±∇ci (x⋆ ) falls i ∈ JG und y = ∇ci (x⋆ ) falls i ∈ JU ∩ JA (x⋆ ) ein, so erhalten wir i ∈ JG =⇒ ±∇ci (x⋆ )T d ≥ 0 also ∇ci (x⋆ )T d = 0,

i ∈ JU ∩ JA (x⋆ ) =⇒ ∇ci (x⋆ )T d ≥ 0.

Daraus folgt d ∈ K(x⋆ ).

Definition 7.9 Die zum Minmierungsproblem mit Nebenbedingungen (7.1) gehörige Lagrange-Funktion lautet X L(x, λ) := f (x) − λi ci (x). i∈JG ∪JU

Die Parameter λ = [λi ]i∈JG ∪JG werden Lagrange-Parameter genannt. Mit Hilfe der Lagrange-Funktion kann man nur ein Minimium von (7.1) charakterisieren. Satz 7.10 (Karush, Kuhn und Tucker) Es sei x⋆ ∈ G eine lokales Minimum des Minimierungsproblems mit Nebenbedingungen (7.1) und x⋆ erfülle die LICQ-Bedingung. Dann gibt es genau einen Vektor von Lagrange-Parametern λ⋆ , so dass die folgenden KKTBedingungen erfüllt sind: ∇x L(x⋆ , λ⋆ ) = 0, ci (x⋆ ) = 0 für ci (x⋆ ) ≥ 0 für λ⋆i ≥ 0 für λ⋆i ci (x⋆ ) = 0 für

alle alle alle alle

i ∈ JG , i ∈ JU , i ∈ JU , i ∈ JG ∪ JU .

62

Kapitel 7. Optimierungsprobleme mit Nebenbedingungen

Beweis. (i.) Angenommen, x⋆ ∈ G ist ein lokales Minimum von (7.1), an dem die LICQBedingung erfüllt ist. Dann gilt nach Satz 7.3 ∇f (x⋆ )T d ≥ 0 für alle dazugehörigen Grenzrichtungen d. Nach Lemma 7.6 sind alle Richtungen, für die die Bedingungen ∇ci (x⋆ )T d = 0 für alle i ∈ JG ,

∇ci (x⋆ )T d ≥ 0 für alle i ∈ JU ∩ JA (x⋆ ). efüllt sind, auch Grenzrichtungen und erfüllen somit ∇f (x⋆ )T d ≥ 0. Mit anderen Worten, für alle d ∈ K(x⋆ ) gilt ∇f (x⋆ )T d ≥ 0. Dies impliziert nach Lemma 7.8 die Existenz von λi ∈ R, so dass X ∇f (x⋆ ) = λi ∇ci (x⋆ ) mit λi ≥ 0 für i ∈ JU ∩ JA (x⋆ ). i∈JA (x⋆ )

Hierin sind die Lagrange-Parameter λi aufgrund der LICQ-Bedingung eindeutig bestimmt. (ii.) Wir definieren ( λi , i ∈ JA (x⋆ ) λ⋆i := 0, sonst und weisen nach, dass damit die KKT-Bedingungen erfüllt sind. Es gilt nach (i.) ∇x L(x⋆ , λ⋆ ) = ∇f (x⋆ ) − = ∇f (x⋆ ) − = 0.

X

i∈JG ∪JU

X

i∈JA (x⋆ )

λi ∇ci (x⋆ ) λi ∇ci (x⋆ )

Die Bedingungen ci (x⋆ ) = 0 für alle i ∈ JG und ci (x⋆ ) ≥ 0 für alle i ∈ JU folgen aus x⋆ ∈ G. Weiter gilt λ⋆i = 0 für i ∈ JU \ JA (x⋆ ) und λ⋆i ≥ 0 für i ∈ JU ∩ JA (x⋆ ), das heißt, es ist λ⋆i ≥ 0 für i ∈ JU . Schließlich gilt λ⋆i = 0 für i ∈ JU \ JA (x⋆ ), sowie ci (x⋆ ) = 0 für alle i ∈ JA (x⋆ ), und somit λ⋆i ci (x⋆ ) = 0 für alle i ∈ JG ∪ JU . Weil ci (x⋆ ) > 0 ist für alle i ∈ JU \ JA , folgt die Eindeutigkeit der λ⋆i aus der letzten KKT-Bedingung.

7.2 Optimalitätsbedingungen zweiter Ordnung Wir wenden uns nun Bedingungen zweiter Ordnung zu, wobei wir voraussetzen, dass f und ci für alle i ∈ JG ∪ JU zweimal stetig differenzierbar seien. Definition 7.11 Zu gegebenen Lagrange-Parametern λ⋆ sei K(x⋆ , λ⋆ ) = {d ∈ K(x⋆ ) : ∇ci (x⋆ )T d = 0 für alle i ∈ JU ∩ IA (x⋆ ) mit λ⋆i > 0}.

7.2. Optimalitätsbedingungen zweiter Ordnung

63

Bemerkung Wegen K(x⋆ , λ⋆ ) ⊂ K(x⋆ ) gilt offensichtlich  ⋆ T  ∇ci (x ) d = 0, i ∈ JG , d ∈ K(x⋆ , λ⋆ ) ⇐⇒ ∇ci (x⋆ )T d = 0, i ∈ JU ∩ JA (x⋆ ) mit λ⋆i > 0,   ∇ci (x⋆ )T d ≥ 0, i ∈ JU ∩ JA (x⋆ ) mit λ⋆i = 0.

Deshalb folgt aus der ersten KKT-Bedingung aus Satz 7.10 die Aussage X λ⋆i ∇ci (x⋆ )T d = 0. d ∈ K(x⋆ , λ⋆ ) =⇒ ∇f (x⋆ )T d = i∈JG ∪JU

△ Satz 7.12 (hinreichende Bedingung 2. Ordnung) Für einen zulässigen Punkt x⋆ ∈ G gebe es einen Vektor von Lagrange-Parametern λ⋆ , so dass die KKT-Bedingungen erfüllt sind. Weiter sei dT ∇2x L(x⋆ , λ⋆ )d > 0 für alle d ∈ K(x⋆ , λ⋆ ) \ {0}.

Dann ist x⋆ Lösung des Minimierungsproblems mit Nebenbedingungen (7.1).

Beweis. Der Satz ist bewiesen, wenn wir zeigen können, dass zu jeder zulässigen Folge {xk } ein M ∈ N existiert, so dass f (xk ) > f (x⋆ ) für alle k ≥ M ist. Zu einer beliebigen zulässigen Folge {xk } sei d eine beliebige Grenzrichtung. Nach Lemma 7.6 gilt d ∈ K(x⋆ ). Nach der Definition der Grenzrichtung gibt es eine Teilfolge {xkℓ }, welche xkℓ − x⋆ = kxkℓ − x⋆ k2 d + O(kxkℓ − x⋆ k2 ) (7.5) erfüllt. Wegen den KKT-Bedingungen gilt für die Lagrange-Funktion X L(xkℓ , λ⋆ ) = f (xkℓ ) − λ⋆i ci (xkℓ ) ≤ f (xkℓ ). i∈JA (x⋆ )

Taylor-Entwicklung ergibt L(xkℓ , λ⋆ ) = L(x⋆ , λ⋆ ) + ∇x L(x⋆ , λ⋆ ) T (xkℓ − x⋆ ) | {z } | {z } =0

=f (x⋆ )

1 + (xkℓ − x⋆ )T ∇2x L(x⋆ , λ⋆ )(xkℓ − x⋆ ) + O(kxkℓ − x⋆ k22 ) 2 1 = f (x⋆ ) + (xkℓ − x⋆ )T ∇2x L(x⋆ , λ⋆ )(xkℓ − x⋆ ) + O(kxkℓ − x⋆ k22 ), 2 wobei wieder die KKT-Bedingungen verwendet wurden. Wir unterscheiden nun zwei Fälle: d 6∈ K(x⋆ , λ⋆ ) und d ∈ K(x⋆ , λ⋆ ). Für d 6∈ K(x⋆ , λ⋆ ) gibt es ein i⋆ ∈ JU ∩ IA (x⋆ ), für welches λ⋆i⋆ ∇ci⋆ (x⋆ )T d > 0

gilt. Eine Taylor-Entwicklung führt auf λ⋆i⋆ ci⋆ (xkℓ ) = λ⋆i⋆ ci⋆ (x⋆ ) +λ⋆i⋆ ∇ci⋆ (x⋆ )T (xkℓ − x⋆ ) + O(kxkℓ − x⋆ k2 ) | {z } =0

(7.5)

=

λ⋆i⋆ kxkℓ

− x⋆ k2 ∇ci⋆ (x⋆ )T d + O(kxkℓ − x⋆ k2 ).

64

Kapitel 7. Optimierungsprobleme mit Nebenbedingungen

Für die Lagrange-Funktion gilt somit X L(xkℓ , λ⋆ ) = f (xkℓ ) − λ⋆i ci (xkℓ ) i∈JA (x⋆ )

≤ f (xkℓ ) − λ⋆i⋆ ci⋆ (xkℓ )

= f (xkℓ ) − λ⋆i⋆ kxkℓ − x⋆ k2 ∇ci⋆ (x⋆ )T d + O(kxkℓ − x⋆ k2 ). Andererseits wurde oben gezeigt, dass L(xkℓ , λ⋆ ) = f (x⋆ ) + O(kxkℓ − x⋆ k2 ). Daraus ergibt sich nun f (xkℓ ) ≥ f (x⋆ ) + λ⋆i⋆ kxkℓ − x⋆ k2 ∇ci⋆ (x⋆ )T d + O(kxkℓ − x⋆ k2 ) und wegen λi⋆ ∇ci⋆ (x⋆ )T d > 0 impliziert dies f (xkℓ ) > f (x⋆ ), falls nur kℓ genügend groß ist. Für d ∈ K(x⋆ , λ⋆ ) erhalten wir direkt f (xkℓ ) ≥ L(xkℓ , λ⋆ ) 1 = f (x⋆ ) + (xkℓ − x⋆ )T ∇2x L(x⋆ , λ⋆ )(xkℓ − x⋆ ) + O(kxkℓ − x⋆ k22 ) 2 1 (7.5) = f (x⋆ ) + kxkℓ − x⋆ k22 dT ∇2x L(x⋆ , λ⋆ )d + O(kxkℓ − x⋆ k22 ), 2 woraus sich wegen dT ∇2x L(x⋆ , λ⋆ )d > 0 wiederum f (xkℓ ) > f (x⋆ ) für genügend große kℓ ergibt. Bemerkung Die notwendige Bedingung zweiter Ordnung für ein Minimum von (7.1) lautet: Ist x⋆ ∈ G ein Minimum von (7.1), an dem die LICQ-Bedingung erfüllt ist, und genügen die zugehörigen Lagrange-Parameter λ⋆ den KKT-Bedingungen, dann gilt dT ∇2x L(x⋆ , λ⋆ )d ≥ 0 für alle d ∈ K(x⋆ , λ⋆ ). △

65

8. Projiziertes Gradientenverfahren 8.1 Konvergenzeigenschaften Ziel dieses Kapitels ist es, das Gradientenverfahren aus Kapitel 2 auf das Minimierungsproblem mit Nebenbedingungen (7.1) zu verallgmeinern. Dazu sei neben der stetigen Differenzierbarkeit der Funktion f : Rn → R vorausgesetzt, dass die zulässige Menge G ⊂ Rn konvex ist. Grundlage des projizierten Gradientenverfahren ist die orthogonale Projektion: Definition 8.1 Es sei G ⊂ Rn eine abgeschlossene, konvexe Menge. Dann ist die orthogonale Projektion PG : Rn → G definiert durch die Bedingung kPG (x) − xk2 = min ky − xk2 . y∈G

Der Punkt PG (x) ∈ G besitzt also die Eigenschaft, den kürzesten Abstand zu einem gegebenen Punkt x ∈ Rn zu besitzen. Die Berechnung des projizierten Punktes PG (x) lässt sich für viele spezielle Nebenbedingungen relativ einfach umsetzen. Beispielsweise gilt dies für affine Nebenbedingungen, die später genauer behandelt werden. Die Grundversion des projizierten Gradientenverfahren ist im folgenden Algorithmus beschrieben: Algorithmus 8.2 (projiziertes Gradientenverfahren) input: Funktion f : Rn → R, konvexe zulässige Menge G ⊂ Rn und Startnäherung x0 ∈ G output: Folge von Iterierten {xk }k∈N ➀ Initialisierung: wähle σ ∈ (0, 1) und setze k := 0 ➁ berechne den Antigradienten dk := −∇f (xk ) und setze αk := 1 ➂ solange

f PG (xk + αk dk ) > f (xk ) − σdTk PG (xk + αk dk ) − xk setze αk := αk /2 ➃ setze xk+1 := PG (xk + αk dk ) ➄ erhöhe k := k + 1 und gehe nach ➁

(8.1)

66

Kapitel 8. Projiziertes Gradientenverfahren

Für den Fall der Minimierung ohne Nebenbedingungen, das heißt G = Rn , stellt obiger Algorithmus genau das bekannte Gradientverfahren 2.4 dar. Insbesondere geht die Bedingung an die Reduktion des Funktionals über in die Armijo-Goldstein-Bedingung f (xk+1 ) ≤ f (xk ) − σαk k∇f (xk )k22 . Ähnlich wie früher kann man zeigen, dass auch für das projizierte Gradientenverfahren ein αk > 0 existiert, für das die Reduktionsbedingung erfüllt ist. Lemma 8.3 Ist die zulässige Menge G ⊂ Rn konvex, so erfüllt die orthogonale Projektion PG die folgenden Eigenschaften:
T
(i.) Es gilt PG (x) − x PG (x) − y ≤ 0 für alle x ∈ Rn , y ∈ G.
T (ii.) Es gilt PG (y) − PG (x) (y − x) ≥ kPG (y) − PG (x)k22 ≥ 0 für alle x, y ∈ Rn , das heißt, PG ist monoton. (iii.) Es gilt kPG (y) − PG (x)k2 ≤ ky − xk2 für alle x, y ∈ Rn , das heißt, PG ist nicht expandierend. b := (1 − t)PG (x) + ty ∈ G für Beweis. (i.) Wegen der Konvexität folgt aus y ∈ G auch y alle t ∈ [0, 1]. Aus kb y − xk22 = kb y − PG (x) + PG (x) − xk22


T
b = kb y − PG (x)k22 + kPG (x) − xk22 − 2 PG (x) − x PG (x) − y

folgt aufgrund der Minimierungseigenschaft von PG , dass
T
b = kb kb y − PG (x)k22 − 2 PG (x) − x PG (x) − y y − xk22 − kPG (x) − xk22 ≥ 0.
b = t PG (x) − y führt auf Einsetzen von PG (x) − y
T
t2 ky − PG (x)k22 − 2t PG (x) − x PG (x) − y ≥ 0, was für t → 0 die gewünschte Aussage liefert. (ii.) Die bereits bewiesene Aussage (i.) impliziert
T
PG (x) − x PG (x) − PG (y) ≤ 0,
T
PG (y) − y PG (y) − PG (x) ≤ 0. Zusammen führt dies auf

PG (y) − y + x − PG (x)


T


PG (y) − PG (x) ≤ 0,

das ist Aussage (ii.). (iii.) Diese Aussage folgt sofort aus Aussage (ii.) durch Anwenden der Cauchy-Schwarzschen Ungleichung. Bemerkung Aus der Monotonieeigenschaft (ii.) folgt wegen xk = PG (xk ) für die neue
Iterierte xk+1 = PG xk − αk ∇f (xk ) des projizierten Gradientenverfahrens, dass
kxk+1 − xk k22 ≤ (xk+1 − xk )T (xk − αk ∇f (xk ) − xk .

8.1. Konvergenzeigenschaften

67

Wir erhalten daher − ∇f (xk )T (xk+1 − xk ) ≥

1 kxk+1 − xk k22 . αk

(8.2)

Dies bedeutet, dass durch die Abstiegsbedingung (8.1) des Algorithmus 8.2 tatsächlich f (xk+1) < f (xk ) erreicht wird. △ Lemma 8.4 Ist die zulässige Menge G ⊂ Rn konvex, so ist für beliebige x ∈ Rn und d ∈ Rn die Funktion kPG (x + αd) − xk2 ϕ(α) := α für alle α > 0 monoton fallend. Beweis. (i.) Für 0 < α < β setzen wir u := PG (x + αd) − x,

v := PG (x + βd) − x

und erhalten unter Verwendung von Lemma 8.3 (i.) uT (u − v) = {PG (x + αd) − (x + αd) + αd}T {PG (x + αd) − PG (x + βd)} ≤ αdT {PG (x + αd) − PG (x + βd)}

und analog vT (v − u) ≤ βdT {PG (x + βd) − PG (x + αd)}.

Zusammen ergibt dies

uT (u − v) vT (u − v) ≤ . α β

(8.3)

(ii.) Weiter erhalten wir mit Lemma 8.3 (ii.) uT (u − v) ≤ αdT {PG (x + αd) − PG (x + βd)} α =− (βd − αd)T {PG (x + βd) − PG (x + αd)} β−α α ≤− kPG (x + βd) − PG (x + αd)k22 β−α ≤ 0. Aus der Cauchy-Schwarzschen Ungleichung ergibt sich uT v(kuk2 + kvk2 ) ≤ kuk2 kvk2 (kuk2 + kvk2 ), woraus kuk2 vT (u − v) = kuk2 (uT v − kvk22 ) ≤ kvk2 (kuk22 − uT v) = kvk2 uT (u − v)

(8.4)

folgt. (iii.) Wir unterscheiden nun zwei Fälle: Für uT (u − v) = 0 gilt PG (x + αd) = PG (x + βd) und somit u = v. Hieraus folgt unmittelbar auch ϕ(α) =

kuk2 kvk2 ≥ = ϕ(β). α β

68

Kapitel 8. Projiziertes Gradientenverfahren

Für den Fall uT (u − v) < 0 folgt aus (8.3) vT (u − v) β ≥ T α u (u − v) und aus (8.4) kvk2 ≤ kuk2

vT (u − v) . uT (u − v)

Kombiniert man diese zwei Ungleichungen, so erhält man wieder ϕ(α) =

kuk2 kvk2 ≥ = ϕ(β). α β

Satz 8.5 Die zulässige Menge G ⊂ Rn sei konvex und die Funktion f : Rn → R sei auf G stetig differenzierbar und nach unten beschränkt. Weiter sei ∇f auf G gleichmäßig stetig. Dann gilt für die Iterierten {xk }k∈N des projizierten Gradientenverfahrens kxkℓ +1 − xk k2 = 0. k→∞ αk lim

Beweis. Wir führen einen Widerspruchsbeweis. Angenommen, es existiert zu jedem ε > 0 eine unendliche Teilfolge {kℓ }ℓ∈N , so dass kxkℓ +1 − xkℓ k2 ≥ ε. αk ℓ Dann gilt insbesondere auch kxkℓ +1 − xkℓ k22 ≥ ε max{εαkℓ , kxkℓ +1 − xkℓ k2 }. αk ℓ

(8.5)

Da die Folge {f (xkℓ )}ℓ∈N monoton fallend und nach unten beschränkt ist, folgt aus der Abstiegsbedingung (8.1) des projizierten Gradientenverfahrens lim ∇f (xkℓ )T (xkℓ +1 − xkℓ ) = 0,

ℓ→∞

was wiederum gemäß (8.2) kxkℓ +1 − xkℓ k22 = 0. ℓ→∞ αk ℓ lim

nach sich zieht. Aufgrund von (8.5) erhalten wir hieraus lim αkℓ = 0 und

k→∞

lim kxkℓ +1 − xkℓ k2 = 0.

k→∞

Für ykℓ +1 := PG (xkℓ + 2αkℓ dkℓ ) gilt aufgrund der algorithmischen Umsetzung des projizierten Gradientenverfahrens f (ykℓ +1 ) > f (xkℓ ) + σ∇f (xkℓ )T (ykℓ +1 − xkℓ ),

8.1. Konvergenzeigenschaften

69

also auch Aus Lemma 8.4 folgt

f (xkℓ ) − f (ykℓ +1 ) < σ∇f (xkℓ )T (xkℓ − ykℓ +1 ).

(8.6)

kykℓ +1 − xkℓ k2 kykℓ +1 − xkℓ k2 ε kxkℓ +1 − xkℓ k22 ≥ kxkℓ +1 − xkℓ k2 ≥ αk ℓ ε = kykℓ +1 − xkℓ k2 . αk ℓ 2αkℓ 2αkℓ 2 Weiter ergibt sich mit Lemma 8.3 (ii.)

(ykℓ +1 − xkℓ +1 )

 T

 (xkℓ +1 − xkℓ )T xkℓ − αkℓ ∇f (xkℓ ) − xkℓ ≥ kxkℓ +1 − xkℓ k22 , {z } | =−αkℓ ∇f (xkℓ )


xkℓ − 2αkℓ ∇f (xkℓ ) − xkℓ − αkℓ ∇f (xkℓ ) ≥ 0. | {z } =−αkℓ ∇f (xkℓ )

Zusammen führt dies auf

(xkℓ − ykℓ +1 )T ∇f (xkℓ ) ≥ (xkℓ − xkℓ +1 )T ∇f (xkℓ ) kxkℓ +1 − xkℓ k22 ≥ αk ℓ ε ≥ kykℓ +1 − xkℓ k2 , 2 woraus sich speziell auch kykℓ +1 −xkℓ k2 → 0 für ℓ → ∞ ergibt. Die gleichmäßige Stetigkeit von ∇f impliziert nun 1 − f (xkℓ ) − f (ykℓ +1 ) = O(kykℓ +1 − xkℓ k2 ) T (xkℓ − ykℓ +1 ) ∇f (xkℓ ) (xkℓ − ykℓ +1 )T ∇f (xkℓ ) 2 O(kykℓ +1 − xkℓ k2 ) ≤ ε kxkℓ − ykℓ +1 k2 ℓ→∞

−→ 0.

Dies steht jedoch im Widerspruch zu der aus (8.6) folgenden Abschätzung f (xkℓ ) − f (ykℓ +1 ) < σ < 1. (xkℓ − ykℓ +1 )T ∇f (xkℓ )

Bemerkung Da αk ≤ 1 ist, folgt aus Satz 8.5 speziell kxk+1 − xk k2 → 0 für k → ∞, das heißt, die Iterierten aus Algorithmus 8.2 konvergieren gegen einen Punkt x⋆ ∈ G. △ Definition 8.6 Die zulässige Menge G ⊂ Rn sei konvex. Für jeden Punkt x ∈ G ist der Tangentialkegel TG (x) definiert als kleinster abgeschlossener Kegel, der die Menge {d = y − x : y ∈ G} enthält.

70

Kapitel 8. Projiziertes Gradientenverfahren

Bemerkung Der Tangentialkegel TG (x⋆ ) ist die Menge aller Grenzrichtungen zulässiger Folgen für das Minimierungsproblem (7.1), skaliert mit einem beliebigen positiven Faktor. Wir erinnern an den Begriff des Linearisierungskegels K(x⋆ ) aus dem vorigen Kapitel. Er besteht aus allen Richtungen d ∈ Rn mit dT ∇ci (x⋆ ) = 0 für i ∈ JG und dT ∇ci (x⋆ ) ≥ 0 für i ∈ JU ∩ JA (x⋆ ). Nach Lemma 7.6 gilt TG (x⋆ ) ⊂ K(x⋆ ) und für den Fall, dass die LICQ-Bedingung erfüllt ist, sogar TG (x⋆ ) = K(x⋆ ). △ Lemma 8.7 Die zulässige Menge G ⊂ Rn
sei konvex. Für jeden Punkt x ∈ G erfüllt die orthogonale Projektion PTG (x) − ∇f (x) der Richtung des steilsten Abstiegs auf den Tangentialkegel TG (x) die folgenden Eigenschaften: (i.) Es gilt



2 ∇f (x)T PTG (x) − ∇f (x) = − PTG (x) − ∇f (x) 2 . (ii.) Es ist


min{∇f (x)T d : d ∈ TG (x), dk2 ≤ 1} = −kPTG (x) − ∇f (x) 2 .

(iii.) Der Punkt x ist genau dann ein stationärer Punkt
des Minimierungsproblems mit Nebenbedingungen (7.1), wenn PTG (x) − ∇f (x) = 0. Beweis. (i.) Nach Definition der Orthogonalprojektion besitzt die Funktion

2
1 g(λ) := λPTG (x) − ∇f (x) + ∇f (x) 2 2

ein Minimum bei λ = 1. Daher gilt


2
g ′ (1) := PTG (x) − ∇f (x) 2 + ∇f (x)T PTG (x) − ∇f (x) = 0.

(ii.) Wegen Aussage (i.) gilt





PT (x) − ∇f (x) + ∇f (x) 2 = k∇f (x)k22 − PT (x) − ∇f (x) 2 . G G 2 2


Für alle d ∈ TG (x) mit kdk2 ≤ PTG (x) −∇f (x) 2 gilt nach Definition der orthogonalen Projektion




PT (x) − ∇f (x) + ∇f (x) 2 ≤ kd + ∇f (x)k22 G 2


2 ≤ PT (x) − ∇f (x) + 2∇f (x)T d + k∇f (x)k2. G

2

2

Zusammen ergibt dies


2 ∇f (x)T d ≥ − PTG (x) − ∇f (x) 2 .


b = d/ PT (x) − ∇f (x) setzt. Das Behauptete erhält man, indem man d G 2

(iii.) Definitionsgemäß ist x ∈ G genau dann ein stationärer Punkt, wenn ∇f (x)T d ≥ 0 ist für alle Grenzrichtungen zulässiger Folgen mit Grenzwert x. Dies ist gleichbedeutend damit, dass ∇f (x)T d ≥ 0 für alle d ∈ TG
(x) ist. Aussage (ii.) impliziert, dass dies genau dann der Fall ist, wenn PTG (x) − ∇f (x) = 0 erfüllt ist.

8.1. Konvergenzeigenschaften

71

Bemerkung Aussage (ii.) des obigen Lemmas kann auch als


min{∇f (x⋆ )T d : d ∈ TG (x⋆ ), dk2 = 1} = −kPTG (x⋆ ) − ∇f (x⋆ ) 2 .

(8.7)

geschrieben werden, da das Minimum für kdk2 = 1 angenommen wird.

△

Satz 8.8 Die zulässige Menge G ⊂ Rn sei konvex und die Funktion f : Rn → R sei auf G stetig differenzierbar und nach unten beschränkt. Weiter sei ∇f auf G gleichmäßig stetig. Dann gilt für die Iterierten {xk }k∈N des projizierten Gradientenverfahrens
lim PTG (xk ) − ∇f (xk ) = 0. k→∞

Beweis. Zu beliebigen ε > 0 gibt es nach Lemma 8.7 (ii.) zu jeder Iterierten xk ein dk ∈ TG (xk ) mit kdk k2 = 1, so dass


∇f (xk )T dk ≤ − PTG (xk ) − ∇f (xk ) 2 + ε (8.8) gilt. Da dk Grenzrichtung einer zulässigen Folge ist, gibt es ein yk ∈ G mit



yk − xk

≤ ε.

− d k

kyk − xk k2 2

Aus Lemma 8.3 (i.) folgt 
T xk+1 − xk − αk ∇f (xk ) (xk+1 − yk+1 ) 

T 
= PG xk − αk ∇f (xk ) − xk − αk ∇f (xk ) PG xk − αk ∇f (xk ) − yk+1 ≤ 0,

was auf beziehungsweise

αk ∇f (xk )T (xk+1 − yk+1 ) ≤ kxk+1 − xk k2 kxk+1 − yk+1 k2 ,

kxk+1 − xk k2 ∇f (xk )T (yk+1 − xk+1 ) ≤ , kxk+1 − yk+1 k2 αk führt. Insgesamt erhalten wir deshalb

T

yk+1 − xk+1 T

− ∇f (xk ) (yk+1 − xk+1 ) − d −∇f (xk ) dk+1 ≤ k∇f (xk )k2 k+1

kyk+1 − xk+1 k2 kyk+1 − xk+1 k2 2 kxk+1 − xk k2 ≤ εk∇f (xk )k2 + . αk −

Die Kombination mit (8.8) ergibt




PT (x ) − ∇f (xk+1) G k+1 2

≤ −∇f (xk+1 )T dk+1 + ε

≤ −∇f (xk )T dk+1 + k∇f (xk+1) − ∇f (xk )k2 kdk+1 k2 +ε | {z } =1

kxk+1 − xk k2 + k∇f (xk+1 ) − ∇f (xk )k2 + ε. ≤ εk∇f (xk )k2 + αk

72

Kapitel 8. Projiziertes Gradientenverfahren

Weil ε > 0 beliebig war, folgt hieraus schließlich


kxk+1 − xk k2 + lim k∇f (xk+1 ) − ∇f (xk )k2 = 0 lim PTG (xk+1 ) − ∇f (xk+1 ) 2 ≤ lim k→∞ k→∞ k→∞ αk

wobei Satz 8.5 und die gleichmäßige Stetigkeit von ∇f zur Anwendung kommt.
In der Regel folgt aus der Stetigkeit von ∇f nicht, dass auch PTG (x) − ∇f (x) stetig ist. Um sicherzustellen, dass die Iterierten des projizierten Gradientenverfahrens tatsächlich gegen einen stationären Punkt konvergieren, benötigen wir daher das folgende Resultat. Satz 8.9 Die zulässige Menge G ⊂ Rn sei konvex und die Funktion f : Rn → R sei auf G stetig differenzierbar. Dann folgt für jede Folge {xk }k∈N ⊂ G mit xk → x⋆ ∈ G





PT (x⋆ ) − ∇f (x⋆ ) ≤ lim inf PT (x ) − ∇f (xk ) . G G k 2 2 k→∞

Beweis. Aus Lemma 8.7 (ii.) folgt für jedes y ∈ G


−∇f (xk )T (y − xk ) ≤ PTG (xk ) − ∇f (xk ) 2 ky − xk k2 , woraus sich für k → ∞ die Ungleichung


−∇f (x⋆ )T (y − x⋆ ) ≤ PTG (xk ) − ∇f (xk ) 2 ky − x⋆ k2

ergibt. Jedes d ∈ TG (x⋆ ) mit kdk2 = 1 ist Grenzrichtung einer zulässigen Folge {yk }k∈N ⊂ G, das heißt, es gilt yk − x⋆ k→∞ kyk − x⋆ k2

d = lim

lim yk = x⋆ .

und

k→∞

Somit erhalten wir


−∇f (x⋆ )T d ≤ lim inf PTG (xk ) − ∇f (xk ) 2 k→∞

und daraus wegen (8.7) die Behauptung:




PT (x⋆ ) − ∇f (x⋆ ) = max{−∇f (x⋆ )T d : d ∈ TG (x⋆ ), dk2 = 1} G 2


≤ lim inf PT (x ) − ∇f (xk ) . k→∞

G

k

2

Bemerkung Die Kombination der Sätze 8.5, 8.8 und 8.9 liefert die folgende Aussage: Ist die zulässige Menge G ⊂ Rn konvex und ist die Funktion f : Rn → R auf G stetig differenzierbar mit gleichmäßig stetigem Gradienten und nach unten beschränkt, dann konvergieren die Iterierten {xk }k∈N des projizierten Gradientenverfahrens 8.2 gegen ein x⋆ ∈ G mit
PTG (x⋆ ) − ∇f (x⋆ ) = 0. Gemäß Lemma 8.7 (iii.) bedeutet dies, dass x⋆ ein stationärer Punkt ist.

△

8.2. Affine Nebenbedingungen

73

8.2 Affine Nebenbedingungen Bei Anwendungsproblemen treten häufig affine Ungleichungsnebenbedingungen auf, weshalb wir diesen Spezialfall eingehender untersuchen wollen. Wir nehmen demnach an, dass das Minimierungsproblem unter Nebenbedingungen von der speziellen Form min f (x) unter den Nebenbedingungen

x∈Rn

aTi x ≤ bi für alle i = 1, 2, . . . , m

(8.9)

mit ai ∈ Rn und bi ∈ R für alle i = 1, 2, . . . , m ist. Dies ist also ein Spezialfall von (7.1), für den JG = ∅, JU = {1, 2, . . . , m} und ci (x) := bi − aTi x gesetzt wird. Die Menge der zulässigen Punkte ist in diesem Fall gegeben durch G = {x ∈ Rn : aTi x ≤ bi für alle i = 1, 2, . . . , m}. Für x ∈ G betrachten wir die Menge der aktiven Indizes  JA (x) := i ∈ {1, 2, . . . , m} : aTi x = bi .

Die LICQ-Bedingung bedeutet hier gerade die lineare Unabhängigkeit der Vektoren ai für i ∈ JA (x). Wir wenden uns der Charakterisierung von stationären Punkten zu, die für den Spezialfall affiner Nebenbedingungen besonders handlich ist.

Lemma 8.10 Für den Punkt x⋆ ∈ G sei die Menge {ai : i ∈ JA (x⋆ )} linear unabhängig. Dann ist x⋆ genau dann ein stationärer Punkt für das Minimierungsproblem unter affinen Nebenbedingungen (8.9), wenn es für jedes i ∈ JA (x⋆ ) eine eindeutige Zahl λ⋆i ≥ 0 gibt, so dass X ∇f (x⋆ ) + λ⋆i ai = 0 (8.10) i∈JA (x⋆ )

gilt.

Beweis. Nach Lemma 8.7 (iii.) ist x⋆ ∈ G genau dann ein stationärer Punkt des Problems
(8.9), wenn für die orthogonale Projektion auf den Tangentialkegel PTG (x⋆ ) −∇f (x⋆ ) = 0 gilt. Dies ist nach Lemma 8.7 (ii.) genau dann der Fall, wenn −∇f (x⋆ )T d ≤ 0 für alle d ∈ TG (x⋆ ). gilt. Da die LICQ-Bedingung erfüllt ist, gilt nun nach Lemma 7.6 TG (x⋆ ) = K(x⋆ ) = {d ∈ Rn : dT ai ≤ 0 für alle i ∈ JA (x⋆ )}. Schließlich ist nach Lemma 7.8 die Aussage −∇f (x⋆ )T d ≤ 0 für alle d ∈ K(x⋆ ) äquivalent zur Existenz und Eindeutigkeit nichtnegativer Zahlen λ⋆i ∈ R, i ∈ JA (x⋆ ), mit X ∇f (x⋆ ) = − λ⋆i ai . i∈JA (x⋆ )

74

Kapitel 8. Projiziertes Gradientenverfahren

Bemerkung Die Notwendigkeit der Bedingung (8.10) ergibt sich auch direkt aus den KKT-Bedingungen (vergleiche Satz 7.10). Lemma 8.10 sagt jedoch zusätzlich aus, dass diese Bedingung bei Vorliegen affiner Nebenbedingungen auch hinreichend für das Vorliegen eines stationären Punktes ist. △ Definition 8.11 Ein stationärer Punkt heißt nicht entartet, wenn die Menge {ai : i ∈ JA (x⋆ )} linear unabhängig ist und für die Lagrange-Parameter in (8.10) λ⋆i > 0 für alle i ∈ JA (x⋆ ) gilt. Satz 8.12 Die Funktion f : Rn → R sei stetig differenzierbar auf der zulässigen Menge G = {x ∈ Rn : aTi x ≤ bi für alle i = 1, 2, . . . , m}. Die Folge {xk }k∈N ∈ G konvergiere gegen den nicht entarteten stationären Punkt x⋆ ∈ G und es gelte
lim PTG (xk ) − ∇f (xk ) = 0. k→∞

Dann gibt es ein L ∈ N mit JA (xk ) = JA (x⋆ ) für alle k ≥ L, das heißt, die Menge der aktiven Indizes ändert sich für genügend große k nicht mehr.

Beweis. Aus xk → x⋆ und xk ∈ G folgt JA (xk ) = JA (x⋆ ) für genügend große k. Wir nehmen nun das Gegenteil der Aussage des Satzes an, nämlich dass eine unendliche Teilfolge {xkℓ }ℓ∈N und ein Index p ∈ JA (x⋆ ) existieren, so dass p 6∈ JA (xkℓ ) für alle ℓ ∈ N. Wir verwenden den orthogonalen Projektor PH auf den Unterraum  H := x ∈ Rn : aTi x = 0 für alle i ∈ JA (x⋆ ) \ {p} .

Für i ∈ JA (x⋆ ) \ {p} gilt insbesondere aTi PH ai = 0, woraus PH ai = 0 folgt. Aufgrund der Eigenschaft der orthogonalen Projektion steht ap − PH ap senkrecht auf  dem Unterraum H, was wiederum bedeutet, dass a − P a im Raum span a : i ∈ p H p i ⋆ JA (x ) \ {p} enthalten ist. Angenommen, es würde PH ap = 0 gelten, dann ließe sich ap als Linearkombination der Vektoren ai , i ∈ JA (x⋆ ) \ {p}, schreiben. Die Menge {ai : i ∈ JA (x⋆ )} ist aber linear unabhängig, da der stationäre Punkt x⋆ nicht entartet ist. Somit muss PH ap 6= 0 gelten. Wegen PH ap ∈ H, das heißt, aTi PH ap = 0 für i ∈ JA (x⋆ ) \ {p} ⊃ JA (xkℓ ) schließen wir insbesondere PH ap ∈ TG (xkℓ ) für alle ℓ ∈ N. Aus Lemma 8.7 (ii.) erhalten wir


∇f (xk )T PH ap ≥ − PTG (xk ) − ∇f (xk ) 2 kPH ap k2 .

Für k → ∞ ergibt sich wegen xk → x⋆ damit

∇f (x⋆ )T PH ap ≥ 0. Andererseits gilt, da der stationäre Punkt x⋆ nicht entartet ist, X ∇f (x⋆ ) + λ⋆i ai = 0 i∈JA (x⋆ )

(8.11)

8.2. Affine Nebenbedingungen

75

mit λ⋆i > 0 für i ∈ JA (x⋆ ). Dies impliziert


T ∇f (x⋆ )T PH ap = PH ∇f (x⋆ ) ap  X T ⋆ =− λi PH ai ap i∈JA (x⋆ )

= −λ⋆p (PH ap )T ap = −λ⋆p kPH ap k22 0 für alle p ∈ kern(AT ) \ {0}, so ist die KKT-Matrix   H A K := AT 0 regulär und folglich (9.3) eindeutig lösbar. Beweis. Sei p ∈ Rn und q ∈ Rm derart, dass      H A p p =0 = K AT 0 q q gilt. Aus der zweiten Zeile folgt AT p = 0 und daher ist     T  H A p (9.4) T T 0= p q = p Hp. AT 0 q

Dies kann aber nach Voraussetzung nur für p = 0 erfüllt sein. Wegen Hp + Aq = 0 folgt schließlich Aq = 0 und aufgrund des vollen Spaltenrangs von A ist auch q = 0. Folglich ist K regulär und das Behauptete bewiesen.

78

Kapitel 9. SQP-Verfahren

Bemerkung Man kann sogar zeigen, dass die Matrix K genau n positive und m negative Eigenwerte besitzt. Die Lösung indefiniter Systeme kann mit geeigneten Erweiterungen der Cholesky-Zerlegung oder aber bei sehr großen Systemen auch iterativ, etwa mit dem MINRES-Verfahren oder dem Bramble-Pasciak-CG, geschehen. △ Der nächste Satz besagt, dass die Lösung von (9.3) die Lösung von (9.2) liefert. Satz 9.2 Es gelten die Voraussetzungen von Satz 9.1. Dann ist die Lösung x⋆ ∈ Rn von (9.3) die eindeutige Lösung von (9.2). Beweis. Sei x 6= x⋆ derart, dass AT x = b und setze p := x⋆ − x. Dann ist AT p = 0 und p 6= 0. Für q(x) := 21 xT Hx + gT x gilt dann offensichtlich mit x = x⋆ − p 1 q(x) = pT Hp − pT Hx⋆ − gT p + q(x⋆ ). 2

(9.5)

Aus (9.3) folgt, dass −Hx⋆ = Aλ + g, so dass −pT Hx⋆ = pT (Aλ + g) = gT p ist. Eingesetzt in (9.5) ergibt sich 1 q(x) = pT Hp + q(x⋆ ) > q(x⋆ ), 2 dies bedeutet, x⋆ ist das eindeutige Minimum von (9.2).

9.2 Bestimmung aktiver Nebenbedingungen Wir wenden uns nun der Lösung der Minimierungsaufgabe (9.1) zu, wobei wir annehmen, dass die Vektoren ai linear unabhängig sind. Indem wir wie bisher für jeden Punkt x ∈ Rn die Menge der aktiven Indizes bezeichnen als JA (x) = {i ∈ JG ∪ JU : aTi x = bi } = JG ∪ {i ∈ JU : aTi x = bi }, können wir (9.1) umschreiben gemäß 1 minn xT Hx + gT x unter den Nebenbedingungen x∈R 2 aTi x = bi für alle i ∈ JA (x), aTi x ≤ bi für alle i ∈ JU \ JA (x).

Für jede lokale Lösung x⋆ ∈ Rn des Minimierungsproblems (9.1) gelten somit nach Satz 7.10 die folgenden KKT-Bedingungen X Hx⋆ + g + λ⋆i ai = 0, i∈JA (x⋆ )

aTi x⋆ aTi x⋆

= bi für i ∈ JA (x⋆ ),

≤ bi für i ∈ JU \ JA (x⋆ ), λi ≥ 0 für i ∈ JU \ JA (x⋆ ).

9.2. Bestimmung aktiver Nebenbedingungen

79

Die Grundidee des folgenden Vorgehens ist die Beobachtung, dass bei Kenntnis der Indizes JA (x⋆ ) der aktiven Nebenbendigungen am Lösungspunkt die Lösung von (9.1) gleichzeitig die Lösung des quadratischen Minimierungsproblem unter Gleichungsnebenbedingungen 1 minn xT Hx + gT x unter den Nebenbedingungen x∈R 2 aTi x = bi für alle i ∈ JA (x⋆ )

ist. Die aktiven Indizes werden über ein iteratives Vorgehen gefunden, ausgehend von einem Sattelpunkt x0 , der zulässiger Punkt für das Problem (9.1) ist, und einer Indexmenge J0 mit JG ⊂ J0 ⊂ JA (x0 ). Ausgehend von einem zulässigem Punkt xk und einer aktuellen Indexmenge Jk mit JG ⊂ Jk ⊂ JA (xk ) besteht ein Schritt dieses iterativen Vorgehens aus den folgenden beiden Halbschritten: 1. Löse das quadratische Minimierungsproblem mit Gleichungsnebenbedingungen an der aktuellen Indexmenge: 1 min (xk + d)T H(xk + d) + gT (xk + d) unter den Nebenbedingungen 2 aTi (xk + d) = bi für alle i ∈ Jk .

d∈Rn

(9.6)

2. Konstruiere aus dem aus dem ersten Halbschritt erhaltenen Punkt xk + dk einen Punkt xk+1 , der in der zulässigen Menge G = {x ∈ Rn : aTi x = bi für alle i ∈ JG und aTi x ≤ bi für alle i ∈ JU } enthalten ist. Die Lösung von (9.6) ist offenbar äquivalent zu 1 1 minn dT Hd + (g + Hxk )T d + xTk Hxk + gT xk unter den Nebenbedingungen d∈R 2 2 T ai d = bi − aTi xk für alle i ∈ Jk . Mit gk := g + Hxk und unter der Beachtung von aTi xk = bi sowie der Tatsache, dass Konstanten in der zu minimierenden Funktion keine Rolle spielen, lässt sich dies auch umschreiben in 1 minn dT Hd + gkT d unter den Nebenbedingungen d∈R 2 (9.7) T ai d = 0 für alle i ∈ Jk .

Dieses quadratische Minimierungsproblem unter Gleichungsnebenbedingungen hat die Form (9.2) und kann entsprechend gelöst werden. Bei der Durchführung des zweiten Halbschritts wird wie folgt vorgegangen. Erfüllt xk +dk auch die Nebenbedingungen aTi (xk + dk ) ≤ bi für alle i ∈ JU \ Jk ,

dann wird xk+1 := xk + dk gesetzt. Falls nicht alle Nebenbedingungen erfüllt sind, bestimmen wir xk+1 := xk + αk dk mit αk ∈ [0, 1) größtmöglich, so dass noch alle Nebenbedingungen aTi (xk + αk dk ) ≤ bi für alle i ∈ JU \ Jk erfüllt sind. Die Bestimmung von αk geschieht durch eine Fallunterscheidung:

80

Kapitel 9. SQP-Verfahren 1. Für alle i ∈ JU \ Jk mit aTi dk ≤ 0 gilt aTi (xk + αk dk ) ≤ aTi xk ≤ bi für alle αk ≥ 0. Somit verursachen diese Nebenbedingungen keine Einschränkung an αk . 2. Für alle i ∈ JU \ Jk mit aTi dk > 0 haben wir jedoch die Einschränkung αk ≤

bi − aTi xk . aTi dk

Insgesamt erhalten wir daher ) bi − aTi xk . αk = min 1, min aTi dk i6∈Jk , aT i dk >0 (

Alle Nebenbedingungen, für die das Minimum angenommen wird, bezeichnen wir als blockierende Nebenbedingungen. Falls αk < 1, legen wir zu xk+1 := xk + αk dk eine neue Indexmenge Jk+1 fest, indem wir eine oder mehrere blockierende Nebenbedingungen hinzufügen. b angelangt sind, der Diese Iteration wird solange durchgeführt, bis wir an einem Punkt x das quadratische Minimierungsproblem mit Gleichungsnebenbedingungen an der aktuellen Indexmenge Jb löst. Dies ist leicht daran zu erkennen, dass das zugehörige quadratische b = 0 besitzt. Minimierungsproblem (9.7) die Lösung d Aufgrund der Minimierungseigenschaft gilt nach Satz 7.10 für das Problem (9.7) die erste KKT-Bedingung X bi ai = 0 Hb x+g+ λ i∈Jb

bi , über deren Vorzeichen aber keine Aussage möglich ist. Setmit Lagrange-Parametern λ zen wir b λi = 0 für i ∈ JU \ Jb, so haben wir auch die erste KKT-Bedingung für das Minimierungsproblem (9.1) erfüllt: X bi ai = 0. λ Hb x+g+ i∈JG ∪JU

Die zweite und dritte KKT-Bedingungen lauten

b = bi für i ∈ JG , aTi x

b ≤ bi für i ∈ JU . aTi x

bi = 0 für alle i ∈ JU \ Jb, weshalb Sie sind nach Konstruktion ebenfalls erfüllt. Ferner ist λ auch die fünfte KKT-Bedingung bi (bi − aT x λ i b) = 0 für i ∈ JG ∪ JU

gilt. bi für Für die vierte KKT-Bedingung ist noch das Vorzeichen der Lagrange-Parameter λ b i ∈ J ∩ JU zu bestimmen. Sind diese sämtlich nichtnegativ, dann ist auch die vierte KKT-Bedingung erfüllt und es gelten alle für das Minimierungsproblem (9.1) notwendigen

9.2. Bestimmung aktiver Nebenbedingungen

81

Bedingungen aus Satz 7.10. Für den Fall, dass H positiv definit ist, folgt insbesondere b tatsächlich eine lokale Lösung des Minimierungsproblems (9.1) ist. sogar, dass x bj mit j ∈ Jb ∩ JU auf, so ist die vierte Tritt dagegen ein negativer Lagrange-Parameter λ b keine Lösung des Minimierungsproblems (9.1) KKT-Bedingung verletzt. Daher kann x sein. In diesem Fall lässt sich der Wert des quadratischen Funktionals verkleinern, indem man die zu j gehörige Nebenbedingung aus der Indexmenge Jb streicht. Dies folgt aus dem Beweis von Lemma 7.8. Man kann zeigen, dass sich auf diese Weise alle Nebenbedingungen b resultiert, mit negativen Lagrange-Parametern entfernen lassen und am Ende ein Punkt x der alle KKT-Bedingungen erfüllt. Beispiel 9.3 Wir illustrieren das beschriebene Vorgehen an folgendem konkreten Minimierungsproblem 2  5 2 unter den Nebenbedingungen min (x1 − 1) + x2 − (x1 ,x2 )∈R2 2 x1 − 2x2 + 2 ≥ 0 −x1 − 2x2 + 6 ≥ 0 −x1 + 2x2 + 2 ≥ 0 x1 ≥ 0 x2 ≥ 0. In der Notation (9.1) haben wir also       2 0 −2 x H= , g= , x= 1 0 2 −5 x2         −1 1 1 −1 a1 = , a2 = , a3 = , a4 = , 2 2 −2 0 b1 = 2,

b2 = 6,

b3 = 2,

b4 = 0,



0 a5 = −1



b5 = 0.

Wir starten das Verfahren mit   2 x0 = 0

und J0 = {3, 5}.

Das Minimierungsproblem (9.6) liefert für k = 0 als Lösung x1 = x0 , da nur dieser Punkt in der zulässigen Menge liegt und die beiden Gleichungsnebenbedingungen erfüllt. Die zugehörigen Lagrange-Parameter erhalten wir aus dem linearen Gleichungssystem " # b 2 + λ3 b3 a3 + λ b5 a5 = 0 = Hx1 + g + λ b3 − λ b5 , −5 − 2λ

b3 = −2 und b dessen Lösung λ λ5 = −1 ist. Wir streichen die Nebenbedingungen zum kleineren Lagrange-Parameter und erhalten J1 = {5}. Das Minimierungsproblem (9.6) für k = 1 liefert als Lösung   1 x2 = 0

82 b5 errechnet sich aus und der Lagrange-Parameter λ  b 0 = Hx2 + g + λ5 a5 =

Kapitel 9. SQP-Verfahren

 0 b5 , −5 − λ

b5 = −5. Wir streichen folglich auch diese Nebenbedingung und landen bei das heißt, λ J2 = ∅. Nun lösen wir das Minimierungsproblem (9.7) für k = 2, das heißt, ein Minimierungsproblem ohne Nebenbedingungen, was auf     3 1 0 und somit x3 = x2 + d2 = d2 = 3/2 5/2 5 führt. Wir addieren die blockierende Nebenbedingung zu unserer leeren Indexmenge J2 hinzu und erhalten J3 = {1}. Die Lösung des Minimierungsproblems (9.7) für k = 3 ergibt     2/5 7/5 d3 = und somit x4 = x3 + d3 = , 1/5 17/10 wobei keine weitere blockierende Nebenbedingung auftritt. Der zugehörige LagrangeParameter folgt aus " # b1 4/5 − λ b1 a1 = 0 = Hx4 + g + λ b1 , −8/5 + 2λ

b1 = 4/5. Folglich erfüllt x4 alle KKT-Bedingungen des Minimierungsproblems (9.1) also λ und ist die gesuchte Lösung. △

Index Algorithmus CG-Verfahren, 39 Gauß-Newton–Verfahren, 49 Gradientenverfahren, 10, 13, 65 Levenberg-Marquardt-Verfahren, 52 modifiziertes Verfahren von Polak und Ribière, 43 Newton-Verfahren, 15, 20, 22 nichtlineares CG-Verfahren, 40, 43 Quasi-Newton-Verfahren, 32 Trust-Region-Verfahren, 25 Verfahren von Fletcher und Reeves, 40 Verfahren von Polak und Ribière, 40 Armijo-Goldstein-Bedingung, 13, 66 Armijo-Schrittweitenregel, 13 BFGS-Verfahren, 32 Cauchy-Punkt, 25 CG-Verfahren, 39 nichtlineares, 40, 43 DFP-Verfahren, 32 Dogleg-Strategie, 25 Energienorm, 40 Folge zulässige, 56 Funktion konvexe, 6 Lagrange-, 61 quadratische, 7 Gauß-Newton-Verfahren, 48 Gradientenverfahren, 10 Grenzrichtung, 56 Hebden-Verfahren, 55 KKT-Bedingungen, 61

Konvexität, 6 gleichmäßige, 6 strikte, 6 Krylov-Raum, 39 Lagrange -Funktion, 61 -Parameter, 52, 61 LICQ-Bedingung, 57 Linearisierungskegel, 60 Minimum globales, 5 lokales, 5, 56 striktes, 5 Nebenbedingung, 56 affine, 73 aktive, 57 blockierende, 80 Gleichungs-, 56 Ungleichungs-, 56 Newton-Verfahren, 15 globalisiertes, 20 inexaktes, 22 Norm Energie-, 40 Projektion orthogonale, 65 Punkt Cauchy-, 25 stationärer, 6, 57 nicht entarteter, 74 zulässiger, 56 Quasi-Newton-Gleichung, 31 Quasi-Newton-Verfahren BFGS-Verfahren, 32 DFP-Verfahren, 32 Limited-Memory-, 36

84

Index Rang-2-Verfahren, 31 symmetrisches Rang-1-Verfahren von Broydon, 32

Satz von Karush, Kuhn und Tucker, 61 von Newton-Kantorovich, 15 SQP-Verfahren, 76 Tangentialkegel, 69 Trust-Region-Verfahren, 25, 51 Vektoren konjugierte, 34, 37 Verfahren der konjugierten Gradienten, 39 des steilsten Abstiegs, 10 Gauß-Newton-, 48 globalisiertes Newton-, 20 Gradienten-, 10 Hebden-, 55 inexaktes Newton-, 22 Levenberg-Marquardt-, 52 modifiziertes Verfahren von Polak und Ribière, 43 Newton-, 15 Quasi-Newton-, 31 BFGS-Verfahren, 32 DFP-Verfahren, 32 Limited-Memory-, 36 Rang-2-Verfahren, 31 symmetrisches Rang-1-Verfahren von Broydon, 32 SQP-, 76 Trust-Region-, 25, 51 von Fletcher und Reeves, 40 von Polak und Ribière, 40 Zielfunktion, 5