Skripts - Universität der Bundeswehr München

Lineare und nichtlineare Optimierung
AXEL DREVES
Institut für Mathematik und Rechneranwendung
Fakultät für Luft- und Raumfahrttechnik
Universität der Bundeswehr München
Werner-Heisenberg-Weg 39
85577 Neubiberg/München
E-Mail: [email protected]
Version: 7. Dezember 2015
Vorwort
Dies ist ein Skript zur Vorlesung Lineare und nichtlineare Optimierung“ im Bachelor”
studiengang Luft- und Raumfahrttechnik an der Universität der Bundeswehr München.
Wesentliche Teile des Skripts sind aus dem gleichnamigen Skript von Prof. Dr. Matthias
Gerdts von der Universität der Bundeswehr München, sowie den Skripten zu den Vorlesungen Optimierungsmethoden“ und Numerische Mathematik“ von Prof. Dr. Christian
”
”
Kanzow von der Universität Würzburg entstanden.
Empfohlene Literatur:
• C. Geiger und C. Kanzow: Numerische Verfahren zur Lösung unrestringierter Optimierungsaufgaben. Springer-Verlag, Berlin-Heidelberg-New York, 1999.
• C. Geiger und C. Kanzow: Theorie und Numerik restringierter Optimierungsaufgaben. Springer-Verlag, Berlin-Heidelberg-New York, 2002.
• M. Gerdts und F. Lempio: Mathematische Optimierungsverfahren des Operations
Research. DeGruyter, Berlin, 2011.
i
Inhaltsverzeichnis
1 Einleitung
1.1 Problemstellung . . . . . .
1.2 Graphische Darstellungen
1.3 Beispiele . . . . . . . . . .
1.4 Klassifikation . . . . . . .
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
2 Unrestringierte nichtlineare Optimierung
2.1 Optimalitätsbedingungen . . . . . . . . . .
2.2 Ein Allgemeines Abstiegsverfahren . . . .
2.3 Schrittweitenstrategien . . . . . . . . . . .
2.4 Gradientenverfahren . . . . . . . . . . . .
2.5 Newton-Verfahren . . . . . . . . . . . . . .
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
1
1
2
4
7
.
.
.
.
.
8
9
13
17
19
21
3 Restringierte nichtlineare Optimierung
30
3.1 Optimalitätsbedingungen . . . . . . . . . . . . . . . . . . . . . . . . . . . . 30
3.2 Lagrange-Newton-Verfahren . . . . . . . . . . . . . . . . . . . . . . . . . . 37
4 Lineare Optimierung
4.1 Lineare Programme . . . . . . . . . . . . . . . .
4.2 Ecken, Basisvektoren und der Fundamentalsatz
4.3 Das primale Simplexverfahren . . . . . . . . . .
4.3.1 Basiswechsel beim Simplexverfahren . . .
4.3.2 Der Algorithmus . . . . . . . . . . . . .
4.3.3 Updateformeln für das Simplextableau .
4.3.4 Phase 1 des Simplexverfahrens . . . . . .
4.3.5 Endlichkeit des Simplexverfahrens . . . .
4.4 Dualität . . . . . . . . . . . . . . . . . . . . . .
ii
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
40
40
43
53
53
57
60
62
65
68
Kapitel 1
Einleitung
1.1
Problemstellung
Seien eine Menge X ⊆ Rn und eine Funktion f : X → R gegeben. Gegenstand der
Vorlesung sind Optimierungsprobleme, d.h. wir wollen die Funktion f auf der Menge X
minimieren oder maximieren. Wir betrachten daher das Minimierungsproblem
min f (x) unter der Nebenbedingung (u.d.N.) x ∈ X,
(1.1)
oder das Maximierungsproblem
max f (x) u.d.N. x ∈ X.
(1.2)
Die Funktion f heißt dabei Zielfunktion und X heißt zulässiger Bereich. Das Problem
eine Funktion f auf der Menge X zu minimieren ist äquivalent zu dem Problem die
Funktion −f auf der Menge zu maximieren. Daher ist es keine Einschränkung, dass wir
uns nur mit Minimierungsproblemen und nicht mit Maximierungsproblemen beschäftigen.
Gilt X = Rn , so liegt ein unrestringiertes Optimierungsproblem vor, für X ( Rn
sprechen wir von einem restringierten Optimierungsproblem.
Definition 1.1.1
(a) x ∈ Rn heißt genau dann zulässig für (1.1), wenn x ∈ X.
(b) x∗ ∈ X heißt globales Minimum von (1.1), wenn f (x∗ ) ≤ f (x) für alle x ∈ X
gilt.
(c) x∗ ∈ X heißt striktes globales Minimum von (1.1), wenn f (x∗ ) < f (x) für alle
x ∈ X mit x 6= x∗ gilt.
(d) x∗ ∈ X heißt lokales Minimum von (1.1), wenn ein ε > 0 existiert, so dass
f (x∗ ) ≤ f (x) für alle x ∈ X ∩ Uε (x∗ ) mit
Uε (x∗ ) = {x ∈ Rn | kx − x∗ k < ε}
gilt.
(e) x∗ ∈ X heißt striktes lokales Minimum von (1.1), wenn ein ε > 0 existiert, so
dass f (x∗ ) < f (x) für alle x ∈ X ∩ Uε (x∗ ) mit x 6= x∗ gilt.
1
(Strikte) globale und lokale Maxima werden analog definiert. Die Begriffe werden in Abbildung 1.1 erläutert.
x1
x2
x3
x4
x5
x6
x7
Abbildung 1.1: Lokale und globale Minima und Maxima einer Funktion: x1 : striktes globales Minimum, x2 : lokales Maximum, x3 : lokales Minimum; (x2 , x3 ): gleichzeitig lokales
Minimum und Maximum, x4 : striktes globales Maximum, x5 : striktes lokales Minimum,
x6 , x7 : lokale Maxima, (x6 , x7 ): gleichzeitig lokales Minimum und Maximum.
1.2
Graphische Darstellungen
Für Raumdimensionen n = 1 oder n = 2 kann man die Zielfunktionen graphisch veranschaulichen.
So erhält man beispielsweise für die
Funktion
0.4
0.3
2
f (x) := sin(x) cos (x),
x ∈ [−4; 4]
0.2
0.1
durch das Matlab-Programm
0
x=-4:0.01:4;
fx=sin(x).*(cos(x).^2);
plot(x,fx)
−0.1
−0.2
−0.3
den nebenstehenden Graphen.
−0.4
−4
2
−3
−2
−1
0
1
2
3
4
Auch für n = 2 kann man den Graphen der Funktion darstellen lassen. So erhalten wir
für die Funktion
f (x1 , x2 ) := sin(3x2 − x21 + 1) + cos(2x22 − 2x1 )
mit Hilfe des Matlab-Programms
[x1,x2]=meshgrid(-2:0.01:2);
fx=sin(3*x2-x1.^2+1)+cos(2*x2.^2-2*x1);
mesh(x1,x2,fx)
den nebenstehenden Graphen auf dem Rechteck
[−2; 2] × [−2; 2].
Die Bestimmung der Lage der lokalen Maxima und Minima ist hier schwierig. Besser
geeignet ist ein Bild der Höhenlinien oder Niveaulinien der Funktion f zum Niveau
c ∈ R, definiert durch
Nf (c) := {x ∈ Rn | f (x) = c},
c ∈ R.
Entlang einer Höhenlinie ist die Funktion f konstant.
Für unser Beispiel ergibt sich mit dem
Matlab-Programm
2
1.5
[x1,x2]=meshgrid(-2:0.01:2);
fx=sin(3*x2-x1.^2+1)+cos(2*x2.^2-2*x1);
contour(x1,x2,fx,20)
1
0.5
0
das nebenstehende Höhenlinienbild, bei
dem sich die Lage der Extremwerte schon
besser sehen lässt.
−0.5
−1
−1.5
−2
−2
3
−1.5
−1
−0.5
0
0.5
1
1.5
2
Als Warnung sei bemerkt, dass Bilder zur Anschauung dienen, aber auch irreführend sein
können.
1.3
Beispiele
Optimierungsprobleme spielen in vielen Gebieten der Mathematik eine wichtige Rolle
und haben Anwendungen bei wirtschaftswissenschaftlichen, naturwissenschaftlichen sowie
technischen Fragestellungen. Physikalische Prozesse werden häufig durch das Prinzip der
Energieminimierung gesteuert, so dass die Optimierung die mathematische Modellierung
von physikalischen Problemen direkt beeinflusst. Es gibt eine Vielzahl verschiedenster
Optimierungsprobleme und wir wollen im Folgenden ein paar wenige konkrete Beispiele
geben.
Beispiel 1.3.1 (Lineares Optimierungsproblem)
Ein Landwirt bewirtschaftet ein Grundstück von 40 Hektar Größe mit Zuckerrüben und
Weizen. Er kann hierzu 2400 Euro und 312 Arbeitstage einsetzen. Pro Hektar betragen
seine Anbaukosten bei Rüben 40 Euro und bei Weizen 120 Euro. Für Rüben benötigt
er 6 Arbeitstage, für Weizen 12 Arbeitstage pro Hektar. Der Reingewinn bei Rüben sei
100 Euro pro Hektar, bei Weizen sei er 250 Euro pro Hektar. Wie muss der Bauer sein
Grundstück bewirtschaften, um einen maximalen Gewinn zu erzielen?
Mathematische Formulierung: Wir bezeichnen mit x1 die Fläche, die mit Rüben bepflanzt
wird und mit x2 die Fläche, die mit Weizen bepflanzt wird. Der zu maximierende Gewinn
lautet
f (x1 , x2 ) = 100x1 + 250x2 =: c> x.
Aus der Aufgabenstellung lassen sich folgende Beschränkungen ableiten:
Grundstücksgröße:
Geld:
Arbeitstage:
keine negativen Flächen:
g1 (x1 , x2 ) := x1 + x2 ≤ 40
g2 (x1 , x2 ) := 40x1 + 120x2 ≤ 2400
g3 (x2 , x2 ) := 6x1 + 12x2 ≤ 312
x1 , x2 ≥ 0
Der zulässige Bereich des Optimierungsproblems ist durch den blauen Bereich in der
folgenden Abbildung gegeben. Die rote Gerade stellt die Höhenlinie der Zielfunktion zum
Niveau 3500, die grüne diejenige zum Niveau 5500 dar.
4
x2
40
30
20
100x 1 +250x 2 =5500
10
40x 1 +120x 2 =2400
10
20
30
100x 1 +250x 2 =3500
40
x 1 +x 2 =40
-10
50
60
70
x1
6x 1 +12x 2 =312
Beobachtungen:
• Die Höhenlinien einer affin linearen Funktion sind Geraden!
• Lineare Optimierungsprobleme lassen sich grafisch lösen, indem die Höhenlinie der
Zielfunktion in Richtung wachsender Niveaus (bei Maximierungsaufgaben) bzw. fallender Niveaus (bei Minimierungsaufgaben) bis an den Rand des zulässigen Bereichs
verschoben“ wird. Konkret:
”
(1) Skizziere den zulässigen Bereich.
(2) Wähle einen beliebigen Punkt x̂, setze ihn in die Zielfunktion ein und berechne
w = c> x̂. Skizziere die durch c> x = w gegebene Gerade.
(3) Bewege die Gerade aus (2) in Richtung −c bei Minimierungsaufgaben bzw. in
Richtung c bei Maximierungsaufgaben solange der Durchschnitt der Geraden
und des zulässigen Bereichs nicht leer ist.
(4) Die extremste Gerade aus (3) ist optimal. Alle zulässigen Punkte auf dieser
Geraden sind optimal. Der optimale Zielfunktionswert kann in der Zeichnung
5
abgelesen werden.
Durch Ablesen ergibt sich hier die grafische Lösung x1 = 30, x2 = 10 mit Zielfunktionswert f (x1 , x2 ) = 5500.
• Das Optimum wird in einer Ecke“ des zulässigen Bereichs angenommen.
”
Beispiel 1.3.2 (Transportproblem)
Ein Transportunternehmer hat m Vorratslager, aus denen n Verbraucher mit einem Produkt beliefert werden können. Die Lieferkosten von Lager i zu Verbraucher j betragen
cij Einheiten pro Produkteinheit. In Lager i sind ai Einheiten des Produktes vorrätig.
Verbraucher j hat einen Bedarf von bj Einheiten. Um die Kunden nicht zu verärgern,
muss der Lieferant den Bedarf der Kunden befriedigen. Andererseits möchte der Lieferant
seine Lieferkosten minimieren.
Lieferung
a1
Verbraucher
Vorratslager
b1
a2
bn
am
Bezeichnet xij die Liefermenge von Lager i zu Verbraucher j, so führt das Problem auf
das folgende Transportproblem, welches ein spezielles lineares Optimierungsproblem ist:
min
m X
n
X
i=1 j=1
cij xij u.d.N.
n
X
j=1
m
X
xij ≤ ai ,
i = 1, . . . , m,
xij ≥ bj ,
j = 1, . . . , n,
i=1
xij ≥ 0,
i = 1, . . . , m, j = 1, . . . , n.
Die erste Nebenbedingung besagt, dass aus Lager i maximal ai Einheiten abtransportiert
werden können. Die zweite Nebenbedingung besagt, dass der Bedarf bj befriedigt werden
muss. Die letzte Nebenbedingung verbietet negative Liefermengen.
Beispiel 1.3.3 (Parameterschätzung bei Differentialgleichungen)
Sei y(t) die Population von Tieren zur Zeit t. Die Population möge eine maximale Größe
6
K > 0 nicht überschreiten können und entwickle sich proportional zur aktuellen Größe
y(t) sowie zum verbleibenden Rest“ K − y(t). Dies liefert die Wachstumsgleichung
”
y 0 (t) = λy(t)(K − y(t)).
Seien yi Messungen der exakten Werte y(ti ) zu gewissen Zeitpunkten ti , i = 1, . . . , m.
Sei y(·; λ, K) die exakte Lösung der Wachstumsgleichung. Um die Parameter der Differentialgleichung möglichst optimal zu schätzen, sucht man eine Lösung des nichtlinearen
Optimierungsproblems
m
1X
(y(ti ; λ, K) − yi )2 .
min
λ,K>0 2
i=1
1.4
Klassifikation
Meistens ist die zulässige Menge des Optimierungsproblems durch Gleichheits– und Ungleichheitsrestriktionen beschrieben, d.h. es gibt Funktionen h : Rn → Rp , g : Rn → Rm ,
so dass
X := {x ∈ Rn | h(x) = 0, g(x) ≤ 0},
wobei g(x) ≤ 0 als Ungleichungen für jede Komponente, also gi (x) ≤ 0 für alle i =
1, . . . , m, zu verstehen ist. Nun kann man das Minimierungsproblem
min f (x) u.d.N. x ∈ X
(1.3)
in eine der 4 Kategorien einsortieren:
p = 0, m = 0 unrestringiert
p = 0, m > 0 nur Ungleichheitsrestriktionen
p > 0, m = 0 nur Gleichheitsrestriktionen
p > 0, m > 0 Gleichheits– und Ungleichheitsrestriktionen
Bei der nichtlinearen Optimierung können die Funktionen f, g, h in (1.3) beliebig sein.
Im Gegensatz dazu steht die lineare Optimierung, bei der man fordert, dass f, g, h
affin lineare Funktionen sind. Weitere häufig betrachtete Spezialfälle der nichtlinearen
Optimierung sind die quadratische Optimierung (f quadratische Funktion, g, h affin
lineare Funktionen), konvexe Optimierung (f und alle gi sind konvexe Funktionen, h
ist eine affin lineare Funktion), oder die linear restringierte Optimierung (g, h sind
affin lineare Funktionen).
7
Kapitel 2
Unrestringierte nichtlineare Optimierung
Wir beschäftigen uns in diesem Abschnitt mit dem Problem das Minimum einer Funktion
f : Rn → R auf ganz Rn zu finden, d.h. wir betrachten das Problem
min f (x) s.t. x ∈ Rn .
Wir wollen uns zunächst einige Bezeichnungen in Erinnerung rufen:
• Der Gradient einer Funktion f : Rn → R an der Stelle x = (x1 , . . . , xn )> ∈ Rn ist
definiert als der Spaltenvektor


∇f (x1 , . . . , xn ) := 
∂f
(x1 , . . . , xn )
∂x1


..
.
.
∂f
(x1 , . . . , xn )
∂xn
Im Spezialfall n = 1 ist der Gradient gerade die erste Ableitung von f . Bekanntlich
zeigt der Gradient einer Funktion f in die Richtung des steilsten Anstiegs von
f . Außerdem steht der Gradient senkrecht auf den Höhenlinien von f .
• Die Hessematrix von f an der Stelle x ist gegeben durch
 2
∂ f (x)
∂ 2 f (x)
· · · ∂x
∂x1 ∂x1
1 ∂xn

..
..
..
Hf (x) = ∇2 f (x) = 
.
.
.

∂ 2 f (x)
∂ 2 f (x)
· · · ∂xn ∂xn
∂xn ∂x1


.

Im Spezialfall n = 1 ist die Hessematrix gerade die zweite Ableitung von f . Die
Hessematrix beschreibt anschaulich die lokale Krümmung einer Funktion.
• Die Richtungsableitung f 0 (x; d) der Funktion f : Rn → R an der Stelle x ∈ Rn
in Richtung d ∈ Rn ist definiert durch
f 0 (x; d) := lim
t↓0
f (x + td) − f (x)
.
t
Für stetig differenzierbare Funktionen gilt
f 0 (x; d) = ∇f (x)> d.
8
• Die mehrdimensionale Taylorsche Formel liefert für eine dreimal stetig differenzierbare Funktion f : Rn → R eine Zwischenstelle ξ auf der Verbindungsstrecke des
Entwicklungspunktes x und eines beliebigen Punktes x + d, so dass
1
f (x + d) = f (x) + ∇f (x)> d + d> ∇2 f (ξ)d.
2
• Vektorungleichungen wollen wir stets komponentenweise interpretieren, d.h. für zwei
Vektoren x, y ∈ Rn bedeutet x ≤ y, dass xi ≤ yi für alle i = 1, . . . , n gilt.
2.1
Optimalitätsbedingungen
Als erstes wollen wir nun Optimalitätsbedingungen für ein lokales Minimum herleiten.
Wir unterscheiden dabei
• notwendige Bedingungen, d.h. Bedingungen, die ein lokales Minimum zwangsläufig
erfüllen muss, und
• hinreichenden Bedingungen, d.h. Bedingungen, mit denen entschieden werden
kann, ob ein lokales Minimum vorliegt oder nicht.
Satz 2.1.1 (Notwendige Bedingung 1. Ordnung)
Sei f : Rn → R stetig differenzierbar. Ist x̂ ein lokales Minimum von f , dann gilt ∇f (x̂) =
0.
Beweis: Beweis durch Widerspruch: Angenommen, es gilt ∇f (x̂) 6= 0. Dann gibt es
einen Vektor d ∈ Rn mit ∇f (x̂)> d < 0. Damit ist
lim
t↓0
f (x̂ + td) − f (x̂)
= f 0 (x̂; d) = ∇f (x̂)> d < 0.
t
Es folgt für alle hinreichend kleinen t > 0 muss f (x̂ + td) − f (x̂) < 0 und damit
f (x̂ + td) < f (x̂)
gelten, ein Widerspruch zu lokalen Minimalität von x̂.
Der Spezialfall n = 1 ist die notwendige Bedingung 1. Ordnung wohlbekannt und besagt,
dass in einem lokalen Minimum f 0 (x̂) = 0 gilt. Da aber auch in einem lokalen Maximum
∇f (x̂) = 0 gilt, ist die Bedingung nicht hinreichend für ein lokales Minimum.
Beispiel 2.1.2
Betrachten wir die Funktionen f1 (x) = x2 , f2 (x) = −x2 und f3 (x) = x3 . Dann gilt
fi0 (0) = 0 für i = 1, 2, 3 und die Funktion f1 hat in x̂ = 0 ein globales Minimum, f2
9
hat dort ein globales Maximum, und f3 hat weder ein lokales Minimum noch ein lokales
Maximum.
Definition 2.1.3 (stationärer Punkt)
Jeder Punkt x mit ∇f (x) = 0 heißt stationärer Punkt von f .
Wir wollen nun noch ein weiteres notwendiges Kriterium herleiten.
Satz 2.1.4 (Notwendige Bedingung 2. Ordnung)
Sei f : Rn → R zweimal stetig differenzierbar. Ist x̂ ein lokales Minimum von f , dann ist
die Hessematrix Hf (x̂) positiv semidefinit, d.h. es gilt für alle d ∈ Rn ist d> Hf (x̂) d ≥ 0.
Beweis: Beweis durch Widerspruch: Angenommen, Hf (x̂) ist nicht positiv semidefinit.
Dann gibt es ein d ∈ Rn mit d> Hf (x̂) d < 0. Verwenden wir nun den Satz von Taylor, so
erhalten wir
1
f (x̂ + td) = f (x̂) + t∇f (x̂)> d + t2 d> Hf (ξt ) d
2
für eine von t abhängige Zwischenstelle ξt ∈ [x̂, x̂+td]. Da x̂ ein lokales Minimum ist, liefert
die notwendige Bedingung 1. Ordnung ∇f (x̂) = 0. Da f zweimal stetig differenzierbar ist,
ist Hf (x) stetig, und somit folgt aus d> Hf (x̂) d < 0, dass auch d> Hf (ξt ) d < 0 für alle
hinreichend kleinen t > 0 gilt. Damit folgt aber
f (x̂ + td) < f (x̂)
für alle hinreichend kleinen t > 0, ein Widerspruch zur Minimalität von x̂.
Die Überprüfung, ob eine Matrix A positiv oder auch negativ (semi-)definit ist, ist im Fall
von symmetrischen Matrizen, d.h. für A> = A mittels der Eigenwerte möglich. Dabei gilt
sämtliche Eigenwerte von A sind positiv ⇔ A ist positiv definit
sämtliche Eigenwerte von A sind negativ ⇔ A ist negativ definit
sämtliche Eigenwerte von A sind ≥ 0 ⇔ A ist positiv semidefinit
sämtliche Eigenwerte von A sind ≤ 0 ⇔ A ist negativ semidefinit
A besitzt positive und negative Eigenwerte ⇔ A ist indefinit
Beispiel 2.1.5
Betrachten wir erneut die Funktionen f1 (x) = x2 , f2 (x) = −x2 und f3 (x) = x3 . Dann
gilt f10 (0) = 0, f100 (0) = 2 > 0. Damit erfüllt f1 die notwendigen Bedingungen 1. und
2. Ordnung. Die Funktion f2 hat die Hessematrix f200 (0) = −2 < 0 und ist damit negativ
definit. Sie verletzt die notwendige Bedingung 2. Ordnung und hat damit in x̂ = 0 kein
lokales Minimum. Die Funktion f3 erfüllt wegen f30 (0) = 0, f300 (0) = 0 beide notwendi10
gen Bedingungen. Da x̂ = 0 aber kein lokales Minimum ist, ist keine der notwendigen
Bedingungen hinreichend.
Wir haben also gesehen, dass alle Punkte, die die notwendigen Bedingungen 1. Ordnung
(∇f (x) = 0) und 2. Ordnung (Hf (x) positiv semidefinit) erfüllen, lediglich Kandidaten
für lokale Minima darstellen. Um nun feststellen zu können, ob tatsächlich ein lokales
Minimum vorliegt, benötigen wir hinreichende Bedingungen.
Satz 2.1.6 (Hinreichende Bedingung 2. Ordnung)
Sei f : Rn → R zweimal stetig differenzierbar und x̂ ein stationärer Punkt mit positiv
definiter Hessematrix Hf (x̂). Dann ist x̂ ein striktes lokales Minimum von f .
Beweis: Die stetige Funktion d> Hf (x̂) d nimmt auf der kompakten Menge {d ∈ Rn |
kdk = 1} ihr Minimum an. Wegen der positiven Definitheit von Hf (x̂) gibt es daher ein
µ > 0 mit
>
d
d
Hf (x̂)
≥µ
kdk
kdk
für alle d ∈ Rn \ {0}. Damit erhält man
d> Hf (x̂) d ≥ µkdk2
für alle d ∈ Rn . Mit dem Satz von Taylor folgt die Existenz eines ξt ∈ [x̂, x̂ + td] mit
1
f (x̂ + td) = f (x̂) + t ∇f (x̂)> d + t2 d> Hf (ξt ) d
| {z }
2
=0
1
1
= f (x̂) + t2 d> Hf (x̂)d + t2 d> (Hf (ξt ) − Hf (x̂))d
2
2
1 2
≥ f (x̂) + t kdk2 (µ − kHf (ξt ) − Hf (x̂)k),
2
wobei die letzte Ungleichung mit Hilfe der Cauchy-Schwarz Ungleichung folgt. Wegen der
Stetigkeit der Hessematrix ist µ − kHf (ξt ) − Hf (x̂)k > 0 für alle hinreichend kleinen t > 0
und somit folgt
f (x̂ + td) > f (x̂)
für alle hinreichend kleinen t > 0 und somit ist x̂ ein striktes lokales Minimum.
Im Spezialfall n = 1 besagt diese Resultat, dass wenn f 0 (x̂) = 0 und f 00 (x̂) > 0 gilt, ein
striktes lokales Minimum vorliegt.
Die hinreichende Bedingung 2. Ordnung ist nicht notwendig, wie das folgende Beispiel
zeigt.
Beispiel 2.1.7
Die Funktion f (x) = x4 hat in x̂ = 0 ein striktes lokales Minimum. x̂ ist zwar ein stati11
onärer Punkt, die Hessematrix Hf (0) = 0 ist aber nicht positiv definit (nur semidefinit)
und erfüllt daher das hinreichende Kriterium 2. Ordnung nicht.
Wir wollen nun noch einige mögliche Fälle im R2 graphisch veranschaulichen:
f (x, y) = x2 + y 2 !
2 0
0 2
Hf (0, 0) =
f (x, y) = −x2 − y 2 !
−2 0
Hf (0, 0) =
0 −2
positiv definit
8
0
6
−2
4
−4
2
−6
0
2
−8
2
1
2
1
0
2
1
0
−1
−1
−2 −2
f (x, y) = x2 − y 2 !
2 0
Hf (0, 0) =
0 −2
1
0
0
−1
negativ definit
−1
−2 −2
indefinit
4
2
0
−2
−4
2
2
1
1
0
0
−1
−1
−2 −2
Abbildung 2.1: Lokales Minimum, lokales Maximum, Sattelpunkt
12
Nun wollen wir noch spezielle quadratische Funktionen untersuchen.
Satz 2.1.8
Eine quadratische Funktion
1
f (x) = x> Qx + c> x + γ
2
mit einer symmetrisch, positiv definiten Matrix Q ∈ Rn×n und c ∈ Rn , γ ∈ R hat genau
ein (lokales = globales) Minimum in x̂ = −Q−1 c.
Beweis: Durch Nachrechnen zeigt man, dass ∇f (x) = Qx + c und Hf (x) = Q gilt. Da
Q als positiv definite Matrix regulär ist, liefert die notwendige Bedingung 1. Ordnung x̂ =
−Q−1 c als einzigen Kandidaten für ein lokales Minimum, und die hinreichende Bedingung
2. Ordnung zeigt, dass x̂ striktes lokales Minimum ist. Weiter gilt für alle x 6= x̂:
1 >
x Qx + c> x + γ
2
1
1
=
(x − x̂)> Q(x − x̂) − x̂> Qx̂ + x> Qx̂ + c> x + γ
2
2
1
1
(x − x̂)> Q(x − x̂) − c> Q−1 c − c> x + c> x + γ
=
{z
} 2
2|
f (x) =
>0
1
> − c> Q−1 c + γ
2
= f (x̂),
also ist x̂ = −Q−1 c auch globales Minimum von f.
2.2
Ein Allgemeines Abstiegsverfahren
In diesem Abschnitt konstruieren wir ein allgemeines Konzept zur Minimierung einer
stetig differenzierbaren Funktion f : Rn → R. Hierzu definieren wir zunächst eine Abstiegsrichtung.
Definition 2.2.1 (Abstiegsrichtung)
Seien f : Rn → R und x ∈ Rn . d ∈ Rn heißt eine Abstiegsrichtung von f in x, falls
es ein t̄ > 0 gibt mit
f (x + td) < f (x) für alle t ∈ (0, t̄].
Die Definition einer Aufstiegsrichtung ist analog. Eine hinreichende Bedingung für eine
Abstiegsrichtung liefert das folgende Lemma.
Lemma 2.2.2
Sei f : Rn → R stetig differenzierbar in x ∈ Rn . d ∈ Rn ist eine Abstiegsrichtung von f
13
in x, wenn ∇f (x)> d < 0 gilt.
Beweis: Nach Definition der Richtungsableitung von f in x in Richtung d ist
lim
t↓0
f (x̂ + td) − f (x̂)
= f 0 (x̂; d) = ∇f (x̂)> d < 0.
t
Hieraus folgt die Existenz eines t̄ > 0 mit f (x + td) − f (x) < 0 für alle 0 < t ≤ t̄, d.h. d
ist Abstiegsrichtung.
Bemerkung 2.2.3
• Die Bedingung ∇f (x)> d < 0 bedeutet geometrisch, dass der Winkel zwischen Gradient und Abstiegsrichtung zwischen 90◦ und 270◦ liegt. Dies lässt sich aus der für
Vektoren a, b ∈ Rn allgemein gültigen Beziehung
cos ∠ (a, b) =
a> b
ha, bi
=
kak · kbk
kak · kbk
ableiten.
• Die Bedingung in Lemma 2.2.2 ist nur hinreichend für eine Abstiegsrichtung aber
nicht notwendig. Betrachte z.B. ein striktes lokales Maximum x̂ mit ∇f (x̂) = 0. Für
alle Richtungen d ∈ Rn gilt dann ∇f (x̂)> d = 0. Andererseits ist in einem strikten
lokalen Maximum jede Richtung d ∈ Rn eine Abstiegsrichtung.
Ist f : Rn → R eine stetig differenzierbare Funktion und x ∈ Rn kein stationärer Punkt
(also ∇f (x) 6= 0), so ist jedes d := −M ∇f (x) mit einer positiv definiten Matrix M ∈ Rn×n
eine Abstiegsrichtung von f in x, da
∇f (x)> d = −∇f (x)> M ∇f (x) < 0.
Algorithmus 2.2.4 (Allgemeines Abstiegsverfahren)
(S.0) Wähle einen Startpunkt x0 ∈ Rn und setze k := 0.
(S.1) Genügt xk einem geeigneten Abbruchkriterium: STOP.
(S.2) Bestimme eine Abstiegsrichtung dk von f in xk
(S.3) Berechne eine Schrittweite tk > 0 mit f (xk + tk dk ) < f (xk ).
(S.4) Setze xk+1 := xk + tk dk , k := k + 1 und gehe zu (S.1).
14
Dieser Algorithmus hat die Wahl der Abstiegsrichtung dk und die Wahl der Schrittweite
tk als Freiheitsgrade. Wie wir gesehen haben, ist die Richtung d = −∇f (xk ) ein Beispiel
einer Abstiegsrichtung, da die Einheitsmatrix M = I positiv definit ist. Wir wollen nun
in einem Beispiel sehen, dass die Wahl der Schrittweite nicht beliebig sein kann.
Beispiel 2.2.5
Wir wenden Algorithmus 2.2.4 auf die Funktion f (x) = x2 , die Abstiegsrichtung dk =
1
−f 0 (xk ) = −2xk , und den Startpunkt x0 = 1 an. Als Schrittweite wählen wir tk := 2k+3
.
Wir erhalten
1
k+1
k
k
k
x
= x − 2tk x = (1 − 2tk )x = 1 − k+2 xk .
2
Mit x0 = 1 folgt xk ∈ (0, 1) für alle k ∈ N und weiter
x0 − x k =
k−1
k−1
k−1
X
X
1
1 i X 1
x
≤
≤ .
(xi − xi+1 ) =
i+2
i+2
2
2
2
i=0
i=0
i=0
Also ist |x0 −xk | = |1−xk | ≤ 12 , und die Folge xk kann nicht gegen den einzigen stationären
Punkt x̂ = 0 konvergieren.
Definition 2.2.6 (Schrittweitenstrategie)
Jede Abbildung T : Rn × Rn → 2R++ (2R++ ist die Potenzmenge, also die Menge aller
Teilmengen, der positiven reellen Zahlen R++ ), (x, d) 7→ T (x, d) heißt Schrittweitenstrategie. Sie heißt wohldefiniert, wenn T (x, d) 6= ∅ für alle x, d mit ∇f (x)> d 6= 0 gilt.
Wie das Beispiel 2.2.5 gezeigt hat, sind nicht alle Schrittweitenstrategien geeignet. Vielmehr müssen sie eine hinreichend große Abnahme des Funktionswerts sicherstellen.
Definition 2.2.7 (Effiziente Schrittweiten)
Seien f : Rn → R stetig differenzierbar, x ∈ Rn und d ∈ Rn eine Abstiegsrichtung von f
in x. Eine Schrittweitenstrategie heißt effizient, wenn es eine von x und d unabhängige
Konstante θ > 0 gibt, so dass
2
∇f (x)> d
f (x + td) ≤ f (x) − θ
kdk
für alle t ∈ T (x, d) gilt. Jedes t ∈ T (x, d) wird dann als effiziente Schrittweite bezeichnet.
Da eine Funktion entlang ihrer Höhenlinien konstant ist, sollte man die Abstiegsrichtung möglichst so wählen, dass man sich nicht nahezu entlang der derzeitigen Höhenlinie
bewegt. Daher definieren wir noch eine Winkelbedingung.
15
Definition 2.2.8 (Winkelbedingung)
Seien f : Rn → R stetig differenzierbar, x ∈ Rn und d ∈ Rn eine Abstiegsrichtung von f
in x. Die Richtung d erfüllt die Winkelbedingung, wenn es eine Konstante C > 0 gibt, so
dass
∇f (xk )> dk
≥C
−
k∇f (xk )k · kdk k
für alle k ∈ N gilt.
Anschaulich besagt die Winkelbedingung, dass der Winkel zwischen −∇f (xk ) und dk
immer zwischen 90◦ und −90◦ liegt und dabei gleichmäßig (also unabhängig von k) von
90◦ und −90◦ wegbleibt, vgl. Abbildung 2.2.
xk
α
α
{x | f (x) = f (xk )
dk
−∇f (xk )
Abbildung 2.2: Winkelbedingung für die Suchrichtungen.
Satz 2.2.9 (Konvergenzsatz für das Allgemeine Abstiegsverfahren)
Seien f : Rn → R stetig differenzierbar und {xk } eine durch den Algorithmus 2.2.4
erzeugte Folge, so dass die Schrittweiten tk > 0 effizient sind und die Abstiegsrichtungen
die Winkelbedingung erfüllen. Dann ist jeder Häufungspunkt der Folge {xk } ein stationärer
Punkt von f.
Beweis: Verwenden wir zuerst die Effizienzbedingung und dann die Winkelbedingung
so gilt
f (xk+1 ) = f (xk + tk dk )
2
∇f (xk )> dk
k
≤ f (x ) − θ
kdk k
≤ f (xk ) − θC 2 k∇f (xk )k2 .
Sei nun x̂ ein Häufungspunkt von {xk }. Dann konvergiert, da f stetig ist, {f (xk )} auf einer
Teilfolge gegen f (x̂), und da die Folge {f (xk )} monoton fallend ist, konvergiert die gesamte
16
Folge {f (xk )} gegen f (x̂). Damit folgt f (xk+1 ) − f (xk ) → 0 und somit k∇f (xk )k2 → 0.
Daher ist ∇f (x̂) = 0, also ist jeder Häufungspunkt ein stationärer Punkt von f.
2.3
Schrittweitenstrategien
In diesem Abschnitt geht es um die Realisierung von Schritt (3) von Algorithmus 2.2.4,
der Bestimmung einer Schrittweite.
Sei f : Rn → R stetig differenzierbar, und seien β, σ ∈ (0, 1). Zu gegebenem x, d ∈ Rn mit
∇f (x)> d < 0 bestimmt man
t := max{β j | j = 0, 1, 2, . . .},
so dass die Armijo-Bedingung
f (x + td) ≤ f (x) + σt∇f (x)> d
(2.1)
erfüllt ist. Diese Schrittweitenstrategie nennt man Armijo-Regel, und dabei besteht
T (x, d) nur aus einem Element.
Definieren wir ϕ(t) := f (x + td), so lautet die Armijo-Bedingung
ϕ(t) ≤ ϕ(0) + σtϕ0 (0).
Anschaulich ergibt sich daher folgendes Bild:
ϕ(t)
ϕ(0) + σ · t · ϕ0 (0)
t=0
t
ϕ(0) + t · ϕ0 (0)
Abbildung 2.3: Schrittweiten, die der Armijo-Bedingung genügen.
Die fett eingezeichneten Intervalle in Abbildung 2.3 erfüllen die Armijo-Bedingung (2.1).
Die Armijo-Regel ist, wie der folgende Satz zeigt wohldefiniert.
17
Satz 2.3.1
Sei f : Rn → R stetig differenzierbar, und seien β, σ ∈ (0, 1). Zu x, d ∈ Rn mit ∇f (x)> d <
0 existiert stets ein endliches j ∈ N mit
f (x + β j d) ≤ f (x) + σβ j ∇f (x)> d.
Beweis: Beweis durch Widerspruch: Angenommen, für alle j ∈ N gilt
f (x + β j d) > f (x) + σβ j ∇f (x)> d.
Dann folgt
f (x + β j d) − f (x)
> σ∇f (x)> d.
βj
Da f stetig differenzierbar ist folgt durch den Grenzübergang j → ∞
∇f (x)> d ≥ σ∇f (x)> d,
also (1 − σ)∇f (x)> d ≥ 0 ein Widerspruch zu σ ∈ (0, 1) und ∇f (x)> d < 0.
Die Realisierung der Armijo-Regel ist einfach und kann mit dem folgenden Algorithmus
erfolgen:
Algorithmus 2.3.2 (Armijo-Regel)
(S.0) Wähle β ∈ (0, 1), σ ∈ (0, 1) und setze t := 1.
(S.1) Falls die Bedingung
ϕ(t) ≤ ϕ(0) + σ · t · ϕ0 (0)
erfüllt ist, setze tk := t und beende das Verfahren.
Andernfalls gehe zu (S.2).
(S.2) Setze t := β · t und gehe zu (S.1).
Bemerkung 2.3.3
• Die Armijo-Schrittweite liegt immer in (0, 1]. Da die Armijo-Bedingung aber auch
durch größere Schrittweiten erfüllt werden kann, gibt es Varianten, wie die skalierte
Armijo-Regel, bei der für ein Skalierungsfaktor s > 0 die Schrittweite mittels
t = max{sβ j | j = 0, 1, 2, . . .}
bestimmt wird.
18
• Die Armijo-Regel ist im Allgemeinen nicht effizient.
• Eine weitere gängige Schrittweitenstrategien ist die Wolfe-Powell-Regel. Hier hat
man für eine Konstante ρ ∈ [σ, 1) zusätzlich zur Armijo-Bedingung (2.1) noch die
Forderung
∇f (x + td)> d ≥ ρ∇f (x)> d.
• Bei der strengen Wolfe-Powell-Regel hat man für eine Konstante ρ ∈ [σ, 1)
zusätzlich zur Armijo-Bedingung (2.1) noch die Forderung
|∇f (x + td)> d| ≤ −ρ∇f (x)> d.
• Die (strenge) Wolfe-Powell-Regel ist wohldefiniert, falls f nach unten beschränkt
ist. Ist ferner ∇f auf der Levelmenge Lf (x0 ) := {x ∈ Rn | f (x) ≤ f (x0 )} Lipschitzstetig, d.h. gibt es eine Konstante L > 0, so dass
k∇f (y) − ∇f (z)k ≤ Lky − zk
für alle y, z ∈ Lf (x0 ),
so ist die (strenge) Wolfe-Powell-Regel auch effizient.
• Die exakte Liniensuche, bei der man t so bestimmt, dass
f (x + td) = min{f (x + αd) | α > 0},
ist in der Regel zu teuer für praktische Zwecke. Sie wird häufiger bei theoretischen
Untersuchungen verwendet.
• Ist der Gradient ∇f auf der Levelmengen Lf (x0 ) Lipschitz-stetig, und ist die Levelmenge kompakt, so ist die exakte Liniensuche wohldefiniert und effizient.
2.4
Gradientenverfahren
Für eine stetig differenzierbare Funktion f : Rn → R und ein x ∈ Rn mit ∇f (x) 6= 0 ist
die Lösung von
min ∇f (x)> d u.d.N. kdk ≤ 1
gegeben durch
∇f (x)
.
k∇f (x)k
Daher zeigt der negative Gradient in Richtung des steilsten Abstiegs von f. Wählen
wir also als Abstiegsrichtung die negative Gradientenrichtung und verwenden die ArmijoRegel zur Schrittweitenbestimmung, so erhalten wir das folgende Gradientenverfahren,
das auch Verfahren des steilsten Abstiegs genannt wird.
Algorithmus 2.4.1 (Gradientenverfahren mit Armijo-Regel)
d=−
19
(S.0) Wähle x0 ∈ Rn , σ, β ∈ (0, 1), ε ≥ 0 und setze k := 0.
(S.1) Ist k∇f (xk )k ≤ ε: STOP.
(S.2) Berechne dk = −∇f (xk ) und bestimme tk := max{β j | j = 0, 1, 2, . . .} mit
f (xk + tk dk ) ≤ f (xk ) + σtk ∇f (xk )> dk .
(S.3) Setze xk+1 := xk + tk dk , k := k + 1 und gehe zu (S.1).
Für das Gradientenverfahren mit Armijo-Regel gilt der folgende Konvergenzsatz:
Satz 2.4.2
Ist f : Rn → R stetig differenzierbar, so ist jeder Häufungspunkt einer durch das Gradientenverfahren mit Armijo-Regel und ε = 0 erzeugten Folge {xk } ein stationärer Punkt
von f.
Beweis: Beweis durch Widerspruch: Sei x̂ ein Häufungspunkt von {xk }. Angenommen,
es gilt ∇f (x̂) 6= 0. Da f stetig ist, konvergiert {f (xk )} auf einer Teilfolge gegen f (x̂), und
da die Folge {f (xk )} monoton fallend ist, konvergiert die gesamte Folge {f (xk )} gegen
f (x̂). Damit folgt f (xk+1 ) − f (xk ) → 0. Mit der Armijo-Bedingung und dk = −∇f (xk )
folgt
f (xk+1 ) − f (xk ) ≤ σtk ∇f (xk )> dk = −σtk k∇f (xk )k2 ≤ 0,
und somit tk k∇f (xk )k2 → 0. Aus der Annahme, dass ∇f (x̂) 6= 0 ist und der Stetigkeit
von ∇f , erhalten wir tk → 0. Setzt man tk = β jk für ein jk ∈ N in Schritt (S.2), so folgt
für hinreichend große k ∈ N
f (xk + β jk −1 dk ) − f (xk ) > σβ jk −1 ∇f (xk )> dk .
Nach dem Mittelwertsatz gibt es ein ξ k ∈ [xk , xk + β jk −1 dk ] mit
f (xk + β jk −1 dk ) − f (xk ) = β jk −1 ∇f (ξ k )> dk .
Somit ergibt sich
β jk −1 ∇f (ξ k )> dk > σβ jk −1 ∇f (xk )> dk
und nach Kürzen von β jk −1 , einsetzen von dk = −∇f (xk ) und Grenzübergang k → ∞
erhalten wir, wegen der Stetigkeit von ∇f
k − ∇f (x̂)k2 ≥ −σk∇f (x̂)k2 .
20
Aus σ ∈ (0, 1) ergibt sich ein Widerspruch zu ∇f (x̂) 6= 0.
Beispiel 2.4.3 (Rosenbrock-Funktion)
Das Gradienten-Verfahren ist im Allgemeinen nur sehr langsam konvergent. Das Beispiel
der sogenannten Rosenbrock-Funktion
f (x1 , x2 ) = 100(x2 − x21 )2 + (1 − x1 )2
verdeutlicht dies. Mit der Parameterwahl β = 0.5, σ := 10−4 und dem Abbruchparameter
ε = 10−4 benötigt das Verfahren 8058 Iterationen und dabei 80045 Funktionsauswertungen um vom Startvektor x0 = (−1.2, 1) beginnend das Abbruchkriterium zu erfüllen.
Die gefundene Lösung x8058 = (0.9999, 0.9998) liegt dann aber in der Nähe des exakten
Minimums x̂ = (1, 1).
2.5
Newton-Verfahren
Um eine Funktion f zu minimieren wollen wir sie in der Nähe des aktuellen Iterationspunktes xk durch eine quadratische Funktion approximieren. Motiviert durch die Taylorentwicklung definieren wir hierzu
1
qk (x) := f (xk ) + ∇f (xk )> (x − xk ) + (x − xk )> Hf (xk )(x − xk ).
2
Nun bestimmen wir die neue Iterierte xk+1 als Minimum der Funktion qk (x) für x ∈ Rn .
Nehmen wir an, dass die Hessematrix Hf (xk ) symmetrisch und positiv definit ist, so liefert
der Satz 2.1.8 ein eindeutiges (globales) Minimum von qk und es gilt
0 = ∇qk (xk+1 ) = ∇f (xk ) + Hf (xk )(xk+1 − xk ).
Da Hf (xk ) als positiv definite Matrix invertierbar ist, erhalten wir die Iterationsvorschrift
xk+1 = xk − Hf (xk )−1 ∇f (xk ).
Die Richtung dk := −Hf (xk )−1 ∇f (xk ) heißt Newton-Richtung, das Gleichungssystem
Hf (xk )dk = −∇f (xk )
heißt Newton-Gleichung. Für positiv definite Matrizen Hf (xk ) ist auch die Inverse
Hf (xk )−1 positiv definit und daher ist die Newton-Richtung für ∇f (xk ) 6= 0 eine Abstiegsrichtung, denn es gilt
∇f (xk )> dk = −∇f (xk )> Hf (xk )−1 ∇f (xk ) < 0.
21
Verzichtet man auf die restriktive Voraussetzung der positiven Definitheit der Hessematrix
in allen Iterationspunkten und fordert lediglich Invertierbarkeit von Hf (x), so ist das
Verfahren noch immer sinnvoll, allerdings im Allgemeinen kein Abstiegsverfahren mehr.
Einen alternativen Zugang zum Newton-Verfahren erhalten wir über nichtlineare Gleichungssysteme. Da ein lokales Minimum x̂ von f notwendigerweise ∇f (x̂) = 0 erfüllt,
liefert das Newton-Verfahren für nichtlineare Gleichungssysteme ebenfalls die Iterationsvorschrift
Hf (xk )dk = −∇f (xk ),
xk+1 = xk + dk ,
k = 0, 1, 2, . . . .
Algorithmus 2.5.1 (Lokales Newton-Verfahren)
(S.0) Wähle einen Startvektor x0 ∈ Rn , ε ≥ 0 und setze k := 0.
(S.1) Ist k∇f (xk )k < ε: STOP.
(S.2) Berechne (falls möglich) dk als Lösung von
Hf (xk )dk = −∇f (xk ).
(S.3) Setze xk+1 := xk + dk , k := k + 1 und gehe zu (S.1).
Beim lokalen Newton-Verfahren gibt es keine Schrittweitenbestimmung, wir wählen stets
tk = 1. Das Verfahren ist unter bestimmten Voraussetzungen sehr schnell konvergent. Um
dies zu messen definieren wir verschiedene Konvergenzgeschwindigkeiten.
Definition 2.5.2 (lineare, superlineare, quadratische Konvergenz)
Eine Folge {xk } ⊂ Rn konvergiert (mindestens)
(a) linear gegen x̂, wenn es eine Konstante c ∈ (0, 1) gibt, so dass
kxk+1 − x̂k ≤ c · kxk − x̂k
für alle hinreichend großen k ∈ N gilt.
(b) superlinear gegen x̂, wenn
kxk+1 − x̂k
= 0.
k→∞ kxk − x̂k
lim
(c) quadratisch gegen x̂, wenn es ein C ≥ 0 gibt, so dass
kxk+1 − x̂k ≤ C · kxk − x̂k2 .
für alle k ∈ N gilt.
22
Definition 2.5.3 (Landau-Symbole)
Für zwei Folgen positiver reeller Zahlen {αk }, {βk } ⊂ R+ gilt:
αk = O(βk )
⇔
αk = o(βk )
⇔
αk
<∞
k→∞ βk
αk
= 0.
lim
k→∞ βk
lim sup
⇔
∃γ > 0 :
αk
≤ γ ∀ k ∈ N;
βk
Mit Hilfe der Landau-Symbole erhalten wir
(a) xk → x̂ superlinear, wenn kxk+1 − x̂k = o(kxk − x̂k).
(b) xk → x̂ quadratisch, wenn kxk+1 − x̂k = O(kxk − x̂k2 ).
Beispiel 2.5.4
(a) Die Folge {xk } mit xk := q k und q ∈ (0, 1) konvergiert linear gegen x̂ = 0, denn für
alle k ∈ N gilt
kxk+1 − x̂k = q k+1 = q · q k = q · kxk − x̂k.
(b) Die Folge {xk } mit xk :=
für alle k ∈ N gilt
1
k
konvergiert sehr langsam (sublinear) gegen x̂ = 0, denn
kxk+1 − x̂k =
1
k
k
· =
·kxk − x̂k.
k + 1 k |k {z
+ 1}
→1
(c) Die Folge {xk } mit xk :=
1
k!
konvergiert superlinear gegen x̂ = 0, denn es gilt
kxk+1 − x̂k
k!
1
=
lim
=
lim
= 0.
k→∞ kxk − x̂k
k→∞ (k + 1)!
k→∞ k + 1
lim
(d) Die rekursiv definierte Folge {xk } mit x0 > a > 0 und xk+1 = 12 xk + xak ist
√
quadratisch konvergent gegen a, denn dies ist die vom Newton-Verfahren erzeugte
Folge für das Gleichungssystem 0 = f (x) = x2 − a.
Die schnelle Konvergenz des lokalen Newton-Verfahrens sichert der folgende Satz.
Satz 2.5.5 (Konvergenzsatz für das lokale Newton-Verfahren)
Sei f : Rn → R eine zweimal stetig differenzierbare Funktion, x̂ ∈ Rn ein stationärer
Punkt von f und die Hessematrix Hf (x̂) sei invertierbar. Dann gibt es ein r > 0 so dass
23
das lokale Newton-Verfahren für alle Startvektoren x0 mit kx0 − x̂k < r wohldefiniert ist.
Die vom lokalen Newton-Verfahren erzeugte Folge {xk } konvergiert superlinear gegen x̂.
Die Konvergenzrate ist quadratisch, wenn Hf lokal Lipschitz-stetig um x̂ ist, d.h. wenn es
ein L > 0 gibt, so dass
kHf (x) − Hf (y)k ≤ Lkx − yk
für alle x, y in einer Umgebung von x̂ gilt.
Zum Beweis dieses Satzes benötigen wir folgendes Lemma aus der Analysis, auf dessen
Beweis wir hier verzichten wollen.
Lemma 2.5.6
Sei f : Rn → R eine zweimal stetig differenzierbare Funktion, x̂ ∈ Rn ein stationärer
Punkt von f und sei die Hessematrix Hf (x̂) invertierbar. Dann existiert ein r1 > 0 und
ein c > 0 mit kHf (x)−1 k ≤ c für alle x ∈ Ur1 (x̂) := {x ∈ Rn | kx − x̂k ≤ r1 }.
Ferner benötigen wir noch das folgende Approximationslemma:
Lemma 2.5.7
Sei f : Rn → R zweimal stetig differenzierbar und {xk } ⊂ Rn eine gegen x̂ ∈ Rn konvergente Folge. Dann gilt
k∇f (xk ) − ∇f (x̂) − Hf (xk )(xk − x̂)k = o(kxk − x̂k).
Ist ferner Hf lokal Lipschitz-stetig, so gilt
k∇f (xk ) − ∇f (x̂) − Hf (xk )(xk − x̂)k = O(kxk − x̂k2 ).
Beweis: Es gilt
k∇f (xk ) − ∇f (x̂) − Hf (xk )(xk − x̂)k
= k∇f (xk ) − ∇f (x̂) − Hf (x̂)(xk − x̂) + Hf (x̂)(xk − x̂) − Hf (xk )(xk − x̂)k
≤ k∇f (xk ) − ∇f (x̂) − Hf (x̂)(xk − x̂)k + kHf (x̂) − Hf (xk )k kxk − x̂k
|
{z
} |
{z
}
o(kxk −x̂k), da ∇f stetig differenzierbar in x̂
k
= o(kx − x̂k).
24
→0, da Hf stetig
Ferner folgt mit dem Mittelwertsatz und der lokalen Lipschitz-Konstante L > 0
k
k
k
k∇f (x ) − ∇f (x̂) − Hf (x )(x − x̂)k =
=
≤
≤
=
Z 1
k
k
k
k
H
(x̂
+
t(x
−
x̂))(x
−
x̂)dt
−
H
(x
)(x
−
x̂)
f
f
0
Z 1
[Hf (x̂ + t(xk − x̂)) − Hf (xk )]dt (xk − x̂)
0
Z 1
kHf (x̂ + t(xk − x̂)) − Hf (xk )kdt kxk − x̂k
Z0 1
Lk(1 − t)(x̂ − xk )kdt kxk − x̂k
0
Z 1
(1 − t)dt kxk − x̂k2
L
0
L k
=
kx − x̂k2
2
= O(kxk − x̂k2 ).
Damit können wir nun den Konvergenzsatz 2.5.5 für das lokale Newton-Verfahren beweisen:
Beweis: Sei x̂ ein stationärer Punkt der zweimal stetig differenzierbaren Funktion
f : Rn → R, und sei Hf (x̂) regulär. Mit Lemma 2.5.6 folgt die Existenz eines r1 > 0 und
eines c > 0, mit
kHf (x)−1 k ≤ c ∀x ∈ Ur1 (x̂).
Mit Lemma 2.5.7 folgt die Existenz eines r2 > 0 mit
k∇f (x) − ∇f (x̂) − Hf (x)(x − x̂)k ≤
1
kx − x̂k ∀x ∈ Ur2 (x̂).
2c
Sei r := min{r1 , r2 }. Für x0 ∈ Ur (x̂) ist dann Hf (x0 ) invertierbar, also ist d0 in Schritt
(S.2) berechenbar und es folgt:
kx1 − x̂k = kx0 + d0 − x̂k
= kx0 − x̂ − Hf (x0 )−1 ∇f (x0 )k
≤ kHf (x0 )−1 k kHf (x0 )(x0 − x̂) − ∇f (x0 )k
≤ ckHf (x0 )(x0 − x̂) − ∇f (x0 ) + ∇f (x̂) k
| {z }
=0
1
≤ c kx0 − x̂k
2c
1 0
=
kx − x̂k
2
25
Dies zeigt, dass x1 ∈ Ur (x̂) liegt, und somit die Wohldefiniertheit von d1 . Nun folgt
induktiv die Wohldefiniertheit des gesamten Verfahrens und die Konvergenz von {xk }
gegen x̂. Ferner folgt induktiv
kxk+1 − x̂k ≤ ckHf (xk )(xk − x̂) − ∇f (xk ) + ∇f (x̂)k.
Mit Lemma 2.5.7 folgt hieraus die superlineare bzw. quadratische Konvergenz.
Bemerkung 2.5.8
Das lokale Newton-Verfahren konvergiert nur gegen einen stationären Punkt, wenn man
genügend nahe bei diesem startet. Die Konvergenzgeschwindigkeit ist dann superlinear
oder quadratisch, was sehr schnell im Vergleich zum linear konvergierenden Gradientenverfahren ist. Man benötigt daher deutlich weniger Iterationen. Die Invertierbarkeitsvoraussetzung an die Hessematrix sichert die lokale Eindeutigkeit des stationären Punktes.
Startet man daher genügend dicht bei einem lokalen Minimum, so konvergiert das Verfahren auch gegen dieses.
Um ein global konvergentes Verfahren zu erhalten, muss man zunächst einmal eine Richtung festlegen, falls die Newton-Gleichung nicht lösbar ist. Ferner wird man einen Schrittweitenstrategie benötigen um globale Konvergenz zeigen zu können. Um diese durchführen
zu können sollte man die Newton-Richtung nicht verwenden, falls sie nicht zu einem hinreichend guten Abstieg führt. Diese Ideen führen auf das folgende Verfahren:
Algorithmus 2.5.9 (Globalisiertes Newton-Verfahren)
(S.0) Wähle einen Startvektor x0 ∈ Rn , ε ≥ 0, ρ > 0, p > 2, β ∈ (0, 1), σ ∈ 0, 12 und setze
k := 0.
(S.1) Ist k∇f (xk )k < ε: STOP.
(S.2) Berechne (falls möglich) dk als Lösung von
Hf (xk )dk = −∇f (xk ).
Ist dieses System nicht lösbar, oder die Bedingung
∇f (xk )> dk ≤ −ρkdk kp
nicht erfüllt, so setze dk := −∇f (xk ).
(S.3) Bestimme tk := max{β j | j = 0, 1, 2, . . .} mit
f (xk + tk dk ) ≤ f (xk ) + tk σ∇f (xk )> dk .
(S.4) Setze xk+1 := xk + tk dk , k := k + 1 und gehe zu (S.1).
26
Für dieses globalisierte Verfahren gilt der folgende Konvergenzsatz, auf dessen Beweis wir
hier verzichten wollen. Man findet ihn in dem Buch Numerische Verfahren zur Lösung
”
unrestringierter Optimierungsaufgaben“ von Geiger und Kanzow (S.92).
Satz 2.5.10 (Konvergenz des globalisierten Newton-Verfahrens)
Sei f : Rn → R zweimal stetig differenzierbar und {xk } eine durch das globalisierte
Newton-Verfahren erzeugte Folge. Für einen Häufungspunkt x̂ der Folge mit positiv definiter Hessematrix Hf (x̂) gelten:
(a) Die gesamte Folge {xk } konvergiert gegen ein striktes lokales Minimum x̂ und von
f.
(b) Für hinreichend großes k ∈ N ist dk stets die Newton-Richtung und die Schrittweite
ist tk = 1.
(c) {xk } konvergiert superlinear gegen x̂.
(d) Ist Hf sogar (lokal) Lipschitz-stetig, so konvergiert {xk } quadratisch gegen x̂.
Aus Teil (b) dieses Satzes folgt, dass das globalisierte Newton-Verfahren irgendwann in
das lokale Newton-Verfahren übergeht und damit die superlineare oder quadratische Konvergenz des lokalen Verfahrens hat.
Beispiel 2.5.11 (Rosenbrock-Funktion)
Wie beim Gradientenverfahren testen wir das globalisierte Newton-Verfahren mit der
Rosenbrock-Funktion aus Beispiel 2.4.3. Als Parameter haben wir ρ = 10−8 , p = 2.1, β =
0.5, σ = 10−4 und als Abbruchparameter ε = 10−12 gewählt. Mit dem Startvektor x0 =
(−1.2, 1) ergibt sich der in der folgenden Tabelle dargestellte Iterationsverlauf. Im Gegensatz zum Gradientenverfahren benötigen wir hier nur 22 Iterationen und können eine viel größere Genauigkeit erzielen. In den letzten drei Iterationsschritten kann man
die lokal sehr schnelle Konvergenz des Verfahrens beobachten. Allerdings verwendet das
Newton-Verfahren im Gegensatz zum Gradientenverfahren auch Informationen über die
Hessematrix und nicht nur über den Gradienten.
27
k
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
tk
k∇f (xk )k
f (xk ) #f
1.0000 4.639426214067 4.731884325267
2
0.1250 28.550080477673 4.087398662072
6
1.0000 11.571520864202 3.228672588622
7
1.0000 30.325894552346 3.213898091447
8
1.0000 3.604102254225 1.942585420622
9
0.2500 9.248418367866 1.600193693647 12
1.0000 4.919800679767 1.178389561005 13
1.0000 8.664340213955 0.922411581923 14
1.0000 1.778814287247 0.597488616790 15
0.5000 5.877803909942 0.452625097830 17
1.0000 2.176370736680 0.280762438188 18
1.0000 9.401289652782 0.211393396424 19
1.0000 0.492062573533 0.089019501316 20
0.5000 3.924249192225 0.051535404875 22
1.0000 1.242107702378 0.019992777968 23
1.0000 2.533067361610 0.007169243633 24
1.0000 0.237581805138 0.001069613679 25
1.0000 0.348272078357 0.000077768464 26
1.0000 0.003874187292 0.000000282467 27
1.0000 0.000118716759 0.000000000009 28
1.0000 0.000000000447 0.000000000000 29
1.0000 0.000000000000 0.000000000000 30
x1
x2
-1.17528090 1.38067416
-0.93298143 0.81121066
-0.78254008 0.58973638
-0.45999712 0.10756339
-0.39304563 0.15000237
-0.20941191 0.00677013
-0.06571902 -0.01632866
0.14204255 -0.02298878
0.23110720 0.04547802
0.37974282 0.11814580
0.47959489 0.22004082
0.65340583 0.39672894
0.70262363 0.49125758
0.80278553 0.63322101
0.86349081 0.74193125
0.94207869 0.88133620
0.96799182 0.93633667
0.99621031 0.99163870
0.99947938 0.99894834
0.99999889 0.99999751
1.00000000 1.00000000
1.00000000 1.00000000
Bemerkung 2.5.12 (Varianten des Newton-Verfahrens)
Beim lokalen Newton-Verfahren muss in jedem Schritt die Hessematrix Hf (xk ) berechnet
werden, und sodann ein lineares Gleichungssystem mit dieser Matrix gelöst werden. Dies
ist bei einigen Beispielen eventuell sehr aufwendig. Daher gibt es einige Varianten des
Newton-Verfahrens, die darauf beruhen die Hessematrix Hf (xk ) durch eine Matrix A(xk )
zu approximieren, die leichter berechenbar ist, und sich leichter Faktorisieren lässt:
• Vereinfachtes Newton-Verfahren: Hier wählt man die konstante Matrix A(xk ) :=
Hf (x0 ). Diese muss man nur einmal berechnen und Faktorisieren, so dass das Lösen
der linearen Gleichungssysteme billig wird. Man kann zeigen, dass das Vereinfachte
Newton-Verfahren lokal noch linear konvergiert.
• Finite Differenzen Approximation: Hier wählt man A(x) = (Aij (x))i,j=1,...,n so,
28
dass
∇fi (x + hk ei ) − ∇f (x)
, i, j = 1, . . . , n,
hk
wobei ei der i-te Einheitsvektor ist, und hk die Gittergröße für die Approximation
im k-ten Iterationsschritt ist, für die limk→∞ hk = 0 gilt.
Aij (x) =
• Quasi-Newton-Verfahren: Hier startet man mit einer symmetrisch positiv definiten Matrix A0 . Dann berechnet man rekursiv Ak+1 aus Ak mittels gewisser update
Formeln, so dass alle Ak symmetrisch positiv definit bleiben und man eine Faktorisierung von Ak+1 relativ günstig aus einer Faktorisierung von Ak berechnen kann.
• Inverse Quasi-Newton-Verfahren: Diese funktionieren ähnlich wie die QuasiNewton-Verfahren, mit dem Unterschied, dass man nicht die Matrix Hf (xk ) approximiert, sondern die Matrix Hf (xk )−1 .
29
Kapitel 3
Restringierte nichtlineare Optimierung
In diesem Kapitel beschäftigen wir uns mit restringierten Problemen, d.h. wir betrachten
für eine Teilmenge X ⊆ Rn das Problem eine Zielfunktion f : Rn → R zu minimieren,
also
min f (x) u.d.N. x ∈ X.
Um allgemeine Aussagen herleiten zu können, benötigen wir eine explizite Darstellung
der Menge X. Daher definieren wir für Funktionen g : Rn → Rm und h : Rn → Rp die
Menge X mittels
X := {x ∈ Rn | g(x) ≤ 0, h(x) = 0}.
Die zulässige Menge wird also durch m ≥ 0 Ungleichheitsrestriktionen und p ≥ 0 Gleichheitsrestriktionen beschrieben. Das Problem hat die Form
min f (x) u.d.N. g(x) ≤ 0, h(x) = 0.
3.1
(3.1)
Optimalitätsbedingungen
Wie in der unrestringierten Optimierung sind wir an notwendigen und hinreichenden
Optimalitätsbedingungen interessiert. Hierzu benötigen wir folgende Definitionen:
Definition 3.1.1 (Lagrange-Funktion)
Die Funktion
p
m
X
X
L(x, λ, µ) := f (x) +
λi gi (x) +
µj hj (x)
i=1
j=1
heißt Lagrange-Funktion für das Problem (3.1).
Definition 3.1.2 (KKT-Bedingungen)
(a) Die Bedingungen
∇f (x) +
m
X
i=1
λi ∇gi (x) +
p
X
µj ∇hj (x) = ∇x L(x, λ, µ) = 0,
(3.2)
j=1
h(x) = 0,
g(x) ≤ 0,
30
λ ≥ 0,
λ> g(x) = 0,
(3.3)
heißen Karush-Kuhn-Tucker-Bedingungen oder kurz KKT-Bedingungen von
(3.1). Die Bedingung (3.2) ist dabei die Stationarität der Lagrange-Funktion,
die Bedingungen in (3.3) heißen Komplementaritätsbedingungen.
(b) Jeder Vektor (x̂, λ̂, µ̂) der den KKT-Bedingungen genügt, heißt KKT-Punkt von
(3.1), λ̂ = (λ̂1 , . . . , λ̂m )> und µ̂ = (µ̂1 , . . . , µ̂p )> heißen Lagrange-Multiplikatoren.
Bemerkung 3.1.3 (Komplementaritätsbedingungen)
Die Komplementaritätsbedingungen in (3.3) sind äquivalent zu
gi (x) ≤ 0,
λi ≥ 0,
λi gi (x) = 0
für alle i = 1, . . . , m.
Ist (x̂, λ̂, µ̂) ein KKT-Punkt und existiert kein Index i ∈ {1, . . . , m} mit λ̂i = 0 und
gi (x̂) = 0, so genügt (x̂, λ̂, µ̂) der strikten Komplementarität.
In Abbildung 3.1 ist die geometrische Bedeutung der KKT-Bedingungen für zwei Ungleichheitsrestriktionen und keine Gleichheitsrestriktionen skizziert.
−∇f (x̂)
∇g2 (x̂)
∇g1 (x̂)
g2 (x) = 0
X
g1 (x) = 0
Abbildung 3.1: Illustration der KKT-Bedingungen: Falls in einem lokalen Minimum x̂
die Gradienten der aktiven Ungleichungsrestriktionen linear unabhängig sind, kann der
negative Gradient der Zielfunktion als nicht-negative Linearkombination der Gradienten
der aktiven Beschränkungen dargestellt werden.
Die KKT-Bedingungen verallgemeinern den unbeschränkten Fall, wo wir ∇f (x̂) = 0 hatten. Sind nur Gleichheitsrestriktionen vorhanden, dann sind die KKT-Bedingungen unter
31
Multiplikatorregel nach Lagrange bekannt und lauten:
∇f (x̂) +
p
X
µj ∇hj (x̂) = 0,
j=1
h(x̂) = 0.
Die KKT-Bedingungen sind im Allgemeinen keine notwendigen Optimalitätsbedingungen. Man benötigt zusätzlich eine Regularitätsbedingung (constraint qualification
(CQ)) an die zulässige Menge X. Für eine dieser benötigen wir den Begriff der Konvexität
einer Funktion:
Definition 3.1.4 (Konvexe Funktion)
Sei X ⊂ Rn nichtleer und konvex. Eine Funktion f : X → R heißt konvexe Funktion
(auf X), wenn
f (λx + (1 − λ)y) ≤ λf (x) + (1 − λ)f (y)
für alle x, y ∈ X und alle λ ∈ (0, 1) gilt.
Für stetig differenzierbare Funktionen gilt die folgende Charakterisierung konvexer Funktionen:
Satz 3.1.5 (Konvexe Funktionen)
Sei X ⊂ Rn nichtleer und konvex und f : X → R eine stetig differenzierbare Funktion.
Dann gilt: f ist genau dann konvex, wenn
f (x) − f (y) ≥ ∇f (y)> (x − y)
für alle x, y ∈ X.
Ist f sogar zweimal stetig differenzierbar, so gilt: f ist genau dann konvex, wenn ∇2 f
positiv semidefinit ist.
Beispiel 3.1.6 (Konvexe Funktionen)
(a) Affin lineare Funktionen f (x) := a> x + b mit a ∈ Rn , b ∈ R sind konvex.
(b) Quadratische Funktionen f (x) := 12 x> Qx + c> x + γ mit symmetrischer Matrix
Q ∈ Rn×n und c ∈ Rn , γ ∈ R sind genau dann konvex, wenn Q positiv semidefinit
ist.
(c) Die Funktion f (x) := exp(x) ist konvex.
Geeignete Regularitätsbedingungen sind nun etwa die folgenden.
Definition 3.1.7 (Regularitätsbedingungen)
32
(a) Ein zulässiger Punkt x von (3.1) erfüllt die Regularitätsbedingung der linearen Unabhängigkeit (LICQ: linear independence constraint qualification),
wenn die Gradienten
∇hj (x), j = 1, . . . , p,
∇gi (x), i ∈ I(x) := {i ∈ {1, . . . , m} | gi (x) = 0}
linear unabhängig sind.
(b) Ein zulässiger Punkt x von (3.1) erfüllt die Regularitätsbedingung von MangasarianFromovitz (MFCQ), wenn die Gradienten
∇hj (x), j = 1, . . . , p
linear unabhängig sind, und ein d ∈ Rn existiert mit
∇hj (x)> d = 0, j = 1, . . . , p
und ∇gi (x)> d < 0, i ∈ I(x)
(c) Die zulässiger Menge X von (3.1) erfüllt die Slater-Bedingung (Slater CQ),
wenn die Funktionen gi , i = 1, . . . , m konvex sind, die Funktionen hj , j = 1, . . . , p
affin linear sind, und ein Punkt x̂ ∈ X existiert mit
hj (x̂) = 0, j = 1, . . . , p,
und
gi (x̂) < 0, i = 1, . . . , m.
Bemerkung 3.1.8
Die ersten beiden Regularitätsbedingungen beziehen sich jeweils auf einen zulässigen Punkt,
während die letzte von einem konkreten Punkt unabhängig ist. MFCQ ist eine schwächer
Bedingung als LICQ, denn gilt LICQ in einem Punkt x, so hat die Matrix
!
∇hj (x) j = 1, . . . , p
A :=
∈ R(p+|I(x)|)×n
∇gi (x)
i ∈ I(x)
vollen Zeilenrang, und mit b := (0> , −1> )> ∈ Rp+|I(x)| hat das lineare Gleichungssystem
A d = b eine Lösung d ∈ Rn , die dann nach Definition MFCQ genügt. Es gilt also
LICQ
⇒
M F CQ,
Allerdings lässt sich LICQ, insbesondere bei linearen Restriktionen, leichter überprüfen.
Die Regularitätsbedingungen sind aber keinesfalls immer erfüllt, wie das folgende Beispiel
verdeutlicht.
Beispiel 3.1.9 (Verletzte CQ)
Betrachten wir das Optimierungsproblem
min x21 + x22
u.d.N. x ∈ X := {x ∈ R2 | x1 x2 = 0}.
33
Dann ist offensichtlich x = (0, 0)> Lösung des Problems und der Gradient der Gleichheitsrestriktion in der Lösung ist (0, 0)> . Daher ist LICQ und MFCQ verletzt. Da die
Gleichheitsrestriktion nicht affin linear ist, kann man hier die Slater-Bedingung nicht verwenden.
Unter den ersten beiden Regularitätsbedingungen sind die KKT-Bedingungen notwendige
Optimalitätskriterien.
Satz 3.1.10 (KKT-Bedingungen unter MFCQ, LICQ)
Sei x̂ ein lokales Minimum von (3.1) das MFCQ oder LICQ erfüllt. Dann gibt es LagrangeMultiplikatoren λ̂ ∈ Rm und µ̂ ∈ Rp , so dass (x̂, λ̂, µ̂) ein KKT-Punkt ist. Falls LICQ gilt,
so sind die Multiplikatoren sogar eindeutig.
Auf die recht technischen Beweise wollen wir hier verzichten. Es sei jedoch erwähnt, dass
das Farkas Lemma, welches eine zum starken Dualitätssatz (siehe Satz 4.4.4) äquivalente
Aussage liefert, das wesentliche Hilfsmittel ist.
Betrachten wir nun einige Beispiele.
Beispiel 3.1.11
Betrachten wir ein Optimierungsproblem mit einer Gleichheitsrestriktion
min x1 + x2
u.d.N. h(x) := x21 − x2 = 0.
Als Lagrange-Funktion erhalten wir
L(x1 , x2 , µ) = x1 + x2 + µ(x21 − x2 ).
!
!
2x1
0
Wegen ∇h(x) =
6=
für alle x1 ∈ R gilt LICQ in jedem zulässigen Punkt.
−1
0
Sei (x1 , x2 ) ein lokales Minimum. Nach Satz 3.1.10 gibt es dann ein µ ∈ R so dass die
KKT-Bedingungen
!
!
!
1
2x1
0
+µ
=
, x21 − x2 = 0
1
−1
0
gelten. Das Gleichungssystem hat die Lösung µ = 1, x1 = − 21 und damit erhält man aus
der verbleibenden Gleichung x2 = 14 . Also ist (x1 , x2 , µ) = − 12 , 14 , 1 ein KKT-Punkt und
somit ist x = − 12 , 14 ein Kandidat (der einzige) für ein lokales Minimum.
Beispiel 3.1.12
Betrachten wir ein Optimierungsproblem mit einer Ungleichheitsrestriktion
min x1
u.d.N. g(x) := (x1 − 1)2 + x22 − 1 ≤ 0.
Als Lagrange-Funktion erhalten wir
L(x1 , x2 , µ) = x1 + λ((x1 − 1)2 + x22 − 1).
34
!
!
2x1 − 2
0
Es gilt ∇g(x) =
6=
für alle x ∈ R2 \ {(1, 0)}. Da g((1, 0)) = −1 < 0 gilt,
2x2
0
ist die Ungleichung in diesem Punkt nicht aktiv. Daher gilt LICQ in jedem zulässigen
Punkt. Sei (x1 , x2 ) ein lokales Minimum. Nach Satz 3.1.10 gibt es dann ein λ ∈ R so dass
die KKT-Bedingungen
!
!
!
1
2x1 − 2
0
+λ
=
,
0
2x2
0
(x1 − 1)2 + x22 − 1 ≤ 0,
λ ≥ 0,
λ((x1 − 1)2 + x22 − 1) = 0
gelten. Nach der ersten Gleichung der Stationaritätsbedingung kann λ = 0 nicht auftreten.
Wegen λ > 0 folgt aus der zweiten Gleichung der Stationaritätsbedingung x2 = 0, und
ferner muss nach den Komplementaritätsbedingungen die Ungleichung aktiv sein, d.h.
(x1 − 1)2 + x22 − 1 = 0. Somit muss x1 = 0 oder x1 = 2 gelten. Hieraus folgt λ = 12 bzw.
λ = − 21 . Letzteres widerspricht λ ≥ 0 und somit erhalten wir den einzigen KKT-Punkt
(x1 , x2 , λ) = 0, 0, 21 .
Beispiel 3.1.13
Betrachten wir ein Optimierungsproblem mit zwei Ungleichheitsrestriktion
min x2
u.d.N. g1 (x) := x21 − x2 ≤ 0, g2 (x) := x41 − x2 ≤ 0
Als Lagrange-Funktion erhalten wir
L(x1 , x2 , λ1 , λ2 ) = x2 + λ1 (x21 − x2 ) + λ2 (x41 − x2 ).
!
!
2x1
4x31
Es gilt ∇g1 (x) =
, ∇g2 (x) =
für alle x ∈ R2 . Da beide Restriktionen gleich−1
−1
zeitig nur in den Punkten (0, 0), (1, 1), (−1, 1) aktiv sind, gilt in allen übrigen Punkten
LICQ. Ferner gilt LICQ auch in (1, 1) und (−1, 1). In (0, 0) hingegen ist LICQ verletzt,
da dort ∇g1 (x) = ∇g2 (x) gilt. Allerdings gilt z.B. für den Vektor d = (0, 1)> , dass
∇gi (x)> d = −1 < 0 für i = 1, 2, und somit ist MFCQ überall erfüllt. Sei (x1 , x2 ) ein lokales Minimum. Nach Satz 3.1.10 gibt es dann λ1 , λ2 ∈ R so dass die KKT-Bedingungen
!
!
!
!
0
2x1
4x31
0
+ λ1
+ λ2
=
,
1
−1
−1
0
x21 − x2 ≤ 0,
λ1 ≥ 0,
λ1 (x21 − x2 ) = 0
x41 − x2 ≤ 0,
λ2 ≥ 0,
λ2 (x41 − x2 ) = 0
gelten. Ist keine der beiden Ungleichungen aktiv, so ist λ1 = 0 = λ2 und die zweite
Stationaritätsgleichung liefert einen Widerspruch. Ist eine der beiden Ungleichungen inaktiv, so liefert die erste Stationaritätsgleichung x1 = 0 und damit folgt aus der aktiven
35
Ungleichung x2 = 0 und somit sind beide Ungleichungen aktiv. Damit genügt es den
Fall, dass beide Ungleichungen aktiv sind zu diskutieren. Hier gibt es genau drei Punkte
(0, 0), (1, 1), (−1, 1). In den Punkten (1, 1) und (−1, 1) folgt aus der ersten Stationaritätsgleichung zusammen mit λi ≥ 0, i = 1, 2, dass λi = 0 sein muss, was aber der zweiten
Stationaritätsgleichung widerspricht. Damit bleibt nur (0, 0) übrig. Hier erhält man keine
eindeutigen Multiplikatoren, sondern alle Elemente der Menge
{(λ, 1 − λ) | 0 ≤ λ ≤ 1}
sind Lagrange Multiplikatoren. Es gibt also keinen eindeutigen KKT-Punkt, nur die xKomponenten sind eindeutig.
Nun wollen wir noch die Klasse der konvexen Optimierungsprobleme genauer untersuchen.
Das Problem (3.1) heißt konvexes Optimierungsproblem, wenn die Funktionen f und
gi , i = 1, . . . , m konvexe, stetig differenzierbare Funktionen und hj (x) := b>
j x − βj , j =
n
1, . . . , p mit bj ∈ R , βj ∈ R, affin lineare Funktionen sind. Wir erhalten also die Form
min f (x) u.d.N.
gi (x) ≤ 0, i = 1, . . . , m,
b>
j x = βj , j = 1, . . . , p.
(3.4)
Bei konvexen Problemen sind die KKT-Bedingungen nicht nur notwendig sondern auch
hinreichend, wenn die Slater-Bedingung erfüllt ist.
Satz 3.1.14 (KKT-Bedingungen unter Slater CQ)
Sei x̂ ein lokales Minimum des konvexen Optimierungsproblems (3.4) und die Slater CQ
erfüllt. Dann gibt es Lagrange-Multiplikatoren λ̂ ∈ Rm und µ̂ ∈ Rp , so dass (x̂, λ̂, µ̂) ein
KKT-Punkt ist. Ferner ist jeder KKT-Punkt des konvexen Problems (3.4) ein Minimum
von (3.4).
Auf den Beweis wollen wir hier verzichten, es sei jedoch erwähnt, dass die hinreichende
Bedingung aus den KKT-Bedingungen mit der Charakterisierung konvexer Funktionen
aus Satz 3.1.5 folgt.
Bemerkung 3.1.15
Die oben betrachteten Beispiele 3.1.11, 3.1.12, 3.1.13 sind allesamt konvexe Optimierungsprobleme. Daher sind die KKT-Bedingungen auch hinreichende Bedingungen und
die gefundenen KKT-Punkte globale Minima der entsprechenden Probleme. Man hätte
auch relativ einfach die Slater CQ überprüfen können, denn es ist leicht einzusehen, dass
die Punkte (0, 0), (1, 0) bzw. (0, 1) diese erfüllen.
Es gibt auch im restringierten Fall notwendige bzw. hinreichende Kriterien zweiter Ordnung. Hierfür benötigen wir die folgende Notation. Sei (x, λ, µ) ein KKT-Punkt von (3.1)
36
und
I0 (x) := {i ∈ {1, . . . , m} | gi (x) = 0, λi = 0}
(strikte Komplementarität verletzt),
I> (x) := {i ∈ {1, . . . , m} | gi (x) = 0, λi > 0}
(strikte Komplementarität erfüllt).
Dann definieren wir die Menge
T (x) := {d ∈ Rn | ∇hj (x)> d = 0, j = 1, . . . , p,
∇gi (x)> d = 0 (i ∈ I> (x)),
∇gi (x)> d ≤ 0 (i ∈ I0 (x))
}.
Damit können wir die folgenden Resultate angeben.
Satz 3.1.16 (Notwendige Optimalitätsbedingung 2. Ordnung)
Sei x̂ ein lokales Minimum von (3.1), so dass LICQ erfüllt ist. Dann gilt
d> ∇2xx L(x̂, λ̂, µ̂) d ≥ 0 ∀d ∈ T (x̂),
wobei λ̂ ∈ Rm , µ̂ ∈ Rp die nach Satz 3.1.10 eindeutig bestimmten Lagrange-Multiplikatoren
zu x̂ sind.
Satz 3.1.17 (Hinreichende Optimalitätsbedingung 2. Ordnung)
Sei (x̂, λ̂, µ̂) ein KKT-Punkt von (3.1) mit
d> ∇2xx L(x̂, λ̂, µ̂) d > 0
∀d ∈ T (x̂) \ {0}.
Dann ist x̂ ein striktes lokales Minimum von (3.1).
3.2
Lagrange-Newton-Verfahren
In diesem Abschnitt wollen wir ein Verfahren zur Lösung restringierter Optimierungsprobleme kennen lernen. Wie im unrestringierten Fall versuchen wir Punkte zu finden,
die notwendigen Optimalitätskriterien genügen, d.h. wir wollen KKT-Punkte finden. Beschränken wir uns zunächst nur auf Gleichheitsrestriktionen, d.h. betrachten wir Probleme
der Form
min f (x) u.d.N. h(x) = 0
(3.5)
mit zweimal stetig differenzierbaren Funktionen f : Rn → R, h : Rn → Rp . Die zugehörigen KKT-Bedingungen führen mit der Lagrange-Funktion
L(x, µ) := f (x) +
p
X
j=1
37
µj hj (x)
auf ein nichtlineares Gleichungssystem
0 = Φ(x, µ) :=
!
∇x L(x, µ)
.
h(x)
Wendet man das Newton-Verfahren auf dieses Gleichungssystem an, so ergibt sich die
sogenannte Lagrange-Newton-Iteration
!
!
xk+1
xk
=
− Φ0 (xk , µk )−1 Φ(xk , µk ), k = 0, 1, 2, . . . .
µk+1
µk
Algorithmus 3.2.1 (Lagrange-Newton-Verfahren)
(S.0) Wähle (x0 , µ0 ) ∈ Rn × Rp , ε > 0 und setze k := 0.
(S.1) Ist kΦ(xk , µk )k ≤ ε: STOP.
(S.2) Berechne (∆xk , ∆µk ) als Lösung des lineare Gleichungssystem
!
!
!
∆x
∇x L(xk , µk )
∇2xx L(xk , µk ) h0 (xk )>
=−
.
h0 (xk )
0
∆µ
h(xk )
(S.3) Setze
(xk+1 , µk+1 ) := (xk , µk ) + (∆xk , ∆µk ),
k := k + 1 und gehe zu (S.1).
Aus dem allgemeinen Konvergenzsatz 2.5.5 des lokalen Newton-Verfahrens folgt speziell
für das Lagrange-Newton-Verfahren:
Satz 3.2.2 (lokale Konvergenz des Lagrange-Newton-Verfahrens)
Sei (x̂, µ̂) ein KKT-Punkt von (3.5), in dem die Matrix
!
∇2xx L(x̂, µ̂) h0 (x̂)>
h0 (x̂)
0
invertierbar ist. Dann gibt es ein r > 0, so dass das Lagrange-Newton-Verfahren für
alle Startwerte (x0 , µ0 ) ∈ Ur (x̂, µ̂) wohldefiniert ist und die Folge {(xk , µk )} konvergiert
superlinear gegen (x̂, µ̂). Sind die zweiten Ableitungen ∇2 f und ∇2 hj , j = 1, . . . , p sogar
Lipschitz-stetig, so ist die Konvergenz quadratisch.
Bemerkung 3.2.3
Ein hinreichendes Kriterium für die Invertierbarkeit der Matrix beim Lagrange-NewtonVerfahren sind die folgenden Bedingungen:
38
(i) Die Gradienten ∇hj (x̂), j = 1, . . . , p sind linear unabhängig (LICQ).
(ii) Es gilt (die hinreichende Optimalitätsbedingung 2. Ordnung)
d> ∇2xx L(x̂, µ̂)d > 0
∀d 6= 0 mit h0 (x̂)d = 0.
Setzen wir Hk := ∇2xx L(xk , µk ) im Gleichungssystem des Lagrange-Newton-Verfahrens,
so erhalten wir
Hk ∆xk + h0 (xk )> (∆µk + µk ) = −f 0 (xk ),
h0 (xk )∆xk = −h(xk ).
Dies sind genau die KKT-Bedingungen des quadratischen Problems
min
∆x
u.d.N.
1
f 0 (xk )∆x + ∆x> Hk ∆x
2
k
0
hj (x ) + hj (xk )∆x = 0
∀j = 1, . . . , p.
Will man nun auch noch Ungleichungen mit einbeziehen, also restringierte Probleme der
Form
min f (x) u.d.N. g(x) ≤ 0, h(x) = 0,
mit zweimal stetig differenzierbaren Funktionen f : Rn → R, g : Rn → Rm , h : Rn → Rp
betrachten, so motiviert obige Interpretation die Betrachtung von
min
∆x
u.d.N.
1
f 0 (xk )∆x + ∆x> Hk ∆x
2
k
0
hj (x ) + hj (xk )∆x = 0
gi (xk ) + gi0 (xk )∆x ≤ 0
∀j = 1, . . . , p,
∀i = 1, . . . , m,
mit Hk := ∇2xx L(xk , λk , µk ). Man löst also in jedem Iterationsschritt ein quadratisches
Problem. Hierdurch entsteht das SQP-Verfahren (sequential quadratic programming), eines der bekanntesten und besten Verfahren zur Lösung restringierter Optimierungsprobleme.
39
Kapitel 4
Lineare Optimierung
4.1
Lineare Programme
Wir wollen nun lineare Optimierungsprobleme genauer untersuchen. Hierzu wiederholen
wir zunächst den Begriff der Linearität.
Definition 4.1.1 (lineare Funktion)
Eine Funktion f : Rn → Rm heißt linear, wenn sie für alle x, y ∈ Rn und alle λ ∈ R die
folgenden beiden Eigenschaften besitzt:
f (x + y) = f (x) + f (y),
f (cx) = λf (x).
(4.1)
(4.2)
Eine Funktion f : Rn → Rm heißt affin linear, wenn es eine lineare Funktion f˜ : Rn →
Rm und einen konstanten Vektor c ∈ Rm gibt, so dass
f (x) = f˜(x) + c
für alle x ∈ Rn gilt.
Lineare und affin lineare Funktionen stehen in engem Zusammenhang zu Matrizen.
Lemma 4.1.2 (Matrixdarstellung)
Jede lineare Funktion f : Rn → Rm kann in der Form f (x) = Ax mit einer Matrix
A ∈ Rm×n dargestellt werden. Jede affin lineare Funktion f : Rn → Rm kann in der Form
f (x) = Ax + b mit einer Matrix A ∈ Rm×n und einem Vektor b ∈ Rm dargestellt werden.
Nun wollen wir eine Definition eines lineare Optimierungsproblems, auch lineares Programm genannt, geben.
Definition 4.1.3 (Lineares Programm (LP))
Gegeben seien die Matrizen A= ∈ Rm= ×n , A≤ ∈ Rm≤ ×n , A≥ ∈ Rm≥ ×n sowie die Vektoren
b= ∈ Rm= , b≤ ∈ Rm≤ , b≥ ∈ Rm≥ , c ∈ Rn . Dann nennt man das Problem
min c> x u.d.N. A= · x = b= ,
A≤ · x ≤ b≤ ,
A≥ · x ≥ b≥
40
ein lineares Programm. Ebenso kann ein lineares Programm auch in Maximierungsform
vorliegen, also
max c> x u.d.N. A= · x = b= ,
A≤ · x ≤ b≤ ,
A ≥ · x ≥ b≥ .
Als nächstes wollen wir uns überlegen, dass man lineare Programme stets auf eine bestimmt Form transformieren kann. Wir definieren diese zunächst.
Definition 4.1.4 (Primale Normalform)
Gegeben seien die Matrix A ∈ Rm×n , sowie die Vektoren b ∈ Rm und c ∈ Rn . Dann nennt
man das Problem
min c> x u.d.N. Ax = b,
x≥0
ein lineares Programm in primaler Normalform.
Um ein beliebiges lineares Programm auf Normalform zu bringen, gehen wir folgendermaßen vor: Liegt das Problem in Maximierungsform vor, so ändern wir das Vorzeichen der
Zielfunktion um die Minimierungsform zu erhalten. Haben wir eine Ungleichung der Form
n
a>
j x ≤ bj , mit aj ∈ R , bj ∈ R, so können wir eine Schlupfvariable sj ≥ 0 einführen, so
dass
a>
j x + s j = bj
gilt. Für Ungleichungen der Form a>
j x ≥ bj erhält man entsprechend eine Schlupfvariable
sj ≥ 0 mit
a>
j x − s j = bj .
Damit lassen sich alle Ungleichungen zu Gleichungen transformieren. Gilt für eine der
Variablen xj nicht xj ≥ 0, so spricht man von einer freien Variablen. Diese kann man
−
in den positiven Anteil x+
j und den negativen Anteil xj so aufsplitten, dass
−
xj = x+
j − xj ,
x+
j ≥ 0,
x−
j ≥ 0
gilt. Damit haben wir ein lineares Programm, durch Einführen von zusätzlichen Variablen,
auf primale Normalform transformiert. Wir wollen dies nun noch einmal an einem Beispiel
erläutern:
Beispiel 4.1.5 (Transformation auf primale Normalform)
41
Betrachte das folgende lineare Optimierungsproblem:
max −x1 − 2x2 u.d.N. x1 ≥ −1,
x1 − 2x2 ≤ 0,
x1 + x2 = 1,
x2 ≥ 0.
Das Problem liegt nicht in primaler Normalform vor und soll darauf transformiert werden.
Durch Vorzeichenänderung der Zielfunktion erhalten wir ein Minimierungsproblem, und
durch Einführung der Schlupfvariablen s1 , s2 ≥ 0 ergibt sich das äquivalente Problem
min x1 + 2x2 u.d.N. x1 − s1 = −1,
x1 − 2x2 + s2 = 0,
x1 + x2 = 1,
x2 ≥ 0,
s1 ≥ 0, s2 ≥ 0.
Die Ungleichung x2 ≥ 0 wurde nicht in eine Gleichung transformiert, da sie bereits die
gewünschte Form hat. Wir erkennen, dass die Variable x1 eine freie Variable ist. Sie wird
ersetzt durch
−
−
x1 = x+
mit x+
1 − x1
1 ≥ 0, x1 ≥ 0.
Damit erhalten wir das Problem
−
+
−
min x+
1 − x1 + 2x2 u.d.N. x1 − x1 − s1 = −1
−
x+
1 − x1 − 2x2 + s2 = 0
−
x+
1 − x1 + x2 = 1,
−
x2 ≥ 0, s1 ≥ 0, s2 ≥ 0, x+
1 ≥ 0, x1 ≥ 0.
In der Matrixschreibweise erkennen wir nun die primale Normalform


x+
1
 − 
 x1 



Minimiere (1, −1, 2, 0, 0) 
 x2 


 s1 
s2


x+

  1−  

−1
1 −1
0 −1 0  x1 
 


= 0 
u.d.N.
x
0 1 
 1 −1 −2
,
2




1 −1
1
0 0  s1 
1
s2
42








x+
1
x−
1
x2
s1
s2




 ≥ 0.



Bemerkung 4.1.6
Die primale Normalform ist in der Regel nur für den Fall m < n interessant. Ist A
invertierbar, also m = n, so ist x bereits durch x = A−1 b eindeutig festgelegt, so dass nichts
zu optimieren ist. Im Fall m > n und Rang(A) = n besitzt das Gleichungssystem Ax = b
höchstens eine Lösung x, so dass ebenfalls nichts zu optimieren ist. Es gelte im folgenden
also stets m < n. Ferner kann man linear abhängige Nebenbedingungen eliminieren und
daher ohne Beschränkung der Allgemeinheit davon ausgehen, dass Rang(A) = m gilt.
4.2
Ecken, Basisvektoren und der Fundamentalsatz
In diesem Abschnitt untersuchen wir die zulässige Menge eines linearen Programms näher
und stellen einen Zusammenhang zwischen gewissen Punkten der Menge und der Lösung
des linearen Programms her. In Beispiel 1.3.1 haben wir gesehen, dass man im zweidimensionalen Fall die zulässige Menge zeichnen kann. Jede lineare Ungleichung führt
dabei auf einen zulässigen Halbraum und der Durchschnitt dieser ist die zulässige Menge.
Der Durchschnitt von endlich vielen abgeschlossenen Halbräumen ist eine polyedrische
Menge (Polyeder). Ein beschränktes Polyeder ist ein Polytop. Wenn man das Polyeder zeichnen kann, gibt es Punkte die wir anschaulich als Ecken bezeichnen. In höheren
Dimensionen als 3 ist eine Zeichnung der zulässigen Menge und eine graphische Lösung
des linearen Programms aber nicht mehr möglich. Zur systematischen Vorgehensweise
definieren wir auch für Polyeder eine Normalform.
Definition 4.2.1 (Polyeder in Normalform)
Eine Menge der Gestalt
P := {x ∈ Rn | Ax = b, x ≥ 0}
mit A ∈ Rm×n und b ∈ Rm heißt Polyeder in Normalform.
Wie wir im vorherigen Abschnitt gesehen haben, lässt sich jedes Polyeder auf Normalform
transformieren.
Abbildung 4.1 skizziert verschiedene Konstellationen, die in der linearen Optimierung
auftreten können.
43
P
P
−c
−c
P
P
−c
−c
Abbildung 4.1: Verschiedene Geometrien des zulässigen Bereichs P (grau), Höhenlinie der
Zielfunktion (rote Gerade), eindeutige optimale Lösung (links), unendlich viele optimale
Lösungen (rechts oben) und keine optimale Lösung (rechts unten).
Die Beispiele erlauben folgende Beobachtungen:
• Der zulässige Bereich P eines linearen Optimierungsproblems kann beschränkt
oder unbeschränkt sein.
• Falls der zulässige Bereich beschränkt ist, dann existiert stets eine Lösung. Dies
folgt aus dem Satz von Weierstrass, welcher besagt, dass eine stetige Funktion auf
einer kompakten Menge ihr Minimum und Maximum annimmt. Falls der zulässige
Bereich unbeschränkt ist, ist die Existenz einer Lösung nicht gesichert.
• Es kann keine, genau eine oder unendlich viele Lösungen geben.
• Existiert eine Lösung, so befindet sich stets ein Eckpunkt unter den Lösungen.
Die letzte Vermutung (sie gilt allgemein) ist die Hauptmotivation für das berühmte Simplexverfahren. Das Simplexverfahren ist eines der populärsten Verfahren zur Lösung linearer Optimierungsprobleme und wurde 1947 von G.B. Dantzig entwickelt. Die wesentliche
Idee des Simplexverfahrens ist, von einem Eckpunkt des zulässigen Bereichs zu einem
benachbarten Eckpunkt des zulässigen Bereichs zu gelangen, dabei den Zielfunktionswert
(möglichst) zu verringern und dies solange zu wiederholen, bis ein optimaler Eckpunkt
44
erreicht wird. Der Fundamentalsatz der linearen Programmierung wird diese Idee rechtfertigen. Zur Konstruktion des Simplexverfahrens müssen allerdings folgende Fragen beantwortet werden:
• Wie können (zulässige) Eckpunkte berechnet werden?
• Wie kann zu einem gegebenen Eckpunkt ein benachbarter Eckpunkt berechnet werden?
• Wie kann festgestellt werden, ob ein Eckpunkt bereits optimal ist?
Zur Beantwortung dieser Fragen muss zunächst der anschauliche Begriff eines Eckpunkts
mathematisch präzisiert werden.
Definition 4.2.2 (Ecke)
Sei P := {x ∈ Rn | Ax = b, x ≥ 0} ein Polyeder in Normalform mit A ∈ Rm×n und
b ∈ Rm . Ein Vektor x ∈ P heißt Ecke von P, wenn aus
x = λx1 + (1 − λ)x2
für x1 , x2 ∈ P und λ ∈ (0, 1) bereits x1 = x2 folgt.
Für eine anschauliche Interpretation benötigen wir den Begriff der Konvexität einer Menge.
Definition 4.2.3 (Konvexität)
Eine Menge M ⊆ Rn heißt konvex, falls mit zwei Punkten stets auch deren gesamte
Verbindungsstrecke zu M gehört, wenn also
x, y ∈ M
⇒
z = λx + (1 − λ)y ∈ M
für alle λ ∈ [0; 1]
gilt. Jeder solche Punkt z heißt Konvexkombination von x und y.
Konvexe und nicht konvexe Menge
Dass jedes Polyeder eine konvexe Menge ist, ist eine leichte Übungsaufgabe. Anschaulich
sind daher die Vektoren Ecken von P, die sich nicht als echte Konvexkombination zweier
verschiedener Punkte aus P darstellen lassen. In der Definition spielt es keine Rolle, dass
45
wir uns auf Polyeder in Normalform beschränkt haben. Sie gilt genauso für beliebige
abgeschlossene und konvexe Mengen. Beispielsweise hat das Einheitsquadrat {(x1 , x2 ) ∈
R2 | 0 ≤ x1 ≤ 1, 0 ≤ x2 ≤ 1} die bekannten 4 Ecken (0, 0), (0, 1), (1, 0), (1, 1). Nach obiger
Definition hat der Kreis unendlich viele Ecken, nämlich alle Punkte auf dem Kreisrand.
Leider ist die Definition 4.2.2 nicht konstruktiv, so dass man mit ihr die Ecken eines
Polyeders berechnen könnte. Daher benötigen wir einen alternativen Zugang, der sich auf
folgenden Satz stützt und für einen Vektor x ∈ Rn die Indexmenge
B> := {i ∈ {1, . . . , n} | xi > 0}
verwendet.
Satz 4.2.4
Sei P := {x ∈ Rn | Ax = b, x ≥ 0} ein Polyeder in Normalform mit A ∈ Rm×n und b ∈
Rm . Dann ist x ∈ P genau dann eine Ecke von P, wenn die zu den positiven Komponenten
von x zugehörigen Spaltenvektoren von A linear unabhängig sind, wenn also die Spalten
ai von A mit i ∈ B> linear unabhängig sind.
Beweis: Sei x ∈ P zunächst eine Ecke von P. Angenommen die Spalten ai , i ∈ B> sind
linear abhängig. Dann gibt es Koeffizienten βi , i ∈ B> mit
X
βi ai = 0
i∈B>
und βi 6= 0 für mindestens ein i ∈ B> . Wegen xi > 0 für alle i ∈ B> existiert ein δ > 0
mit xi + δβi ≥ 0 und xi − δβi ≥ 0 für alle i ∈ B> . Setzt man nun
(
(
x
+
δβ
,
falls
i
∈
B
,
xi − δβi , falls i ∈ B> ,
i
i
>
x1i :=
, und x2i :=
0,
sonst
0,
sonst
so gilt x1 ≥ 0, x2 ≥ 0 und
1
Ax =
n
X
ai x1i =
i=1
sowie
2
Ax =
n
X
i=1
X
ai (xi + δβi ) = b + δ
i∈B>
ai x2i =
X
X
βi ai = b,
i∈B>
ai (xi − δβi ) = b − δ
i∈B>
X
βi ai = b.
i∈B>
Damit ist x1 , x2 ∈ P und x = 21 x1 + 12 x2 . Wegen x1 6= x2 steht dies im Widerspruch dazu,
dass x eine Ecke von P ist.
Seien nun umgekehrt die Spaltenvektoren ai , i ∈ B> linear unabhängig. Ferner sei
x = λx1 + (1 − λ)x2
46
für x1 , x2 ∈ P und λ ∈ (0, 1). Aus x1 ≥ 0, x2 ≥ 0 und xi = 0, i 6∈ B> folgt mit λ ∈ (0, 1)
sofort
x1i = 0 und x2i = 0 für alle i 6∈ B> .
Daher gilt
0 = b − b = Ax1 − Ax2 = A(x1 − x2 ) =
X
(x1i − x2i )ai .
i∈B>
Die lineare Unabhängigkeit der ai , i ∈ B> impliziert
gilt x1 = x2 , d.h. x ist eine Ecke von P.
x1i − x2i
= 0 für alle i ∈ B> und somit
Nun wollen wir den Begriff eines Basisvektors einführen.
Definition 4.2.5 (Basisvektor)
Sei P := {x ∈ Rn | Ax = b, x ≥ 0} ein Polyeder in Normalform mit A ∈ Rm×n und
b ∈ Rm . Dann heißt x ∈ P Basisvektor von P , falls die Spaltenvektoren ai , i ∈ B> =
{i ∈ {1, . . . , n} | xi > 0} linear unabhängig sind.
Nach dieser Definition ist ein Basisvektor stets in P, so dass man auch von einem zulässigen Basisvektor spricht. Ferner ist, aus historischen Gründen, auch die Bezeichnung
zulässige Basislösung üblich. In der Literatur finden sich auch leicht andere Definitionen, die wegen unserer Voraussetzung Rang(A) = m zu unserer äquivalent sind.
Wegen m < n kann B> für einen Basisvektor höchstens m verschiedene Indizes haben.
Sind es weniger als m, so kann man, wegen Rang(A) = m, die Menge B> zu einer Menge
B ergänzen, so dass B> ⊆ B und die Spalten ai , i ∈ B linear unabhängig sind.
Schließlich bekommen wir mit Satz 4.2.4, dass die Basisvektoren genau die Ecken eines
Polyeders sind.
Satz 4.2.6
Sei P := {x ∈ Rn | Ax = b, x ≥ 0} ein Polyeder in Normalform mit A ∈ Rm×n und
b ∈ Rm . Dann ist x genau dann eine Ecke von P, wenn x ein Basisvektor ist.
Der Grund für unser Interesse an Ecken des Polyeders liegt in dem folgenden zentralen
Resultat.
Satz 4.2.7 (Fundamentalsatz der linearen Programmierung)
Gegeben sei das lineare Programm in primaler Normalform
min c> x u.d.N. Ax = b, x ≥ 0
mit c ∈ Rn , A ∈ Rm×n , b ∈ Rm . Der zulässige Bereich P = {x ∈ Rn : Ax = b, x ≥ 0} sei
nichtleer. Dann gelten folgende Aussagen:
(a) P hat mindestens eine Ecke, und die Gesamtzahl der Ecken ist endlich.
(b) Existiert eine Lösung, so ist auch mindestens ein Eckpunkt von P Lösung.
47
(c) Ist P beschränkt, dann existiert eine Lösung und x ∈ P ist optimal genau dann,
wenn x eine Konvexkombination von optimalen Eckpunkten ist. Hierin heißt x Konvexkombination von Punkten x1 , . . . , xN , falls x sich darstellen lässt als
x=
N
X
λi x
i
N
X
mit
i=1
λi = 1, λi ≥ 0, i = 1, . . . , N.
i=1
Beweis:
(a) Wir zeigen zunächst, dass P eine Ecke hat. Ist der Nullvektor in P, so ist dieser
offensichtlich Basisvektor und somit nach Satz 4.2.6 Ecke. Andernfalls gibt es einen
Vektor x ∈ P mit einer minimalen Anzahl positiver Komponenten. Für diesen
ist B> = {i | xi > 0} 6= ∅ und wir wollen zeigen, dass die Spalten ai , i ∈ B>
linear unabhängig sind. Angenommen, sie sind linear abhängig. Dann existieren
Koeffizienten βi , i ∈ B> mit
X
β i ai = 0
i∈B>
und βi 6= 0 für mindestens ein i ∈ B> . O.B.d.A. ist βi < 0 für eine i ∈ B> (sonst
multipliziere die Gleichung mit -1). Definieren wir
xi δ := min − i ∈ B> , βi < 0 ,
βi
und setzen
(
yi :=
xi + δβi , i ∈ B> ,
0,
sonst,
so gilt y ≥ 0 und yi = 0 für mindestens ein i ∈ B> . Ferner ist
Ay =
n
X
i=1
ai yi =
X
ai (xi + δβi ) = b + δ
i∈B>
X
βi ai = b,
i∈B>
also gilt y ∈ P. Da y aber weniger positive Komponenten als x hat, haben wir einen
Widerspruch. Also sind die Spalten ai , i ∈ B> linear unabhängig. Damit ist x ein
Basisvektor und nach Satz 4.2.6 eine Ecke.
Da es nur endlich viele Möglichkeiten gibt aus den n Spalten von A jeweils m linear
unabhängige auszuwählen, und ein Basisvektor dann durch die lineare Gleichung
Ax = b eindeutig festgelegt ist, kann P nur endlich viele Basisvektoren haben. Also
ist die Gesamtzahl der Ecken endlich.
(b) Nach Voraussetzung existiert eine Lösung. Daher ist die Zielfunktion nach unten
beschränkt. Dann ist
f ∗ := inf{c> x | x ∈ P } > −∞.
48
Wir betrachten nun das modifizierte LP
min c> x u.d.N. x ∈ P̄ := {x ∈ Rn | Ax = b, c> x = f ∗ , x ≥ 0}.
Das LP hat nach Voraussetzung eine Lösung und somit ist P̄ 6= ∅. Ist c linear
abhängig von den Zeilenvektoren von A, so gilt P = P̄ , und der nach Teil (a)
existierende Eckpunkt ist Lösung. Ist c linear unabhängig von den Zeilenvektoren
von A, so hat P̄ als Polyeder in Normalform nach Teil (a) eine Ecke x∗ . Wir zeigen,
dass x∗ auch Ecke von P ist. Dazu wählen wir Punkte x1 , x2 ∈ P und λ ∈ (0, 1), so
dass
x∗ = λx1 + (1 − λ)x2 .
Wegen x1 , x2 ∈ P ist f ∗ ≤ c> x1 und f ∗ ≤ c> x2 . Andererseits ist x∗ ∈ P̄ , also
f ∗ = c> x∗ und somit folgt
f ∗ = c> x 1
und f ∗ = c> x2 ,
also x1 , x2 ∈ P̄ . Da x∗ Ecke von P̄ ist, folgt x1 = x2 , und somit ist x∗ Ecke von P̄ .
Ferner gilt
c> x∗ = f ∗ ≤ c> x für alle x ∈ P,
und somit ist x∗ Lösung.
(c) Ist P beschränkt, so folgt mit der Abgeschlossenheit, dass P kompakt ist. Nach dem
Satz von Weierstraß hat die stetige Zielfunktion c> x somit auf P ein Minimum, d.h.
es existiert eine Lösung.
Da ein beschränktes Polyeder die konvexe Hülle, d.h. die Menge aller Konvexkombinationen, seiner Ecken ist, lässt sich jedes x ∈ P als Konvexkombination von
P
i
Ecken schreiben, d.h. es gibt Ecken x1 , . . . , xN mit x = N
i=1 λi x , λi ≥ 0 für alle
PN
i = 1, . . . , N und i=1 λi = 1. Dabei gilt
>
c x=
N
X
λi c> xi .
(4.3)
i=1
Ist x ∈ P optimal, so gilt c> x ≤ c> xi für alle i = 1, . . . , N. Wegen (4.3) muss dann
aber c> x = c> xi für alle i = 1, . . . , N gelten und somit ist x Konvexkombination
von optimalen Ecken.
Andersherum ist jede Konvexkombination von optimalen Ecken wegen (4.3) ebenfalls optimal.
Der Satz zeigt, dass es ausreicht die Eckpunkte des zulässigen Bereichs zu untersuchen um
eine Lösung des linearen Programms in Normalform zu erhalten. Allerdings erhält man
49
nicht alle Lösungen. Er nutzt wesentlich, dass das Problem in Normalform vorliegt, da ansonsten keine Ecken existieren müssen. Beispielsweise hat das Polyeder {x ∈ R2 | x1 ≥ 0}
keine Ecke. Ein naiver Ansatz zur Lösung eines linearen Programms in Normalform besteht nun darin alle Ecken des zulässigen Bereichs zu berechnen. Da Ecken auch Basisvektoren sind, und wir, wie oben gesehen, wegen Rang(A) = m zu jedem Basisvektor eine
Indexmenge
B mit m linear unabhängigen Spalten von A bekommen, kann es höchstens
!
n
n!
Ecken geben. Da diese Zahl mit größer werdenden n, m zwar endlich
= m!(n−m)!
m
ist, aber sehr schnell (exponentiell) wächst, ist diese Vorgehensweise nicht effizient. Zudem
bekommen wir alleine durch Betrachtung der Ecken keinen Hinweis, ob das Problem überhaupt lösbar ist, oder ob die Zielfunktion im zulässigen Bereich nach unten unbeschränkt
ist.
Wir wollen nun einige wichtige Bezeichnungen einführen.
Definition 4.2.8 (Notationen)
(a) Sei x ein Basisvektor. Jedes System {ai | i ∈ B} von m linear unabhängigen Spalten
von A, welches jene Spalten ai mit xi > 0 enthält, heißt Basis von x.
(b) Sei {ai | i ∈ B} eine Basis von x. Die Indexmenge B heißt Basisindexmenge, die
Indexmenge N := {1, . . . , n} \ B heißt Nichtbasisindexmenge, die Matrix AB :=
(ai )i∈B heißt Basismatrix, die Matrix AN := (ai )i∈N heißt Nichtbasismatrix,
der Vektor xB := (xi )i∈B heißt Basisvariable und der Vektor xN := (xi )i∈N heißt
Nichtbasisvariable.
Wir verdeutlichen die Begriffe im folgenden Beispiel.
Beispiel 4.2.9
Gegeben seien die Nebenbedingungen
x1 + 4x2 ≤ 24,
3x1 + x2 ≤ 21,
x1 + x2 ≤ 9,
x1 ≥ 0, x2 ≥ 0.
Durch Einführung von Schlupfvariable x3 ≥ 0, x4 ≥ 0, x5 ≥ 0 erhalten wir die Nebenbe50
dingungen in primaler Normalform Ax = b, x ≥ 0 mit



1 4 1 0 0


A =  3 1 0 1 0 ,
1 1 0 0 1


24


b =  21  ,
9



x=



x1
x2
x3
x4
x5




.



Die Projektion des zulässigen Bereichs P = {x ∈ R5 | Ax = b, x ≥ 0} auf die (x1 , x2 )Ebene ergibt sich wie folgt:
x2
3x1 + x2 = 21
(0, 6)
(4, 5)
x1 + 4x2 = 24
(6, 3)
(0, 0)
x1 + x2 = 9
x1
(7, 0)
Betrachte x = (6, 3, 6, 0, 0)> ∈ P . 
Die ersten
sind positiv und die
  drei
 Komponenten
 
1
4
1
   
 
zugehörigen Spalten von A, nämlich  3 ,  1  und  0 , sind linear unabhängig.
1
1
0
x ist also Basisvektor und nach Satz 4.2.4 Ecke. Wir erhalten die Basis
      

4
1 
 1

     
 3 , 1 , 0 




1
1
0
mit Basisindexmenge B = {1, 2, 3}, Nichtbasisindexmenge N = {4, 5}, Basisvariable xB =
(6, 3, 6)> , Nichtbasisvariable xN = (0, 0)> und Basismatrix bzw. Nichtbasismatrix




1 4 1
0 0




AB :=  3 1 0  ,
AN :=  1 0  .
1 1 0
0 1
51
Mit den eingeführten Bezeichnungen kann man Ecken eines Polyeders mit Hilfe von linearen Gleichungssystemen bestimmen.
Algorithmus 4.2.10 (Berechnung einer Ecke)
(S.1) Wähle m linear unabhängige Spalten ai , i ∈ B, mit B ⊆ {1, . . . , n} von A und setze
N = {1, . . . , n} \ B.
(S.2) Setze xN = 0 und löse das lineare Gleichungssystem
AB xB = b.
(S.3) Falls xB ≥ 0 gilt, ist x Basisvektor: STOP.
Falls es einen Index j ∈ B mit xj < 0 gibt, dann ist x nicht zulässig. Gehe zu
(S.1) und wiederhole die Prozedur mit einer anderen Wahl von linear unabhängigen
Spalten von A.
Beispiel 4.2.11
Betrachte P = {x ∈ R4 : Ax = b, x ≥ 0} mit
!
!
2 3 1 0
6
A=
,
b=
.
1 0 0 1
2
!
4
Es gibt
= 6 mögliche Kombinationen von zwei Spalten von A:
2
!
>
2
2
−1
1
B1 = {1, 2} ⇒ AB1 b =
,
x = 2, , 0, 0 ,
2
3
3
!
2
B2 = {1, 3} ⇒ A−1
,
x2 = (2, 0, 2, 0)> ,
B2 b =
2
!
3
B3 = {1, 4} ⇒ A−1
,
x3 = (3, 0, 0, −1)> ,
B3 b =
−1
B4 = {2, 3} ⇒
B5 = {2, 4} ⇒
B6 = {3, 4} ⇒
AB4 ist singulär,
!
2
A−1
,
B5 b =
2
!
6
A−1
,
B6 b =
2
x5 = (0, 2, 0, 2)> ,
x6 = (0, 0, 6, 2)> .
x3 ist nicht zulässig. Die Eckpunkte des zulässigen Bereichs sind also gegeben durch
x1 , x2 , x5 , x6 .
52
4.3
Das primale Simplexverfahren
Wir haben in Satz 4.2.4 gesehen, dass ein Eckpunkt algebraisch über linear unabhängige
Spalten von A charakterisiert werden kann. Nach dem Fundamentalsatz genügt es, die
Eckpunkte des zulässigen Bereichs zu berechnen, um mindestens eine optimale Lösung –
falls eine solche existiert – zu erhalten. Darauf basiert das Simplexverfahren.
Der schematische Ablauf des Simplexverfahrens lautet wie folgt:
Phase 0: Transformiere das vorliegende lineare Optimierungsproblem auf primale Normalform.
Phase 1: Bestimme einen Basisvektor (Ecke) x mit Basisindexmenge B, Nichtbasisindexmenge N , Basisvariable xB ≥ 0 und Nichtbasisvariable xN = 0.
Phase 2: Bestimme solange benachbarte Basisvektoren (Ecken) x+ mit Basisindexmenge
B + und Nichtbasisindexmenge N + bis entweder ein Optimum erreicht ist, oder bis
entschieden werden kann, dass das Problem keine Lösung besitzt.
Im Folgenden widmen wir uns
• der effizienten Berechnung von Ecken,
• der Berechnung benachbarter Eckpunkte,
• Optimalitätskriterien für Eckpunkte.
4.3.1 Basiswechsel beim Simplexverfahren
In diesem Abschnitt wird auf die Konstruktion eines benachbarten Basisvektors x+ in Phase 2 des Simplexverfahrens eingegangen, wenn ein Basisvektor x gegeben ist. Ausgehend
von Basis- bzw. Nichtbasisindexmengen B und N in Iteration j des Simplexverfahrens
werden neue Indexmengen
B + = (B \ {p}) ∪ {q},
N + = (N \ {q}) ∪ {p}
konstruiert, wobei p ∈ B und q ∈ N geeignet gewählte Indizes sind. Da hier lediglich der
Basisindex p und der Nichtbasisindex q gegeneinander ausgetauscht werden, spricht man
auch von einem Basiswechsel. Die Indizes p und q werden dabei nicht beliebig gewählt,
sondern so, dass die folgenden beiden Mindestanforderungen erfüllt sind:
(i) Zulässigkeit:
Mit x ∈ P muss auch x+ ∈ P gelten, d.h
Ax = b, x ≥ 0
⇒
53
Ax+ = b, x+ ≥ 0
(ii) Abstiegseigenschaft:
Der Zielfunktionswert fällt monoton, d.h.
c> x+ ≤ c> x.
Es sei eine Basis mit Basisindexmenge B gegeben. Da wir ohne Beschränkung der Allgemeinheit vorausgesetzt haben, dass Rang(A) = m gilt, ist AB invertierbar. Somit kann
die Nebenbedingung Ax = b nach den Basisvariablen xB aufgelöst werden:
Ax = AB xB + AN xN = b
−1
N
xB = A−1
B b − AB AN xN =: βB − ΓB xN .
} | {z }
| {z
⇒
=βB
(4.4)
=ΓN
B
Einsetzen in die Zielfunktion liefert die Darstellung
>
> N
>
>
>
c> x = c>
B xB + cN xN = cB βB − (cB ΓB − cN ) xN =: d − ζN xN .
| {z } |
{z
}
(4.5)
>
=ζN
=d
Hierin besitzt der Vektor βB die Komponenten βi , i ∈ B, der Vektor ζN die Komponenten
ζi , i ∈ N , und die Matrix ΓN
B die Komponenten γij , i ∈ B, j ∈ N .
Für einen Basisvektor x gilt xN = 0 in (4.4) und (4.5) und somit
xB = βB ≥ 0,
xN = 0,
c> x = d.
(4.6)
Wir wollen nun die Basis abändern, um zu einer benachbarten Basis zu gelangen. Dazu
wird die Kantenhalbgerade
z(t) = x + ts,
t ≥ 0,
(4.7)
mit noch zu bestimmender Suchrichtung s und Schrittweite t ≥ 0 betrachtet, vgl. Abbildung 4.2.
x+
s
z(t)
x
Abbildung 4.2: Idee des Basiswechsels beim Simplexverfahren: Es wird entlang einer Kante
gesucht, so dass die Zielfunktion monoton fällt.
54
Die Suchrichtung s wird speziell so gewählt, dass lediglich eine Nichtbasisvariable xq mit
geeignetem q ∈ N abgeändert wird, während die übrigen Nichtbasisvariablen mit xj = 0,
j ∈ N , j 6= q, unverändert bleiben. Für einen geeigneten Index q ∈ N gelte also
sj := 0 für j ∈ N, j 6= q,
sq := 1,
(4.8)
und somit
zq (t) = t ≥ 0,
zj (t) = 0 für j ∈ N, j 6= q.
(4.9)
Natürlich sollen die Werte z(t) darüber hinaus auch zulässig sein, d.h. wie in (4.4) muss
gelten
b = Az(t) = |{z}
Ax +tAs ⇒ 0 = As = AB sB + AN sN .
(4.10)
=b
Auflösen ergibt
sB = −ΓN
B sN ,
(4.11)
so dass die Suchrichtung s mit (4.8) und (4.11) vollständig festgelegt ist. Die Zielfunktionswerte berechnen sich zu
c> z(t)
(4.7)
c> x + tc> s
(4.6)
=
>
d + tc>
B sB + tcN sN
N
>
d − t c>
B ΓB − cN sN
=
d − tζN> sN
=
=
(4.11)
(4.8)
=
d − tζq .
(4.12)
Aus dieser Darstellung lässt sich sofort ablesen, dass der Zielfunktionswert entlang der
Richtung s für t ≥ 0 abnimmt, wenn ζq > 0 gilt.
Gilt andererseits ζj ≤ 0 für alle j ∈ N , so ist entlang s kein Abstieg in der Zielfunktion
möglich. Es bleibt zu klären, ob es dann möglicherweise andere Punkte gibt, die nicht auf
der Kantenhalbgeraden z liegen, und einen kleineren Zielfunktionswert annehmen. Dies
ist nicht der Fall, denn für beliebiges x̂ ∈ P gilt nach (4.5) und mit ζj ≤ 0 und x̂N ≥ 0
>
>
> N
>
>
>
c> x̂ = c>
B x̂B + cN x̂N = cB βB − (cB ΓB − cN )x̂N = d − ζN x̂N ≥ d = c x.
Damit ist der aktuelle Basisvektor x bereits optimal. Wir fassen zusammen:
Zwischenresultat zur Wahl der Pivotspalte q:
Um die Abstiegseigenschaft in (ii) zu erreichen, muss ein Index q ∈ N mit ζq > 0 gewählt
werden. Gilt ζj ≤ 0 für alle j ∈ N , so ist der aktuelle Basisvektor x optimal.
Nun wird versucht, die Zulässigkeitsbedingung in (i) zu erfüllen. Es muss gelten
(4.11)
zB (t) = xB + tsB = βB − tΓN
B sN ≥ 0
55
bzw. in Komponentenschreibweise
zi (t)
=
βi − t
X
γij sj
j∈N
=
βi − tγiq sq −t
|{z}
=1
(4.8)
=
βi − γiq t ≥ 0,
X
j∈N,j6=q
γij sj
|{z}
=0
i ∈ B.
Die Bedingungen βi −γiq t ≥ 0, i ∈ B, schränken die Schrittweite t ≥ 0 ein. Wir diskutieren
zwei Fälle:
(a) Fall 1: Es gelte γiq ≤ 0 für alle i ∈ B.
Wegen βi ≥ 0 gilt dann zi (t) = βi − γiq t ≥ 0 für alle t ≥ 0 und alle i ∈ B. Mit
anderen Worten: z(t) bleibt für alle t ≥ 0 zulässig.
Gilt nun zusätzlich ζq > 0, so ist die Zielfunktion wegen (4.12) für t → ∞ nicht
nach unten beschränkt. Das lineare Optimierungsproblem besitzt keine Lösung!
Zwischenresultat: Unlösbarkeit des linearen Optimierungsproblems
Gilt für ein q ∈ N ζq > 0 und γiq ≤ 0 für alle i ∈ B, so besitzt das lineare
Optimierungsproblem keine Lösung. Die Zielfunktion ist nach unten unbeschränkt.
(b) Fall 2: Es gelte γiq > 0 für mindestens ein i ∈ B.
Aus βi − γiq t ≥ 0 folgt dann t ≤ γβiqi . Diese Forderung schränkt die Schrittweite t
ein.
Zwischenresultat zur Wahl der Pivotzeile p:
Die Zulässigkeit in (i) wird durch die Wahl eines Index p ∈ B mit
βp
βi tmin :=
:= min
γiq > 0, i ∈ B
γpq
γiq
erreicht.
Es bleibt zu zeigen, dass der neue Vektor z(tmin ) tatsächlich ein Basisvektor ist. Nach
Konstruktion und der Definition von tmin gilt zi (tmin ) = 0 für alle i ∈ N + = (N \ {q}) ∪
{p}. Zum Nachweis der linearen Unabhängigkeit der Vektoren ai , i ∈ B + = (B \ {p})∪{q}
sei
X
α i ai + α q aq = 0
i∈B\{p}
56
für Zahlen αi , i ∈ B + . Wir müssen zeigen, dass αi = 0, i ∈ B + gilt. Aus (4.10) ergibt sich
zusammen mit (4.8)
X
aq = −AB sB = −
s i ai
i∈B
Setzen wir dies ein, so erhalten wir
X
X
0 =
αi ai − αq
s i ai
i∈B
i∈B\{p}
=
X
(αi − αq si )ai − αq sp ap .
i∈B\{p}
Da x ein Basisvektor ist, sind die Vektoren ai , i ∈ B linear unabhängig und somit folgt
αi − αq si = 0 für alle i ∈ B \ {p} und αq sp = 0.
Nach Definition gilt tmin = γβpqp = − xspp und somit sp < 0. Wir erhalten also αq = 0 und
damit αi = 0 für alle i ∈ B \ {p}. Dies zeigt die lineare Unabhängigkeit der ai , i ∈ B +
und somit ist z(tmin ) ein Basisvektor.
4.3.2 Der Algorithmus
Mit diesen Betrachtungen ist das Simplexverfahren definiert und wir fassen alles in einem
Algorithmus zusammen:
Algorithmus 4.3.1 (Primales Simplexverfahren)
(S.0) Phase 1:
Bestimme einen Basisvektor (Eckpunkt) x mit Basisindexmenge B, Nichtbasisindexmenge N , Basismatrix AB , Nichtbasismatrix AN , Basisvariable xB ≥ 0 und
Nichtbasisvariable xN = 0 der primalen Normalform mit Rang(A) = m.
−1
−1
>
(S.1) Berechne ΓN
B := AB AN , βB := AB b und die negativen reduzierten Kosten ζN :=
>
N
c>
B ΓB − cN .
(S.2) Überprüfung der Optimalität:
Falls ζj ≤ 0 für alle j ∈ N gilt, so beende das Verfahren. Der aktuelle Basisvektor
xB = βB , xN = 0 ist optimal. Der Zielfunktionswert beträgt d = c>
B βB .
(S.3) Überprüfung der Unbeschränktheit:
Gibt es einen Index q mit ζq > 0 und γiq ≤ 0 für alle i ∈ B, so besitzt das lineare
Optimierungsproblem keine Lösung. STOP.
(S.4) Bestimmung des Pivotelements:
Wähle einen Index q mit ζq > 0. q legt die Pivotspalte fest. Wähle einen Index p mit
βp
βi = min
γiq > 0, i ∈ B .
γpq
γiq p legt die Pivotzeile fest.
57
(S.5) Führe Basiswechsel durch:
Setze B := (B \ {p}) ∪ {q} und N := (N \ {q}) ∪ {p} und gehe zu (S.1).
Das Simplexverfahren wird häufig in einer kompakten Darstellung – dem sogenannten
Simplextableau – durchgeführt. Die Beziehungen
xB = (AB )−1 b − ΓN
B xN ,
> N
>
c> x = c>
B xB − cB ΓB − cN xN ,
werden in einem Tableau nach folgendem Schema zusammengefasst:
xN
xB
−1
ΓN
B = (AB ) AN
βB = (AB )−1 b
N
>
ζN> = c>
B ΓB − cN
c>
B xB
Da xN = 0 in einem Basisvektor gilt, kann der aktuelle Wert der Basisvariablen xB direkt
im Tableau abgelesen werden: xB = (AB )−1 b. Der zugehörige Zielfunktionswert lautet
−1
>
c> x = c>
B xB = cB (AB ) b.
xB
Hat man sich für eine Pivotspalte q ∈ N entschieden, so werden die Quotienten q
ΓB
benötigt, welche aus Gründen der Übersichtlichkeit in einer weiteren Spalte des Tableaus
aufgelistet werden können:
xN
xB
ΓN
B
βB
N
>
c>
B ΓB − cN
c>
B xB
xB
ΓqB
q
Hierin bezeichnet ΓqB die Spalte von ΓN
B , die zum Index q ∈ N gehört. Die Division xB /ΓB
ist komponentenweise zu verstehen.
Beispiel 4.3.2 (vgl. Beispiel 4.2.9)
Betrachte das lineare Optimierungsproblem
min c> x u.d.N. Ax = b, x ≥ 0
mit den Daten (x3 , x4 , x5 sind Schlupfvariable):



−2







 −5 

1 4 1 0 0
24







 , A =  3 1 0 1 0  , b =  21  , x = 
c=
0






0
1
1
0
0
1
9



0
58
x1
x2
x3
x4
x5




.



Ein Basisvektor ist gegeben durch die Basisindexmenge B = {3, 4, 5} mit N = {1, 2} und




 
!
1 0 0
1 4
0
−2




 
AB =  0 1 0  , AN =  3 1  , cB =  0  , cN =
.
−5
0 0 1
1 1
0
und

 
 
 

γ31 γ32
1 4
x3
24

 
 
 

ΓN
B =  γ41 γ42  =  3 1  ,  x4  =  21  ,
γ51 γ52
1 1
x5
9
ζ1
ζ2
!
=
2
5
!
.
Es ergibt sich das folgende Starttableau für das Simplexverfahren:
x3
x4
x5
x1
1
3
1
2
x2
4
1
1
5
24 6
21 21
9 9
0
Das Tableau ist noch nicht optimal und das Problem ist nicht unbeschränkt. Wir wählen
den Nichtbasisindex q = 2 ∈ N mit ζq = 5 > 0 als Pivotspalte (wir hätten auch q = 1
wählen können) und berechnen die Quotienten
x3
24
x4
21
x5
9
=
= 6,
=
= 21,
= = 9.
γ32
4
γ42
1
γ52
1
Diese Werte stehen in der rechten Spalte des Starttableaus. Der erste Quotient ist minimal
und definiert die Pivotzeile p = 3 ∈ B .
Basiswechsel liefert
sowie die Matrizen

4

AB =  1
1
die Basisindexmenge B = {2, 4, 5}, die Nichtbasismenge N = {1, 3},





0 0
1 1
−5





1 0  , AN =  3 0  , cB =  0  , cN =
0 1
1 0
0
−2
0
!
.
und

 
γ21 γ23

 
ΓN
B =  γ41 γ43  = 
γ51 γ53
1
4
11
4
3
4
 

x2
6
 
 

− 14  ,  x4  =  15  ,
− 14
x5
3
1
4
 
Damit folgt Tableau 1:
x2
x4
x5
x1
x3
1
4
11
4
3
4
3
4
1
4
− 14
− 14
− 54
59
6
15
24
3
4
-30
60
11
ζ1
ζ3
!
=
3
4
− 54
!
.
Das Tableau ist nicht optimal und das Problem ist nicht unbeschränkt. Wir wählen die
Pivotspalte q = 1 und die Pivotzeile p = 5 und führen einen Basiswechsel durch. Dieser
liefert die Basisindexmenge B = {2, 4, 1}, die Nichtbasismenge N = {5, 3}, sowie die
Matrizen






!
0 1
−5
4 0 1
0






AB =  1 1 3  , AN =  0 0  , cB =  0  , cN =
.
0
1 0 1
1 0
−2
und
 
  
 
1
x
5
γ25 γ23
− 31
2
3
  11 2  
  

N
,
=  4 ,
ΓB =  γ45 γ43  =  − 3
x
3   4 
4
1
−3
x1
4
γ15 γ13
3

ζ5
ζ3
!
=
−1
−1
!
.
Damit folgt Tableau 2:
x2
x4
x1
x5
− 13
− 11
3
x3
4
3
− 13
-1
-1
1
3
2
3
5
4
4
-33
Tableau 2 ist optimal, da ζN < 0 gilt. Der zugehörige optimale Basisvektor lautet x2 = 5,
x4 = 4, x1 = 4, x5 = 0 und x3 = 0. Der optimale Zielfunktionswert beträgt c> x = −33.
Bemerkung 4.3.3
Anstatt die Inverse A−1
B in Schritt (S.1) von Algorithmus 4.3.1 explizit zu berechnen,
werden die linearen Gleichungssysteme AB ΓN
B = AN und AB xB = b gelöst. Hierbei
können effiziente numerische Verfahren (LR-Zerlegung, QR-Zerlegung, ggf. Verfahren zur
Lösung von dünnbesetzten Gleichungssystemen) eingesetzt werden. Ebenso gibt es UpdateTechniken zur effizienten Berechnung der erwähnten Zerlegungen beim Basiswechsel. Man
spricht hier auch vom revidierten Simplexverfahren.
4.3.3 Updateformeln für das Simplextableau
Die Neuberechnung des Tableauinhalts bei einem Basiswechsel mit Hilfe der Inversen
Basismatrix A−1
B ist für Handrechnungen mitunter sehr mühsam. Es lassen sich jedoch
effiziente Updateformeln herleiten. Hierzu nehmen wir an, dass ein Simplextableau zur
Basisindexmenge B vorliegt. Wir wollen nun die Basis B abändern, um zu einer benachbarten Basis zu gelangen. Dazu seien nach den Regeln des Simplexverfahrens bereits die
Pivotspalte q ∈ N und die Pivotzeile p ∈ B gewählt worden.
60
Betrachte dazu die Zeile p im Tableau:
X
xp = βp −
X
γpj xj = βp − γpq xq −
j∈N
γpj xj .
j∈N,j6=q
Für das Pivotelement gilt γpq 6= 0 und die Gleichung kann nach xq aufgelöst werden:
xq
!
1
=
γpq
X
β p − xp −
γpj xj
=
j∈N,j6=q
X γpj
1
βp
xj
−
xp −
γpq
γpq
γ
pq
j∈N,j6=q |{z}
|{z} |{z}
=βq+
= βq+ −
X
+
=γqp
+
=γqj
+
γqj
xj ,
j∈N +
wobei N + = (N \ {q}) ∪ {p} wieder die neue Nichtbasisindexmenge bezeichnet.
Ersetzen von xq in den restlichen Gleichungen
xi = βi −
X
X
γij xj = βi − γiq xq −
j∈N
γij xj ,
i 6= p, i ∈ B,
j∈N,j6=q
liefert für i 6= p, i ∈ B, die Formeln
X βp
γiq
γpj
xi =
βi − γiq
− −
xp −
γij − γiq
xj
γpq
γpq
γ
pq
j∈N,j6=q |
{z
} | {z }
{z
}
|
=βi+
= βi+ −
+
=γip
X
+
=γij
γij+ xj .
j∈N +
Einsetzen in die Zielfunktion liefert
c> x = d −
X
ζj xj = d − ζq xq −
j∈N
X
ζj xj
j∈N,j6=q
X ζq
γpj
βp
− −
xp −
ζj − ζq
xj
= d − ζq
γpq
γpq
γ
pq
j∈N,j6=q |
| {z } | {z }
{z
}
=d+
= d+ −
=ζp+
X
=ζj+
ζj+ xj
j∈N +
Für die neue Basisindexmenge B + = (B \ {p}) ∪ {q} und die neue Nichtbasisindexmenge
N + = (N \ {q}) ∪ {p} ergibt sich somit die neue Basisdarstellung
xB + = β + − Γ + xN + ,
c> x+ = d+ − (ζ + )> xN + ,
61
welche wieder in einem Tableau dargestellt werden kann. Beachtet man die oben hergeleiteten Updateformeln, so berechnet sich der Tableauinhalt für die neue Basis zur Indexmenge B + aus dem alten Tableauinhalt (mit Daten γij , βi , ζj ) nach folgendem Schema:
xj , j ∈ N \ {q}
xi , i ∈ B \ {p}
γij −
γiq γpj
γpq
γpj
γpq
xq
ζj −
xp
γ
iq
βi −
− γpq
βp
γpq
1
γpq
ζq γpj
γpq
γiq βp
γpq
− γζpqq
d−
ζq βp
γpq
Bemerkung 4.3.4 (Aufwand des Simplexverfahrens)
Zwar ist die Anzahl der Ecken endlich, aber V. Klee und G.J. Minty haben an einem
Beispiel – dem sogenannten Klee und Minty Würfel – gezeigt, dass das Simplexverfahren
im schlimmsten Fall exponentielle Komplexität besitzt und somit nicht zur Klasse der polynomialen Algorithmen gehört. Dies war ein enttäuschendes Resultat aus theoretischer
Sicht und trug maßgeblich zur Entwicklung von Innere Punkte-Verfahren bei, für die unter geeigneten Voraussetzungen eine polynomiale Komplexität gezeigt werden kann. Die
Tatsache, dass das Simplexverfahren – selbst für große Probleme – immer noch eines der
am häufigsten benutzten Verfahren im Operations Research ist, zeigt, dass die worst-case
Komplexität in der Praxis in der Regel nicht beobachtet wird.
4.3.4 Phase 1 des Simplexverfahrens
Zur Durchführung des Simplexverfahrens wird ein Basisvektor benötigt. Für lineare Optimierungsprobleme der Form
min c> x u.d.N. Ax ≤ b, x ≥ 0
mit b ≥ 0 kann eine Basisvektor sofort abgelesen werden. Die Einführung von Schlupfvariablen s ≥ 0 führt auf
Ax + Is = b, x ≥ 0, s ≥ 0
bzw. in Matrixschreibweise
A I
x
s
!
= b,
x
s
!
≥ 0.
Da die Einheitsmatrix I invertierbar ist, ist ein Basisvektor durch
s = b ≥ 0,
62
x=0
gegeben.
Häufig liegt das Problem allerdings nicht in der oben beschriebenen Form vor, so dass ein
Basisvektor im Allgemeinen nicht sofort abgelesen werden kann, und die Bestimmung ist
eine nichttriviale Aufgabe. Es gibt jedoch einen Trick, bei dem ein Basisvektor mit Hilfe
des Simplexverfahrens berechnet werden kann. Hierzu werden künstliche Schlupfvariable
s ≥ 0 eingeführt, die dann zu Null gemacht werden sollen. Diese Idee führt auf das
Hilfsproblem
>
min e s =
m
X
si
u.d.N. Ax + Is = b, x ≥ 0, s ≥ 0,
(4.13)
i=1
wobei e = (1, . . . , 1)> ∈ Rm ist. Hierbei gelte o.B.d.A. b ≥ 0, was durch Multiplikation
entsprechender Gleichungen in Ax = b mit −1 stets erreicht werden kann.
In diesem Fall kann für das Hilfsproblem (4.13) sofort ein Basisvektor angegeben werden, nämlich s = b, x = 0 mit Basismatrix I und es kann dann das Simplexverfahren
angewendet werden.
Für das Hilfsproblem gilt:
(i) Hat das Hilfsproblem die Simplexlösung s = 0, so ist das zugehörige x ein Basisvektor für das Ausgangsproblem.
(ii) Besitzt das Hilfsproblem eine optimale Lösung s 6= 0 und somit e> s > 0, so ist dies
gleichbedeutend damit, dass das Ausgangsproblem keinen zulässigen Punkt besitzt.
Wir illustrieren die Vorgehensweise am folgenden Beispiel.
Beispiel 4.3.5 (Bestimmung eines Basisvektors)
Betrachte Ax = b, x ≥ 0 mit
A=
1 1 2
−1 1 0
!
, b=
1
1
!


x1


, x =  x2  ,
x3
sowie die zu minimierende Zielfunktion
c> x = x1 + 2x2 + x3 .
Zur Bestimmung eines Basisvektors lösen wir zunächst das folgende Hilfsproblem mit dem
Simplexverfahren:
min s1 + s2
u.d.N. Ax + Is = b,
x ≥ 0,
s = (s1 , s2 )> ≥ 0.
Das Simplexverfahren liefert die folgenden Tableaus, wobei x4 := s1 und x5 := s2 gesetzt
wurden.
63
Starttableau:
x4
x5
x1
1
-1
0
x2
1
1
2
x3
2
0
2
1 1
1 1
2
x5
-1
1
-2
x3
2
0
2
0 0
1 –
0
Pivotzeile und -spalte: p = 5 , q = 2 .
Tableau 1:
x4
x2
x1
2
-1
2
Pivotzeile und -spalte: p = 4 , q = 1 .
Tableau 2:
x4
x1
x2
1
2
1
2
-1
x5
− 12
1
2
-1
x3
1
1
0
0
1
0
Dieses Tableau ist optimal mit Zielfunktionswert 0 und s1 = s2 = 0(= x4 = x5 ). Daher
erfüllt x = (x1 , x2 , x3 )> = (0, 1, 0)> die Restriktionen Ax = b und x ≥ 0. x ist Basisvektor
mit Basisindexmenge B = {1, 2} und Nichtbasisindexmenge N = {3}.
Nun soll die Zielfunktion
c> x = x1 + 2x2 + x3
minimiert werden. Das Starttableau ist durch Tableau 2 gegeben, wenn die zu s1 und s2
(bzw. x4 und x5 ) gehörenden Spalten gelöscht werden. Da nun eine andere Zielfunktion
als im Hilfsproblem betrachtet wird, muss die untere Zeile des Tableaus neu berechnet
werden:
!
0
= 2,
c>
B xB = (1, 2)
1
!
1
N
>
c>
− (1) = 2.
B ΓB − cN = (1, 2)
1
Starttableau:
x1
x2
x3
1
1
2
64
0 0
1 1
2
Pivotzeile und -spalte: p = 1 , q = 3 .
Tableau 1:
x3
x2
x1
1
-1
-2
0
1
2
Dieses Tableau ist optimal mit Basisvektor x3 = 0, x2 = 1, x1 = 0 und Zielfunktionswert
2.
Bemerkung 4.3.6
Am Ende der Phase 1 kann es vorkommen, dass einige Komponenten von s Basisvariable
sind. In diesem Fall definieren die Spalten von A, die zu den Basiskomponenten xi (i ∈ B)
gehören, keine volle Basismatrix. Es müssen zusätzliche linear unabhängige Spalten aj von
A mit j ∈ N bestimmt werden, um eine volle Basismatrix zu erhalten. Dies kann durch
zusätzliche Simplexschritte erreicht werden mit dem Ziel, die Komponenten von s in die
Nichtbasis zu schieben.
4.3.5 Endlichkeit des Simplexverfahrens
Es gibt noch zwei Freiheitsgrade in Schritt (S.4) des Algorithmus 4.3.1:
• Welcher Index q mit ζq > 0 soll gewählt werden? Das Beispiel 4.3.2 zeigt, dass es
im allgemeinen mehrere Möglichkeiten gibt.
• Welcher Index p mit
βp
= min
γpq
βi γiq > 0, i ∈ B
γiq soll gewählt werden? Im Allgemeinen gibt es auch hier mehrere Möglichkeiten.
Eine übliche Wahl geht auf Dantzig zurück und wählt diejenige Pivotspalte k mit
ζk = max{ζj | ζj > 0, j ∈ N }.
(4.14)
Von dieser Wahl erhofft man sich eine möglichst starke Abnahme der Zielfunktion beim
Basiswechsel. Diese Wahl kann jedoch zu Zyklen im Simplexverfahren führen, so dass das
Verfahren nicht immer terminiert. Dieser Fall tritt allerdings nur in sogenannten entarteten Basisvektoren auf.
Beispiel 4.3.7 (Schleifenbildung beim Simplexverfahren)
Das folgende Beispiel wurde von K.T. Marshall und J.W. Suurballe zum Nachweis von
65
Zyklen im Simplexverfahren verwendet:
min c> x u.d.N. Ax = b, x ≥ 0
mit den Daten (x5 , x6 , x7 sind Schlupfvariable)



−10



 57 






 



9
0.5
−5.5
−2.5
9
1
0
0
0







 

c =  24  , A =  0.5 −1.5 −0.5 1 0 1 0  , b =  0  , x = 






0
1
0
0
0
0
0
1
1






0 


0
x1
x2
x3
x4
x5
x6
x7







.





Anwendung des Simplexverfahrens, wobei die Pivotspalte stets gemäß (4.14) gewählt wird,
führt auf die folgenden Tableaus. Dabei entspricht das letzte Tableau bis auf Permutationen der Spalten dem ersten Tableau. Folglich führt eine wiederholte Wahl derselben
Pivotelemente zu einem Zyklus.
Starttableau:
x5
x6
x7
x1
0.5
0.5
1
10
x2
x3
-5.5 -2.5
-1.5 -0.5
0
0
-57
-9
Pivotzeile und -spalte: p = 5 , q = 1 .
Tableau 1:
x5 x2 x3
x1 2 -11 -5
x6 -1
4
2
x7 -2 11 5
-20 53 41
Pivotzeile und -spalte: p = 6 , q = 2 .
Tableau 2:
x5
x6
x1 -0.75 2.75
x2 -0.25 0.25
x7 0.75 -2.75
-6.75 -13.25
66
x4
9
1
0
-24
x4
18
-8
-18
-204
x3
0.5
0.5
-0.5
14.5
0 0
0 0
1 1
0
0
0
1
0
x4
-4
-2
4
-98
–
0
1
11
0 0
0 0
1 –
0
Pivotzeile und -spalte: p = 1 , q = 3 .
Tableau 3:
x5
x3 -1.5
x2 0.5
x7
0
15
x6
5.5
-2.5
0
-93
x1
2
-1
1
-29
Pivotzeile und -spalte: p = 2 , q = 4 .
Tableau 4:
x5
x6
x3 0.5
-4.5
x4 0.25 -1.25
x7
0
0
10.5 -70.5
x1
-2
-0.5
1
-20
Pivotzeile und -spalte: p = 3 , q = 5 .
Tableau 5:
x3
x5
2
x4 -0.5
x7
0
-21
x6
-9
1
0
24
Pivotzeile und -spalte: p = 4 , q = 6 .
Tableau 6:
x3
x5 -2.5
x6 -0.5
x7
0
-9
x4 x1
9 0.5
1 0.5
0
1
-24 10
x1
-4
0.5
1
22
x4
-8
2
0
18
0 –
0 0
1 –
0
x2
4
0.5
0
-9
x2
8
-1.5
0
-93
x2
-5.5
-1.5
0
-57
0 0
0 0
1 –
0
0
0
1
0
–
0
–
0
0
1
0
Es ist interessant zu beobachten, dass alle Tableaus denselben (entartete) Basisvektor x =
(0, 0, 0, 0, 0, 0, 1)> beschreiben, obwohl sich die Tableaus selbst unterscheiden. Geometrisch
bleibt das Verfahren also in derselben Ecke des zulässigen Bereichs hängen.
Wir wollen einige Überlegungen anstellen, wann die Schleifenbildung eintreten kann.
Zunächst nehmen wir an, dass bei der Durchführung des Simplexverfahrens nur Pivotzeilen p mit βp > 0 gewählt werden (dies war beim vorangehenden Beispiel nicht der
Fall). Die Updateformeln aus Abschnitt 4.3.3 zeigen, dass die Zielfunktion in diesem Fall
strikt fällt, d.h. es gilt d+ < d. Da das Simplexverfahren per Konstruktion ausschließlich
Basisvektoren berechnet, von denen es nur endlich viele gibt und von denen mindestens
67
eine optimal ist, endet das Verfahren in diesem Fall nach endlich vielen Schritten (falls
überhaupt eine Lösung existiert).
Eine Schleife kann also nur dann auftreten, wenn eine Pivotzeile p mit βp = 0 gewählt
wird, da dann d+ = d aus den Updateformeln aus Abschnitt 4.3.3 folgt. Der Zielfunktionswert ändert sich somit nicht beim Übergang von dem Basisvektor x zum benachbarten
Basisvektor x+ . Des weiteren ändert sich gemäß Abschnitt 4.3.3 auch β + nicht, d.h. es
gilt xB + = βB++ = βB = xB (vgl. Beispiel 4.3.7). Das Verfahren stagniert also beim Basisvektor x = x+ . Beachte, dass sich die Darstellung des Basisvektors, also das Tableau,
sehr wohl ändert, da ein Basiswechsel stattfindet. In diesem Fall spricht man von einem
entarteten Basisvektor.
Definition 4.3.8 (Entartung)
Ein Basisvektor x mit xi = 0 für mindestens ein i ∈ B heißt entartet. Andernfalls heißt
er nicht entartet.
Als Fazit halten wir fest, dass eine Schleifenbildung nur auftreten kann (aber nicht muss),
wenn entartete Basisvektoren auftreten.
Es gibt eine einfache Regel, die die Schleifenbildung verhindert (der Beweis ist allerdings
sehr kompliziert, so dass wir ihn hier nicht führen wollen).
Blandsche Regel
Wähle für die Pivotspalte q ∈ N und die Pivotzeile p ∈ B in Schritt (S.4) des Algorithmus
4.3.1 stets den Kandidaten mit dem kleinsten Index. Dann treten keine Schleifen auf und
das Simplexverfahren endet nach endlich vielen Schritten.
Bemerkung 4.3.9
In der Praxis tritt die Schleifenbildung in der Regel nicht auf, so dass dort häufig auf die
Blandsche Regel verzichtet wird.
4.4
Dualität
In diesem Abschnitt kümmern wir uns um duale lineare Programme. Wir beginnen mit
der Definition.
Definition 4.4.1 (Duales lineares Programm)
Zu einem gegebenem linearen Programm in primaler Normalform
min c> x u.d.N. Ax = b, x ≥ 0
(4.15)
mit A ∈ Rm×n , b ∈ Rm , c ∈ Rn bezeichnen wir mit
max b> y
u.d.N. A> y ≤ c.
68
(4.16)
das zugehörige duale lineare Programm.
Der Nutzen dualer Programme ist vielfältig. Die dualen Programme liefern beispielsweise
untere Schranken für den Zielfunktionswert des primalen Problems. Ferner lassen sich mit
Hilfe dualer Probleme durch Störung der Problemdaten hervorgerufene Änderungen des
Minimalwertes des primalen Problems abschätzen. Manchmal sind die dualen Probleme
auch einfacher zu lösen als die primalen.
Man kann nachrechnen (Übungsaufgabe), dass das duale Problem des dualen Problems
wieder das primale Problem ist. Einen weiteren Zusammenhang zwischen den Problemen
liefert der folgende Satz
Satz 4.4.2 (Schwacher Dualitätssatz)
Für alle primal zulässigen Punkte x (mit Ax = b, x ≥ 0) und alle dual zulässigen Punkte
y (mit A> y ≤ c) gilt
b> y ≤ c> x.
Beweis: Die Behauptung folgt wegen b = Ax, x ≥ 0 und A> y ≤ c aus
b> y = (Ax)> y = x> A> y ≤ x> c = c> x.
Der schwache Dualitätssatz liefert eine Motivation für das duale Problem. Dual zulässige
Werte liefern untere Schranken für den optimalen Zielfunktionswert des primalen Problems. Umgekehrt liefern primal zulässige Werte obere Schranken für den optimalen Zielfunktionswert des dualen Problems. Diese Eigenschaft kann z.B. beim Branch & BoundVerfahren für ganzzahlige Probleme zur Bestimmung unterer Schranken ausgenutzt werden.
Wenn b> y < c> x gilt, so spricht man von einer Dualitätslücke. Hat man diese nicht, so
gilt:
Korollar 4.4.3
Sei x primal und y dual zulässig und gelte c> x = b> y. Dann löst x das primale und y das
duale lineare Programm.
Beweis: Sei x̃ ∈ Rn ein beliebiger für das primale LP zulässiger Punkt. Dann gilt mit
dem schwachen Dualitätssatz
c> x = b> y ≤ c> x̃,
also löst x das primale LP.
69
Sei ỹ ∈ Rn ein beliebiger für das duale LP zulässiger Punkt. Dann gilt mit dem schwachen
Dualitätssatz
b> ỹ ≤ c> x = b> y,
also löst y das duale LP.
Satz 4.4.4 (Starker Dualitätssatz)
Das primale lineare Programm (4.15) besitzt genau dann eine Lösung x, wenn das duale
lineare Programm (4.16) eine Lösung y besitzt. Ferner gilt dann c> x = b> y.
Beweis: Das primalen LP habe eine Lösung. Dann lässt sich mit dem Simplexverfahren
auch ein Basisvektor x mit Basisindexmenge B und Nichtbasisindexmenge N finden, der
Lösung ist. Es gilt also
xB ≥ 0,
xN = 0,
−1
>
c>
B AB AN − cN =: ζN ≤ 0.
Wir definieren nun
−1
y := (A>
B ) cB
und zeigen, dass y das duale LP löst. Mit dieser Definition und ζN ≤ 0 gilt
A>
B y = cB ,
>
> −1
A>
N y = AN (AB ) cB = ζN + cN ≤ cN ,
also ist A> y ≤ c und y ist zulässig für das duale Problem. Ferner folgt aus der primalen
Zulässigkeit AB xB = b
> −1
>
>
c> x = c>
B xB = cB AB b = y b = b y,
und somit liefert Korollar 4.4.3, dass y Lösung des dualen LP ist.
Umgekehrt besitze das duale LP eine Lösung. Man kann das duale LP auf primale Normalform bringen. Die gleiche Argumentation wie oben liefert dann eine Lösung des dazu
dualen Problems, was wiederum das primale LP ist (Übung) und vervollständigt den
Beweis.
Nun wollen wir noch zeigen, dass für Lösungen der linearen Programme Komplementaritätsbedingungen gelten.
Korollar 4.4.5 (Komplementaritätsbedingungen)
Seien x eine Lösung des primalen und y eine Lösung des dualen linearen Programms.
Dann gelten die Komplementaritätsbedingungen:
xi (c − A> y)i = 0 für alle i = 1, . . . , N,
70
d.h. es gilt für alle i = 1, . . . , m
xi > 0 ⇒ (c − A> y)i = 0,
(c − A> y)i > 0 ⇒ xi = 0.
Beweis: Nach dem starken Dualitätssatz gilt
>
>
>
>
>
>
>
0 = c x − b y = x c − x A y = x (c − A y) =
m
X
xi (c − A> y)i .
i=1
Wegen der Zulässigkeit von x und y sind alle Summanden nicht negativ, und daher muss
xi (c − A> y)i = 0
für alle i = 1, . . . , m gelten.
Schließlich bleibt der folgende Satz zu erwähnen, auf dessen Beweis wir hier verzichten
wollen. Er sichert die Existenz von Lösungen für lineare Programme und unterstreicht die
Bedeutung der dualen Programme.
Satz 4.4.6 (Existenzsatz)
Haben sowohl das primale als auch das duale lineare Programm einen zulässigen Punkt,
so haben beide Probleme eine Lösung.
Betrachten wir nun eine Anwendung, nämlich Matrixspiele in gemischten Strategien. Gegeben seien 2 Spieler und eine Auszahlungsmatrix A ∈ Rm×n . Spieler 1 wählt eine Strategie
in
(
)
m
X
X = x ∈ Rm xi ≥ 0, i = 1, . . . , m,
xi = 1 ,
i=1
wobei xi die Wahrscheinlichkeit für die Wahl der reinen Strategie i darstellt. Analog wählt
Spieler 2 eine Strategie in
(
)
n
X
Y = y ∈ Rn yi ≥ 0, i = 1, . . . , n,
yi = 1 .
i=1
Die Auszahlung an Spieler 1 ist dann durch x> Ay gegeben. Da wir Nullsummenspiele
betrachten wollen ist die Auszahlung an Spieler 2 entsprechend −x> Ay. Ziel beider Spieler
ist es ihren Gewinn zu maximieren. Ein zulässiges Strategiepaar (x̄, ȳ) ∈ X × Y heißt
Nash-Gleichgewicht in gemischten Strategien, wenn
x̄> Aȳ ≥ x> Aȳ
für alle x ∈ X
und
71
− x̄> Aȳ ≥ −x̄> Ay
für alle y ∈ Y
gilt. J. Nash hat gezeigt, dass es immer ein Nash-Gleichgewicht in gemischten Strategien gibt. Weiter kann man folgenden Zusammenhang zu primalen und dualen linearen
Programmen herstellen, den wir hier nicht beweisen wollen.
Satz 4.4.7 (Nash-Gleichgewichte in gemischten Strategien)
Gegeben sei ein Matrixspiel mit Auszahlungsmatrix A ∈ Rm×n und X, Y wie oben. Dann
gilt: (x̄, ȳ) ist genau dann ein Nash-Gleichgewicht in gemischten Strategien, wenn folgende
Bedingungen erfüllt sind:
(a) x̄ (zusammen mit einem v̄ ∈ R) löst das primale lineare Programm
max v
v,x
u.d.N. A> x ≥ v1n ,
m
X
xi = 1,
x ≥ 0.
i=1
(b) ȳ (zusammen mit einem w̄ ∈ R) löst das duale lineare Programm
min w
w,y
u.d.N. Ay ≤ w1m ,
n
X
yi = 1,
y ≥ 0.
i=1
Beispiel 4.4.8 (Stein-Schere-Papier Spiel)

0
1 −1


Das Spiel Stein-Schere-Papier lässt sich durch die Auszahlungsmatrix A = −1 0
1
1 −1 0
beschreiben. Löst man die linearen Programme aus Satz 4.4.7 so erhält man für beide Spie
ler jeweils die gemischte Strategie 13 , 31 , 13 als Nash-Gleichgewicht. Es lohnt sich also für
keinen Spieler eine der drei Strategien Stein, Schere oder Papier zu bevorzugen.

72