Date post: | 30-Dec-2015 |
Category: |
Documents |
Upload: | aretha-rios |
View: | 28 times |
Download: | 0 times |
WIR
TS
CH
AF
TS
INF
OR
MA
TIK
WestfälischeWilhelms-Universität Münster
WIRTSCHAFTSINFORMATIK
Lösung nichtlinearer Lösung nichtlinearer GleichungssystemeGleichungssysteme
Präsentation in Rahmen des Seminars
„Parallele Programmierung“
SS 2003
Referent: Roman Achziger
Datum: 11.06.2003
2
WIRTSCHAFTSINFORMATIK
Einleitung: Einleitung: MotivationMotivation
Anwendungsgebiete Gebiete der Technik, Wirtschafts- und Naturwissenschaften
Lösung nichtlinearen Abhängigkeiten Approximation durch lineare Modelle
Einsatz numerischen Algorithmen
Die Standardform nichtlinearer Gleichungssysteme:
= oder kurz F(x)
)...,,(.
.)...,,(
)...,,(
21
212
211
nm
n
n
xxxf
xxxf
xxxf
0.
.0
0
3
WIRTSCHAFTSINFORMATIK
Grundlagen der IterationsverfahrenGrundlagen der Iterationsverfahren
Fixpunktiteration: Betrachtung des skalaren Falls Transformation der vorliegende Gleichung f(x)=0 in eine
äquivalente Fixpunktform (x) = x
Definition: Ein Element x* R heißt Fixpunkt der Abbildung : R R falls
(x*)=x* gilt.
Wähle die Abbildung so, dass x* die Abbildung (x*) = x* und f(x*) = 0 erfüllt
Ausgehend vom Startwert x(0) gilt die Rekursion x(k+1) = (x(k)),
Abbildung verkleinert sukzessive den Abstand zwischen x(k+1) und x*
0x
4
WIRTSCHAFTSINFORMATIK
Grundlagen der IterationsverfahrenGrundlagen der Iterationsverfahren
Spinnweb-Diagramm für (x)= exp(-x2)
Die Folge oszilliert um den Fixpunkt
-0.2 0.2 0.4 0.6 0.8 1
0.2
0.4
0.6
0.8
1y=x
x(0)=0
(x) = exp(-x2)
5
WIRTSCHAFTSINFORMATIK
Grundlagen der IterationsverfahrenGrundlagen der Iterationsverfahren
Die Wahl einer gute Iterationsfunktion ist entscheidend
Definition:Eine auf einer Teilmenge A R definierte Abbildung : A A heißt Kontraktion, wenn eine Lipschitz-Konstante 0<L<1 existiert,
so dass für alle Paare x,y A gilt: | (x) - (y)| L |x-y|.
L = |´(Z)|, mit einem Zwischenwert Z zwischen x und y
Die Abbildung konvergiert, wenn dieses Maximum der ersten Ableitung kleiner als 1
AZmax
6
WIRTSCHAFTSINFORMATIK
Grundlagen der Iterationsverfahren: Grundlagen der Iterationsverfahren: BeispielBeispiel
Gesucht ist die Nullstelle von f(x) = x3+x-1,5
Einsetzen: f(0,8) = -0,188 und f(0,9) = 0,129
Da die Funktion stetig ist, existiert eine Nullstelle x* mit 0,8<x*<0,9
1) 1(x)= -x3+1,5 ; | 1´(x)|=3x2
| 1´(0,8)| <| 1´(x*)|<| 1´(0,9)|<=> 1,92 <| 1´(x*)|< 2,43
abstoßender Fixpunkt: x(0)=0,9; x(1)=0,77; x(2)=1,04; x(3)=0,37
2) 2(x)=(-x3+x+1,5)/2 ; | 2´(x)| =| (-3x2+1)/2|
| 2´(0,8) |<| 2´(x*)|<| 2´(0,9)| d.h. 0,46<| 2´(x*)|< 0,715
Die Folge konvergiert gegen x* für Startwerte aus [0,8; 0,9]
x(0)=0,9; x(1)=0,8355; x(2)=0,8761; x(3)=0,8518; x(4)=0,8667;
x(5)=0,8577; x(5)=0,8633
.
7
WIRTSCHAFTSINFORMATIK
Grundlagen der IterationsverfahrenGrundlagen der Iterationsverfahren
Banachsche Fixpunktsatz:
Es sei : A A eine Kontraktion über der abgeschlossenen Menge AR mit der Kontraktionskonstanten L< 1. Dann :
1) existiert ein eindeutiger Fixpunkt x* = (x*) in A;
2) konvergiert die durch x(k+1) = (x(k)) definierte Folge für jeden Startwert x(0) A gegen den Fixpunkt;
3) gelten die Abschätzungen:
| x-x |1
L |x-x|
1 |x*-x | (0)(1)
k1)-(k(k)(k)
LL
L
“ priori-„a “ posteriori-a „
8
WIRTSCHAFTSINFORMATIK
Grundlagen der IterationsverfahrenGrundlagen der Iterationsverfahren: : EinzugsbereichEinzugsbereich
Definition:
Die Menge E R heißt Einzugsbereich der Lösung x*, wenn x(k)= x* für alle x(0) E gilt.
globale Konvergenz: Ist der gesamte Definitionsbereich von Einzugsbereich
Iterationsfolge konvergiert für alle Startwerte
lokale Konvergenz: Iterationsfolge konvergiert nur auf einer Umgebung um x*
k
lim
9
WIRTSCHAFTSINFORMATIK
Grundlagen der IterationsverfahrenGrundlagen der Iterationsverfahren: : KonvergenzgeschwindigkeitKonvergenzgeschwindigkeit
Der Fehler (k) := x(k) -x* stellt auf einer kleinen Umgebung von x* das Konvergenzverhalten der Folge dar
Lineare Konvergenz: Es gibt eine Konstante q (0; 1), Konvergenzordnung p=1, so dass
für alle k=0, 1, 2, . . . gilt.
der Fehler (k+1) nimmt linear mit dem Faktor q ab
Quadratische Konvergenz:
Es gibt eine Konvergenzfaktor K > 0 und eine Konvergenzordnung
p=2, so dass für alle k=0, 1, 2, . . . gilt.
Fehler (k+1) ist proportional zum Quadrat von (k) mit dem Faktor K
1)()1( ** xxqxx kk
2)()1( ** xxKxx kk
10
WIRTSCHAFTSINFORMATIK
Grundlagen der IterationsverfahrenGrundlagen der Iterationsverfahren: : AbbruchkriterienAbbruchkriterien
Kompromiss zwischen der Genauigkeit und der Effizienz
Abbruch der Verfahren beim Erreichen einer geforderten Genauigkeit
Fehler-Kriterium: Prozess stoppt, wenn || x(k+1)-x(k)|| < abs
Residuums-Kriterium: Prozess stoppt, wenn ||F(x)|| < f
Abbruch der Verfahren bei Nicht-Konvergenz Kriterien: ||x(k)||>xmax oder || F(x(k))|| > fmax oder
||F(x(k))|| > ||F(x(k-1))||> ||F(x(k-2))||
sicherste Methode: Beschränkung der maximalen Anzahl der
Iterationsschritte Verhinderung einer unendlichen Folge
11
WIRTSCHAFTSINFORMATIK
Lösung nichtlinearer Gleichungen: Lösung nichtlinearer Gleichungen: BisektionsverfahrenBisektionsverfahren
Nullstellensatz für stetige Funktionen Sei eine stetige Funktion f : [a,b] R mit f(a)*f(b) < 0 gegeben.
Dann existiert ein x* (a,b), so dass f(x*) = 0 gilt.
Algorithmus: Startintervall: I0 =[x(k-1),x(k)] mit x(k-1) < x(k) und f(x(k-1))f(x(k))< 0
while (Abbruchkriterium nicht erfüllt)
{ (x(k+1) = x(k-1) + x(k)) /2
if (f(x(k-1))f(x(k+1)) 0 ) x(k)= x(k+1) else x(k-1)= x(k+1)
} return x(k)
Lineare Konvergenz mit dem Konvergenzfaktor ½
Global konvergent
12
WIRTSCHAFTSINFORMATIK
Lösung nichtlinearer Gleichungen: Lösung nichtlinearer Gleichungen: Bisektionsverfahren Bisektionsverfahren
Parallelen Realisierung
Unterteilung des Intervalls in in p+1 Teilintervalle
parallelen Berechnung der p Funktionswerte
Dasjenige mit einen Vorzeichenwechsel wird aus den p+1 Teilintervallen ausgewählt
wiederhole bis Abbruchkriterium erfüllt
Nach k-Schritt wird das Anfangsintervall auf das
( p+1)-k –fache reduziert speed-up von Sp= O(log2 p)
13
WIRTSCHAFTSINFORMATIK
Lösung nichtlinearer Gleichungen: Lösung nichtlinearer Gleichungen: Linearisiertes Newton-VerfahrenLinearisiertes Newton-Verfahren
Ersetzung der Funktion durch eine lineare Modellfunktion
Geometrisch: Der Schnittpunkt der Tangente mit der x-Achse gilt als Näherungswert für die Nullstelle x*
Iterationsvorschrift: x(k+1)= x(k) - , k =0,1,2,...,
Fixpunktiterationen: (x(k))= x(k)-
Newton-Iteration: x(k+1) = (x(k)), k = 0, 1, 2,....,
Eigenschaften: quadratische Konvergenz lokale Konvergenz selbstkorrigierend
)´(
)()(
)(
k
k
xf
xf
)´(
)(
xf
xf
14
WIRTSCHAFTSINFORMATIK
Lösung nichtlinearer Gleichungssysteme: Lösung nichtlinearer Gleichungssysteme: FixpunktiterationFixpunktiteration
Übertragung des skalaren Falls der Fixpunktiteration auf ein
n-dimensionales Problem
vorliegende Problem:
sind : Rn Rn und F: Rn Rn gegeben, finde x* Rn, so dass
(x*) = x*. Für einen Fixpunkt x* von gilt F(x*) = 0
Ausgehend von x(0) Rn konvergiert die Iterationsvorschrift
x(k+1) = ( x(k)), für k = 0, 1, 2..... gegen x*
15
WIRTSCHAFTSINFORMATIK
Lösung nichtlinearer Gleichungssysteme:Lösung nichtlinearer Gleichungssysteme:FixpunktiterationFixpunktiteration
parallele Realisierung vergleichbar zur der parallelen Realisierung eines Iterations-
verfahrens zur Lösung linearer Gleichungssysteme
Aufteilung der Komponenten i(x(k)) der Kontraktionsabbildung
auf n Prozessoren blockweise Aufteilung Fixpunktabbildung auf die Prozessoren
jeder einzelne Prozessor
berechnet n/p Spalten des nächsten Iterationsvektors
Verteilung des Iterationsvektor per Multi-Broadcastoperation
P1 P2 P3 P4
12345678
1 2 3 4 5 6 7 8
16
WIRTSCHAFTSINFORMATIK
Lösung nichtlinearer Gleichungssysteme: Lösung nichtlinearer Gleichungssysteme: Newton-VerfahrenNewton-Verfahren
Linearisierung der Funktion F: Rn Rn
Iterationsvorschrift: x(k+1)= x(k)-[DF(x(k))]-1 F(x(k)) für k=0, 1, 2, . . ..
Fixpunktproblem: (x) = x(k)-[DF(x(k))]-1 F(x(k))
Fixpunktiteration x(k+1) = (x(k))
DF(x)=
Eigenschaften: quadratische Konvergenz
lokale Konvergenz
)(....)(
.
.
.
.
.
.
)(....)(
1
1
1
1
xx
Fx
x
F
xx
Fx
x
F
n
nn
n
17
WIRTSCHAFTSINFORMATIK
Lösung nichtlinearer Gleichungssysteme: Lösung nichtlinearer Gleichungssysteme: Newton-VerfahrenNewton-Verfahren
Einführung eines Korrekturterms (k):1) x(k+1) = x(k) + (k)
2) DF(x(k))*(k) = - F(x(k))
Überführung des Problem des Lösens eines Systems nichtlinearer
Gleichungen in das Problem des Lösens eines linearen Gleichungs-
systems
Rechenaufwand pro Iterationsschritt: n2 +n Komponenten-Funktionsauswertung
2n3/3+O(n2) arithmetischen Operationen
n2 + O(n) Speicheroperationen
18
WIRTSCHAFTSINFORMATIK
Lösung nichtlinearer Gleichungssysteme: Lösung nichtlinearer Gleichungssysteme: Modifiziertes Newton-VerfahrenModifiziertes Newton-Verfahren
Zyklisches Aktualisieren der Jacobimatrix LU-Faktorisierung der Jakobimatrix nur alle n Iterationsschritte
Differenzenapproximation der Jacobimatrix Bestimmung einer Näherung der Jacobimatrix durch
n-dimensionale Differenzen
Berechne für jeden Eintrag aus DF(x(k)): (x(k))
Möglichkeit der Wahl für rj und ej: ||F(x(k) + rjej)- F(x(k))|| F(x(k))
j
kijj
ki
r
xFerxF )()( )()(
j
i
x
F
19
WIRTSCHAFTSINFORMATIK
Lösung nichtlinearer Gleichungssysteme:Lösung nichtlinearer Gleichungssysteme: Parallele Implementierung des Newton-VerfahrensParallele Implementierung des Newton-Verfahrens
Verwendung des Ergebnisvektors des k-ten Iterationsschrittes als Datum bei der Berechnung des nächsten Iterationsschrittes
die Iterationen müssen nacheinender durchgeführt werden
Teilaufgaben pro Schritt:
Auswertung der Funktion F
Berechnung der Jacobimatrix
Durchführung der Gauß-Elimination
Berechnung des nächsten Approximationsvektors
Kontrolle der Konvergenz des Verfahrens mit Bestimmung der Norm des Korrekturterms
20
WIRTSCHAFTSINFORMATIK
Lösung nichtlinearer Gleichungssysteme: Lösung nichtlinearer Gleichungssysteme: Parallele Implementierung des Newton-VerfahrensParallele Implementierung des Newton-Verfahrens
Parallelisierung der Teilaufgaben des Iterationsschrittes unter Beachtung der Datenabhängigkeiten
Datenabhängigkeit in der Newton-Iteration:
x(0)
funktion jacobian gauss norm
x(k)
w(k)DF(x(k))x(k) F(x(k))update
F(x(k))x(final)
w(k)
x(k+1)
x(k+1)
x(k)
||w||
21
WIRTSCHAFTSINFORMATIK
Lösung nichtlinearer Gleichungssysteme: Lösung nichtlinearer Gleichungssysteme: Parallele Implementierung des Newton-VerfahrensParallele Implementierung des Newton-Verfahrens
Bei der parallelen Ausführung auf eine gemeinsame Datenverteilung achten
Gauß-Elimination verursacht den meisten Rechenaufwand
Datenverteilung bzgl. Gauß-Elimination optimieren
Durchführung der Gauß-Elimination auf mehreren Prozessoren:
zeilenzyklische Berechnung
gesamtzyklische Berechnung
22
WIRTSCHAFTSINFORMATIK
Lösung nichtlinearer Gleichungssysteme:Lösung nichtlinearer Gleichungssysteme:Parallele zeilenzyklische ImplementierungParallele zeilenzyklische Implementierung
Verteilung der Jacobimatrix zyklisch auf die vorhandenen Prozessoren
Berechnung der Jacobimatrix in der entsprechenden
Verteilung
Approximation der Jacobimatrix durch:
(x(k))
Berechnung der Jacobimatrix: Ein Prozessor Pq speichert und berechnet alle Zeilen i mit q = ((i-
1) mod p) +1 maximal Zeilen
Vektor x(k) wird repliziert gehalten keine Kommunikation
Gesamtaufwand: T1(n,p)= (n+1) Teval(F)
pn /
pn /
j
kijj
ki
r
xFerxF )()( )()(
j
i
x
F
23
WIRTSCHAFTSINFORMATIK
Lösung nichtlinearer Gleichungssysteme:Lösung nichtlinearer Gleichungssysteme:Parallele zeilenzyklische ImplementierungParallele zeilenzyklische Implementierung
Durchführung der Gauß-Elimination Zeilenzyklische Ablage der Jacobimatrix
Zeilenzyklische Ablage der rechten Seite des zu lösenden Gleichungssystems
zeilenzyklischen Variante der Gauß-Elimination ohne
Umverteilung der Daten
Ein Prozessor speichert die
benötigten Zeilen der Jacobimatrix und
den entsprechenden Funktionswert ab
a11 a1n b1. . . .
an1 ann bn. . . .
24
WIRTSCHAFTSINFORMATIK
Lösung nichtlinearer Gleichungssysteme:Lösung nichtlinearer Gleichungssysteme:Parallele zeilenzyklische ImplementierungParallele zeilenzyklische Implementierung
Berechnung des neuen Iterationsvektors: Replizierte Ablage des Iterationsvektors x(k) der Newton-Iteration
die Prozessoren führen die Berechnung ohne
Kommunikation durch
Jeder Prozessor p berechnet: Werte der gespeicherten Funktion
Zugehörige Zeilen der Jacobimatrix
Zugehörige Werte des Vektors (k)
Zugehörige Werte des Vektors x(k+1)
Zugehörige Komponenten der Norm ||(k)||
25
WIRTSCHAFTSINFORMATIK
Lösung nichtlinearer Gleichungssysteme:Lösung nichtlinearer Gleichungssysteme:Parallele zeilenzyklische ImplementierungParallele zeilenzyklische Implementierung
double *newton_zyklisch (double *F( )) { double *x_neu, *x_alt, *temp, *y, f[MAX_SIZE], fg[MIN_SIZE] ; int i, j ; x_neu = x_buf1; x_old= x_buf2; initialize (x_neu); do { temp=x_neu ; x_neu= x_alt ; x_alt=temp ; for (i=me; i<n; i+=p) {
fg[ i ] = -F(i, x_alt);
for (j=0; i<n; i++) {
x_alt[ j ] += r[ j ]; Berechnung der Jacobimatrix f [ j ] = F(i, x_alt); x_alt[ j ] - = r[ j ];
DF[ i ][ j ] = ( fg[ i ] + f [ j ] ) /r [ j ]; }
}
y = gauss_zyklisch (DF, fg); zeilenzyklische Gauß-Elimination
for ( j=0 ; j<n; j++ ) Berechnung des neuen Iterationsvektors x_neu[ j ] = x_alt[ j ] + y[ j ]; }
while (norm(x_alt, x_neu ) > (1-L) /L ); Überprüfung des Abbruchkriteriums return x_neu ; }
26
WIRTSCHAFTSINFORMATIK
Lösung nichtlinearer Gleichungssysteme:Lösung nichtlinearer Gleichungssysteme:Parallele gesamtzyklische ImplementierungParallele gesamtzyklische Implementierung
Berechnung der Jacobimatrix: Blockzyklische schachbrettartige Datenverteilung der Matrix
Bei gleichen Anzahl der Prozessoren in jeder Dimension berechnet ein
Prozessor etwa Einträge
Der benötigte Iterationsvektor wird repliziert abgelegt
lokales Arbeiten mit der Funktion F
P1 P2 P1 P212345678
1 2 3 4 5 6 7 8
P3 P4 P3 P4
P1 P2 P1 P2
P3 P4 P3 P4p
n
p
n
p
n 2
27
WIRTSCHAFTSINFORMATIK
Lösung nichtlinearer Gleichungssysteme:Lösung nichtlinearer Gleichungssysteme:Parallele gesamtzyklische ImplementierungParallele gesamtzyklische Implementierung
Vermeidung von Kommunikation bei der Berechnung der Jacobimatrix zwischen den Prozessoren
Fi(x(k)) mindestens Mal errechnen, da sich an der
Berechnung einer Zeile der Jacobimatrix Prozessoren
beteiligen
Es ergibt sich der Gesamtaufwand von:
T2(n,p) = n/p( )Teval(F)
Durchführung der Gauß-Elimination: Allokierung der Jacobimatrix DF vergleichbar mit der
zeilenzyklischen Variante
Eine allokierte Zeile enthält nur noch +1 Einträge
pn
p
p
ypn /
28
WIRTSCHAFTSINFORMATIK
Lösung nichtlinearer Gleichungssysteme :Lösung nichtlinearer Gleichungssysteme :Vergleich der beiden VariantenVergleich der beiden Varianten
Zeilenzyklische Implementierung: geringere Anzahl der Funktionsauswertungen
Gesamtzyklische Implementierung: besseres Kommunikationsverhalten der Gauß-Elimination
Wahl abhängig davon, welcher der beiden Aspekte überwiegt
Auswertungsaufwand der Funktion F
Zeilenzyklische Implementierung
Kommunikationsverhalten der Zielmaschine
Gesamtzyklische Implementierung
29
WIRTSCHAFTSINFORMATIK
ZusammenfassungZusammenfassung
Vielzahl von unterschiedlichen Iterationsverfahren Beurteilung anhand des Effizienzvergleich
Anzahl der Iterationsschritte
Rechenaufwand pro Iterationsschritt Bisektionsverfahren:
für skalaren Fall
einfaches Verfahren mit globaler, lineare Konvergenz Newton-Verfahren:
für skalaren und auch mehrdimensionalen Fall
quadratische, lokale Konvergenz Parallele Realisierung:
Verteilung des Berechnungsaufwands auf mehrere Prozessoren
30
WIRTSCHAFTSINFORMATIK
Herzlichen Dank
für Ihre Aufmerksamkeit!
Noch Fragen?