Seminar zur Approximationstheorie im Wintersemester 2009 ... · uber diese Methode, wobei dort mehr...

32
Seminar zur Approximationstheorie im Wintersemester 2009/2010 Monge-Amp` ere-Gleichung Numerische Verfahren zur L¨ osung der Monge-Amp` ere-Gleichung, Teil II Andreas Platen 29.01.2010 1

Transcript of Seminar zur Approximationstheorie im Wintersemester 2009 ... · uber diese Methode, wobei dort mehr...

Page 1: Seminar zur Approximationstheorie im Wintersemester 2009 ... · uber diese Methode, wobei dort mehr auf eine Methode zur L osung eines Minimierungs- problems mit Nebenbedingung mit

Seminar zur Approximationstheorieim Wintersemester 2009/2010

Monge-Ampere-Gleichung

Numerische Verfahren zur Losung derMonge-Ampere-Gleichung, Teil II

Andreas Platen

29.01.2010

1

Page 2: Seminar zur Approximationstheorie im Wintersemester 2009 ... · uber diese Methode, wobei dort mehr auf eine Methode zur L osung eines Minimierungs- problems mit Nebenbedingung mit

Inhaltsverzeichnis

Inhaltsverzeichnis

1 Einleitung 3

2 Monge-Ampere-Gleichung 3

3 Dritte Methode: Erweiterte Lagrange-Methode 43.1 Einfuhrung der erweiterten Lagrange-Methode . . . . . . . . . . . . . . . 4

3.1.1 Lagrange-Funktion . . . . . . . . . . . . . . . . . . . . . . . . . . . 53.1.2 Strafterm-Methode . . . . . . . . . . . . . . . . . . . . . . . . . . . 83.1.3 Erweiterte Lagrange-Funktion . . . . . . . . . . . . . . . . . . . . . 10

3.2 Sattelpunktproblem . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 103.2.1 Sattelpunkt-Formulierung der Monge-Ampere-Gleichung . . . . . . 113.2.2 Losungsverfahren fur das Sattelpunktproblem . . . . . . . . . . . . 143.2.3 Konvergenz des Algorithmus . . . . . . . . . . . . . . . . . . . . . 153.2.4 Losen von Teilproblem (21) . . . . . . . . . . . . . . . . . . . . . . 153.2.5 Losen von Teilproblem (22) . . . . . . . . . . . . . . . . . . . . . . 15

3.3 Finite-Elemente-Approximation . . . . . . . . . . . . . . . . . . . . . . . . 183.3.1 Finite-Elemente-Methode . . . . . . . . . . . . . . . . . . . . . . . 183.3.2 Gemischte Finite-Elemente-Approximation von (1) . . . . . . . . . 203.3.3 Iteratives Losen von Problem (36) . . . . . . . . . . . . . . . . . . 223.3.4 Losen von Teilproblem (39) . . . . . . . . . . . . . . . . . . . . . . 233.3.5 Losen von Teilproblem (40) . . . . . . . . . . . . . . . . . . . . . . 26

4 Numerische Experimente 264.1 Beispiele . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 264.2 Vergleich . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 30

5 Ausblick 30

2

Page 3: Seminar zur Approximationstheorie im Wintersemester 2009 ... · uber diese Methode, wobei dort mehr auf eine Methode zur L osung eines Minimierungs- problems mit Nebenbedingung mit

1 Einleitung

In dieser Ausarbeitung geht es um numerische Verfahren zur Losung der Monge-Ampere-Gleichung mit Dirichlet-Randbedingung im Zweidimensionalen. Es wird dabei direkt anYasemin Hafizogullaris Ausarbeitung [16] angeknupft. Dort wurden bereits in Kapitel 1und 2 die Grundlagen geklart und in Kapitel 3 und 4 zwei verschiedene numerische Ver-fahren vorgestellt. Das erste dieser beiden Verfahren ist ein Finite-Differenzen-Verfahrenund beim zweiten Verfahren wurde das Losen der nichtlinearen Monge-Ampere-Glei-chung auf ein iteratives Losen von Poisson-Gleichungen, welche linear sind, reduziert.Beide Verfahren stammen von Benamou, Froese und Oberman [1].

In dieser Ausarbeitung wird das von Dean und Glowinski [5, 6, 7] vorgeschlagene Ver-fahren vorgestellt. Sie formulierten das Dirichlet-Problem der Monge-Ampere-Gleichungmit Hilfe eines erweiterten Lagrange-Funktionals in ein Sattelpunktproblem um, welchesanschließend mit einem Uzawa-Douglas-Rachford-Algorithmus gelost werden kann. Fureine Implementierung wahlen sie zur Diskretisierung der Problemstellung eine gemischteFinite-Elemente-Methode.

In Kapitel 2 wird zunachst noch einmal die Problemstellung wiederholt. Anschließendfolgt in Kapitel 3 die Vorstellung des Verfahrens von Dean und Glowinski [5, 6, 7].Im letzten Kapitel werden dann alle drei numerischen Verfahren, das heißt die beidenVerfahren von Benamou, Froese und Oberman [1] und das hier vorgestellte, miteinanderverglichen.

2 Monge-Ampere-Gleichung

Sei Ω ⊂ R2 ein Gebiet mit Rand ∂Ω. Die Monge-Ampere-Gleichung im Zweidimensio-nalen mit Dirichlet-Randbedingung ist durch

detD2u = uxxuyy − u2xy = f in Ω,

u = v auf ∂Ω

(1)

definiert, wobei D2u die Hessematrix von u : Ω → R bezeichne und v ∈ C0(∂Ω) undf > 0 gegebene Funktionen seien.

Bemerkung 1. (Zusammenhang zwischen Monge-Ampere- und Poisson-Gleichung)Das Poisson-Dirichlet-Problem auf Ω lautet

−∆u = SpurD2u = f in Ω,u = v auf ∂Ω.

(2)

Seien λ+, λ− : Ω → R mit λ+ ≥ λ− die beiden Eigenwerte der Hessematrix D2u : Ω →R2×2 in den entsprechenden Punkten aus Ω. Fur das Dirichlet-Problem der Monge-

3

Page 4: Seminar zur Approximationstheorie im Wintersemester 2009 ... · uber diese Methode, wobei dort mehr auf eine Methode zur L osung eines Minimierungs- problems mit Nebenbedingung mit

3 Dritte Methode: Erweiterte Lagrange-Methode

Ampere-Gleichung (1) kann man nun auch

λ+ · λ− = f in Ω,u = v auf ∂Ω

(3)

und fur die Poisson-Gleichung

−(λ+ + λ−) = f in Ω,u = v auf ∂Ω

(4)

schreiben. Da die Monge-Ampere-Gleichung mit f > 0 betrachtet wird, haben die Eigen-werte das gleiche Vorzeichen. Damit zeigt die Monge-Ampere-Gleichung (3) eine Ver-bindung zum geometrischen und die Poisson-Gleichung (4) zum arithmetischen Mittelder Eigenwerte der Hessematrix D2u.

Im folgenden Abschnitt wird sich zeigen, dass man mit geeigneten iterativen Methodenund gemischten Finite-Elemente-Verfahren das Dirichlet-Problem der Monge-Ampere-Gleichung (1) auf das Losen einer Folge von diskreten Varianten des Poisson-Dirichlet-Problems (2) und (4) reduzieren kann.

3 Dritte Methode: Erweiterte Lagrange-Methode

Die in diesem Abschnitt vorgestellte Methode zur Losung der Monge-Ampere-Gleichungmit Dirichlet-Randbedingung (1) stammt von Dean und Glowinski [5, 6, 7].

Zunachst wird zu dem Problem (1) ein Sattelpunktproblem formuliert, wobei dies mitHilfe eines erweiterten Lagrange-Funktionals geschieht. Aus diesem Grund wird zunachstdie erweiterte Lagrange-Funktion in Abschnitt 3.1 in seiner ursprunglichen Form ein-gefuhrt. Danach folgt in Abschnitt 3.2 die Formulierung des Sattelpunktproblems fur(1) und ein Algorithmus, mit dem man dieses Losen kann. Der Abschnitt 3.3 befasstsich schließlich mit einer Diskretisierung der Problemstellung und der Losungsmethodeaus Abschnitt 3.2, um dieses Verfahren auf einem Rechner umsetzen zu konnen.

3.1 Einfuhrung der erweiterten Lagrange-Methode

Die folgende Motivation der erweiterten Lagrange-Methode stammt von Hestenes [17],Abschnitt 2 und 3. Powell [20] veroffentlichte im gleichen Jahr ebenfalls einen Artikeluber diese Methode, wobei dort mehr auf eine Methode zur Losung eines Minimierungs-problems mit Nebenbedingung mit Hilfe der erweiterten Lagrange-Methode und umKonvergenzaussagen als auf eine Herleitung eingegangen wird.

4

Page 5: Seminar zur Approximationstheorie im Wintersemester 2009 ... · uber diese Methode, wobei dort mehr auf eine Methode zur L osung eines Minimierungs- problems mit Nebenbedingung mit

3.1 Einfuhrung der erweiterten Lagrange-Methode

Sei zu gegebenen n ∈ N und Funktionen J, g ∈ C2(Rn) das folgende Minimierungs-problem gegeben:

Finde x ∈ Rn,mit J(x) = minx∈Rn : g(x)=0 J(x)

unter der Nebenbedingung g(x) = 0.

(5)

Die Funktion g kann auch vektorwertig sein, wobei man dies dann als mehrere Neben-bedingungen auffassen kann. Im folgenden wird jedoch der Einfachheit halber lediglicheine einzelne Nebenbedingung verwendet.

Mit Hilfe der erweiterten Lagrange-Funktion kann man das Problem (5) in eine Pro-blemstellung ohne Nebenbedingung umschreiben. Um dies zu motivieren, werden in Ab-schnitt 3.1.1 und 3.1.2 zwei weitere Verfahren vorgestellt, welche schließlich gemeinsamdie Methode mit der erweiterten Lagrange-Funktion liefern, siehe dazu Abschnitt 3.1.3.

3.1.1 Lagrange-Funktion

Das Ziel ist es, die Problemstellung (5) mit der geforderten Nebenbedingung g(x) = 0 zulosen. Bei Minimierungsproblemen ohne Nebenbedingung kann auf bekannte Verfahren,wie dem Newton-Verfahren, siehe zum Beispiel [8], Abschnitt 5.5, zuruckgegriffen werden.

Die Lagrange-Funktion ist dabei eine Moglichkeit von verschiedenen, vergleiche dazuzum Beispiel den nachsten Unterabschnitt, die Nebenbedingung g(x) = 0 von (5) mitder Funktion J als neue Funktion L zu verknupfen und eine zugehorige Problemstellungaufzustellen. Dieses neue Problem soll dabei keine Nebenbedingung mehr enthalten unddie Funktion L so gewahlt werden, dass die ursprungliche Nebenbedingung g(x) = 0

”automatisch“ durch das Losen des neuen Problems erfullt wird.Als erstes folgt ein Hilfssatz, der als Motivation fur die Lagrange-Funktion dient.

Satz 1. (Lagrange-Multiplikator, vergleiche [21], Proposition 1.19)Seien n ∈ N, J, g ∈ C2(Rn) und x ∈ Rn eine Losung von (5). Falls ∇g(x) 6= 0 gilt, dannexistiert ein Lagrange-Multiplikator λ ∈ R, so dass fur

L(x, λ) := J(x) + λg(x) (6)

gilt, dass

∇L(x, λ) = 0. (7)

L wird auch Lagrange-Funktion des Problems (5) genannt.

Beweis. Da ∇g(x) 6= 0 ist, existiert ein 0 6= h ∈ Rn, so dass c := (∇g(x))Th 6= 0 gilt.Definiere zu einem beliebigen aber festen y ∈ Rn die Funktionen

Φ(ε, t) := J(x+ εy + th) und

Γ(ε, t) := g(x+ εy + th).

5

Page 6: Seminar zur Approximationstheorie im Wintersemester 2009 ... · uber diese Methode, wobei dort mehr auf eine Methode zur L osung eines Minimierungs- problems mit Nebenbedingung mit

3 Dritte Methode: Erweiterte Lagrange-Methode

Es gilt Γ(0, 0) = g(x) = 0 und Γt(0, 0) = (∇g(x))Th = c 6= 0. Mit dem Satz uber impliziteFunktionen folgt, dass ein 0 < ε0 ∈ R und eine Funktion τ = τ(ε) ∈ C1((−ε0, ε0)) mitτ(0) = 0 existieren, so dass Γ(ε, τ(ε)) = 0 fur alle ε ∈ (−ε0, ε0) gilt. Differentiation nachε liefert

0 =d

∣∣∣ε=0

Γ(ε, τ(ε)) = Γε(0, 0) + Γt(0, 0)τε(0) = (∇g(x))T y + cτε(0)

und damit

τε(0) = −1

c(∇g(x))T y,

da c 6= 0 ist. Nach Voraussetzung ist x ein Minimum von J unter der Nebenbedingungg(x) = 0. Da samtliche x+ εy+ τ(ε)h fur ε ∈ (−ε0, ε0) die Nebenbedingung erfullen undτ sowie J und g stetig differenzierbar sind, gilt Φ(0, 0) ≤ Φ(ε, τ(ε)) fur alle ε ∈ (−ε0, ε0)und damit

0 =d

∣∣∣ε=0

Φ(ε, τ(ε)) = Φε(0, 0) + Φt(0, 0)τε(0)

= (∇J(x))T y +

[−1

c(∇J(x))Th

]︸ ︷︷ ︸

=:λ

(∇g(x))T y. (8)

Der Wert λ = − (∇J(x))T h(∇g(x))T h

ist unabhangig von y, so dass (8) fur alle Richtungsableitungen

in den x-Koordinaten von L(x, λ) = J(x) + λg(x) gilt. Nach Voraussetzung gilt bereits∂L∂λ (x, λ) = g(x) = 0, woraus schließlich ∇L(x, λ) = 0 folgt.

Bemerkung 2. In Satz 1 kann

λ := −(∇J(x))T∇g(x)

|∇g(x)|2

gewahlt werden, da im Beweis h := ∇g(x) die Bedingung (∇g(x))Th 6= 0 erfullt.

Satz 1 bietet jedoch nur eine notwendige Bedingung an die Losung von (5), da furein Maximierungsproblem wie (5), wobei lediglich das

”min“ durch

”max“ ersetzt wird,

Satz 1 vollig analog bewiesen werden kann. Das heißt also, falls (x, λ) ∈ Rn × R Glei-chung (7) erfullt, dann kann dies eine Losung von (5) sein, aber auch vom analogenMaximierungsproblem.

Eine Moglichkeit, um an eine hinreichende Bedingung zu gelangen, ist das Aufstel-len eines Sattelpunktproblems. Dazu folgt zunachst die Definition eines Sattelpunktes,welche speziell an diese Situation angepasst ist.

6

Page 7: Seminar zur Approximationstheorie im Wintersemester 2009 ... · uber diese Methode, wobei dort mehr auf eine Methode zur L osung eines Minimierungs- problems mit Nebenbedingung mit

3.1 Einfuhrung der erweiterten Lagrange-Methode

Definition 1. (Sattelpunkt, vergleiche [9], Teil Eins, Kapitel 3, Definition 3.2)Seien m,n ∈ N und die Funktion L : Rn×Rm → R gegeben. Ein Punkt (x, λ) ∈ Rn×Rmmit der Eigenschaft

L(x, λ) ≤ L(x, λ) ≤ L(x, λ) fur alle (x, λ) ∈ Rn × Rm (9)

wird Sattelpunkt von L genannt.

Nun folgt das Resultat, welches eine Problemstellung ohne Nebenbedingung liefert,dessen Losung das Minimierungsproblem (5) lost.

Satz 2. (vergleiche [9], Teil Eins, Kapitel 3, Proposition 3.1)Sei das Minimierungsproblem (5) gegeben. Dann gilt fur jeden Sattelpunkt (x, λ) ∈ Rn×Rvon

L(x, λ) := J(x) + λg(x), (10)

dass x eine Losung von (5) ist.

Beweis. Sei (x, λ) ein Sattelpunkt von L. Dann gilt

L(x, λ) ≤ L(x, λ) fur alle λ ∈ R ⇔ λg(x) ≤ λg(x) fur alle λ ∈ R⇔ g(x) = 0.

Die zweite Ungleichung des Sattelpunkts liefert unter Berucksichtigung von g(x) = 0,dass

L(x, λ) ≤ L(x, λ) fur alle x ∈ Rn ⇔ J(x) ≤ J(x) + λg(x) fur alle x ∈ Rn

Insbesondere folgt dann J(x) ≤ J(x) fur alle x ∈ Rn mit g(x) = 0, so dass x schließlicheine Losung von (5) ist.

Satz 2 zeigt also, dass man anstelle einer Minimierung unter einer Nebenbedingungauch versuchen kann ein geeignetes Sattelpunktproblem ohne Nebenbedingung zu losen.Fur die Ruckrichtung, das heißt fur den Nachweis der Existenz eines λ ∈ R, so dass fureine Losung x von (5) der Punkt (x, λ) ein Sattelpunkt von (10) ist, verlangt jedoch einigeVoraussetzungen, vergleiche dazu zum Beispiel [9], Teil Eins, Kapitel 3, Proposition 3.1.

Sei nun speziell das Problem (5) mit der symmetrisch positiv definiten Matrix A ∈Rn×n, den Vektoren b, c ∈ Rn und den Funktionen

J(x) :=1

2xTAx− bTx,

g(x) := cTx

7

Page 8: Seminar zur Approximationstheorie im Wintersemester 2009 ... · uber diese Methode, wobei dort mehr auf eine Methode zur L osung eines Minimierungs- problems mit Nebenbedingung mit

3 Dritte Methode: Erweiterte Lagrange-Methode

gegeben. Das zugehorige Sattelpunktproblem aus Satz 2 kann dann beispielsweise mitdem sogenannten Uzawa-Algorithmus, vergleiche [10], Gleichungen (2.1) bis (2.3) furr = 0, und [2], Kapitel 4, §5, berechnet werden. Es handelt sich dabei um ein iterativesVerfahren zur Losung solcher Sattelpunktprobleme, welches von besonderem Interesseist, da Dean und Glowinski eine Abwandlung dieses Algorithmus fur das Losen derMonge-Ampere-Gleichung verwenden. Es zeigt sich jedoch, dass man die Konvergenz-Eigenschaften verbessern kann, indem man die Lagrange-Funktion durch einen weiterenTerm erganzt. Um diese Erganzung zu motivieren, wird zunachst noch eine weitereMethode zur Losung von (5) eingefuhrt.

3.1.2 Strafterm-Methode

Neben der Umformulierung von (5) in ein Sattelpunktproblem der Lagrange-Funktionaus Satz 2, ist es moglich, die Nebenbedingung als sogenannten Strafterm zur Minimie-rungsbedingung zu addieren. Man wahlt dazu einen Parameter r > 0 und minimiertnun

Lr(x) := J(x) +r

2g2(x). (11)

Auf diese Weise entsteht wieder ein Optimierungsproblem ohne Nebenbedingung. DerUnterschied von (11) zur Lagrange-Funktion (10) ist, dass der Zusatzterm nun nichtnegativ werden kann, so dass hier tatsachlich wieder wie in der ursprunglichen Problem-stellung eine Minimierungsaufgabe zu losen ist.

Bemerkung 3. Die Funktion g2 in (11) ist der sogenannte quadratische Strafterm. Mankann auch eine beliebige Funktion g ∈ C1(R) nehmen, so dass g g ≥ 0 auf ganz Rn und

g(g(x)) = 0 ⇔ g(x) = 0

gilt. Die Ableitung ist dann gerade durch ∇(g g)(x) = g′(g(x))∇g(x) gegeben.Speziell die Funktion g(g(x)) := g2(x) liefert damit eine sehr einfach zu bestimmende

Ableitung, weshalb diese haufig in der Literatur und auch in dieser Ausarbeitung ver-wendet wird.

Es zeigt sich aber, dass Losungen von (11), das heißt Minimierer, fur beliebiges r > 0im Allgemeinen keine Losung von (5) sind. Um dies zu zeigen, folgt ein kleines Gegen-beispiel.

Beispiel 1. (Losung von (11) ist im Allgemeinen keine Losung von (5))Seien J(x) := x2, g(x) := x − 1 und r = 2. Dann ist x = 1 die einzige Losung von (5),da es als einziges die Bedingung g(x) = 0 erfullt. Aber es gilt

Lr(x) = x2 + (x− 1)2 = 2x2 − 2x+ 1 = 2(x− 1

2

)2+

1

2

8

Page 9: Seminar zur Approximationstheorie im Wintersemester 2009 ... · uber diese Methode, wobei dort mehr auf eine Methode zur L osung eines Minimierungs- problems mit Nebenbedingung mit

3.1 Einfuhrung der erweiterten Lagrange-Methode

mit Minimum an der Stelle 12 und damit

Lr

(1

2

)=

1

2< 1 = Lr(x),

g(1

2

)=

1

26= 0.

Das bedeutet, dass das Minimum von Lr die Nebenbedingung g(x) = 0 nicht erfullt unddamit keine Losung von (5) sein kann.

Mit Hilfe der Funktion Lr aus (11) kann man jedoch trotzdem an die Losung von (5)gelangen: Seien fur r > 0 der Punkt xr ∈ Rn ein Minimierer von Lr und x eine Losungvon (5). Dann gilt

Lr(xr) = J(xr) +r

2g2(xr) ≤ Lr(x) = J(x) +

r

2g2(x)︸ ︷︷ ︸

=0

= J(x) (12)

fur alle r > 0. Auf der linken Seite der Ungleichung muss g2(xr) jedoch nicht Null sein,da J(xr) < J(x) gelten und damit die Ungleichung weiterhin gultig sein kann, sieheBeispiel 1. Um mit Hilfe von Lr an x zu gelangen, kann der Grenzwert fur r gegen ∞betrachtet werden, falls g(xr) 6= 0 fur alle r > 0 ist. Der Grund dafur ist folgender: In(12) hangt die rechte Seite nicht mehr von r ab, woraus schließlich

limr→∞

Lr(xr) = limr→∞

[J(xr) +

r

2g2(xr)

]≤ J(x) <∞ (13)

folgt. Dann musslimr→∞

g2(xr) = 0

gelten, da sonst limr→∞r2g

2(xr) =∞ und damit limr→∞ J(xr) + r2g

2(xr) =∞ gilt, wasein Widerspruch zu (13) ist. Falls nun limr→∞ xr existiert, setze

x := limr→∞

xr.

Da g und J stetig sind, minimiert x ∈ Rn die Funktion J mit Nebenbedingung g(x) = 0und ist damit eine Losung von (5). Es muss jedoch nicht x = x gelten, es sei denn dieLosung von (5) ist eindeutig.

Beim Ausfuhren einer solchen Methode auf einem Rechner treten jedoch fur sehr großer aufgrund der Maschinengenauigkeit starke Rundungsfehler im Term r

2g2 auf. Dadurch

lasst sich nur ein sehr ungenaues Ergebnis der xr fur große r und schließlich auch fur xerwarten.

Bemerkung 4. (Verbindung zum Lagrange-Multiplikator)Da xr jeweils ein Minimum von Lr ohne Nebenbedingung und Lr stetig differenzierbarist, gilt

0 = ∇Lr(xr) = ∇J(xr) + rg(xr)∇g(xr).

Wenn nun ∇g(x) 6= 0 ist, dann konvergiert λr := rg(xr) gegen einen Lagrange-Multi-plikator λ aus Satz 1, da dann die geforderte Eigenschaft in (7) erfullt sind.

9

Page 10: Seminar zur Approximationstheorie im Wintersemester 2009 ... · uber diese Methode, wobei dort mehr auf eine Methode zur L osung eines Minimierungs- problems mit Nebenbedingung mit

3 Dritte Methode: Erweiterte Lagrange-Methode

3.1.3 Erweiterte Lagrange-Funktion

Die ursprungliche Motivation, die zur sogenannten erweiterten Lagrange-Funktion ge-fuhrt hat, war es, eine Methode zur Losung von (5) zu entwickeln, die das Problem dereinfach anzuwendenden Strafterm-Methode fur große r beseitigt.

Hestenes [17] und Powell [20] verfolgten dabei den gleichen Ansatz und erhielten alserstes das folgende Resultat:

Satz 3. (erweiterte Lagrange-Funktion, siehe [17], Theorem 2.1)Sei x ∈ Rn eine Losung von (5). Dann existiert ein Multiplikator λ ∈ R und eineKonstante r > 0, so dass x zusatzlich ein Minimum von

Lr(x) = J(x) + λg(x) +1

2r(g(x))2 (14)

ist. Umgekehrt gilt, dass jedes Minimum x von Lr mit g(x) = 0 eine Losung von (5) ist.

Lr ist die erweiterte Lagrange-Funktion zum Problem (5). Man sieht, dass dies eineMischung aus Strafterm und Lagrange-Funktion ist, vergleiche (6) und (11).

Bemerkung 5. Betrachtet man die Funktion

Lr(x, λ) := J(x) + λg(x) +1

2r(g(x))2, (15)

dann gilt Satz 2 vollig analog mit der Funktion Lr anstelle von L in (10), da der Zu-satzterm 1

2r(g(x))2 unabhangig von λ ist und im Sattelpunkt (x, λ), wobei x gleichzeitigProblem (5) lost, verschwindet und damit an dieser Stelle sein Minimum annimmt.

Auf den ersten Blick hat die erweiterte Lagrange-Funktion sehr ahnliche Eigenschaftenzur Lagrange-Funktion und es ist daher nicht klar, aus welchem Grund man diese neueFunktion wahlen sollte. Jedoch bezieht sich der in Abschnitt 3.1.1 erwahnte Uzawa-Algorithmus in [10], Gleichungen (2.1) bis (2.3), gerade auf Lr und damit auch aufL0 = L. Es zeigt sich, dass sich die Konvergenz-Eigenschaften des Uzawa-Algorithmusfur r > 0 verbessern, siehe dazu [10], Abschnitt 2.3. Aufgrund dieser Verbesserung, wirdim spateren Verlauf eine erweiterte Lagrange-Funktion mit zugehorigem Sattelpunktpro-blem verwendet.

3.2 Sattelpunktproblem

Es wird nun in Abschnitt 3.2.1 ein Sattelpunktproblem fur das Dirichlet-Problem derMonge-Ampere-Gleichung (1) aufgestellt. Zur Losung dieser neuen Problemstellung wirdin Abschnitt 3.2.2 ein Losungsverfahren vorgestellt, wobei dort das Losen zweier Teil-probleme gefordert wird, womit sich Abschnitt 3.2.4 und 3.2.5 befassen. Zur Konvergenzdes Verfahrens folgen in Abschnitt 3.2.3 einige Bemerkungen.

10

Page 11: Seminar zur Approximationstheorie im Wintersemester 2009 ... · uber diese Methode, wobei dort mehr auf eine Methode zur L osung eines Minimierungs- problems mit Nebenbedingung mit

3.2 Sattelpunktproblem

3.2.1 Sattelpunkt-Formulierung der Monge-Ampere-Gleichung

Es werden nun starke Losungen des Dirichlet-Problems der Monge-Ampere-Gleichung(1) betrachtet, das heißt die Losungen u sind in H2(Ω). Um einen Uberblick uber dieim Verlaufe dieses Abschnittes benotigten Mengen zu erhalten, werden diese bereits zuBeginn durch

Q := q ∈ (L2(Ω))2×2 : q = (qij)1≤i,j≤2, q21 = q12,Qf := q ∈ Q : det q = f fast uberall in Ω,Vv := u : Ω→ R : u ∈ H2(Ω), u = v auf ∂Ω,

Ef,v := u ∈ Vv : detD2u = f fast uberall in Ω,Ef,v := u,q ∈ Vv ×Qf : q = D2u fast uberall in Ω

(16)

definiert, wobei D2u im Sinne der schwachen Ableitungen zu verstehen ist.

Bemerkung 6. a) Wenn v ∈ H32 (∂Ω) ist, dann ist der Raum Vv nicht leer.

b) Falls f ∈ L1(Ω), dann ist Qf 6= ∅.

Beweis. a) Diese Aussage folgt aus der Umkehrung des Spursatzes fur Sobolev-Raume,siehe [19], Theorem 4.2.4.

b) Wahlt man zu f ∈ L1(Ω) eine Funktion q : Ω→ R2×2 mit

q1,1(x) := q2,2(x) :=

√f(x), falls f(x) ≥ 0,

0, falls f(x) < 0

und

q1,2(x) := q2,1(x) :=

0, falls f(x) ≥ 0,√|f(x)|, falls f(x) < 0,

dann folgt fur i, j ∈ 1, 2, dass

‖qi,j‖2L2(Ω) =

∫Ω|qi,j(x)|2 dx ≤

∫Ω|f(x)|dx

f∈L1(Ω)< ∞

und damit q ∈ Q ist. Zudem gilt

det q(x) =

f(x), falls f(x) ≥ 0,

−|f(x)|, falls f(x) < 0

= f(x)

fur alle x ∈ Ω und damit q ∈ Qf , das heißt Qf 6= ∅.

11

Page 12: Seminar zur Approximationstheorie im Wintersemester 2009 ... · uber diese Methode, wobei dort mehr auf eine Methode zur L osung eines Minimierungs- problems mit Nebenbedingung mit

3 Dritte Methode: Erweiterte Lagrange-Methode

Sei nun Efv 6= ∅, das heißt (1) hat Losungen in Vv, wobei an dieser Stelle keine Aussageuber die Eindeutigkeit der Losung gemacht wird. Es macht also Sinn, folgendes Problemaus der Variationsrechnung zu betrachten:

Finde u ∈ Ef,v,mit J(u) ≤ J(u)

fur alle u ∈ Ef,v,

(17)

wobei

J(u) :=1

2

∫Ω|∆u|2 dx =

1

2‖∆u‖2L2(Ω).

Es wird dabei der Laplace-Operator gewahlt, da dies spater auf das iterative Losenvon Poisson-Gleichungen fuhrt, welche lineare partielle Differentialgleichungen sind unddamit im Vergleich zur Monge-Ampere-Gleichung leichter zu losen sind.

Die folgenden Ideen wurde unter anderem durch Fortin und Glowinski [11], eine An-wendung zusammen mit Marrocco [12] und Glowinski [13], Kapitel 6, motiviert.

Das Problem (17) wird zunachst umgeschrieben, so dass man in einem weiteren Schrittdie Nebenbedingung detD2u = f von u trennen kann. Aus (17) wird dann mit q = D2uund q = D2u die aquivalente Problemstellung:

Finde u, q ∈ Ef,v,mit j(u, q) ≤ j(u,q)

fur alle u,q ∈ Ef,v,

(18)

wobei

j(u,q) := J(u) =1

2

∫Ω|∆u|2 dx.

Ziel ist es nun, das Variationsproblem (18), welches die Nebenbedingung

D2u− q = 0 fast uberall in Ω

beinhaltet, zu losen. Dazu wird mit Hilfe des erweiterten Lagrange-Funktionals, ahnlichwie in (15), die Nebenbedingung in das zu minimierende Funktional geschrieben und derLagrange-Multiplikator variabel gehalten, so dass man wieder ein Sattelpunktproblemaufstellen kann, um die Idee von Satz 2 auszunutzen. Zu einem positiven Parameter rund u,q,µ ∈ (Vv ×Qf )×Q sei das erweiterte Lagrange-Funktional durch

Lr(u,q;µ) :=1

2

∫Ω|∆u|2 dx+

r

2

∫Ω‖D2u− q‖2F dx+

∫Ωµ : (D2u− q) dx (19)

mit µ : q =∑2

i,j=1 µijqij und der Frobeniusnorm ‖D2u‖2F =∑2

i,j=1 |uxixj |2 definiert. Der

Term D2u−q kann dabei als Storung von D2u aufgefasst werden, da fur die Losung von(18) gerade D2u = q gelten soll. Es wird also die gestorte Hessematrix von u betrachtet,

12

Page 13: Seminar zur Approximationstheorie im Wintersemester 2009 ... · uber diese Methode, wobei dort mehr auf eine Methode zur L osung eines Minimierungs- problems mit Nebenbedingung mit

3.2 Sattelpunktproblem

welches ein typischer Ansatz in der konvexen Analysis und der Variationsrechnung zurLosung von Minimierungsproblemen ist, vergleiche dazu zum Beispiel [9], Teil Eins,Kapitel 3, Abschnitt 1–3.

Folgendes Resultat sagt aus, dass ein Sattelpunkt des erweiterten Lagrange-Funktio-nals eine gesuchte Losung ist, vergleiche dazu auch Satz 2.

Satz 4. (Sattelpunktproblem, vergleiche [11], Theorem 2.1)Seien Lr wie in (19), die Mengen Vv, Qf und Q wie in (16) definiert und folgendesSattelpunktproblem gegeben:

Finde u, q, µ ∈ (Vv ×Qf )×Q,mit Lr(u, q;µ) ≤ Lr(u, q; µ) ≤ Lr(u,q; µ)

fur alle u,q,µ ∈ (Vv ×Qf )×Q.

(20)

Falls u, q, µ eine Losung von (20) ist, das heißt u, q, µ ist ein Sattelpunkt vonLr, dann lost u, q Problem (18) und schließlich u das Dirichlet-Problem der Monge-Ampere-Gleichung (1).

Beweis. Aus der ersten Ungleichung Lr(u, q;µ) ≤ Lr(u, q; µ) folgt sofort, dass∫Ωµ : (D2u− q) dx ≤

∫Ωµ : (D2u− q) dx fur alle µ ∈ Q

und damit

0 ≤∫

Ω(µ− µ) : (D2u− q) dx fur alle µ ∈ Q.

Sei S ∈ Q, so dass jeder Eintrag von S gerade das Vorzeichen des entsprechendenEintrags von (D2u− q) sei, das heißt Si,j = sgn(D2u− q)i,j , fur i, j ∈ 1, 2. Setzt manµ = µ− S ∈ Q, dann ergibt sich

0 ≤ −∫

Ω

2∑i,j=1

|(D2u− q)i,j |dx︸ ︷︷ ︸≥0

.

Es muss also D2u = q fast uberall in Ω und damit u, q ∈ Ef,v gelten.Die zweite Ungleichung Lr(u, q; µ) ≤ Lr(u,q; µ) liefert schließlich

j(u, q) =1

2

∫Ω|∆u|2 dx ≤ 1

2

∫Ω|∆u|2 dx = j(u,q) fur alle u,q ∈ Ef,v,

da die Zusatzterme wegen u, q, u,q ∈ Ef,v, das heißt D2u = q und D2u = q fastuberall in Ω, verschwinden. Damit ist u, q ∈ Ef,v eine Losung von (18) und u lost dannnach Definition von Ef,v das Dirichlet-Problem der Monge-Ampere-Gleichung (1).

Mit (20) ist also das Sattelpunktproblem gegeben, mit dem eine Losung von (1) be-stimmt werden soll. Im nachsten Abschnitt wird ein Algorithmus vorgestellt, um dasSattelpunktproblem (20) iterativ zu losen.

13

Page 14: Seminar zur Approximationstheorie im Wintersemester 2009 ... · uber diese Methode, wobei dort mehr auf eine Methode zur L osung eines Minimierungs- problems mit Nebenbedingung mit

3 Dritte Methode: Erweiterte Lagrange-Methode

3.2.2 Losungsverfahren fur das Sattelpunktproblem

Fur Sattelpunktprobleme von der Form (20) existieren bereits verschiedene iterativeLosungsmethoden, vergleiche dazu [11], Abschnitt 3, [13], Kapitel 6, Abschnitt 3, und[14], Algorithmus (4.20)-(4.23). Dean und Glowinski [5, 6], Abschnitt 3 respektive 4,wahlten zur Losung von (20) auf Grund seiner Einfachheit den sogenannten AlgorithmusALG2 aus der oben genannten Literatur. Es handelt sich dabei um einen Uzawa-Douglas-Rachford-Algorithmus, welcher speziell an die Problemstellung (20) angepasst wurde,siehe dazu Algorithmus 1.

Algorithmus 1 Uzawa-Douglas-Rachford-Algorithmus/ALG2 zur Losung des Sattel-punktproblems (20), siehe [6], Abschnitt 4.1

1: Wahle Startvektor u(−1),µ(0) ∈ Vv ×Q.2: for n = 0 to nmax do . nmax ist maximale Anzahl an Iterationen3: Lose Minimierungsproblem

Finde q(n) ∈ Qf ,

mit Lr(u(n−1),q(n);µ(n)) ≤ Lr(u(n−1),q;µ(n))fur alle q ∈ Qf .

(21)

4: Lose Minimierungsproblem

Finde u(n) ∈ Vv,mit Lr(u(n),q(n);µ(n)) ≤ Lr(u,q(n);µ(n))

fur alle u ∈ Vv.

(22)

5:

µ(n+1) ← µ(n) + r(D2u(n) − q(n)) (23)

6: end for7: return u(nmax),q(nmax),µ(nmax)

Bemerkung 7. Dean und Glowinski [6], Remark 4.1, empfehlen fur den Algorithmus1 als Startvektor µ(0) = 0 und u(−1) definiert durch die Losung des Poisson-Dirichlet-Problems

−∆u(−1) =√f in Ω, (24)

u(−1) = v auf ∂Ω.

Dies kann dadurch motiviert werden, dass man der Einfachheit halber annimmt, dassin Bemerkung 1 die Eigenwerte der Hessematrix von u gleich seien, das heißt λ+ =

14

Page 15: Seminar zur Approximationstheorie im Wintersemester 2009 ... · uber diese Methode, wobei dort mehr auf eine Methode zur L osung eines Minimierungs- problems mit Nebenbedingung mit

3.2 Sattelpunktproblem

λ− =: λ gilt. Dann gilt fur die Monge-Ampere-Gleichung gerade λ ∈ ±√f und damit

2λ = λ+ +λ− = ∆u ∈ ±2√f. Es macht also Sinn als rechte Seite in (24) die Funktion

2√f zu wahlen. Laut Dean und Glowinski [6], Remark 4.1, verbessert dies jedoch nicht

die Konvergenz des Algorithmus.

Das Vorzeichen der Wurzel hangt davon ab, ob die Losung konvex oder konkav seinsoll. Dean und Glowinski experimentierten unter anderem auch mit Beispielen, die einekonkave Losung verlangen, was das negative Vorzeichen erklart.

3.2.3 Konvergenz des Algorithmus

Angenommen man ersetzt in (21) die Menge Qf durch eine abgeschlossene, konvexe undnichtleere Teilmenge K von Q und die Wahl von µ(n+1) in (23) durch

µ(n+1) := µ(n) + ρ(D2u(n) − q(n))

mit 0 < ρ < r 1+√

52 . Falls Lr aus (19) ein Sattelpunkt in (Vv×K)×Q besitzt, dann wurde

von Fortin und Glowinski [11], Theorem 5.1, und Glowinski [13], Kapitel 6, Theorem 5.1,bewiesen, dass u(n),q(n),µ(n)n∈N aus Algorithmus 1 in (Vv ×K) ×Q gegen einenSattelpunkt von Lr konvergiert.

Die Menge Qf ⊂ Q ist jedoch nicht konvex, da fur zwei ungleiche p,q ∈ Qf gilt, dass

det

(1

2p +

1

2q

)im Allgemeinen

6= 1

2det(p) +

1

2det(q) = f,

woraus folgt, dass Qf im Allgemeinen nicht konvex ist.

Aufgrund der Erfahrung mit nichtlinearen Problemen der finiten Elastizitat, auf denenAlgorithmen wie Algorithmus 1 angewandt wurden, siehe zum Beispiel [14], Kapitel 7,erwarten Dean und Glowinski [6], siehe Abschnitt 4.2, fur r > 0 auch hier die Konvergenz.

3.2.4 Losen von Teilproblem (21)

Die Idee zur Losung von Teilproblem (21) ist es, nicht das Integral zu minimieren,sondern den Integranden. Dies geschieht punktweise an endlich vielen Punkten in Ω.Die Wahl dieser Punkte hangt mit der Diskretisierung des Sattelpunktproblems (20)zusammen, was erst in Abschnitt 3.3 behandelt wird. Aus diesem Grund folgen dieDetails erst in Abschnitt 3.3.4.

3.2.5 Losen von Teilproblem (22)

Das Teilproblem (22) kann mit Mitteln der Variationsrechnung umgeschrieben werden.Sei dazu u(n) ∈ Vv eine Losung von (22) und ϕ ∈ H2(Ω)∩H1

0 (Ω). Sei weiter ε > 0. Nach

15

Page 16: Seminar zur Approximationstheorie im Wintersemester 2009 ... · uber diese Methode, wobei dort mehr auf eine Methode zur L osung eines Minimierungs- problems mit Nebenbedingung mit

3 Dritte Methode: Erweiterte Lagrange-Methode

Definition ist klar, dass u(n) + εϕ ∈ Vv ist, da ϕ|∂Ω ≡ 0 ist. Weil u(n) eine Losung von(22) ist, gilt

Lr(u(n),q(n);µ(n)) ≤ Lr(u(n) + εϕ,q(n);µ(n)) fur alle ϕ ∈ H2(Ω) ∩H10 (Ω)︸ ︷︷ ︸

=:V0

, ε > 0.

Sei

Φ(ε) := Lr(u(n) + εϕ,q(n);µ(n))

=

∫Ω

[1

2|∆u(n) + ε∆ϕ|2 +

r

2‖D2u(n) + εD2ϕ− q(n)‖2F

+ µ(n) : (D2u(n) + εD2ϕ− q(n))

]dx.

Da nach Voraussetzung fur ε = 0 ein Minimum vorliegt und Φ ∈ C1(R) ist, muss

0 = Φ′(0)

=

∫Ω

[∆u(n)∆ϕ+ r(D2u(n) − q(n)) : D2ϕ+ µ(n) : D2ϕ

]dx

=

∫Ω

[∆u(n)∆ϕ+ rD2u(n) : D2ϕ

]dx+

∫Ω

(µ(n) − rq(n)) : D2ϕdx︸ ︷︷ ︸=:−Ln(ϕ)

gelten, wobei das Funktional Ln(·) linear und stetig auf V0 ist. Daher lautet die zulosende Aufgabe:

Finde u(n) ∈ Vv,mit

∫Ω ∆u(n)∆ϕdx+ r

∫ΩD

2u(n) : D2ϕdx = Ln(ϕ)fur alle ϕ ∈ V0,

(25)

Dieses Problem kann mit dem Algorithmus der konjugierten Gradienten gelost werden,wobei auf den Raumen Vv und V0 mit dem Skalarprodukt (u,w) 7→

∫Ω ∆u∆w und

zugehoriger induzierter Norm gearbeitet wird. Als Resultat erhalt man Algorithmus 2,was im folgenden noch erlautert wird.

Das Verfahren der konjugierten Gradienten (CG-Verfahren) ist ursprunglich fur dasLosen linearer Gleichungssysteme Ax = b zu gegebener symmetrisch positiv definiterMatrix A ∈ Rn×n, Vektor b ∈ Rn und der gesuchten Losung x ∈ Rn vorgesehen, sie-he dazu [18], Abschnitt 3. Das Ergebnis der Rechenoperation Ax wird nun durch dasProblem (27) beschrieben. Die rechte Seite b entspricht dem Funktional Ln(·), so dassdas Residuum Ax − b durch Problem (26) ersetzt wird. Zu beachten ist nun, dass derStartvektor un,0 bereits die Randbedingung un,0|∂Ω ≡ v erfullen soll. Wenn man nun denKorrekturschritt (28) durchfuhrt, muss die Randbedingung weiterhin erfullt bleiben, das

16

Page 17: Seminar zur Approximationstheorie im Wintersemester 2009 ... · uber diese Methode, wobei dort mehr auf eine Methode zur L osung eines Minimierungs- problems mit Nebenbedingung mit

3.2 Sattelpunktproblem

Algorithmus 2 Algorithmus zum Losen von Problem (25), siehe [6], Abschnitt 4.4

1: Wahle Startvektor un,0 ∈ Vv. . zum Beispiel un,0 := u(n−1)

2: Lose das bi-harmonische Problem

Finde gn,0 ∈ V0,mit

∫Ω ∆gn,0∆ϕdx =

∫Ω ∆un,0∆ϕdx+ r

∫ΩD

2un,0 : D2ϕdx− Ln(ϕ)fur alle ϕ ∈ V0.

(26)

3: wn,0 ← gn,0.4: for k = 0 to kmax do . kmax ist maximale Iterationsanzahl5: Lose das bi-harmonische Problem

Finde gn,k ∈ V0,mit

∫Ω ∆gn,k∆ϕdx =

∫Ω ∆wn,k∆ϕdx+ r

∫ΩD

2wn,k : D2ϕdxfur alle ϕ ∈ V0.

(27)

6: Berechne . vergleiche CG-Verfahren, siehe [18], Abschnitt 3

ρn,k ←∫

Ω |∆gn,k|2 dx∫

Ω ∆gn,k∆wn,k dx,

un,k+1 ← un,k − ρn,kwn,k, (28)

gn,k+1 ← gn,k − ρn,kgn,k.

7: if∫Ω |∆g

n,k+1|2 dx∫Ω |∆gn,0|2 dx

≤ ε then8:

u(n+1) ← un,k+1,

break.

9: else10:

γn,k ←∫

Ω |∆gn,k+1|2 dx∫

Ω |∆gn,k|2 dx,

wn,k+1 ← gn,k+1 + γn,kwn,k.

11: end if12: end for13: return u(n+1)

17

Page 18: Seminar zur Approximationstheorie im Wintersemester 2009 ... · uber diese Methode, wobei dort mehr auf eine Methode zur L osung eines Minimierungs- problems mit Nebenbedingung mit

3 Dritte Methode: Erweiterte Lagrange-Methode

heißt die Randwerte von w mussen Null sein. Aus diesem Grund wird in Problem (26)und (27) im Raum V0 anstelle von Vv gesucht. Nach diesen Veranderung der ursprung-lichen Problemstellung Ax = b entspricht Algorithmus 2 gerade dem Verfahren derkonjugierten Gradienten, vergleiche dazu [18], Gleichungen (3:1a)–(3:1f).

Um ein solchen Algorithmus implementieren zu konnen, bedarf es vor allem zumLosen von (26) und (27) einer geeigneten Diskretisierung. Eine Moglichkeit dazu wirdim folgenden Abschnitt beschrieben.

3.3 Finite-Elemente-Approximation

Laut Satz 4 kann das Sattelpunktproblem (20) anstelle der Monge-Ampere-Gleichung (1)mit der Dirichlet-Randbedingung gelost werden. Das Losungsverfahren aus Abschnitt 3.2zur Losung dieses Sattelpunktproblems kann in dieser Form jedoch nicht implementiertwerden, da zum Beispiel im Algorithmus 2 der unendlich-dimensionale FunktionenraumH2(Ω) ∩H1

0 (Ω) verwendet wird und Integrale berechnet werden mussen.

Als geeignete Diskretisierung des Problems (1) verwenden Dean und Glowinski [6], Ka-pitel 5, eine Methode basierend auf einem Finite-Elemente-Verfahren. Daher werden inAbschnitt 3.3.1 zunachst einige allgemeine Bemerkungen zur Finite-Elemente-Methodegemacht. Anschließend folgt in Abschnitt 3.3.2 das diskretisierte Dirichlet-Problem derMonge-Ampere-Gleichung (1) und daraufhin in Abschnitt 3.3.3 nach dem Schema vonAbschnitt 3.2.1 ein zugehoriges Sattelpunktproblem. Zusatzlich wird in Abschnitt 3.3.3ein iteratives Losungsverfahren vorgestellt, welches ein diskretisiertes Analogon zu demin Abschnitt 3.2.2 darstellt. Am Ende mussen wieder, wie in Abschnitt 3.2.4 und 3.2.5,zwei Teilprobleme gelost werden. Dies geschieht in Abschnitt 3.3.4 und 3.3.5.

3.3.1 Finite-Elemente-Methode

In ihrer einfachsten Form ist die Finite-Elemente-Methode, zur Losung eines Problemsauf einem geeigneten unendlich-dimensionalen Funktionenraum V auf dem Gebiet Ω ⊂R2, durch die folgenden drei Grundaspekte charakterisiert, vergleiche dazu [3], Einleitungvon Kapitel 2 und Abschnitt 2.1, und [4], Abschnitt 12.4:

a) Die abgeschlossene Menge Ω wird in endlich viele abgeschlossene Teilmengen Tifur i ∈ 1, 2, ...,M, wie zum Beispiel in Dreiecke, zerlegt. Sei hTi := diam(Ti)der Durchmesser der Teilmenge Ti fur i = 1, ...,M . Dann bezeichnet man mith := maxi∈1,2,...,M hTi die Menge Th := T1, T2, ..., TM als Triangulierung von Ω.

b) Es wird ein endlich-dimensionaler Finite-Elemente-Raum Vh ⊂ V gewahlt, zumBeispiel Vh := v ∈ C0(Ω) : v|K ∈ Pk fur alle K ∈ Th, wobei Pk der Raum derPolynome in mehreren Variablen vom Totalgrad kleiner oder gleich k sei. Anschlie-ßend kann die diskrete Problemstellung aufgestellt werden. Die Funktionen in Vhwerden als Finite-Elemente bezeichnet.

18

Page 19: Seminar zur Approximationstheorie im Wintersemester 2009 ... · uber diese Methode, wobei dort mehr auf eine Methode zur L osung eines Minimierungs- problems mit Nebenbedingung mit

3.3 Finite-Elemente-Approximation

c) Es sollte eine Basis von Vh existieren, dessen Elemente einen moglichst kleinenTrager besitzen, um eine Lokalitatseigenschaft zu erhalten. Das heißt, dass jedesElement aus Vh als endlich viele Koeffizienten aus R, welche eindeutig mit dementsprechenden Basiselement identifiziert sind, aufgefasst werden kann. Eine Ande-rung oder Storung eines Koeffizienten wirkt sich dann nur auf den

”kleinen“ Trager

des zugehorigen Basiselements aus, so dass die ubrigen Bereiche unverandert blei-ben.Mit Hilfe dieser Basis kann dann ein entsprechendes diskretisiertes Losungsverfah-ren entwickelt werden, welches schließlich implementiert werden kann.

In der Literatur werden haufig Dreiecke fur die Zerlegung von Ω verwendet, da sich jedespolygonale Gebiet in Dreiecke zerlegen lasst und eine solche Zerlegung im Allgemeinenleichter zu bestimmen ist, als zum Beispiel die Zerlegung in Vierecke, falls dies moglichist. Dean und Glowinski [5, 6, 7] wahlten ebenfalls eine Zerlegung in Dreiecke. Dazufolgen nun zwei wichtige Definitionen.

Definition 2. (Zulassige und uniforme Triangulierung, siehe [2], Kapitel 2, Definition5.1)Sei Ω ⊂ R2 ein Gebiet.

a) Eine Triangulierung Th = T1, ..., TM von Ω in Dreiecke heißt zulassig, wennfolgende Eigenschaften erfullt sind:

(i) Es gilt Ω =⋃Mi=1 Ti.

(ii) Besteht Ti∩Tj aus genau einem Punkt, so ist dieser ein Eckpunkt sowohl vonTi als auch von Tj.

(iii) Besteht Ti ∩Tj fur i 6= j aus mehr als einem Punkt, so ist Ti ∩Tj eine Kantesowohl von Ti als auch von Tj.

b) Eine Familie von Zerlegungen Th von Ω heißt uniform, wenn es eine Zahl κ > 0gibt, so dass jedes Element T von Th einen Kreis mit Radius RT ≥ h

κ enthalt.

Ein Beispiel einer unzulassigen Triangulierung ist in Abbildung 1 (a) zu sehen, da sodie Bedingung (iii) in Definition 2 verletzt ist.

Die Definition einer uniformen Familie von Triangulierungen ist unter anderem eineArt Winkelbedingung. Sie sagt aus, dass die Winkel der Dreiecke nicht beliebig spitzwerden durfen, siehe dazu Abbildung 1 (b), denn je spitzer ein Winkel im Dreieck, destokleiner ist sein Innenkreisradius, welcher durch h

κ nach unten beschrankt ist. Ebenso wirddie Große bzw. der Durchmesser der Dreiecke von dem Dreieck mit dem großten Durch-messer nach unten beschrankt, da bei schrumpfendem Durchmesser der Innenkreisradiusebenfalls kleiner wird.

19

Page 20: Seminar zur Approximationstheorie im Wintersemester 2009 ... · uber diese Methode, wobei dort mehr auf eine Methode zur L osung eines Minimierungs- problems mit Nebenbedingung mit

3 Dritte Methode: Erweiterte Lagrange-Methode

(a) Hangender Knoten (b) Der Innenkreisradius ist bei glei-chem Durchmesser abhangig vonden Winkeln.

Abbildung 1: Unzulassige Triangulierung und Schaubild fur uniforme Familie von Tri-angulierungen.

3.3.2 Gemischte Finite-Elemente-Approximation von (1)

Es wird nun eine gemischte Finite-Elemente-Methode verwendet, vergleiche dazu [3],Kapitel 7 und [15], Kapitel 4, Abschnitt 4.5, und Appendix 4, Abschnitt 3.

Sei Ω ⊂ R2 der Einfachheit halber wieder ein polygonales Gebiet, so dass es in Dreieckezerlegt werden kann, und Th eine zulassige Triangulierung von Ω in Dreiecke. Der Raumder Polynome in zwei Variablen vom maximalen Totalgrad eins sei durch

P1 := P : R2 → R, (x, y) 7→ a1x+ a2y + a3 : a1, a2, a3 ∈ R

definiert. Von Th ausgehend werden die Raume L2(Ω), H1(Ω) und H2(Ω) durch denendlich dimensionalen Raum

Vh := w ∈ C0(Ω) : w|T ∈ P1 fur alle T ∈ Th

approximiert. Da im Algorithmus 2 auch H10 (Ω) gebraucht wird, sei

V0,h := Vh ∩H10 (Ω) = w ∈ Vh : w|∂Ω ≡ 0.

Zu einer genugend oft schwach differenzierbaren Funktion u seien die Ableitungen durch

Di(u) :=∂u

∂xiund D2

ij(u) :=∂2u

∂xi∂xj

20

Page 21: Seminar zur Approximationstheorie im Wintersemester 2009 ... · uber diese Methode, wobei dort mehr auf eine Methode zur L osung eines Minimierungs- problems mit Nebenbedingung mit

3.3 Finite-Elemente-Approximation

fur i, j ∈ 1, 2 definiert, wobei x1 := x und x2 := y jeweils die entsprechenden Koordi-naten seien. Aus der Greenschen Formel, siehe [3], Gleichung (1.2.4), folgt, dass∫

ΩD2ii(u)ϕdx = −

∫ΩDi(u)Di(ϕ) dx fur alle ϕ ∈ H1

0 (Ω) und alle i ∈ 1, 2 (29)

und∫ΩD2

12(u)ϕdx = −1

2

∫Ω

[D1(u)D2(ϕ) +D2(u)D1(ϕ)] dx fur alle ϕ ∈ H10 (Ω) (30)

gilt. Sei nun u ∈ Vh. Die Funktion u ist dann Lipschitz-stetig und damit auch in H1(Ω).Es muss jedoch nicht notwendigerweise u ∈ H2(Ω) gelten. Mit Hilfe von (29) und (30)wird das diskrete Analogon D2

hij des Differentialoperators D2ij so definiert, dass D2

hij(u) ∈V0,h fur i, j ∈ 1, 2 und∫

ΩD2hii(u)ϕdx = −

∫ΩDi(u)Di(ϕ) dx fur alle ϕ ∈ V0,h und alle i ∈ 1, 2 (31)

bzw.∫ΩD2h12(u)ϕdx = −1

2

∫Ω

[D1(u)D2(ϕ) +D2(u)D1(ϕ)] dx fur alle ϕ ∈ V0,h (32)

gilt. Aus [2], Kapitel 2, Satz von Lax-Milgram und Satz 2.2, folgt, dass die FunktionenD2hij(u) eindeutig durch (31) und (32) festgelegt sind. Um die diskreten zweiten Ablei-

tungen auf der linken Seite der Integrale in (31) und (32) zu berechnen, wird das Integralmit der Trapezregel approximiert. Dabei wird wie folgt vorgegangen:

a) Zunachst werden die Menge der Knoten Σh von Th, das heißt Σh ist die Mengealler Eckpunkte von den Dreiecken in Th, und Σ0,h := P ∈ Σh : P /∈ ∂Ωdefiniert. Die Anzahl der Elemente beider Mengen wird durch Nh := #Σh undN0,h := #Σ0,h definiert. Dann gilt dimVh = Nh und dimV0,h = N0,h. Seien nun

Σ0,h = PkN0,h

k=1 und Σh = Σ0,h ∪ PkNhk=N0,h+1.

b) Zu einem Punkt Pk ∈ Σh ist die Funktion wk ∈ Vh durch

wk(Pi) =

1, falls i = k,

0, falls i ∈ 1, ..., Nh und i 6= k

eindeutig definiert, da so in jedem Dreieck T drei Funktionswerte festgelegt werdenund damit das Polynom wk|T ∈ P1 eindeutig bestimmt ist. Die Stetigkeit an denKanten der Dreiecke folgt aus der Linearitat der Polynome, so dass tatsachlichwk ∈ Vh gilt.

21

Page 22: Seminar zur Approximationstheorie im Wintersemester 2009 ... · uber diese Methode, wobei dort mehr auf eine Methode zur L osung eines Minimierungs- problems mit Nebenbedingung mit

3 Dritte Methode: Erweiterte Lagrange-Methode

Die Mengen Bh := wkNhk=1 und B0,h := wk

N0,h

k=1 sind dann jeweils Basen von Vhrespektive V0,h, da die Dimension der Basen mit dem zugehorigen Raum uberein-stimmt und die Elemente in den Mengen linear unabhangig sind, siehe auch [2],Kapitel 2, Abschnitt

”Dreieckelemente mit vollstandigen Polynomen“.

c) Da nun eine Basis fur V0,h bekannt ist und die Integrale in (31) und (32) linear in ϕsind, reicht es diese Gleichheiten lediglich fur die Basiselemente in B0,h zu prufen.Diese Basisfunktionen sind jeweils an allen bis auf einen Knoten gleich Null. Umdies auszunutzen, wird die Trapezregel verwendet, da hier lediglich Funktionswertean den Knoten benotigt werden, welche noch dazu die Werte sind, an denen spatereine diskrete Variante der Monge-Ampere-Gleichung (1) erfullt sein soll.Sei dazu Ak der Flacheninhalt des Gebiets, welches die Vereinigung der Elementeaus Th ist, die den Eckpunkt Pk haben. Wendet man nun auf (31) und (32) dieTrapezregel mit ϕ := wk an, so erhalt man fur i ∈ 1, 2 gerade

D2hii(u)(Pk) = − 3

Ak

∫ΩDi(u)Di(wk) dx, (33)

D2h12(u)(Pk) = D2

h21(u)(Pk) = − 3

2Ak

∫Ω

[D1(u)D2(wk) +D2(u)D1(wk)] dx (34)

fur alle k ∈ 1, 2, ..., N0,h. Das Losen der Integrale auf der rechten Seite von (33)und (34) ist dann weitestgehend nur noch die Berechnung der Flacheninhalte derDreiecke aus Th, da die Ableitungen von u und wk stuckweise konstant sind.

An dieser Stelle kann nun das Dirichlet-Problem der Monge-Ampere-Gleichung (1) dis-kretisiert werden. Angenommen die Funktion v, welche die Randbedingung vorgibt, iststetig und fh sei eine Naherung von f aus (1), welche ebenfalls stetig ist. Dann werdendie Mengen Q, Qf und Vv aus (16) durch

Qh := q ∈ (V0,h)2×2 : q = (qij)1≤i,j≤2, q21 = q12,Qf,h := q ∈ Qh : det q(Pk) = fh(Pk) fur alle k ∈ 1, 2, ..., N0,h,Vv,h := u ∈ Vh : u(P ) = v(P ) fur alle P ∈ Σh ∩ ∂Ω.

(35)

approximiert. Die Diskretisierung von (1) lautet schließlich:

Finde uh ∈ Vv,h,mit D2

h11(uh)(Pk)D2h22(uh)(Pk)− (D2

h12(uh)(Pk))2 = fh(Pk)

fur alle k ∈ 1, 2, ..., N0,h

(36)

3.3.3 Iteratives Losen von Problem (36)

Motiviert durch Abschnitt 3.2 wird nun eine diskrete Variante des Algorithmus 1 vorge-stellt. Als erstes muss dazu das Sattelpunktproblem (20) und damit auch das erweiterte

22

Page 23: Seminar zur Approximationstheorie im Wintersemester 2009 ... · uber diese Methode, wobei dort mehr auf eine Methode zur L osung eines Minimierungs- problems mit Nebenbedingung mit

3.3 Finite-Elemente-Approximation

Lagrange-Funktional Lr in (19) approximiert werden. Das Funktional Lr wird dabei furu,q,µ ∈ (Vh ×Qf,h)×Qh mit

∆hu := D2h11(u) +D2

h22(u), Dh(u) := (D2hij(u))1≤i,j≤2,

(w, w)h :=1

3

N0,h∑k=1

Akw(Pk)w(Pk), ‖w‖h := (w,w)12h ,

((s, t))h :=1

3

N0,h∑k=1

Aks(Pk) : t(Pk), ‖s‖h := ((s, s))12h

fur alle w, w ∈ V0,h und alle s, t ∈ Qh durch

Lr,h(u,q;µ) :=1

2‖∆hu‖2h +

r

2‖D2

h(u)− q‖2h + ((µ,D2h(u)− q))h (37)

definiert. Anschließend kann das diskrete Sattelpunktproblem wie folgt definiert werden:

Finde uh, qh, µh ∈ (Vv,h ×Qf,h)×Qh,mit Lr,h(uh, qh;µ) ≤ Lr,h(uh, qh; µh) ≤ Lr,h(u,q; µh)

fur alle u,q,µ ∈ (Vv,h ×Qf,h)×Qh.

(38)

Zur Losung dieses Sattelpunktproblems wird eine diskrete Variante des Algorithmus 1benutzt, siehe Algorithmus 3.

Bemerkung 8. Dean und Glowinski [6], Gleichung (5.29), wahlten als Startwert furden Algorithmus 3 eine diskrete Variante der Wahl aus Bemerkung 7. Es wird µ0

h = 0gesetzt und u−1

h ∈ Vv,h als Losung von∫Ω

(∇u−1h )T∇ϕ dx =

1

3

N0,h∑k=1

Ak√fh(Pk)ϕ(Pk) fur alle ϕ ∈ V0,h

definiert.

3.3.4 Losen von Teilproblem (39)

Fur das Teilproblem (39) muss das Funktional Lr,h(u(n−1)h ,q;µ

(n)h ) aus (37) in der Va-

riable q ∈ Qf,h minimiert werden. Da der erste Summand in Lr,h nicht von q abhangt,muss lediglich der Ausdruck

r

2‖D2

h(u(n−1)h )− q‖2h + ((µ

(n)h ,D2

h(u(n−1)h )− q))h

=Ak3

N0,h∑k=1

(r2‖(D2

h(u(n−1)h )− q)(Pk)‖2F + µ

(n)h (Pk) : (D2

h(u(n−1)h )− q)(Pk)

)

23

Page 24: Seminar zur Approximationstheorie im Wintersemester 2009 ... · uber diese Methode, wobei dort mehr auf eine Methode zur L osung eines Minimierungs- problems mit Nebenbedingung mit

3 Dritte Methode: Erweiterte Lagrange-Methode

Algorithmus 3 Diskretisierung von Algorithmus 1 zur Losung des Sattelpunktproblems(38), siehe [6], Abschnitt 5.3

1: Wahle Startvektor u(−1)h ,µ

(0)h ∈ Vv,h ×Qh.

2: for n = 0 to nmax do . nmax ist maximale Anzahl an Iterationen3: Lose Minimierungsproblem

Finde q(n)h ∈ Qf,h,

mit Lr,h(u(n−1)h ,q

(n)h ;µ

(n)h ) ≤ Lr,h(u

(n−1)h ,q;µ

(n)h )

fur alle q ∈ Qf,h.

(39)

4: Lose Minimierungsproblem

Finde u(n)h ∈ Vv,h,

mit Lr,h(u(n)h ,q

(n)h ;µ

(n)h ) ≤ Lr,h(u,q

(n)h ;µ

(n)h )

fur alle u ∈ Vv,h.

(40)

5:

µ(n+1)h ← µ

(n)h + r(D2

hu(n)h − q

(n)h ).

6: end for7: return u(nmax)

h ,q(nmax)h ,µ(nmax)

h

in q minimiert werden. Dies wiederum ist ein System von N0,h ungekoppelten Minimie-rungsproblemen, denn es reicht, wenn jeder einzelne Summand einzeln minimiert wird.Problem (39) ist damit aquivalent zum Minimieren von

r

2‖(D2u(n−1) − q)(Pk)‖2F + µ(n)(Pk) : (D2u(n−1) − q)(Pk) =

=r

2

2∑i,j=1

(D2u(n−1)(Pk)− q(Pk))2i,j +

2∑i,j=1

µ(n)i,j (Pk)(D

2u(n−1)(Pk)− q(Pk))i,j

=2∑

i,j=1

r

2qi,j(Pk)

2 −2∑

i,j=1

[2r

2(D2u(n−1))i,j(Pk) + µ

(n)i,j (Pk)

]qi,j(Pk)

+r

2|D2u(n−1)(Pk)|2 + µ(n)(Pk) : D2u(n−1)(Pk) (41)

24

Page 25: Seminar zur Approximationstheorie im Wintersemester 2009 ... · uber diese Methode, wobei dort mehr auf eine Methode zur L osung eines Minimierungs- problems mit Nebenbedingung mit

3.3 Finite-Elemente-Approximation

fur alle k ∈ 1, 2, ..., N0,h. Seien nun

b(n)1 (Pk) := r(D2u(n−1)(Pk))1,1 + µ

(n)1,1 (Pk),

b(n)2 (Pk) := r(D2u(n−1)(Pk))2,2 + µ

(n)2,2 (Pk),

b(n)3 (Pk) := r(D2u(n−1)(Pk))1,2 + µ

(n)1,2 (Pk) + r(D2u(n−1)(Pk))2,1 + µ

(n)2,1 (x),

zk := (z1,k, z2,k, z3,k)T ∈ R3,

z1,k := q1,1(Pk),

z2,k := q2,2(Pk),

z3,k := q1,2(Pk) = q2,1(Pk),

wobei letztere Gleichheit wegen der Symmetrie von q gilt. Da die beiden Summandenin (41) nicht mehr von q abhangen, muss also punktweise

2∑i,j=1

r

2qi,j(Pk)

2 −2∑

i,j=1

[2r

2(D2u(n−1))i,j(Pk) + µ

(n)i,j (Pk)

]qi,j(Pk)

=r

2(z2

1,k + z22,k + 2z2

3,k)− b(n)1 (Pk)z1,k − b

(n)2 (Pk)z2,k − b

(n)3 (Pk)z3,k

unter der Nebenbedingung

z1,kz2,k − z23,k = f(Pk)

minimiert werden. Um die Nebenbedingung zu beseitigen, wird gemaß Abschnitt 3.1.1,siehe Satz 1, die Lagrange-Funktion

L(zk, λ) :=r

2(z2

1,k + z22,k + 2z2

3,k)− b(n)1 (Pk)z1,k − b

(n)2 (Pk)z2,k − b

(n)3 (Pk)z3,k

− λ(z1,kz2,k − z23,k − f(Pk))

betrachtet. Fur einen stationaren Punkt zk, λk ∈ R3 × R von L(·, ·) gilt dann

∇L(zk, λk) =

rz1,k − λkz2,k − b

(n)1,k(Pk)

rz2,k − λkz1,k − b(n)2,k(Pk)

2rz3,k + 2λkz3,k − b(n)3,k(Pk)

z1,kz2,k − z23,k − f(Pk)

= 0. (42)

Dieses nichtlineare Gleichungssystem mit den Unbekannten z1,k, z2,k, z3,k und λk kannman beispielsweise mit dem Newton-Verfahren oder einem Quasi-Newton-Verfahren, ver-gleiche dazu [8], Kapitel 5 respektive 6, losen.

25

Page 26: Seminar zur Approximationstheorie im Wintersemester 2009 ... · uber diese Methode, wobei dort mehr auf eine Methode zur L osung eines Minimierungs- problems mit Nebenbedingung mit

4 Numerische Experimente

3.3.5 Losen von Teilproblem (40)

An dieser Stelle wird ahnlich wie in Abschnitt 3.2.5 vorgegangen. Das Teilproblem (40)ist aquivalent zu dem endlich-dimensionalen linearen Variationsproblem

Finde unh ∈ Vv,h,mit (∆hu

nh∆hϕ)h + r((D2

h(unh),D2h(ϕ))) = Ln,h(ϕ)

fur alle ϕ ∈ V0,h,

(43)

wobei Ln,h(·) in V0,h linear ist. Auch hier wird der Algorithmus des Verfahrens der konju-gierten Gradienten verwendet, wobei das Skalarprodukt (u,w) 7→ (∆hu,∆hw)h benutztwird. Es ergibt sich damit fur Problem (43) eine diskrete Variante des Algorithmus 2,siehe dazu Algorithmus 4.

Mittels des Algorithmus 3 wurde das Dirichlet-Problem der Monge-Ampere-Gleichung(1) auf das Losen folgender zwei Problemstellungen reduziert:

a) Das Losen einer Folge von nichtlinearen Gleichungssystemen, siehe Abschnitt 3.3.4bzw. Abschnitt 3.2.4 und

b) das Losen einer Folge von diskreten Poisson-Dirichlet-Problemen, siehe (44) und(45).

4 Numerische Experimente

Benamou, Froese und Oberman [1], Abschnitt 3, haben ihre beiden Methoden, siehedazu auch [16], Methode 1 und Methode 2, miteinander verglichen. Sie haben dabeiunter anderem die gleichen Beispiele gewahlt, wie schon zuvor Dean und Glowinski[5, 6, 7]. Diese Beispiele werden daher in Abschnitt 4.1 kurz vorgestellt. Anschließendwerden die drei Methoden in Abschnitt 4.2 miteinander verglichen, siehe dazu auch [1],Abschnitt 3.1 und 3.2. Methode 1 und 2 bezeichnen dabei die Verfahren in [16], Kapitel3 respektive 4, und Methode 3 die erweiterte Lagrange-Methode aus Kapitel 3.

4.1 Beispiele

Als Startwerte wahlten Benamou, Froese und Oberman [1] stets die Funktion u, welche

∆u =√

2f

lost. Dean und Glowinski wahlten stets r = 1 und als Startfunktion u die Losung von

−∆u =√f

mit vorgegebener Dirichlet-Randbedingung, siehe Bemerkung 7. Zudem verwendeten sieeine uniforme Familie von Triangulierung.

26

Page 27: Seminar zur Approximationstheorie im Wintersemester 2009 ... · uber diese Methode, wobei dort mehr auf eine Methode zur L osung eines Minimierungs- problems mit Nebenbedingung mit

4.1 Beispiele

Algorithmus 4 Algorithmus zum Losen von Problem (43), siehe [6], Abschnitt 5.5

1: Wahle Startvektor un,0h ∈ Vv,h. . zum Beispiel un,0h := u(n−1)h

2: Lose das diskrete bi-harmonische Problem

Finde gn,0h ∈ V0,h,

mit (∆hgn,0h ,∆hϕ)h = (∆hu

n,0h ,∆hϕ)h + r((D2

h(un,0h ),D2h(ϕ)))h − Lnh(ϕ)

fur alle ϕ ∈ V0,h.

(44)

3: wn,0h ← gn,0h

4: for k = 0 to kmax do . kmax ist maximale Iterationsanzahl5: Lose das bi-harmonische Problem

Finde gn,kh ∈ V0,h,

mit (∆hgn,kh ,∆hϕ)h = (∆hw

n,kh ,∆hϕ)h + r((D2

h(wn,kh ),D2h(ϕ)))h

fur alle ϕ ∈ V0,h.

(45)

6: Berechne . vergleiche CG-Verfahren, siehe [18], Abschnitt 3

ρn,k ←‖∆hg

n,kh ‖

2h

(∆hgn,kh ,∆hw

n,kh )h

,

un,k+1h ← un,k − ρn,kwn,kh ,

gn,k+1h ← gn,k − ρn,kgn,kh .

7: if‖∆hg

n,k+1h ‖2h

‖∆hgn,0h ‖

2h

≤ ε then

8:

u(n+1)h ← un,k+1

h ,

break.

9: else10:

γn,k ←‖∆hg

n,k+1h ‖2h

‖∆hgn,kh ‖2h

,

wn,k+1h ← gn,k+1

h + γn,kwn,kh .

11: end if12: end for13: return u

(n+1)h

27

Page 28: Seminar zur Approximationstheorie im Wintersemester 2009 ... · uber diese Methode, wobei dort mehr auf eine Methode zur L osung eines Minimierungs- problems mit Nebenbedingung mit

4 Numerische Experimente

Beispiel 1: klassische und starke Losung Es wurde das Gebiet Ω = (0, 1) × (0, 1)gewahlt. Die rechte Seite der Monge-Ampere-Gleichung und die Dirichlet-Randbedingungaus (1) wurden mit

f(x, y) := (1 + x2 + y2)ex2+y2

und v(x, y) := e12

(x2+y2)

definiert, siehe dazu [1], Abschnitt 4.1, und [6],”first test problem“ in Abschnitt 6. Dabei

ist v sowohl klassische als auch starke Losung des Problems.Hier konvergieren alle drei Verfahren gegen die exakte Losung v, wobei Methode 2

deutlich schneller arbeitet als Methode 1.

Beispiel 2: klassische und starke Losung Hier wurde ebenfalls das Gebiet Ω = (0, 1)×(0, 1) gewahlt. Die Funktionen sind durch

f(x, y) :=1√

x2 + y2und v(x, y) :=

2√

2

3(x2 + y2)

34

gegeben, siehe dazu [1], Abschnitt 4.2, und [6],”third test problem“ in Abschnitt 6. Auch

hier ist v gleichzeitig eine klassische und starke Losung von (1). Im Nullpunkt weist feine Singularitat auf, was bedeutet, dass die rechte Seite der Monge-Ampere-Gleichungunbeschrankt ist.

Hier ist die gleiche Situation wie in Beispiel 1 zu beobachten, das heißt alle Verfahrenkonvergieren und Methode 2 ist viel schneller als Methode 1.

Beispiel 3: keine glatte Losung In [1] wurde Ω = (−1, 1) × (−1, 1) und in [6] Ω =(0, 1)× (0, 1) gewahlt. Die Monge-Ampere-Gleichung und die Dirichlet-Randbedingungwurden durch

f(x, y) := 1 und

v(x, y) := 1 in [1], Abschnitt 5.3,

v(x, y) := 0 in [6],”fourth test problem“ in Abschnitt 6,

festgelegt. Fur dieses Beispiel existiert trotz glatter Daten f und v keine exakte Losung.Hier zeigt sich bei allen drei Verfahren, dass die berechnete Losung im Inneren konvex

bzw. konkav ist, wie es auch gewunscht ist, jedoch in der Nahe des Randes von dieserForm abweicht, da keine exakte Losung existiert. Bei diesem Beispiel arbeitet dieses malMethode 1 deutlich schneller als Methode 2.

Beispiel 4: klassische aber keine starke Losung In diesem Beispiel verwendeten beideQuellen wieder Ω = (0, 1) × (0, 1). Sie verwendeten die gleiche rechte Seite der Monge-Ampere-Gleichung, aber unterschiedliche Randbedingungen. Diese Funktionen wurden

28

Page 29: Seminar zur Approximationstheorie im Wintersemester 2009 ... · uber diese Methode, wobei dort mehr auf eine Methode zur L osung eines Minimierungs- problems mit Nebenbedingung mit

4.1 Beispiele

durch

f(x, y) :=2

(2− x2 − y2)2und

v(x, y) := −√

2− x2 − y2 in [1], Abschnitt 5.4,

v(x, y) :=√

2− x2 − y2 in [6],”second test problem“ in Abschnitt 6,

definiert. Die klassische Losung v ist dabei nur noch in H1(Ω) und damit keine Losungim starken Sinne.

Methode 3 divergiert bei diesem Beispiel fur jedes r > 0. Die anderen beiden konver-gieren gegen v, wobei die zweite Methode schneller ist.

(a) Exakte Losung von Beispiel 1 (b) Starke Losung von Beispiel 2

(c) Berechnete Losung von Beispiel 3 und Null-randwerten mit Methode 3

(d)”Losung“ von Beispiel 4 und der zweitenWahl von v

Abbildung 2: Die Bilder wurden aus [6], siehe Abbildung 6.1, 6.2, 6.4 und 6.5, entnom-men.

29

Page 30: Seminar zur Approximationstheorie im Wintersemester 2009 ... · uber diese Methode, wobei dort mehr auf eine Methode zur L osung eines Minimierungs- problems mit Nebenbedingung mit

Literatur

4.2 Vergleich

Methode 1 und 2, siehe [16], Kapitel 3 und 4, konvergieren bei allen Beispielen ausAbschnitt 4.1 gegen eine Losung von (1). Bei Methode 3 sieht es ahnlich aus, nur dassdies bei Beispiel 4 fur beliebiges r > 0 divergiert. Methode 2 und 3 benotigen grobgesagt die gleiche Anzahl an Iterationen bis die berechnete Losung nahe an der exaktenLosung ist. Dean und Glowinski [5, 6, 7] haben keine Rechenzeiten angegeben. Es wirdjedoch vermutet, dass Methode 2 schneller arbeitet, da bei dieser in jeder Iterationlediglich eine Poisson-Gleichung gelost werden muss, wohingegen in Methode 3 nochweitere Berechnungen erforderlich sind.

Methode 2 arbeitet bei glatten Losungen schneller als Methode 1. Bei Losungen die

”viele“ Singularitaten enthalt ist Methode 1 schneller. Wahrend die erste Methode immer

ungefahr die gleiche Anzahl an Iterationen unabhangig von der Glattheit der Losungbenotigt, nimmt die Rechenzeit fur Methode 2 bei Abnahme der Glattheit zu.

Weitere Details und Ergebnisse zu den Beispielen sind in [1], Abschnitt 4 und 5 und6, und [6], Abschnitt 4, zu finden.

5 Ausblick

Das Ziel der Ausarbeitung von Yasemin Hafizogullari [16] und dieser war es, drei nume-rische Loser fur das Dirichlet-Problem der Monge-Ampere-Gleichung vorzustellen. Dieersten beiden Methoden stammen von Benamou, Froese und Oberman [1] und die drittevon Dean und Glowinski [5, 6, 7].

Es hat sich dabei gezeigt, dass die ersten beiden Methoden im Vergleich zur drittenleichter zu implementieren sind. Es stellte sich zudem heraus, dass die zweite der drittenMethode vorzuziehen ist, da Methode 2 vermutlich schneller arbeitet, ungefahr die gleicheAnzahl an Iterationen benotigt und sogar in Beispiel 4 konvergiert, bei welchem Methode3 divergiert.

Benamou, Froese und Oberman [1], Abschnitt 8, empfehlen die zweite Methode, wennbekannt ist, dass die Losungen in H2(Ω) oder in C2(Ω) sind. Falls die Glattheit derLosung unbekannt ist, sollte Methode 1 gewahlt werden. Da jedoch bisher noch nichtsuber die Konvergenz dieser Verfahren bekannt ist, kann man im nachsten Schritt versu-chen, dies genauer zu untersuchen.

Literatur

[1] Benamou, Jean-David, Brittany D. Froese und Adam M. Oberman: Twonumerical methods for the elliptic Monge-Ampere equation. Preprint, Simon FraserUniversity, 2009. http://www.divbyzero.ca/froese/w/images/4/40/MA.pdf.

[2] Braess, Dietrich: Finite Elemente. Springer, Heidelberg, 2. Auflage, 1997.

30

Page 31: Seminar zur Approximationstheorie im Wintersemester 2009 ... · uber diese Methode, wobei dort mehr auf eine Methode zur L osung eines Minimierungs- problems mit Nebenbedingung mit

Literatur

[3] Ciarlet, Philippe G.: The finite element method for elliptic problems, Band 4 derReihe Studies in Mathematics and its Applications. North-Holland, Amsterdam,1978.

[4] Dahmen, Wolfgang und Arnold Reusken: Numerik fur Ingenieure und Na-turwissenschaftler. Springer, Heidelberg, 2008.

[5] Dean, Edward J. und Roland Glowinski: Numerical solution of the two-dimensional elliptic Monge-Ampere equation with Dirichlet boundary conditions:An augmented Lagrangian approach. Comptes Rendus. Mathematique. Academiedes Sciences, Paris, 336(9):779–784, 2003.

[6] Dean, Edward J. und Roland Glowinski: An augmented Lagrangian approachto the numerical solution of the Dirichlet problem for the elliptic Monge-Ampereequation in two dimensions. Electronic Transactions on Numerical Analysis, 22:71–96, 2006.

[7] Dean, Edward J. und Roland Glowinski: Numerical methods for fully nonli-near elliptic equations of the Monge-Ampere type. Computer Methods in AppliedMechanics and Engineering, 195(13-16):1344–1386, 2006.

[8] Dennis, John E. und Robert B. Schnabel: Numerical Methods for Unconstrai-ned Optimization and Nonlinear Equations. Prentice-Hall series in computationalmathematics. Prentice-Hall, Englewood Cliffs, 1983.

[9] Ekeland, Ivan und Roger Temam: Convex analysis and variational problems,Band 1 der Reihe Studies in Mathematics and its Applications. North-Holland,Amsterdam, 1976.

[10] Fortin, Michel und Roland Glowinski: Augmented Lagrangian methods inquadratic programming. In: Fortin, Michel und Roland Glowinski (Heraus-geber): Augmented Lagrangian methods: Applications to the numerical solution ofboundary-value problems, Band 15 der Reihe Studies in Mathematics and its App-lications, Kapitel 1, Seiten 1–46. North-Holland, Amsterdam, 1983.

[11] Fortin, Michel und Roland Glowinski: On decomposition-coordination me-thods using an augmented Lagrangian. In: Fortin, Michel und Roland Glowin-ski (Herausgeber): Augmented Lagrangian methods: Applications to the numericalsolution of boundary-value problems, Band 15 der Reihe Studies in Mathematics andits Applications, Kapitel 3, Seiten 97–146. North-Holland, Amsterdam, 1983.

[12] Fortin, Michel, Roland Glowinski und Americo Marrocco: Application tothe solution of strongly nonlinear second-order boundary value problems. In: For-tin, Michel und Roland Glowinski (Herausgeber): Augmented Lagrangian me-thods: Applications to the numerical solution of boundary-value problems, Band 15

31

Page 32: Seminar zur Approximationstheorie im Wintersemester 2009 ... · uber diese Methode, wobei dort mehr auf eine Methode zur L osung eines Minimierungs- problems mit Nebenbedingung mit

Literatur

der Reihe Studies in Mathematics and its Applications, Kapitel 5, Seiten 171–216.North-Holland, Amsterdam, 1983.

[13] Glowinski, Roland: Numerical methods for nonlinear variational problems.Scientific Computation. Springer, Heidelberg, 2008.

[14] Glowinski, Roland und Patrick Le Tallec: Augmented Lagrangian andoperator-splitting methods in nonlinear mechanics, Band 9 der Reihe Studies in Ap-plied and Numerical Mathematics. Society for Industrial and Applied Mathematics,Philadelphia, 1989.

[15] Glowinski, Roland, Jacques-Louis Lions und Raymond Tremolieres: Nu-merical Analysis of variational inequalities, Band 8 der Reihe Studies in Mathema-tics and its Applications. North-Holland, Amsterdam, 1981.

[16] Hafizogullari, Yasemin: Numerische Verfahren zur Losung der Monge-Ampere-Gleichung, Teil I. Seminar zur Approximationstheorie am Institut fur Geometrieund Praktische Mathematik der Rheinisch-Westfalischen Technischen HochschuleAachen, Wintersemester 2009/2010.

[17] Hestenes, Magnus R.: Multiplier and Gradient Methods. Journal of OptimizationTheory and Applications, 4(5):303–320, 1969.

[18] Hestenes, Magnus R. und Eduard Stiefel: Methods of Conjugate Gradients forSolving Linear Systems. Journal of Research of the National Bureau of Standards,49(6):409–436, 1952.

[19] Krantz, Steven G. und Harold R. Parks: The geometry of domains in space.Birkhauser Advanced Texts / Basler Lehrbucher. Birkhauser, Boston, 1999.

[20] Powell, Michael J. D.: A method for nonlinear constraints in minimizationproblems. In: Fletcher, Roger (Herausgeber): Optimization, Kapitel 19, Seiten283–298. Academic Press, New York, 1969.

[21] von der Mosel, Heiko: Variationsrechnung I. Skript zur Vorlesung am Institutfur Mathematik der Rheinisch-Westfalischen Technischen Hochschule Aachen, Win-tersemester 2009/2010. http://www.instmath.rwth-aachen.de/~heiko/lehre/

variWS09-10/VarI.pdf.

32