WIRTSCHAFTSINFORMATIK Westfälische Wilhelms-Universität Münster WIRTSCHAFTS INFORMATIK Lösung...

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


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


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



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



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


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



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


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


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


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


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


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


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


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


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


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


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


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


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


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


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


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



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



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



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


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


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


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


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


Herzlichen Dank

für Ihre Aufmerksamkeit!

Noch Fragen?

WIRTSCHAFTSINFORMATIK Westfälische Wilhelms-Universität Münster WIRTSCHAFTS INFORMATIK Lösung...

Documents

Transcript of WIRTSCHAFTSINFORMATIK Westfälische Wilhelms-Universität Münster WIRTSCHAFTS INFORMATIK Lösung...