Ansätze zur Reduzierung der Anzahl der ausgeführten Operationen

Durch einige einfache Veränderungen kann der Programmcode in Abb.2.4so angepaßt werden, daß die Anzahl der Operationen deutlich vermindert wird.

Abspaltung von Schleifeniterationen

Eine wichtige Technik in diesem Zusammenhang ist das Abspalten von Schleifeniterationen. Das bedeutet, daß man die ersten oder die letzten Iterationen einer Schleife aufrollt und vom Schleifenkörper trennt. Die Grenzen der Schleife müssen danach entsprechend angepaßt werden.

Beispiel2.1. Von der Schleife for(i:=1;i≤N;i++) work(i);

werden die ersten drei und die letzten zwei Iterationen abgespalten:

work(1);

work(2);

work(3);

for(i:=4;i≤N−2;i++) work(i);

work(N−1);

work(N);

Mit Hilfe dieser Technik ist es möglich, bei der Berechnung von Vektorsummen die Initialisierung der Summationsvektoren mit 0 zu vermeiden. Stattdessen wird die erste Iteration der Summationsschleife ab-gespalten, so daß das erste Glied der Summe zur Initialisierung verwendet werden kann. Auf diese Weise spart man für jede Vektorkomponente eine Zuweisungsoperation (z. B. das Initialisieren eines Registers mit 0, was in der Regel einen Prozessortakt benötigt). Dies kann bereits zu einem meßbaren absoluten Laufzeit-gewinn führen, wenn man bedenkt, daß die Dimension des Differentialgleichungssystems unter Umstän-den sehr groß ein kann. Bei der in Abb.2.4dargestellten Implementierung ist der erzielbare Laufzeitgewinn jedoch noch wesentlich größer. Da hier die Initialisierung mit Hilfe einer separaten Schleife geschieht, die über die Vektorkomponenten iteriert, ist mit jeder Zuweisung auch ein Schreibzugriff auf das Speichersy-stem verbunden. Insbesondere dann, wenn die Dimension des DifferentialgleichungssySpeichersy-stems sehr groß ist, so daß bei der Iteration über einen Vektor der Cacheinhalt mehrfach überschrieben wird, könnten durch die Speicherzugriffe große Wartezeiten entstehen. Kann man diese Speicherzugriffe also vermeiden, führt dies oft zu einer signifikanten Verbesserung der Laufzeit. Hinzu kommt im Fall der Implementierung aus Abb.2.4noch, daß die Einsparung der Sprung- und Vergleichsoperationen, die zur Realisierung der Schlei-fe nötig waren, zu einer weiteren Laufzeitverbesserung führt.

Multiplikationen mit der Schrittweite

An vier Stellen des Berechnungskernels in Abb. 2.4finden Multiplikationen mit der Schrittweiteh statt:

Zur Berechnung der Argumentvektorenwl (Zeile 7), zur Berechnung von∆η(Zeile 14), zur Berechnung des Fehlervektorse(Zeile 19) und zur Berechnung der Skalierungsfaktoren (Zeile 21). Der zweite Faktor ist dabei immer entweder ein einzelner Stufenvektor oder eine zuvor berechnete gewichtete Summe von Stu-fenvektoren. Durch Ausnutzung der Distributivität läßt sich die Anzahl der Multiplikationen reduzieren, indem man die modifizierten Stufenvektoren

˜

v_l=hκv_l für l=1, . . . ,s (2.4)

2.3 Implementierung eingebetteter Runge-Kutta-Verfahren 39

1: // Berechnev₁, . . . ,vs nach(2.2a)

2: for(l:=1;l≤s;l++)

3: {

4: for(j:=1;j≤n;j++)w_l[j]:=0;

5: for(i:=1;i<l;i++)

6: for(j:=1;j≤n;j++)wl[j]+=alivi[j];

7: for(j:=1;j≤n;j++)wl[j]:=η[j] +hwl[j];

8: for(j:=1;j≤n;j++)vl[j]:= fj(t+clh,wl);

9: }

10: // Berechne∆η=η_κ+1−η_κ, vgl.(2.2b)

11: for(j:=1;j≤n;j++)∆η[j]:=0;

12: for(l:=1;l≤s;l++)

13: for(j:=1;j≤n;j++)∆η[j]+=b_lv_l[j];

14: for(j:=1;j≤n;j++)∆η[j]:=_h∆η[j];

15: // Berechne Fehlervektore=ηˆκ+1−ηκ+1, vgl.(2.2b)und(2.2c)

16: for(j:=1;j≤n;j++)e[j]:=0;

17: for(l:=1;l≤s;l++)

18: for(j:=1;j≤n;j++)e[j]+= ˜b_lv_l[j]; //b˜_l =b^ˆ_l−b_l

19: for(j:=1;j≤n;j++)e[j]:=he[j];

20: // Berechne Skalierungsfaktoren

21: for(j:=1;j≤n;j++)s[j]:=|η[j]|+|h·v1[j]|;

Abb. 2.4:Berechnungskernel eines eingebetteten Runge-Kutta-Verfahrens, realisiert als direkte vektororientierte Um-setzung der Berechnungsvorschrift (2.2).

einführt. Damit ergibt sich die modifizierte Berechnungsvorschrift

˜

v_l=hκf tκ+c_lhκ,η_κ+

l−1

∑

i=1

a_liv˜i

!

für l=1, . . . ,s , (2.5a)

η_κ+1=η_κ+

∑

s l=1

blv˜l , (2.5b)

ˆ

η_κ+1=ηˆ_κ+

∑

s l=1

bˆlv˜l , (2.5c)

mit den Skalierungsfaktoren

s_j=η_κ,j +v˜_1,j

für j=1, . . . ,n . (2.6)

Die Multiplikation mitherfolgt jetzt nur noch an einer Stelle, da direkt nach der Funktionsauswertung deren Ergebnis mit der Schrittweite multipliziert wird. Dadurch lassen sich nicht nur die eigentlichen Mul-tiplikationsoperationen einsparen, sondern auch einige Schleifen, die nur zum Zweck der Multiplikation mit h implementiert wurden. Somit ist auch eine Verbesserung des Lokalitätsverhaltens möglich. Da die Reihenfolge der mathematischen Operationen verändert wurde, kann sich dieses Vorgehen allerdings auch auf den Rundungsfehler auswirken.

Reduzierung der Dimension von Feldzugriffen

Mehrdimensionale Feldzugriffe sind oft mit höheren Kosten verbunden als Speicherzugriffe über einen skalaren Zeiger, da die Adresse des Feldzugriffs erst berechnet werden muß. Wird das mehrdimensionale Feld als Block gespeichert, wie dies z. B. in FORTRAN der Fall ist, muß anhand der Indizes und der Größe des Typs der Datenelemente ein Offset ermittelt und zur Anfangsadresse des Blocks hinzuaddiert werden, um die Adresse des gewünschten Feldelements zu bestimmen. In C werden mehrdimensionale Feldzugriffe in der Regel über Zeigerfelder, d. h. durch indirekte Adressierung realisiert. Das bedeutet, daß für jede Dimension ein eindimensionaler Feldzugriff erfolgt. Es muß also für jede Dimension durch Multiplikation von Index und Größe des Datentyps ein Offset bestimmt und zur Anfangsadresse des betreffenden Feldes hinzuaddiert werden, damit anschließend durch eine Speicherzugriff die Anfangsadresse des Feldes der nächstniedrigeren Dimension geladen werden kann. Es besteht aber auch in C die Möglichkeit, z. B. mit Hilfe von Präprozessormakros, Feldzugriffe ähnlich wie in FORTRAN zu realisieren.

In jedem Fall ist ein mehrdimensionaler Feldzugriff in beiden Programmiersprachen relativ teuer, da entweder eine Reihe arithmetischer Operationen zur Bestimmung des Offsets oder mehrere Speicherzu-griffe erforderlich sind. Dies kann sich insbesondere dann negativ auf die Laufzeit auswirken, wenn mehr-dimensionale Feldzugriffe innerhalb von Schleifen verwendet werden. Man kann daher in diesem Fall oft eine Laufzeitverbesserung erreichen, indem man die Dimension der Feldzugriffe reduziert. Dies erreicht man, indem man die Anfangsadresse des Feldes der Dimension, über die iteriert wird, in einer Zeigerva-riablen speichert.

Beispiel2.2. Die Reduktion der Dimension des mehrdimensionalen Zugriffs auf das FeldAim Programm S:=0;

for(i:=_1;i≤N;i++) for(j:=1;j≤N;j++)

for(j:=1;j≤N;j++) S:=S+A[i][j][k]; führt zu dem modifizierten Programm

2.3 Implementierung eingebetteter Runge-Kutta-Verfahren 41 S:=0;

for(i:=1;i≤N;i++) {

AI:= A[i];

for(j:=1;j≤N;j++) {

AI J:= AI[j];

for(j:=_1;j≤ N;j++) S:=S+AI J[k]; }

}

Wurde das Feld A als Block der Größe N³ allokiert, könnte man das Beispielprogramm auch wie folgt realisieren:

S:=0;

P:=A[0][0];

for(i:=1;i≤N³;i++)S:=S+P[i];

Der erste im Beispiel gezeigte Transformationsschritt kann von vielen modernen Compilern, die die Optimierungstechnik des Verschiebens von schleifeninvariantem Code beherrschen, automatisch durch-geführt werden. Die zweite vorgeschlagene Transformation kann ein Compiler jedoch nicht automatisch durchführen, da er davon ausgehen muß, daß die Elemente des dreidimensionalen FeldesAinN² separa-ten Feldern der GrößeNgespeichert werden.

Verschmelzen von Schleifen

Da für die Realisierung einer Schleife pro Iteration eine Sprung- und eine Vergleichsoperation erforderlich sind, kann man die Anzahl der ausgeführten Operationen reduzieren, indem man Schleifen mit gleichen Grenzen zu einer einzelnen Schleife verschmilzt. Dies hat allerdings auch zur Folge, daß sich das Loka-litätsverhalten verändert. Denn durch das Verschmelzen der Schleifen ensteht eine neue Schleife, deren Schleifenkörper im Vergleich zu den einzelnen Schleifen mehr Instruktionen enthält und in der Regel auch Speicherzugriffe auf eine größere Menge von Daten durchführt.

Im Dokument Effiziente Implementierung eingebetteter Runge-Kutta-Verfahren durch Ausnutzung der Speicherzugriffslokalität (Seite 54-57)