Parallele Lösungsverfahren

Die von Skeel(1989) vorgeschlagenen Techniken lassen sich auch auf die verallgemeinerten WR-Ver-fahren übertragen (Nevanlinna 1990; Lubich 1992). Darüber hinaus ist eine Verbesserung der Konver-genz durch Multisplitting, d. h. durch eine durch mehrfache Unterteilung erzeugte Überlappung, möglich (Jeltsch u. Pohl 1991;Pohl 1992;Frommer u. Pohl 1995).

Der besondere Vorteil der WR-Verfahren besteht darin, daß eine Entkoppelung von Teilsystemen mög-lich ist, so daß Teilsysteme mit unterschiedmög-lichen Eigenschaften unabhängig voneinander behandelt wer-den können. Dadurch ist es möglich, für jedes Teilsystem ein passendes Integrationsverfahren auszuwäh-len, und die Teilsysteme können mit unterschiedlichen Schrittweiten integriert werden. Neben der An-passung der zeitlichen Diskretisierung ermöglichen Multigrid-Waveform-Relaxationsverfahren auch eine Anpassung des räumlichen Gitters bei der Lösung semi-diskretisierter partieller Differentialgleichungen (abgekürzt: PDEs für englisch:partial differential equations) (siehe Lubich u. Ostermann 1987; Vandewal-le 1993;van Lent u. Vandewalle 2002). Als Basisverfahren können z. B. lineare Mehrschrittverfahren oder Runge-Kutta-Verfahren eingesetzt werden. Für einen ausführlichen Überblick sieheBurrage(1995).

1.2 Parallele Lösungsverfahren

Bei einer parallelen Implementierung von Lösungsverfahren für gewöhnliche Differentialgleichungssyste-me lassen sich verschiedene Parallelitätsarten unterscheiden. NachGear (1987, 1988) unterscheidet man aus Sicht der Lösungsverfahren folgende Parallelitätsarten:

I Parallelität bezüglich des Systems: Verschiedene Komponenten des gewöhnlichen Differentialglei-chungssystems werden parallel behandelt.

I Parallelität bezüglich der Methode: Unabhängige Berechnungen innerhalb eines Zeitschrittes wer-den parallel ausgeführt.

I Parallelität bezüglich der Zeitschritte: Mehrere Zeitschritte werden parallel ausgeführt.

Aus Sicht der Algorithmen, die diese Lösungsverfahren realisieren, unterscheidet man:

I Instruktionsparallelität: Mehrerer Maschineninstruktionen werden durch unterschiedliche Funkti-onseinheiten eines einzelnen Prozessors gleichzeitig ausgeführt.

I Datenparallelität: Durch einen oder mehrere Prozessoren werden gleiche Instruktionen auf un-terschiedliche Daten parallel angewendet.

I Taskparallelität: Unabhängige Instruktionsfolgen (Tasks) werden durch einen oder mehrere Prozessoren parallel ausgeführt, wobei Tasks selbst ebenfalls parallelisierbar sein können.

Die parallele Implementierung von Verfahren zur Lösung von gewöhnlichen Differentialgleichungen wurde lange Zeit vernachlässigt, da die bekannten Verfahren als inhärent sequentiell galten. In der jünge-ren Vergangenheit wurden jedoch von Mathematikern mehrere Verfahjünge-ren zur Lösung gewöhnlicher Diffe-rentialgleichungssysteme vorgeschlagen, die gut für eine parallele Implementierung geeignet sind, da sie ein größeres Potential an Parallelität besitzen als herkömmliche Methoden.

Die verschiedenen Verfahren wurden meist von Mathematikern mit dem Ziel vorgeschlagen, neue merische Verfahren zu entwerfen, die voneinander unabhängige Berechnungen aufweisen und deren nu-merische Eigenschaften mit den nunu-merischen Eigenschaften der bekannten Standardverfahren, wie sie z. B.

vonHairer u. a. (2000) undHairer u. Wanner (2002) beschrieben werden, vergleichbar sind. Die entspre-chenden Veröffentlichungen enthalten meist eine genaue Analyse der numerischen Eigenschaften (Konver-genz, Stabilität) der Verfahren; üblicherweise wird aber ein paralleler Algorithmus, der eine genaue Ab-bildung der durchzuführenden Berechnungsschritte auf parallel arbeitende Berechnungseinheiten angibt, nur sehr vage beschrieben oder fehlt ganz. Oft wird nur ein zusätzlicher Grad an Parallelität identifiziert, Wechselwirkungen mit anderen Parallelitätsquellen werden aber nicht diskutiert.

Im folgenden soll ein kurzer Überblick über die wichtigsten parallelen Lösungsverfahren gegeben wer-den. Eine detaillierte Einführung in dieses Thema gibtBurrage(1995).

1.2.1 Runge-Kutta-Verfahren

Die Möglichkeiten zur parallelen Ausführung von Runge-Kutta-Verfahren werden durch die Berechnungs-vorschrift (1.9) bestimmt. Da, wie bei allen Ein- und Mehrschrittverfahren, Abhängigkeiten zwischen auf-einanderfolgenden Zeitschritten bestehen, ist eine Parallelisierung bezüglich der Zeitschritte nur in Aus-nahmenfällen durch Überlappung von wenigen direkt benachbarten Zeitschritten möglich. Eine Anpas-sung der Schrittweite nach jedem Zeitschritt erschwert dies zusätzlich. Eine Parallelisierung bezüglich des Differentialgleichungssystems und bezüglich der Methode ist jedoch prinzipiell möglich.

Explizite Runge-Kutta-Verfahren

Bei expliziten Runge-Kutta-Verfahren hängen die Stufen (1.9a) sequentiell voneinander ab, so daß eine par-allele Berechnung der Stufen im allgemeinen Fall nicht möglich ist. Die Möglichkeit einer parpar-allelen Aus-führung besteht jedoch für spezielle Verfahren, deren Verfahrensmatrix besondere Eigenschaften aufweist.

Einer der ersten Artikel über parallele Lösungsverfahren für gewöhnliche Differentialgleichungssyste-me vonMiranker u. Liniger(1967) enthält im Anhang die Beschreibung eines parallelen expliziten Runge-Kutta-Verfahrens der Ordnung 3. Hier wird eine parallele Berechnung der Stufen dadurch ermöglicht, daß das Verfahren eingebettete Lösungen der Ordnung 1 und 2 enthält und zusätzlich eine konstante Schrittwei-te vorausgesetzt wird. Die eingebetSchrittwei-teSchrittwei-ten Lösungen und die zugehörigen Stufen werden auf unSchrittwei-terschiedli- unterschiedli-chen Prozessoren berechnet, wobei die Berechnung einer Lösung der Ordnungp der Berechnung der Lö-sung der Ordnungp+1 um einen Zeitschritt vorausläuft. Zur Berechnung der Lösung der Ordnungp+1 muß deshalb nur noch die Stufep+1 berechnet werden, da die Stufeni=1, . . . ,pbereits auf Prozessorizur Bestimmung der Lösung der Ordnungiberechnet wurden.Khalaf u. Hutchinson(1992) zeigen jedoch, daß die Konstruktion derartiger Verfahren für eine Ordnung größer als 3 nicht durchführbar ist. Dies bedeutet insbesondere, daß diese Art der Parallelisierung auf maximal drei Prozessoren beschränkt ist. Darüber hin-aus besitzt das vonMiranker u. Liniger(1967) vorgeschlagene Verfahren schlechte Stabilitätseigenschaften, wie die Autoren selbst einräumen.

Eine detaillierte Analyse des Parallelitätspotentials von Runge-Kutta-Verfahren bezüglich der parallelen Berechnung der Stufen wurde vonIserles u. Nørsett(1990) undJackson u. Nørsett(1995) durchgeführt. Sie fanden heraus, daß eine Parallelisierung expliziter Runge-Kutta-Verfahren nur begrenzt möglich ist:

Satz 1.1. Für ein explizites Runge-Kutta-Verfahren mitσsequentiellen Stufen unterliegt die Ordnung p der Bedin-gung

p≤σ für jede beliebige Anzahl verfügbarer Prozessoren.

Unter der Anzahl der sequentiellen Stufenσversteht man dabei die Länge des längsten Pfades im Pro-duktionsgraphendes Verfahrens, der die Abhängigkeiten zwischen den Stufen beschreibt. Für einen Beweis siehe die Originalarbeiten vonIserles u. Nørsett(1990) undJackson u. Nørsett(1995) oderHairer u. a.(2000).

Verfahren mitp=σbezeichnet man alsp-optimal. Für die Parallelisierung einesp-optimalens-stufigen ex-pliziten Runge-Kutta-Verfahrens bedeutet Satz1.1, daß ein solches Verfahren maximals−pStufen besitzt, die für eine parallele Ausführung in Betracht kommen. Die Verteilung der parallel ausführbaren Stufen auf verschiedene Prozessoren kann mittels der vonIserles u. Nørsett(1990) vorgeschlagenenDigraphmethode durchgeführt werden. Diesen Ansatz verfolgen auchPetcu(1999) undPetcu u. Dragan(2000), deren Soft-wareEpODE(für englisch:expert system for ordinary differential equations) nach Vorgabe eines Runge-Kutta-Verfahrens eine automatische Verteilung der Stufenvektorberechnungen durchführt.

Eine Untersuchung des Parallelitätspotentials eingebetteter Runge-Kutta-Verfahren bezüglich des Dif-ferentialgleichungssystems wurde vonRauber u. Rünger(1999b) durchgeführt. Dazu wurden Implemen-tierungen für gemeinsamen und für verteilten Speicher betrachtet und eine theoretische Abschätzung der Kommunikations- und Berechnungskosten durchgeführt. Laufzeitexperimente wurden auf zwei Parallel-rechnern mit verteiltem Speicher, einer Intel Paragon und einer IBM SP2, durchgeführt. Gute Speedup-Werte konnten dabei jedoch nur beobachtet werden, wenn die Zeit für die Auswertung der rechten Sei-te des DifferentialgleichungssysSei-tems ausreichend hoch war, um die KosSei-ten der aufgrund des verSei-teilSei-ten Adreßraumes auf beiden Parallelrechnern sehr teuren Kommunikationsoperationen auszugleichen.

1.2 Parallele Lösungsverfahren 13 Implizite Runge-Kutta-Verfahren

Für implizite Runge-Kutta-Verfahren ergibt die Berechnungsvorschrift (1.9a) eine implizite Beziehung der Stufenvektoren, so daß für ein allgemeiness-stufiges Verfahren und ein Differentialgleichungssystem der Dimensionnin jedem Zeitschritt die Lösung eines nichtlinearen Gleichungssystems der Dimensions·n erforderlich ist. Im Falle einer vollbesetzten Verfahrensmatrix kann Parallelität daher im wesentlichen nur bei der Lösung des nichtlinearen Gleichungssystems ausgenutzt werden.

Erfüllt jedoch die Struktur der Verfahrensmatrix Abestimmte Voraussetzungen, ist eine parallele Be-rechnung der Stufen prinzipiell möglich. Einfachstes Beispiel sind die sogenannten streng diagonal-impli-ziten Runge-Kutta-Verfahren, die durch ali = 0 fürl 6= i gekennzeichnet sind. Da in diesem Fall jeder Stufenvektor nur von sich selbst abhängig ist, können alle Stufenvektoren gleichzeitig berechnet werden.

Das bedeutet insbesondere auch, daß anstelle des nichtlinearen Gleichungssystems der Dimensions·njetzt sunabhängige nichtlineare Gleichungssysteme der Dimensionngelöst werden müssen, was gleichzeitig zu einer Reduzierung des Berechnungsaufwandes führt.Jackson u. Nørsett(1995) zeigen jedoch, daß derartige Verfahren für allgemeine Probleme nur eine maximale Ordnung von 2 erreichen können.

Bessere Ordnungseigenschaften bieten dagegen die parallelen mehrfach-impliziten Runge-Kutta-Ver-fahren (PMIRK-VerRunge-Kutta-Ver-fahren), die ebenfalls eine Entkoppelung insunabhängige nichtlineare Gleichungssy-steme der Dimensionnermöglichen (vgl.Burrage 1995). Dies wird durch die Voraussetzung möglich, daß die VerfahrensmatrixAdurch eine Ähnlichkeitstransformation in eine DiagonalmatrixDüberführt werden kann, d. h., es gilt:

A=T⁻¹DT , D=diag(λ₁, . . . ,λs) .

Vergleicht man jedoch die für eine parallele Implementierung eines PMIRK-Verfahrens erforderliche Zahl von Operationen mit der erforderlichen Anzahl von Operationen für ein effizientes sequentielles SIRK-Verfahren, so resultieren aus der Parallelität keine nennenswerten Vorteile.

Iterierte Runge-Kutta-Verfahren

Eine Verfahrensklasse mit vielversprechendem Parallelitätspotential sind diePrädiktor-Korrektor-Verfahren.

Prädiktor-Korrektor-Verfahren versuchen die Vorteile der einfachen Berechnungsstruktur expliziter Ver-fahren mit der verbesserten Stabilität und Genauigkeit impliziter VerVer-fahren zu kombinieren. Wie bereits in Abschnitt1.1.2im Zusammenhang mit linearen Mehrschrittverfahren angesprochen, geht man dazu von einem impliziten Basisverfahren aus und verwendet einen iterativen Ansatz zur Auflösung der nichtlinea-ren Gleichungssysteme. Als Startwert der Iteration verwendet man einen Näherungswert der Lösung, der mit Hilfe eines expliziten Verfahrens ermittelt wird. Dieses Verfahren wird daher alsPrädiktorverfahren be-zeichnet. Anschließend wird die Näherungslösung in mehreren Iterationsschritten mit Hilfe des impliziten Basisverfahrens „korrigiert“, weshalb das Basisverfahren auchKorrektorverfahrengenannt wird.

Handelt es sich bei dem Basisverfahren um ein Runge-Kutta-Verfahren, so spricht man von einem ite-rierten Runge-Kutta-Verfahren. Das einfachste und bekannteste Verfahren dieser Art, das daher auch oft als dasiterierte Runge-Kutta-Verfahren bezeichnet wird, verwendet den sogenannten „trivialen“ Prädiktor

Y⁽⁰⁾_l =η_κ für l=1, . . . ,s (1.26)

als Startwert für eine Fixpunktiteration mit einer festen Anzahlmvon Iterationen:

Y^(j)_l =η_κ+hκ

∑

a_lif(tκ+c_lhκ,Y^(j−1)_i ) für l=1, . . . ,s und j=1, . . . ,m . (1.27) Die neue Approximation für die Lösungsfunktion an der Stelletκ+1erhält man anschließend aus

η_κ+1=η_κ+h_κ

∑

b_lf(t_κ+c_lh_κ,Y^(m)_l ) . (1.28) Dieses iterierte Verfahren geht zurück auf die Arbeiten vonNørsett u. Simonsen(1989) undvan der Houwen u. Sommeijer(1990a). Es läßt sich als explizites Verfahren mits(m+1)Stufen auffassen, dessen

Verfahrensmatrix gegeben ist durch:

0 0

c A 0

c 0 A 0

... ... . .. ...

c 0 . . . 0 A 0 0 . . . 0 0 b^T

(1.29)

Bezeichnet man die Ordnung des Basisverfahrens mitp, so beträgt die Ordnung eines solchen Verfahrens p^∗ = min{p,m+1}. Die Anzahl der Korrekturenmwird dementsprechend in der Regel auf m = p−1 gesetzt, um eine höchstmögliche Ordnung bei gleichzeitig möglichst geringer Iterationszahl zu erreichen.

Der Vorteil eines solchen Verfahrens besteht darin, daß innerhalb jedes Iterationsschrittesjdie Vektoren Y^(j)₁ , . . . ,Y^(j)s unabhängig voneinander sind und daher parallel berechnet werden können. Man erhält somit ein Lösungsverfahren der Ordnungp^∗= pmit lediglichpsequentiellen Stufen. Dies bedeutet insbesonde-re, daß das Verfahren optimal in bezug auf Satz1.1ist. Im Vergleich dazu wächst die für die Berechnung eines Zeitschrittes eines sequentiellen expliziten Verfahrens erforderliche Stufenzahl schneller als linear im Verhältnis zur Ordnung.

Eine ausführliche Untersuchung des Parallelitätspotentials iterierter Runge-Kutta-Verfahren, die den trivialen Prädiktor sowie eine feste Anzahl von Iterationen eines impliziten Runge-Kutta-Korrektors ver-wenden, wurde vonRauber u. Rünger(1996,1999b) durchgeführt. Dazu wurden verschiedene Implemen-tierungen für Parallelrechner mit verteiltem Speicher erstellt und untersucht, die sowohl Parallelität be-züglich der Methode als auch Parallelität bebe-züglich des Differentialgleichungssystems ausnutzen, jedoch jeweils unterschiedliche Datenverteilungen verwenden. Für diese Implementierungen wurde eine theoreti-sche Analyse der Berechnungs- und der Kommunikationszeit durchgeführt, die eine sehr genaue Vorhersa-ge der Laufzeit ermöglicht. Eine Untersuchung des realen Laufzeitverhaltens erfolgte auf Parallelrechnern vom Typ Intel iPSC/860, Intel Paragon und IBM SP2. Zufriedenstellende Speedup-Werte konnten jedoch nur gemessen werden, wenn die Funktionsauswertungszeit des zu lösenden Differentialgleichungssystems sehr hoch war.

Neben den eben besprochenen Verfahren existieren eine Reihe weiterer iterierter RK-Verfahren. (Paral-lele) diagonal-implizite iterierte RK-Verfahren (DIIRK- bzw. PDIRK-Verfahren) entstehen durch eine Varia-tion des Korrektorschrittes (1.27), die zu einem impliziten Verfahren führt, das auch zur Lösung von steifen Differentialgleichungen geeignet ist, das aber ebenfalls einen höheren Grad an Parallelität als andere im-plizite Verfahren besitzt (van der Houwen u. a. 1992). Da diese Veränderung von (1.27) dazu führt, daß in jedem Iterationsschrittsimplizite, in der Regel nichtlineare Gleichungssysteme gelöst werden müssen, er-gibt sich jedoch eine wesentlich höhere Berechnungszeit. Eine detaillierte Analyse des Laufzeitverhaltens wurde vonRauber u. Rünger(1999a,c) durchgeführt.

Zur Verbesserung der Konvergenzeigenschaften der DIIRK-Verfahren wurden weiterhin parallele trian-gular-implizite iterierte Verfahren (PTIRK-Verfahren) vorgeschlagen (van der Houwen u. de Swart 1997).

Diese Verfahren besitzen einen etwas geringeren Grad an Parallelität und erfordern einen höheren Berech-nungsaufwand pro Iterationsschritt, sie benötigen jedoch deutlich weniger Iterationsschritte um eine vorge-gebene Genauigkeit zu erreichen und ermöglichen dadurch gleichzeitig eine Überlappung der Berechnung mehrerer Zeitschritte. Die Idee der PTIRK-Verfahren wurde vonBurrage u. Suhartanto(1997a) für nicht-steife und ebenfalls vonBurrage u. Suhartanto (1997b) für steife Probleme auf Mehrschritt-Runge-Kutta-Verfahren des Radau-Typs verallgemeinert.

1.2.2 Extrapolationsverfahren

Wie in Abschnitt1.1.1beschrieben, berechnen Extrapolationsverfahren für die gleiche Argumentstelle ver-schiedene Approximationen der gesuchten Lösungsfunktion, indem sie – ausgehend von der errechneten Lösung für die vorangehende Argumentstelle – verschiedene Schrittweiten verwenden. Durch Kombina-tion dieser ApproximaKombina-tionen entsteht eine ApproximaKombina-tion höherer Konvergenzordnung. Wegen der Un-abhängigkeit der Berechnungen der Approximationen für eine Argumentstelle mit verschiedenen Schritt-weiten sind Extrapolationsverfahren mehrfach für eine parallele Ausführung diskutiert worden, u. a. von