UNIVERSITÄT PASSAU

Fakultät für Mathematik und Informatik

Automatische Methoden zur Parallelisierung im Polyedermodell

Diplomarbeit

Referent: Verfasser: Wieninger

Aufgabensteller: Prof. Christian Lengauer PhD. Lehrstuhl für Informatik Universität Passau

Betreuer: Martin Griebl

Passau, den 14. März 1995

Vorab möchte ich mich bei Herrn Prof. Christian Lengauer PhD. bedanken, der stets für alle Fragen offen war und zur Lösung etwaiger Probleme beitrug. Mein besonderer Dank gilt meinem Betreuer, Herrn Martin Griebl, für seine unermüdliche Hilfsbereitschaft und konstruktiven Vorschläge zu allen Bereichen dieser Arbeit. Auch Herrn Jean-François Collard bin ich für seine hilfreichen Erklärungen zu großem Dank verpflichtet. Außerdem danke ich allen Mitgliedern des LooPo-Teams: Nils Ellmenreich, Peter Faber, Harald Keimer, Klaus Köbler, Frank Schuler, Sabine Wetzel und Alex Wüst, für viele ergiebige Diskussionen und Vorschläge während der Implementierung der Methoden, die in dieser Arbeit behandelt werden. Auch allen anderen Angehörigen des Lehrstuhls sei für ihre stete Hilfsbereitschaft gedankt.

Zusammenfassung

Mit der Einführung paralleler Rechnerarchitekturen wurde auch eine neue Form von Compilern notwendig, die sequentielle Programme parallelisieren, um die erhöhte Rechenleistung dieser Systeme zu nutzen. Wesentliche Bestandteile dieser Parallelisierung sind die Berechnung einer neuen Ausführungsreihenfolge für die Abarbeitung der einzelnen Befehle und die Verteilung der einzelnen Befehlsausführungen auf ein Mehrprozessorsystem. Diese Arbeit behandelt implementierbare Methoden für beide Bereiche, die einerseits effizient arbeiten und andererseits geringe Anforderungen an die zu parallelisierenden Programme und deren Abhängigkeiten stellen. Für den ersten Bereich werden zwei Verfahren mit sehr unterschiedlichen Ansätzen vorgestellt. Für den Bereich der Verteilung wird zunächst eine Methode behandelt, die eine lineare Anordnung der Prozessoren voraussetzt. Dieses Verfahren wird dann auf mehrdimensionale Prozessoranordnungen verallgemeinert.

Inhaltsverzeichnis

1 Einführung 2

2 Notation 2

3 Grundlagen 3

3.1 Formale Darstellung von Eingabe-Programmen und Abhängig-

keiten . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 3 3.2 Grundidee der schleifenbasierten Parallelisierung . . . . . . . . 14

4 Berechnung von Schedules 19 4.1 Die Hyperebenen-Methode . . . . . . . . . . . . . . . . . . . . 19 4.2 Die Methode von Feautrier . . . . . . . . . . . . . . . . . . . . 30

4.2.1 Herleitung des Optimierungsproblems . . . . . . . . . . 30 4.2.2 Bewertung der Methode . . . . . . . . . . . . . . . . . 42 4.2.3 Anpassung an PIP . . . . . . . . . . . . . . . . . . . . 47

5 Berechnung von Allokationen 52

5.1 Eindimensionale Allokation . . . . . . . . . . . . . . . . . . . . 52 5.1.1 Vorstellung des Verfahrens . . . . . . . . . . . . . . . . 54 5.1.2 Einführung einer Gewichtsfunktion . . . . . . . . . . . 62 5.1.3 Bewertung der Methode . . . . . . . . . . . . . . . . . 63

5.2 Mehrdimensionale Allokation . . . . . . . . . . . . . . . . . . . 64 5.2.1 Vorstellung des Verfahrens . . . . . . . . . . . . . . . . 65 5.2.2 Bewertung der Methode . . . . . . . . . . . . . . . . . 68

6 Erzeugung des parallelen Programms 69

7 Schlußbemerkungen 72

1 Einführung

Mit der Einführung paralleler Rechnerarchitekturen wurde auch eine neue Form von Compilern notwendig, die sequentielle Programme parallelisieren, um die erhöhte Rechenleistung dieser Systeme zu nutzen. Der Prozeß dieser Parallelisierung besteht zunächst aus einer Abhängigkeitsanalyse, gefolgt von der Berechnung einer neuen Ausführungshalbordnung, die durch einen sogenannten Schedule repräsentiert wird. Für die Verteilung der Berechnungen auf ein Mehrprozessorsystem wird eine sogenannte Allokation konstruiert. Abschließend wird anhand von Schedule und Allokation das parallele Programm aufgebaut.

Es gibt mittlerweile eine große Anzahl von Verfahren für jede dieser Aufgaben. Die meisten dieser Verfahren besitzen jedoch ernste Einschränkungen bzgl. der zu parallelisierenden Programme und der darin auftretenden Abhängigkeiten.

Gegenstand der vorliegenden Arbeit ist es, für die Bereiche Schedule und Allokation Methoden vorzustellen, die einerseits effizient arbeiten und andererseits geringe Anforderungen an die zu parallelisierenden Programme stellen. Diese Methoden wurden im Rahmen des Projekts LooPo^{1} an der Universität Passau als Bestandteil dieser Arbeit implementiert.

Wir stellen nach der Einführung einiger Notationen in Abschnitt 3 die formalen Grundlagen zur Darstellung von Programmen und den darin auftretenden Abhängigkeiten vor. Die beiden folgenden Abschnitte bilden den Schwerpunkt dieser Arbeit. In Abschnitt 4 sollen zwei Verfahren zur Berechnung von Schedules erläutert werden [5][1]. In Abschnitt 5 gehen wir auf eine Methode zur Berechnung von Allokationen ein [3]. Im letzten Abschnitt wird schließlich die Konstruktion des parallelen Programms bei Vorgabe eines Schedules und einer Allokation angesprochen [3].

2 Notation

Zur Darstellung von Quantoren in formalen Ausdrücken benutzen wir die Notation von Dijkstra und Scholten [11]. Die Quantifizierung über einer Variablen x hat damit die allgemeine Form:

1WWW-Referenz http://brahms.fmi.uni-passau.de/cl/loopo/index.html

( Qx : R.x : P.x)

Q repräsentiert den Quantifikator, R stellt den Wertebereich von x dar, und P ist ein Ausdruck in Abhängigkeit von x. Für formale Herleitungen in Beweisen benutzen wir folgende Darstellungsweise:

formaler Ausdruck 1 op { Kommentar zur Erklärung der Gültigkeit von op ∈ { ⇒, ⇐, ⇔ } }

formaler Ausdruck 2

Zur Darstellung der lexikografischen Ordnung auf einem k-dimensionalen ℝ-Vektorraum benutzen wir den Operator ≺_{l} mit folgender Festlegung:

(∀ i⃗, j⃗ : i⃗, j⃗ ∈ ℝ^{k} : i⃗ ≺_{l} j⃗ ⇔
(∃ l : 1≤l ≤k : i_{l} < j_{l} ∧ (∀ m : 1≤m < l : i_{m} = j_{m})))

Wir verallgemeinern den Operator ≺_{l} durch den Operator ≺ mittels:

(∀ i⃗, j⃗ : i⃗, j⃗ ∈ ℝ^{k} : i⃗ ≺ j⃗ ⇔ (∃ l : 1≤l ≤k : i⃗ ≺_{l} j⃗))

Die Operatoren ≻_{l} und ≻ werden analog festgelegt.

Um formal auszudrücken, daß eine Zuweisung S in einem Eingabe-Pro-gramm syntaktisch vor einer Zuweisung T steht, schreiben wir S <_{L} T.

3 Grundlagen

3.1 Formale Darstellung von Eingabe-Programmen und Abhängigkeiten

Zunächst wollen wir in Anlehnung an die Darstellung in [1] und [5] festlegen, welche Programme als Eingabe für die vorzustellenden Verfahren akzeptiert werden sollen. Die erforderlichen Einschränkungen an diese Programme werden dabei so gering wie möglich gehalten und bei Bedarf später verschärft. Wir betrachten Fragmente von imperativen Programmen mit folgenden Eigenschaften.

• for-Schleifen und Zuweisungen sind die einzigen erlaubten Statements eines Programms und dürfen in beliebiger Kombination angeordnet werden. Die Schrittweite der for-Schleifen beträgt konstant +1. Es werden also auch sogenannte nicht perfekt verschachtelte Programme akzeptiert. Ein Programm heißt perfekt verschachtelt, wenn es die folgende Form besitzt:

for i_{1} := l_{1} to u_{1} do

for i_{2} := l_{2}(i_{1}) to u_{2}(i_{1}) do

for i_{3} := l_{3}(i_{1}, i_{2}) to u_{3}(i_{1}, i_{2}) do
⋮
for i_{n} := l_{n}(i_{1}, …, i_{n−1}) to u_{n}(i_{1}, …, i_{n−1}) do

{Schleifenkörper B} end ⋮

end end

end

Der Schleifenkörper B besteht ausschließlich aus Zuweisungen. Eine Iteration eines Schleifenkörpers bezeichnen wir mit B(i⃗), wobei i⃗ = (i_{1}, …, i_{n}) ist.

• Die einzigen Datenstrukturen sind Skalare und Arrays beliebiger Dimension.

• Schleifengrenzen und Indizes von Array-Zugriffen sind affin lineare Ausdrücke in den Indizes umgebender Schleifen und in Strukturparametern. Mit Strukturparametern sollen Konstanten bezeichnet werden, für die ein Bezeichner vergeben wurde. Sie beschreiben oft die Problemgröße (z. B. Invertierung einer n × n-Matrix mit n als Strukturparameter). Für die Schleifengrenzen gibt es zusätzlich die Möglichkeit, ein Maximum bzw. Minimum affiner Ausdrücke als untere bzw. obere Schranke festzulegen. Die Koeffizienten der Ausdrücke sind rational.

Ein Beispiel für ein Programm, das die Einschränkungen dieser Notation erfüllt, ist folgendes Fragment: c[0] := 1; {S_{1}} for i := 0 to min(3, m/4) do

for j := max(2 ∗ i, 2 ∗ n) to 20 do a[i, j − 2 ∗ i] := a[i, 0]; {S_{2}} b[i] := a[j, m]; {S_{3}}

end

c[i] := a[i, n] + c[i]; {S_{4}} end

Programme, die diesen Einschränkungen genügen, werden im folgenden als ein Schleifenprogramm^{2} L bezeichnet. Die Menge der Zuweisungen eines Schleifenprogramms bezeichnen wir mit Ω_{L}.

Ein Schleifenprogramm besitzt einen statischen Kontrollfluß. Dadurch ist es möglich, jede Ausführung eines Befehls zur Compilezeit durch Angabe der Zuweisung und des Indexvektors zu kennzeichnen.

Definition 3.1.1 (Indexvektor und Indexraum) Gegeben sei ein Schleifenprogramm L mit Zuweisungsmenge Ω_{L}. Für eine Zuweisung S ∈ Ω_{L} seien i_{1}, …, i_{m} die umgebenden m Schleifenindizes und l_{j} bzw. u_{j} die untere bzw. obere Grenze eines Index i_{j}, 1 ≤ j ≤ m. Dann bezeichnet

D_{S} = { (i_{1}, …, i_{m}) : (i_{1}, …, i_{m}) ∈ ℤ^{m} ∧ (∀ j : 1≤j ≤m : l_{j} ≤i_{j} ≤u_{j}) }

den Indexraum von S. Ein Vektor i⃗ ∈ D_{S} heißt Indexvektor von S. Mit

D_{Ω_{L}} = {D_{S} : S ∈ Ω_{L}}

wird die Menge aller Indexräume eines Schleifenprogramms bezeichnet. p_{S} gibt die Dimension eines Indexraums D_{S} an.

Geometrisch betrachtet stellt ein Indexraum einen Polyeder {x : Ax ≥ b} im ℤ^{m} dar, wobei das Ungleichungssystem Ax ≥ b wie in der Definition angegeben durch die Schleifengrenzen festgelegt ist.

2Der in der Theorie geläufige Begriff Schleifensatz bezeichnet Eingabe-Programme, in denen stets eine umgebende Schleife existiert. Da in unserer Festlegung aber auch Zuweisungen ohne umgebende Schleifen und auch aufeinanderfolgende Schleifen erlaubt sind, verwenden wir den Begriff Schleifenprogramm.

Definition 3.1.2 (Instanz einer Zuweisung) Sei L Schleifenprogramm, S ∈ Ω_{L} eine Zuweisung in L und i⃗ ∈ D_{S} ein Indexvektor. Dann heißt das Paar (S, i⃗) Instanz der Zuweisung S mit Indexvektor i⃗.

Damit sind sämtliche Konstrukte definiert, die wir zur formalen Darstellung des Eingabe-Programms benötigen. Als nächsten Punkt wollen wir nun eine Darstellung für auftretende Abhängigkeiten festlegen [5].

Definition 3.1.3 (Abhängigkeit) Gegeben seien ein Schleifenprogramm L, zwei nicht notwendig verschiedene Zuweisungen S, T ∈ Ω_{L} und zwei Indexvektoren i⃗ ∈ D_{S}, j⃗ ∈ D_{T}. Dann heißt die Instanz (S, i⃗) abhängig von der Instanz (T, j⃗), i. Z. (S, i⃗)δ(T, j⃗), falls gilt:

1. (S, i⃗) und (T, j⃗) greifen beide auf dieselbe Speicherzelle M zu, wobei mindestens einer dieser Zugriffe schreibend ist, und

2. (S, i⃗) wird in L vor (T, j⃗) ausgeführt.

Sind zwei Instanzen (S, i⃗) und (T, j⃗) abhängig so heißen auch die Zuweisungen S und T abhängig, i. Z. SδT.

Mit dieser Darstellung von Abhängigkeiten wird ein breites Spektrum abgedeckt. Es werden auch indirekte Abhängigkeiten erfaßt. Für direkte Abhängigkeiten fordert man zusätzlich, daß zwischen den beiden Zugriffen durch S und T kein anderer Zugriff auf die Speicherzelle erfolgt. Eine Spezialfall von Abhängigkeiten, der für Methoden zur Parallelisierung von Programmen im Polyedermodell häufig auftaucht, sind uniforme Abhängigkeiten:

Definition 3.1.4 (uniforme Abhängigkeit) Für ein Schleifenprogramm L mit der Zuweisungsmenge Ω_{L} und S, T ∈ Ω_{L} heißt eine direkte Abhängigkeit SδT uniform, falls D_{S} = D_{T} und ein lexikografisch positiver Vektor d⃗ ∈ ℤ^{p_{S}} existiert, sodaß die Differenz zwischen den Indexvektoren der abhängigen Instanzen von S und T konstant ist. Formal:

(∃ d⃗ : d⃗ ∈ ℤ^{p_{S}} \ {0} ∧ d⃗ ≻ 0 : (∀ i⃗, j⃗ : i⃗, j⃗ ∈ D_{S} ∧ (S, i⃗)δ(T, j⃗) : d⃗ = j⃗ − i⃗))

Der Vektor d⃗ heißt dann der Distanzvektor der Abhängigkeit.

Im Zusammenhang mit uniformen Abhängigkeiten in einem perfekt verschachtelten Schleifenprogramm spricht man auch vom Level einer Abhängigkeit:

Definition 3.1.5 (Level uniformer Abhängigkeiten) Gegeben sei ein perfekt verschachteltes Schleifenprogramm L = (L_{1}, …, L_{m}) mit m Schleifen und der Distanzvektor d⃗ = (d_{1}, …, d_{m}) einer uniformen Abhängigkeit in L. Da ein Distanzvektor stets ungleich dem Nullvektor ist, existiert ein von links gelesen erstes Element d_{l} ≠ 0. Formal:

(∃ l : l ∈ ℕ ∧ 1 ≤ l ≤ m : l = (min k : k ∈ ℕ ∧ 1≤k ≤m ∧ d_{k} ≠ 0 : k))

Dieses l bezeichnet man als den Level von d⃗, i. Z. level(d⃗). Es gilt also

level(d⃗) = l ⇔ d⃗ ≻_{l} 0

Definition 3.1.6 (Abhängigkeiten tragende Schleife) Gegeben sei ein perfekt verschachteltes Schleifenprogramm L = (L_{1}, …, L_{m}) mit m Schleifen und der Distanzvektor d⃗ = (d_{1}, …, d_{m}) einer uniformen Abhängigkeit in L. Ist l der Level von d⃗, so sagt man, daß die Schleife L_{l} diese Abhängigkeit trägt, i. Z. carries_{d}(L_{l}), ansonsten gilt ≠ carries_{d}(L_{l}). Formal:

(∀ l : 1≤l ≤m : carries_{d}(L_{l}) = ( level(d⃗) = l))

Anmerkung 3.1.7 Uniforme Abhängigkeiten sind Abhängigkeiten zwischen den Iterationen eines Schleifenprogramms. Um von den Zuweisungen auf die Iterationen zu abstrahieren benutzen wir die Notation B(i⃗)δB(j⃗) für zwei abhängige Iterationen B(i⃗) und B(j⃗).

In Definition 3.1.3 wird unter 2. implizit von einer Abbildung

Θ : (Ω_{L}, D_{S∈Ω_{L}}) → ℚ

ausgegangen, die jeder Instanz jeder Zuweisung, also jedem Paar (S, i⃗) ∈ Ω_{L} × D_{S∈Ω_{L}} einen Ausführungszeitpunkt zuweist. Eine solche Abbildung bezeichnet man als Schedule, falls sie die Scheduleeigenschaft besitzt.

Definition 3.1.8 (Scheduleeigenschaft) Eine Abbildung

Θ : (Ω_{L}, D_{S∈Ω_{L}}) → ℚ

besitzt für einen Schleifenprogramm L genau dann die Scheduleeigenschaft, wenn für jedes abhängige Paar von Zuweisungen S, T ∈ Ω_{L} gilt:

(∀ i⃗, j⃗ : i⃗ ∈ D_{S} ∧ j⃗ ∈ D_{T} ∧ (S, i⃗)δ(T, j⃗) : Θ(S, i⃗) < Θ(T, j⃗))

Mit

(∀ S : S ∈ Ω_{L} : Θ_{S}(i⃗) := Θ(S, i⃗))

wird der Schedule der Zuweisung S definiert.

Man beachte, daß Schedules von unterschiedlichen Zuweisungen i. a. unterschiedlich sind.

Anmerkung 3.1.9 Für eine uniforme Abhängigkeit mit Distanzvektor d⃗ kann eine affine Abbildung Θ der Form

Θ(S, i⃗) = t⃗i⃗ + c

dadurch auf die Scheduleeigenschaft bzgl. dieser Abhängigkeit getestet werden, daß man prüft, ob t⃗d⃗ > 0 gilt. Dies folgt unmittelbar aus der Definition uniformer Abhängigkeiten und der Festlegung eines Schedules.

Im Zusammenhang mit einem Schedule definiert man als eine weitere Abbildung die Allokation:

Definition 3.1.10 (Allokation) Wird für ein Schleifenprogramm L durch die Abbildung

Π : (Ω_{L}, D_{S∈Ω_{L}}) → ℚ

jeder Instanz einer Zuweisung ein Prozessor zugewiesen, auf dem diese Instanz ausgeführt wird, so heißt die Abbildung Π Allokation für L. Die Allokation einer Zuweisung wird mittels

(∀ S : S ∈ Ω_{L} : Π_{S}(i⃗) := Π(S, i⃗))

festgelegt.

In vielen Verfahren wird die Allokation in Abhängigkeit vom Schedule berechnet. Man fordert dabei, daß zwei Instanzen, die entsprechend dem Schedule zur selben Zeit berechnet werden, nicht auf demselben Prozessor berechnet werden. Formal:

(∀ i⃗, j⃗ : i⃗ ∈ D_{S}, j⃗ ∈ D_{T} : Θ(S, i⃗) = Θ(T, j⃗) ⇒ Π(S, i⃗) ≠ Π(T, j⃗))

Diese Bedingung vereinfacht die im späteren Verlauf der Parallelisierung nötige Programmtransformation. Geht man davon aus, daß ein Prozessor zu einer Zeit nur eine Zuweisung abarbeitet und interpretiert man die Werte des Schedules als Systemzeit, so hat diese Bedingung auch einen technischen Ursprung. Repräsentieren diese Werte jedoch eine virtuelle Zeit, so ist diese Anforderung an die Allokation nicht notwendig, da die verschiedenen Instanzen dann in beliebiger Reihenfolge auf demselben Prozessor ausgeführt werden können.

Wir haben Schedule und Allokation als Abbildungen mit Bild in ℚ festgelegt. Es ist aber auch möglich, diese Abbildungen mit einem höherdimensionalen Bild zu berechnen. Dies kann vor allem bei der Allokation im Hinblick auf die Architektur des Parallelrechners hilfreich sein. Hat dieser beispielsweise ein zweidimensionales Prozessorfeld, so benötigt man auch eine zweidimensionale Allokation.

Abhängigkeiten definieren ein zeitliches Verhältnis zwischen Zuweisungen bzw. Instanzen von Zuweisungen. Diese Relation veranschaulicht der allgemeine Abhängigkeitsgraph.

Definition 3.1.11 (Allgemeiner Abhängigkeitsgraph) Ein allgemeiner Abhängigkeitsgraph G ist ein Tupel (V, E, D, R). Dabei gilt:

• V ist die Knotenmenge von G. Jede Zuweisung S ∈ Ω_{L} ist ein Knoten in G. Ihr wird mittels der Abbildung D : V → D_{Ω_{L}} ein Indexraum zugeordnet.

• E ist die Kantenmenge von G. Zwischen zwei Knoten S, T ∈ V existiert eine Kante e = (S, T), falls es abhängige Instanzen dieser Zuweisungen gibt. Formal:

E = {(S, T) : (∃ i⃗, j⃗ : i⃗ ∈ D_{S} ∧ j⃗ ∈ D_{T} : (S, i⃗)δ(T, j⃗))}

Mit σ(e) bezeichnet man die Zuweisung, von der eine Abhängigkeit ausgeht, mit δ(e) das Ziel der Abhängigkeit. Einer Kante e wird durch die Abbildung R : E → ℤ^{p_{S}+p_{T}} ein Polyeder R_{e} zugeordnet, wobei

R_{e} = {(i⃗, j⃗) : i⃗ ∈ D_{S} ∧ j⃗ ∈ D_{T} ∧ (S, i⃗)δ(T, j⃗)} (1)

R_{e} enthält die Paare von Indexvektoren von Instanzen, zwischen denen eine Abhängigkeit existiert.

Beispiel 3.1.12 Gegeben sei folgendes Programm:

for i := 0 to n do

for j := 0 to n do

a[i, j] := a[i − 2, j + 1] {S_{1}} end s[i] := a[i, n] {S_{2}}

end

Da es sich hier um ein nicht perfekt verschachteltes Programm handelt, sind die Indexräume der Zuweisungen unterschiedlich. Für die Beschreibung des Indexraums D_{S_{1}} ergeben sich die Ungleichungen

0 ≤ i, j ≤ n,

für D_{S_{2}}

0 ≤ i ≤ n

Es existieren zwei Abhängigkeiten. Zum einen entsteht durch das Schreiben von a[i, j] in Zuweisung S_{1} und das spätere Lesen von a[i − 2, j + 1] in Zuweisung S_{1} eine Abhängigkeit mit Kante e_{1} = (S_{1}, S_{1}) im allgemeinen Abhängigkeitsgraphen. Für die Instanzen gilt:

(∀ i_{1}, j_{1}, i_{2}, j_{2} : (i_{1}, j_{1}) ∈ D_{S_{1}} ∧ (i_{2}, j_{2}) ∈ D_{S_{1}} ∧ i_{2} = i_{1} − 2 ∧ j_{2} = j_{1} + 1 : (S_{1}, (i_{1}, j_{1}))δ(S_{1}, (i_{2}, j_{2})))

e_{1}  e_{2} e_{3} S_{1} S^{2}

Abbildung 1: Der allgemeine Abhängigkeitsgraph aus Bespiel 3.1.12.

Damit sind auch die Indexpunkte des entsprechenden Polyeders R_{e_{1}} ∈ ℤ^{4} festgelegt. Mit (1) beschreiben folgende Ungleichungen R_{e_{1}}:

2 ≤ i_{1} ≤ n
0 ≤ j_{1} ≤ n − 1
i_{1} − 2 ≤ i_{2} ≤ i_{1} − 2
j_{1} + 1 ≤ j_{2} ≤ j_{1} + 1

Eine weitere Abhängigkeit mit Kante e_{2} = (S_{1}, S_{2}) im allgemeinen Abhängigkeitsgraphen gibt es zwischen den Zuweisungen S_{1} und S_{2}. Die letzte Speicherzelle, die in der inneren Schleife beschrieben wird, ist a[i, n]. Diese wird im Anschluß von S_{2} gelesen. Die Indexpunkte abhängiger Instanzen sind wie folgt festgelegt:

(∀ i_{1}, j_{1}, i_{2} : (i_{1}, j_{1}) ∈ D_{S_{1}} ∧i_{2} ∈ D_{S_{2}} ∧i_{2} = i_{1}∧j_{1} = n : (S_{1}, (i_{1}, j_{1}))δ(S_{2}, i_{2}))

Der allgemeine Abhängigkeitsgraph ist in Abbildung 1 dargestellt. Wie zuvor ergeben sich mit (1) für die Beschreibung des Polyeders R_{e_{2}} ∈ ℤ^{3} folgende Ungleichungen:

0 ≤ i_{1} ≤ n
n ≤ j_{1} ≤ n
i_{1} ≤ i_{2} ≤ i_{1}

[figure: 0; n; n; n; i_{1}; i_{2}; j_{1}; D_{S_{1}}; D_{S_{2}}]

Abbildung 2: Der Polyeder R_{e_{2}} aus Beispiel 3.1.12.

R_{e_{2}} hat also die in Abbildung 2 gezeigte graphische Darstellung (Punkte in der Graphik repräsentieren Paare von Indexpunkten abhängiger Instanzen).

Beim allgemeinen Abhängigkeitsgraphen wird also in den assoziierten Polyedern R_{e} festgehalten, welche Instanzen voneinander abhängen. Diese Darstellung hat aber zum einen Nachteile in der Implementierung, zum anderen sind aus den Ungleichungen, die den Polyeder beschreiben, die Indexvektoren der zugehörigen Instanzen nur schwer zu erkennen, solange keine graphische Darstellung existiert. Abhilfe schafft eine weitere Möglichkeit zur Repräsentation von Abhängigkeiten:

Mit den anfangs genannten Einschränkungen an die Eingabe-Programme ist es mittels Datenflußanalyse möglich, für jede Abhängigkeit zwischen den Zuweisungen S, T ∈ Ω_{L} mit Kante e = (S, T) im allgemeinen Abhängigkeitsgraphen eine affine Funktion h_{e} : D_{δ(e)} → D_{σ(e)} zu finden, die für jede Instanz im Ziel einer Abhängigkeit die Instanz einer Zuweisung angibt, die Quelle der Abhängigkeit ist [4]. Formal:

(∃ h : h ∈ D_{δ(e)} → D_{σ(e)} : (∀ i⃗, j⃗ : i⃗ ∈ D_{σ(e)} ∧ j⃗ ∈ D_{δ(e)}∧
(S, i⃗)δ(T, j⃗) : h_{e}(j⃗) = i⃗))

for i := 0 to n do s := i; {S_{1}} for j := 0 to n do

a[i, j] := s; {S_{2}} end

end

Abbildung 3: Kein eindeutiges Ziel der Abhängigkeit S_{1}δS_{2}

Diese Funktion h_{e} bezeichnen wir im folgenden als h-Transformation.

Anmerkung 3.1.13 Es stellt sich die Frage, ob die h-Transformation nicht auch als Funktion im Raum D_{σ(e)} → D_{δ(e)} definiert werden kann. Das ist aber nicht möglich, da die Quelle einer Abhängigkeit zwar immer eindeutig ist, nicht aber das Ziel einer Abhängigkeit, wie das Fragment in Abbildung 3 für S_{1}δS_{2} zeigt. Bei der zweiten Wahl wäre hier eine Darstellung der Abhängigkeit durch eine Funktion nicht möglich.

Anmerkung 3.1.14 Bei einer uniformen Abhängigkeit ist die h-Transfor-mation eine Translation um den Distanzvektor.

Definiert man nun den Polyeder P_{e}, der die Indexpunkte von Instanzen im Ziel der Abhängigkeit beschreibt, mit

P_{e} = {j⃗: j⃗ ∈ D_{δ(e)} ∧ (∃i⃗: i⃗ ∈ D_{σ(e)} : i⃗ = h_{e}(j⃗))}

so gilt:

(∀i⃗, j⃗: i⃗ ∈ D_{σ(e)} ∧ j⃗ ∈ D_{δ(e)} : (i⃗, j⃗) ∈ R_{e})
⇔ (∀ i⃗, j⃗ : i⃗ ∈ D_{σ(e)} ∧ j⃗ ∈ D_{δ(e)} : i⃗ = h_{e}(j⃗) ∧ j⃗ ∈ P_{e})

Durch die h-Transformation h_{e} und den Polyeder P_{e} ist es folglich möglich, die Punkte in R_{e} und damit die Abhängigkeit vollständig zu beschreiben. Wir werden diese Darstellung im folgenden benutzen.

Beispiel 3.1.15 Die h-Transformationen aus Beispiel 3.1.12 ergeben sich wie folgt:

h_{e_{1}}(i, j) = (i − 2, j + 1)
h_{e_{2}}(i) = (i, n)

Die zugehörigen Polyeder sind:

P_{e_{1}} = {i, j : 2 ≤ i ≤ n ∧ 0 ≤ j ≤ n − 1}
P_{e_{2}} = {i : 0 ≤ i ≤ n}

3.2 Grundidee der schleifenbasierten Parallelisierung

Bei der Parallelisierung eines Eingabe-Programms unterscheidet man zwischen zwei Ansätzen. Einerseits existieren Verfahren, welche die Instanzen von Zuweisungen eines Programms unter Berücksichtigung der Abhängigkeiten neu ordnen. Wir bezeichnen diesen Ansatz als statement-basiert. Ein anderes und oberflächlicheres Prinzip ist es, eine Schleife mit ihrem Schleifenkörper so umzuformen, daß diese Schleife parallel ausgeführt werden kann. In diesem Abschnitt wollen wir die notwendigen Grundlagen für den letzten, sogenannten schleifenbasierten Ansatz vorstellen. Wir suchen hier nach einem Kriterium für die Ersetzung von sequentiellen Schleifen durch parallele Schleifen. Die Schleifenindizes und ihre Grenzen bleiben dabei also erhalten. Wir werden im nächsten Abschnitt ein schleifenbasiertes Verfahren kennenlernen, das auf diesen Grundlagen aufbaut und zusätzlich eine Transformation der Schleifenindizes vornimmt.

Wir beschränken uns bei den folgenden Untersuchungen auf ein perfekt verschachteltes Schleifenprogramm mit ausschließlich uniformen Abhängigkeiten. Zur Notation paralleler Schleifensprogramme verwenden wir zur Kennzeichnung von sequentiellen Schleifen das Konstrukt seqfor, zur Kennzeichnung paralleler Schleifen das Konstrukt parfor. Für die formale Darstellung der Ersetzung definieren wir das Prädikat parfor mit folgender Festlegung: Wird eine Schleife L_{i} eines Schleifensprogramms L = (L_{1}, …, L_{m}) durch eine parallele parfor-Schleife ersetzt, so gilt parfor(L_{i}).

Allgemein gilt, daß sich durch die Parallelisierung die Ausführungsordnung eines Eingabe-Programms ändert. Im sequentiellen Eingabe-Programm ist diese Ordnung eine Totalordnung, dadurch festgelegt, daß eine Iteration B(i⃗) vor einer Iteration B(j⃗) ausgeführt wird, wenn i⃗ ≺ j⃗ gilt. Bei der Erzeugung des parallelen Programms wird ausgenutzt, daß lediglich abhängige Iterationen nacheinander ausgeführt werden müssen, um die Semantik zu erhalten. Die Ausführungsordnung ist dann eine partielle Ordnung. Wir suchen also nach einer neuen Ausführungsordnung für das gegebene Schleifenprogramm und bezeichnen diese als gültig, wenn für zwei abhängige Iterationen B(i⃗) und B(j⃗) des transformierten Schleifenprogramms gilt, daß B(i⃗) vor B(j⃗) ausgeführt wird.

Um die Semantik eines parallelen Schleifenprogramms bestehend aus seqfor - und parfor-Schleifen definieren zu können, benötigen wir folgende Festlegung:

Definition 3.2.1 (parallele Ausführungsordnung) Für ein paralleles Schleifenprogramm L′ wird die Ausführungsordnung wie folgt definiert:

Eine Iteration B(i⃗) von L′ wird vor einer anderen Iteration B(j⃗) ausgeführt, wenn es eine seqfor-Schleife L′_{r} gibt und i⃗ ≺_{r} j⃗ gilt.

Definition 3.2.2 (Äquivalenz von Schleifensätzen) Sei L ein sequentielles Schleifenprogramm und L′ ein aus L erzeugtes Schleifenprogramm mit neuer Ausführungsordnung, so heißen L und L′ äquivalent, i. Z. L ≡ L′, wenn die Ausführungsordnung von L′ gültig ist.

Folgender Satz gibt uns das gesuchte Kriterium [5]:

Lemma 3.2.3 (Parallelisierbarkeit von Schleifen) Gegeben sei ein sequentielles Schleifenprogramm L = (L_{1}, …, L_{m}). Sei L′ = (L′_{1}, …, L′_{m}) das Schleifenprogramm, das aus L durch Ersetzung einiger Schleifen durch parfor-Schleifen entsteht. Dann ist L äquivalent zu L′, genau dann wenn keine der ersetzten Schleifen eine Abhängigkeit in L trägt. Formal: Ist L′ = {L′_{1}, …, L′_{m}} die Menge der Schleifen von L′, L′_{repl} ⊆ L′ die Menge der ersetzten Schleifen und D die Menge der Distanzvektoren aller Abhängigkeiten, so gilt:

L ≡ L′ ⇔
(∀ i : 1≤i≤m ∧ L′_{i} ∈ L′_{repl} : (∀ d⃗ : d⃗ ∈ D : ¬carries_{L_{i}}(d⃗)))

Beweis:

“⇐”: Voraussetzung: Keine der ersetzten Schleifen trägt eine Abhängigkeit in L. Um die Äquivalenz zu beweisen, muß für zwei abhängige Iterationen in L gezeigt werden, daß diese auch in L′ nacheinander ausgeführt werden. Sei dazu B(j⃗) eine Iteration in L, die von einer Iteration B(i⃗) in L mit Distanzvektor d⃗ abhängt. Dann gilt:

B(j⃗)δB(i⃗) ⇒ { sequentielle Ausführungsordnung in L }

i⃗ ≺ j⃗ ⇔ { lexikographische Ordnung }

(∃ r : 1≤r ≤m : i⃗ ≺_{r} j⃗) ⇔ { d⃗ Distanzvektor der Abhängigkeit }

(∃ r : 1≤r ≤m : d⃗ ≻_{r} 0) ⇔ { Def. 3.1.5 }

(∃ r : 1≤r ≤m : level(d⃗) = r) ⇔ { Def. 3.1.6 }

(∃ r : 1≤r ≤m : carriesd⃗(L_{r})) ⇒ { Voraussetzung }

L_{r} ∉ L′_{repl} ⇒ { Def. 3.2.1 }

B(i⃗) wird in L′ vor B(j⃗) ausgeführt ⇔ { Definition gültiger Ausführungsordnung }

Ausführungsordnung gültig ⇔ { Def. 3.2.2 }

L ≡ L′

“⇒”: Zu zeigen ist:

(L ≡ L′) ⇒ (∀ r : 1≤r ≤m : (∀ d⃗ : d⃗ ∈ D : parfor(L_{r}) ⇒ ¬carriesd⃗(L_{r})))

Wir zeigen die inverse Implikation:

(∃ r : 1≤r ≤m : (∃ d⃗ : d⃗ ∈ D : parfor(L_{r}) ∧ carriesd⃗(L_{r}))) ⇒ (L ≢ L′)

Voraussetzung: L_{r} wird in L′ durch eine parfor-Schleife L′_{r} ersetzt und trägt eine Abhängigkeit mit Distanzvektor d⃗:

(∃ r : 1≤r ≤m : (∃ d⃗ : d⃗ ∈ D : carriesd⃗(L_{r})) ∧ parfor(L_{r})) ⇒ { Def. 3.1.6, Def. 3.1.5 }

(∃ r : 1≤r ≤m : (∃ d⃗ : d⃗ ∈ D : d⃗ ≻_{r} 0) ∧ parfor(L_{r})) ⇒ { uniforme Abhängigkeit }

(∃ r : 1≤r ≤m : (∃ i⃗, j⃗ : : i⃗ ≻_{r} j⃗) ∧ parfor(L_{r})) ⇒ { Def. 3.2.1 }

B(i⃗) wird in L′ nicht notwendig vor B(j⃗) ausgeführt ⇔ { Def. gültige Ausführungsordnung }

Ausführungsordnung ungültig ⇔ { Def. 3.2.2 }

L ≢ L′

Folgendes Beispiel soll die Anwendung von Lemma 3.2.3 veranschaulichen:

Beispiel 3.2.4 Gegeben sei folgendes Programm:

for i := 0 to n do

for j := 0 to n do

for k := 0 to n do a[i, j, k] := a[i − 1, j + 2, k − 3]; b[i, j, k] := b[i, j − 5, k + 2]

end end

end

In diesem Schleifenprogramm existieren zwei uniforme Abhängigkeiten mit den Distanzvektoren d⃗_{1} = (1, −2, 3) und d⃗_{2} = (0, 5, −2). Der Level von d⃗_{1} ist 1, der von d⃗_{2} ist 2. Somit tragen die Schleifen L_{1} und L_{2} eine Abhängigkeit, nicht jedoch L_{3}. Mit Lemma 3.2.3 ergibt sich also folgendes äquivalente Schleifenprogramm:

seqfor i := 0 to n do seqfor j := 0 to n do

parfor k := 0 to n do a[i, j, k] := a[i − 1, j + 2, k − 3]; b[i, j, k] := b[i, j − 5, k + 2]

end end

end

Wir sind nun in der Lage anhand einer Untersuchung der Abhängigkeiten über die Ersetzung von sequentiellen Schleifen durch parallele zu entscheiden. Wir werden im nächsten Abschnitt ein Verfahren zur Berechnung von Schedules kennenlernen, das diese Tatsache nach einer Transformation des Programmcodes ausnutzt.

4 Berechnung von Schedules

In diesem Abschnitt werden Methoden zur Bestimmung von Schedules vorgestellt.

Wir werden zunächst ein einfaches, schleifenbasiertes Verfahren kennenlernen, das sich nur auf die auftretenden Abhängigkeiten stützt und eine sogenannte Transformationsmatrix berechnet.

Als zweites Verfahren wird eine Methode eingeführt, die anhand einer Untersuchung der Indexräume stückweise affine Schedules berechnet. Diese Schedules sind unter noch zu beschreibenden Einschränkungen zeit-optimal.

4.1 Die Hyperebenen-Methode

Bei diesem Verfahren handelt es sich um eine einfache schleifenbasierte Parallelisierungsmethode. Als Ergebnis liefert es eine Transformationsmatrix, in der ein Schedule und eine Allokation enthalten ist. Die Hyperebenen-Methode betrachtet zwar neben den Abhängigkeiten auch den Indexraum, der vorzustellende Algorithmus von Lamport[8] wird diesen aber aus Gründen der Komplexität vernachlässigen. Aus diesem Grund ist der vom Algorithmus implizit berechnete Schedule i. a. nicht zeit-optimal. Ein weiterer Nachteil sind die starken Einschränkungen an die Eingabe-Programme:

• Das Eingabe-Programm ist perfekt verschachtelt.

• Alle Abhängigkeiten sind uniform.

Die Idee des Verfahrens basiert auf einer Partitionierung des Indexraums in Teilmengen von Instanzen, die parallel ausgeführt werden können. Geometrisch betrachtet beschränkt sich das Verfahren dabei auf Teilmengen, die die Form von Hyperebenen haben.

Definition 4.1.1 (Hyperebene) Für einen Zeilenvektor u⃗ ∈ ℝ^{m} \ {0} und eine Konstante c ∈ ℝ heißt eine Menge H ⊂ ℝ^{m} mit

H = {i⃗: i⃗ ∈ ℝ^{m} ∧ u⃗i⃗ = c}

eine Hyperebene mit Richtung u⃗. Durch zwei Konstanten p, q ∈ ℤ wird ein System von Hyperebenen erzeugt mit

H_{0} = {i⃗: i⃗ ∈ ℝ^{m} ∧ u⃗i⃗ = p}
H_{1} = {i⃗: i⃗ ∈ ℝ^{m} ∧ u⃗i⃗ = p + 1}
⋮
H_{q−p} = {i⃗: i⃗ ∈ ℝ^{m} ∧ u⃗i⃗ = p + q}

Betrachtet man das perfekt verschachtelte Eingabe-Programm, so gilt hier für die Werte von Index i_{1} der äußeren Schleife mit unterer Schleifengrenze l_{1} und oberer Schleifengrenze u_{1}:

i_{1} = l_{1}, i_{1} = l_{1} + 1, …, i_{1} = u_{1}

oder

( 1 0 … 0 ) ( i_{1} ; ⋮ ; i_{m} ) = l_{1}, …, ( 1 0 … 0 ) ( i_{1} ; ⋮ ; i_{m} ) = u_{1}

Also definiert ein Schleifenprogramm ein System von Hyperebenen mit Richtung (1, 0, …, 0).

Beispiel 4.1.2 Gegeben sei das Schleifenprogramm

for i := 1 to 4 do
for j := 0 to 5 − i do

for k := 0 to 2 do {Schleifenkörper}

end end

end

mit beliebigem Schleifenkörper. Der Indexraum aller Zuweisungen des Schleifenkörpers ist in der Form von Hyperebenen in Abbildung 4 dargestellt.

Unser Ziel ist es einen neuen Richtungsvektor zu finden, sodaß sich in den damit verbundenen neuen Hyperebenen ausschließlich unabhängige Instanzen von Zuweisungen befinden. Wie zuvor werden diese Hyperebenen

[figure: i; j; 0]

Abbildung 4: Die Hyperebenendarstellung des Indexraums in Beispiel 4.1.2.

durch einen einzigen Schleifenindex aufgezählt, alle anderen Schleifenindizes spannen jeweils eine Dimension einer Hyperebene auf. Da wir die Hyperebenen so konstruieren wollen, daß sie keine abhängigen Instanzen beinhalten, trägt außer des erstgenannten kein Schleifenindex eine Abhängigkeit. Nach der Grundidee des schleifenbasierten Ansatzes können somit alle Schleifen mit den letztgenannten Schleifenindizes parallel ausgeführt werden.

Mathematisch wird diese Transformation durch eine invertierbare Matrix dargestellt, welche die Schleifenindizes des Eingabe-Programms so umformt, daß die genannte Eigenschaft bzgl. der Hyperebenen erfüllt ist. Wir werden sehen, daß es immer möglich ist, einen solche sogenannte Transformationsmatrix zu finden.

Bevor wir den Hauptsatz dieses Abschnitts formulieren können, ist noch theoretische Vorarbeit nötig. Zunächst soll geklärt werden, welche Eigenschaft eine Transformationsmatrix auszeichnet. Um die Semantik des Eingabe-Pro-gramms nicht zu zerstören, muß die durch den Schedule bzw. die Abhängigkeiten festgelegte zeitliche Ordnung auch im transformierten Programm gelten. Eine Transformationsmatrix, die diese Eigenschaft erfüllt, heißt gültig bezüglich eines Schleifenprogramms.

Definition 4.1.3 (gültige Transformation, valid_{L}(T)) Gegeben sei ein perfekt verschachteltes Schleifenprogramm L mit Zuweisungsmenge Ω_{L} und Indexraum D. Eine Transformationsmatrix T heißt gültig bezüglich des Schleifenprogramms L, genau dann wenn die durch die Abhängigkeiten festgelegte zeitliche Ordnung auch nach der Anwendung der Transformation auf die Indexvektoren gilt. Formal:

(∀ S_{1}, S_{2} : S_{1}, S_{2} ∈ Ω_{L} : (∀ i⃗, j⃗ : i⃗, j⃗ ∈ D ∧ (S_{1}, i⃗)δ(S_{2}, j⃗) : valid_{L}(T) ⇔ T i⃗ ≺ T j⃗))

Lemma 4.1.4 (Charakterisierung gültiger Transformationen) Eine Transformation T eines perfekt verschachtelten Schleifenprogramms L mit ausschließlich uniformen Abhängigkeiten ist genau dann gültig, wenn für jeden Distanzvektor d⃗ einer uniformen Abhängigkeit T d⃗ ≻ 0 gilt. Formal:

Sei D die Menge der Distanzvektoren aller uniformen Abhängigkeiten. Dann gilt:

(∀ d⃗ : d⃗ ∈ D : valid_{L}(T) ⇔ T d⃗ ≻ 0))

Beweis:

valid_{L}(T) ⇔ { Def. 4.1.3 }

(∀ S_{1}, S_{2} : S_{1}, S_{2} ∈ Ω_{L} : (∀ i⃗, j⃗ : i⃗, j⃗ ∈ D ∧ (S_{1}, i⃗)δ(S_{2}, j⃗) : T i⃗ ≺ T j⃗)) ⇔ { Arithmetische Umformung }

(∀ S_{1}, S_{2} : S_{1}, S_{2} ∈ Ω_{L} : (∀ i⃗, j⃗ : i⃗, j⃗ ∈ D ∧ (S_{1}, i⃗)δ(S_{2}, j⃗) : T (j⃗ − i⃗) ≻ 0)) ⇔ { Sei d⃗ = j⃗ − i⃗ der Distanzvektor der Abhängigkeit }

(∀ d⃗ : d⃗ ∈ D : T d⃗ ≻ 0)

Lemma 4.1.5 (transformierte Distanzvektoren) Ist T eine gültige Transformation eines Schleifenprogramms L in ein Schleifenprogramm L′ und D die Menge der Distanzvektoren der Abhängigkeiten in L, so gilt für die Menge D′ der Distanzvektoren der Abhängigkeiten in L′:

D′ = {T d⃗: d⃗ ∈ D}

Beweis:

d⃗ ∈ D
⇔ { Def. 3.1.4 }
(∃ i⃗, j⃗ : : d⃗ = j⃗ − i⃗ ∧ B(i⃗)δB(j⃗))
⇔ { Definiere k⃗ = T i⃗ und l⃗ = T j⃗ und mit valid_{L}(T) = true }
B′(k⃗)δB′(l⃗)
⇒ { Def. 3.1.4 }
(∃ d⃗′ : d⃗′ ∈ D′ : d⃗′ = l⃗− k⃗)
⇔ { k⃗ = T i⃗ und l⃗ = T j⃗ }
(∃ d⃗′ : d⃗′ ∈ D′ : d⃗′ = T j⃗ − T i⃗)
⇔ { d⃗ = j⃗ − i⃗ }
(∃ d⃗′ : d⃗′ ∈ D′ : d⃗′ = T d⃗)

Da wegen der Gültigkeit der Transformation eine Abhängigkeit zwischen zwei Iterationen in L′ genau dann besteht, wenn die gleiche Abhängigkeit in L existiert, sind damit alle Abhängigkeiten in L′ erfaßt.

Eine Transformationsmatrix, deren Determinante ±1 ist, erleichtert den Prozess der Programmkonstruktion, wie wir im letzten Abschnitt sehen werden. Eine Matrix mit dieser Eigenschaft heißt unimodular.

Wir sind nun in der Lage, den Hauptsatz dieses Abschnitts zu formulieren, der im Beweis bereits den Algorithmus zur Berechnung der Transformationsmatrix enthält [5]:

Satz 4.1.6 (Existenz einer Transformation) Für ein gegebenes perfekt verschachteltes Schleifenprogramm L = (L_{1}, …, L_{m}) mit ausschließlich uniformen Abhängigkeiten existiert eine gültige unimodulare Transformation T in ein Schleifenprogramm L′ = (L′_{1}, …, L′_{m}), sodaß die Schleifen L′_{2}, …, L′_{m} parallel ausgeführt werden können.

Beweis: Wir unterscheiden zunächst zwei Fälle:

D = ∅:

D = ∅ ⇒ { Lemma 3.2.3 }

L_{1}, …, L_{m} können parallel ausgeführt werden

Also erfüllt T := I_{m} die geforderten Eigenschaften.

D ≠ ∅: Im Beweis benötigen wir folgende Definition einer Partitionierung von D:

(∀ r : 1≤r ≤m : D_{r} = {d⃗ : d⃗ ∈ D ∧ level(d⃗) = r})

Da D die Vereinigung der Mengen D_{r} darstellt und nicht die leere Menge ist, folgt daß zumindest eine der Mengen D_{r} ungleich der leeren Menge ist. Wir nehmen nun an, daß eine Matrix T existiert, für die gilt

(valid_{L}(T) ∧ L′_{2}, …, L′_{m}sind parfor-Schleifen )

und zeigen durch Äquivalenzumformungen, daß T konstruierbar ist:

(valid_{L}(T) ∧ L′_{2}, …, L′_{m}sind parfor-Schleifen ) ⇔ { Def. 4.1.3 }

(∀ d⃗ : d⃗ ∈ D : T d⃗ ≻ 0) ∧ L′_{2}, …, L′_{m} sind parfor-Schleifen ⇔ { Lemma 3.2.3 }

(∀ d⃗ : d⃗ ∈ D : T d⃗ ≻ 0 ∧ (∀ i : 2≤i≤m : ¬carries_{d}(L_{i}))) ⇔ { Lemma 4.1.5 und Def. 3.1.6 }

(∀ d⃗ : d⃗ ∈ D : T d⃗ ≻ 0 ∧ level(T d⃗) = 1) ⇔ { Def. 3.1.5 }

(∀ d⃗ : d⃗ ∈ D : T d⃗ ≻_{1} 0) ⇔ { Definiere t⃗ als erste Zeile in T }

(∀ d⃗ : d⃗ ∈ D : t⃗d⃗≥1) ⇔ { Matrixprodukt }

(∀ d⃗ : d⃗ ∈ D : (∑ i : 1≤i≤m : t_{i}d_{i})≥1) ⇔ { Definition von D_{r} und level }

(∀ r : 1≤r ≤m : (∀ d⃗ : d⃗ ∈ D_{r} : (∑ i : r ≤i≤m : t_{i}d_{i})≥1)) ⇔ { Auflösen nach t_{r}, wobei d_{r} >0 }

(∀ r : 1≤r ≤m : (∀ d⃗ : d⃗ ∈ D_{r} : t_{r} ≥(1− ; (∑ i : r + 1≤i≤m : t_{i}d_{i}))/d_{r}))  (1)

Ist D_{r} nicht leer, so existiert damit eine untere Schranke für t_{r}. Im anderen Fall ist t_{r} unbeschränkt und kann einen beliebigen Wert annehmen.

Die Aussage in (1) kann als Konstruktionsvorschrift verwendet werden, indem man sukzessive für r = m, …, 1 die Werte der Einträge t_{r} jeweils gleich der entsprechenden unteren Schranke setzt.

Somit besitzt also jede Matrix, deren erste Zeile das Ungleichungssystem (1) erfüllt, aufgrund der Äquivalenz im Beweis die in der Annahme geforderten Eigenschaften. Sei also t⃗ = (t_{1}, …, t_{m}) ein solcher Vektor. Es bleibt noch zu zeigen, daß t⃗ zu einer unimodularen Matrix ergänzt werden kann. Mit g = gcd(t_{1}, …, t_{m}) erfüllt auch der Vektor t⃗′ = (t_{1}/g, …, t_{m}/g) das Ungleichungssystem. Mittels Methoden der Zahlentheorie [5] ist es nun möglich den Vektor t⃗′ mittels Echelon-Reduktion zu einer unimodularen Matrix T zu ergänzen, die damit sämtliche geforderten Eigenschaften besitzt.

Anmerkung 4.1.7 Für die berechnete Matrix T gilt laut Beweis:

(∀ d⃗ : : T d⃗ ≻_{1} 0)

Mit dem ersten Zeilenvektor t⃗ gilt damit auch t⃗d⃗ > 0. Mit der Anmerkung 3.1.9 ist also die Abbildung

Θ(i⃗) = t⃗i⃗

ein Schedule für das gegebene Schleifenprogramm.

Der Beweis des Satzes 4.1.6 zeigt, daß es sogar beliebig viele solcher Transformationen gibt, denn für die Einträge der erste Zeile der Matrix existieren lediglich untere Schranken. Unser Ziel ist es, ein laufzeitoptimales Programm zu erhalten. Die Laufzeit eines Programms, dessen innere Schleifen alle parallel ausgeführt werden, wird (abgesehen vom Kommunikationsaufwand) ausschließlich durch die Grenzen der äußeren sequentiellen Schleife bestimmt. Allgemein gilt für die Anzahl von Iterationen bei gegebenem Schedule Θ für das parallele Programm

(max Θ(i⃗) : i⃗ ∈ D_{Ω_{L}} : Θ(i⃗)) − (min Θ(i⃗) : i⃗ ∈ D_{Ω_{L}} : Θ(i⃗)) + 1

Mit Anmerkung 4.1.7 kann damit die Anzahl von Iterationen eines transformierten Programms berechnet werden. Sei dazu t der erste Zeilenvektor einer Transformationsmatrix wie in Satz 4.1.6. Dann ist diese Anzahl gleich folgendem Ausdruck:

(max t⃗i⃗ : i⃗ ∈ D_{Ω_{L}} : t⃗i⃗)) − (min t⃗i⃗ : i⃗ ∈ D_{Ω_{L}} : t⃗i⃗) + 1

Zusammenfassend ist also folgendes Optimierungsproblem zu lösen:

Optimierungsproblem der Hyperebenen-Methode

Gegeben: Ein perfekt verschachteltes Schleifenprogramm L = (L_{1}, …, L_{m}) mit einer Menge D von Distanzvektoren uniformer Abhängigkeiten.

Gesucht: Ein Vektor t⃗ = (t_{1}, …, t_{m}) mit folgenden Eigenschaften:

1. gcd(t_{1}, …, t_{m}) = 1

2. (∀ d⃗ : d⃗ ∈ D ∧ d⃗ = (d_{1}, …, d_{m}) : (∑ k : 1≤k ≤m : t_{k}d_{k})≥1)

3. Der Wert des Ausdrucks

(max t⃗i⃗ : i⃗ ∈ D_{Ω_{L}} : t⃗i⃗) − (min t⃗i⃗ : i⃗ ∈ D_{Ω_{L}} : t⃗i⃗) + 1

ist minimal.

Die Lösung dieses Minimierungsproblems gestaltet sich allerdings sehr aufwendig. Dies gilt vor allem dann, wenn die Grenzen innerer Schleifenindizes von äußeren Schleifenindizes abhängen.

Deshalb stellen wir eine Methode vor, welche die dritte Forderung außer acht läßt [8]. Sie beruht auf der Heuristik, daß bei der Wahl der jeweils unteren Schranken des Ungleichungssystems (1) im Beweis des Satzes 4.1.6 meistens eine Transformation für ein laufzeitoptimales Programm entsteht. Das Verfahren besteht aus den beiden Schritten, die auch im Beweis von Lemma 4.1.6 durchgeführt werden. Zunächst wird ein Vektor t⃗ mit den oben beschriebenen Eigenschaften konstruiert, der dann zu einer unimodularen Matrix ergänzt wird. Zur Konstruktion der Matrix wird jedoch nicht, wie im Beweis von Satz 4.1.6, die Echelon-Reduktion benutzt:

Algorithmus 4.1.8 (Hyperebenen-Methode von Lamport) Eingabe: Eine Menge D ≠ ∅ von Distanzvektoren uniformer Abhängigkeiten

Ausgabe: Eine gültige Transformation T Methode:

1. Konstruktion des ersten Zeilenvektors t⃗ von T:

{ Konstruiere die Mengen D_{r} } for r := 1 to m do

D_{r} := {d⃗: d⃗ ∈ D ∧ level(d) = r} end; { Berechne t_{r} aus den Werten t_{r+1}, …, t_{m} } for r := m to 1 step −1 do

if D_{r} = ∅ then

t_{r} := 0 else

t_{r} := ⌈max_{d⃗∈D_{r}}{(1 − t_{r+1}d_{r+1} − … − t_{m}d_{m})/d_{r}}⌉ end

end

2. Ergänzung zu einer unimodularen Matrix T = (t_{ij}), wobei 1 ≤ i, j ≤m:

r := m; { Suche das erste t_{r} ≠ 0 in der Sequenz (t_{m}, …, t_{1}) } while t_{r} = 0 do

r := r − 1 end; k:=r; { Setze T′ gleich der Einheitsmatrix I_{m−1} } T′ := I_{m−1}; füge einen Nullvektor als k-te Spalte in T′ ein; { Ergänze T′ um t⃗ als erste Zeile }

T := ( t⃗ ; T′ )

Terminierung und Korrektheit:

Terminierung: Die Terminierung der beiden Schleifen im ersten Schritt des Algorithmus ist offensichtlich, da im Schleifenkörper weder die Laufvariable noch die Schleifengrenzen verändert werden. Um die Terminierung der while-Schleife im zweiten Schritt zu sehen, muß die Existenz eines Eintrags t_{r} ≠ 0 gezeigt werden. Da wir davon ausgehen, daß die Menge D nicht leer ist, existiert eine Menge D_{k} ≠ ∅, die diese Eigenschaft als erste in der Sequenz D_{m}, …, D_{1} besitzt. Also folgt nach der Berechnungsvorschrift im ersten Schritt, daß t_{k} = 1 gilt. Somit terminiert auch die while-Schleife.

Korrektheit: Es gilt die beiden ersten Eigenschaften des Optimierungsproblems zu zeigen. Da die Berechnung der Einträge des Vektors t⃗ der Konstruktion im Beweis von Satz 4.1.6 folgt, erfüllt t⃗ die zweite geforderte Eigenschaft. Für den ersten Eintrag t_{r} ≠ 0 in der Sequenz (t_{m}, …, t_{1}) gilt nach der Berechnungsvorschrift t_{r} = 1. Damit ist

gcd(t_{1}, …, t_{m}) = 1,

wodurch die erste geforderte Eigenschaft erfüllt wird.

Wird in der Matrix T die k-te Spalte und k-te Zeile gestrichen, so ergibt sich die Einheitsmatrix I_{m−1}. Die k-te Spalte von T hat die Form (1, 0, …, 0). Damit ergibt sich nach dem Entwicklungssatz von Laplace für die Determinante von T:

|T | = ±1

Somit ist T unimodular.

Abschließend wollen wir die Hyperebenen-Methode für ein Beispiel betrachten, das auch zeigt, daß der Algorithmus, wie bereits erwähnt, nicht immer optimale Ergebnisse liefert:

Beispiel 4.1.9 Gegeben sei folgendes Schleifenprogramm:

for i := 0 to 10 do for j := 0 to 5 do

a[i, j] := a[i − 1, j − 1]; {S_{1}} b[i, j] := a[i − 1, j − 2] {S_{2}}

end end

[figure: i; j; 5; 10; 0]

Abbildung 5: Der Indexraum des Schleifenprogramms aus Beispiel 4.1.9.

Das Schleifenprogramm enthält zwei uniforme Abhängigkeiten mit den Distanzvektoren d⃗_{1} = ( 1 1 ) und d⃗_{2} = ( 1 2 ). Der Indexraum des Schleifenprogramms und die auftretenden Abhängigkeiten sind in Abbildung 5 dargestellt. Für den Level der Abhängigkeiten gilt:

level(d⃗_{1}) = 1
level(d⃗_{2}) = 1

Daraus ergeben sich die Mengen D_{r} wie folgt:

D_{1} = {d⃗_{1}, d⃗_{2}}
D_{2} = ∅

Mit diesen Mengen gilt dann für den ersten Zeilenvektor der Transformationsmatrix T: t⃗ = ( 1 0 ). Die Ergänzung zu einer unimodularen Matrix liefert schließlich als Ergebnis:

T = ( \begin{matrix}1 & 0 \\ 0 & 1 \end{matrix} )

Somit wird also keine Transformation vorgenommen. Lediglich die innere Schleife wird parallel ausgeführt. Die Laufzeit des parallelen Programms beträgt 11 Zeiteinheiten. Berücksichtigt man jedoch auch die Minimierungsbedingung, gilt es den Ausdruck

10t_{1} + 5t_{2} mit den Nebenbedingungen, die durch die erste und zweite Eigenschaft festgelegt sind, zu minimieren. Man erhält dann für den Vektor t⃗:

t⃗ = ( 0 1 )

Dies führt zu der Transformationmatrix:

T = ( \begin{matrix}0 & 1 \\ 1 & 0 \end{matrix} )

Eine derartige Transformation tauscht die innere mit der äußeren Schleife, weshalb die Laufzeit des parallelen Programms nun lediglich 6 Zeiteinheiten beträgt.

4.2 Die Methode von Feautrier

Die zweite Methode zur Parallelisierung verfolgt einen statementbasierten Ansatz. Außer den anfangs genannten Einschränkungen fordert sie keine zusätzlichen Eigenschaften der Eingabe-Programme. Im Vergleich zur Hyperebenen-Methode entfällt also die Beschränkung auf die Uniformität der Abhängigkeiten und die perfekte Verschachtelung.

Anders als bei der Hyperebenen-Methode wird hier für jede Zuweisung ein eigenständiger Schedule berechnet. Wie bereits erwähnt, berücksichtigt dieses Verfahren auch den Indexraum einer Zuweisung. Da dieser von Strukturparametern abhängig sein kann, gilt dies auch für den berechneten Schedule.

Wir leiten zunächst ein Optimierungsproblem her, dessen Lösung den Schedule einer Zuweisung darstellt. Darauf soll eine Möglichkeit zur Lösung des Optimierungsproblems vorgestellt werden.

4.2.1 Herleitung des Optimierungsproblems

Die Grundidee des Verfahrens besteht darin, zunächst eine Menge von Schedules zu berechnen und darauf nach dem Kriterium der kürzesten Laufzeit das Minimum zu suchen. Da die Vorstellung dieser Methode die Implementierung eines Algorithmus ermöglichen soll, werden wir auch auf die Probleme dieser Minimierung eingehen.

Die Grundlage für die Konstruktion der Menge von Schedules leitet sich aus der Definition des Schedules selbst ab. Wie bereits angegeben, besitzt ein Schedule Θ für zwei abhängige Zuweisungen S, T ∈ Ω_{L} eines Schleifenprogramms L folgende Eigenschaft:

(∀ i⃗, j⃗ : i⃗ ∈ D_{S} ∧ j⃗ ∈ D_{T} ∧ (S, i⃗)δ(T, j⃗) : Θ(S, i⃗) < Θ(T, j⃗))

Der Einfachheit halber nehmen wir für die Abarbeitungszeit einer jeden Instanz eine Zeiteinheit an. Damit gilt für den Schedule:

(∀ i⃗, j⃗ : i⃗ ∈ D_{S} ∧ j⃗ ∈ D_{T} ∧ (S, i⃗)δ(T, j⃗) : Θ(S, i⃗)≤Θ(T, j⃗) + 1) (2)

Die Grundform des Verfahrens von Feautrier beschränkt sich auf eindimensionale affine Schedules mit der allgemeinen Form:

Θ(S, i⃗) = α⃗_{S}i⃗ + β⃗_{S}n⃗ + γ_{S}

Dabei bezeichnet n⃗ den Vektor der Strukturparameter, α⃗_{S} und β⃗_{S} die entsprechenden zu bestimmenden Koeffizientenvektoren und γ_{S} eine ebenfalls noch unbekannte Konstante.

Zur Darstellung des Indexraums D_{S} einer Zuweisung S ∈ Ω_{L} verwenden wir folgende Notation für die zugehörigen Ungleichungen:

D_{S} = { i⃗ : (∀ k : , 1≤k ≤m_{S} : a_{S}^{k} ( i⃗ ; n⃗ ) + b_{S}^{k} ≥ 0) }

Die Anzahl der Ungleichungen zur Beschreibung von D_{S} wird von m_{S} angegeben. Der Vektor a_{S}^{k} ist der Koeffizientenvektor der Schleifeninidizes und Sturkturparameter, b_{S}^{k} entspricht dem konstanten Anteil. Analog notieren wir die Polyeder P_{e}, welche die Indexpunkte im Ziel einer Abhängigkeit beschreiben, durch

P_{e} = { j⃗ : (∀ k : 1≤k ≤m_{e} : c_{e}^{k} ( j⃗ ; n⃗ ) + d_{e}^{k} ≥ 0) }

Hier wird die Anzahl der Ungleichungen durch m_{e} wiedergegeben.

Mit folgendem Lemma der linearen Programmierung läßt sich ein erster Zusammenhang zwischen Schedule und Indexraum herstellen:

Lemma 4.2.1 (Affine Form von Farkas’ Lemma) Sei D ein nicht leerer Polyeder, der durch die n Ungleichungen

{ a⃗_{i}i⃗ + b_{i} ≥0 : 1≤i≤n }

beschrieben ist. Dann ist eine affine Abbildung ψ : D → ℚ genau dann in ihrem gesamten Bild nicht negativ, wenn ψ eine positive Linearkombination in diesen Ungleichungen ist, also ψ folgende Form hat:

ψ(x) = λ_{0} + (∑ i : 1≤i≤n ∧ λ_{i} ≥0 : λ_{i}(a⃗_{i}i⃗ + b_{i}))

Die Koeffizienten λ_{i} werden als Farkas-Multiplikatoren bezeichnet.

Wenden wir diese Aussage auf den zu bestimmenden affinen Schedule Θ an, so folgt, daß dieser bei Anwendung auf einen Indexvektor i⃗ ∈ D_{S} positive Werte besitzt, genau dann wenn Farkas-Multiplikatoren μ_{S}^{k} mit folgender Eigenschaft existieren:

Θ(S, i⃗) = μ_{S}^{0} + (∑ k : 1≤k ≤m_{S} : μ_{S}^{k}(a_{S}^{k} ( i⃗ ; n⃗ ) + b_{S}^{k})) (3)

Eine derartige Darstellung des Schedules bezeichnen wir als den Prototypen eines Schedules.

Anmerkung 4.2.2 Der Indexraum einer Zuweisung ist neben den Schleifenindizes auch von der Größe der Strukturparameter abhängig. Aus diesem Grund werden Strukturparameter wie Schleifenindizes ohne Grenzen behandelt. Wir fordern allerdings für Strukturparameter, daß sie stets positiv sind. Dies wird durch die zusätzlichen Ungleichungen

n⃗≥0

zur Beschreibung eines Indexraums ausgedrückt. Die Ungleichungen zur Beschreibung der Polyeder P_{e} sind ebenfalls um diese Ungleichungen zu ergänzen. Eine Einschränkung stellt dies nicht dar, weil ein negatives Vorzeichen stets in den Ausdruck, in dem der Strukturparameter steht, verlagert werden kann. Wir benötigen diese Forderung vor allem für Zuweisungen ohne umgebende Schleifen, also Zuweisungen mit 0-dimensionalem Indexraum. Mit den zusätzlichen Ungleichungen für eine solche Zuweisung ergibt sich damit für ein Schleifenprogramm mit r Strukturparametern n_{1}, …, n_{r} folgender Prototyp:

Θ(S, 0) = μ_{S}^{0} + μ_{S}^{1}n_{1} + … + μ_{S}^{r}n_{r}

Der einzige Indexpunkt soll mit null bezeichnet werden. Bei Zuweisungen mit ein- oder mehrdimensionalem Indexraum sind diese Ungleichungen meist nicht notwendig, da Strukturparameter in Schleifengrenzen häufig wie in

for i := 0 to n

auftreten und damit implizit größer oder gleich null gesetzt werden. Allgemein sollten zur Beschreibung des Indexraums sämtliche redundanten Ungleichungen eliminiert werden, um damit die Anzahl der Farkas-Multiplikatoren zu reduzieren.

Jeder Abhängigkeit mit Kante e im allgemeinen Abhängigkeitsgraphen des Eingabe-Programms läßt sich mittels der Ungleichung (2) eine Verzögerung Δ_{e} zuordnen:

(∀ i⃗, j⃗ : (i⃗, j⃗) ∈ R_{e} : Δ_{e}(i⃗, j⃗) = Θ(δ(e), j⃗) − Θ(σ(e), i⃗) − 1)

Nutzt man zur Darstellung der Abhängigkeiten die anfangs eingeführte h-Transformation und den Polyeder P_{e}, ergibt sich für die Verzögerung:

(∀ j⃗ : j⃗ ∈ P_{e} : Δ_{e}(j⃗) = Θ(δ(e), j⃗) − Θ(σ(e), h(j⃗)) − 1)

Aufgrund der Definition eines Schedules ist Δ_{e} stets positiv und wegen der Affinität des Schedules und der h-Transformation selbst affin. Somit erfüllt Δ_{e} die Anforderungen zur Anwendung des obigen Lemmas, wodurch der Zusammenhang von Schedule und Abhängigkeiten hergestellt wird. Wieder existieren Farkas-Multiplikatoren λ_{e}^{k}, sodaß gilt:

(4)  (∀ j⃗ : j⃗ ∈ P_{e} : Θ(δ(e), j⃗) − Θ(σ(e), h(j⃗)) − 1 = ; λ_{e}^{0} + (∑ k : 1≤k ≤m_{e} : λ_{e}^{k}(c_{e}^{k} ( j⃗ ; n⃗ ) + d_{e}^{k})))

Benutzt man zur Darstellung der Schedules in Gleichung (4) die oben definierten Prototypen und führt nach Einsetzen der Prototypen einen Koeffizientenvergleich für alle Schleifenindizes, Strukturparameter und den konstanten Anteil durch, erhält man für jede Abhängigkeit e ein Gleichungssystem in den Unbekannten μ_{S}^{k} und λ_{e}^{k}. Dieses Gleichungssystem ist erfüllt, falls die eingesetzten Prototypen Schedules sind. Somit schränkt jedes solche Gleichungssystem die Werte der Farkas-Multiplikatoren in den Prototypen derart ein, daß die Prototypen bezüglich e die Scheduleeigenschaft besitzen.

Anmerkung 4.2.3 Bei uniformen Abhängigkeiten ist die Anwendung von Farkas’ Lemma auf die Verzögerung nicht notwendig, wenn diese Abhängigkeit zwischen derselben Zuweisung auftritt. Das liegt daran, daß die zugehörige Verzögerung stets konstant ist und nicht von den Indexvektoren abhängig ist. Hier genügt es also, nur für den konstanten Anteil einen Koeffizientenvergleich durchzuführen. Auf der rechten Seite der erzeugten Gleichung steht lediglich der Farkas-Multiplikator λ_{e}^{0} als Repräsentant für die konstante Verzögerung.

Beispiel 4.2.4 Wir betrachten folgendes Schleifenprogramm:

for i := 0 to 5 do

for j := i to i + 9 do

a[i, j] := a[i − 1, j]+a[i − 2, j + 2] {S_{1}} end b[i] := a[i, i]; {S_{2}}

end

Für die Indexräume D_{S_{1}} und D_{S_{2}} gelten die Ungleichungen:

D_{S_{1}} = {(i, j) : 0≤i≤5 ∧ i≤j ≤i + 9}
D_{S_{2}} = {i : 0≤i≤5}

Es existieren drei Abhängigkeiten: e_{1} = S_{1}δS_{1}, e_{2} = S_{1}δS_{1}, die jeweils uniform sind, und e_{3} = S_{1}δS_{2} mit den folgenden h-Transformationen und Polyedern P_{e_{i}}:

\begin{matrix}e_{1} : & h_{e_{1}}(i, j) = (i − 1, j) & P_{e_{1}} = \{(i, j) : 1≤i≤5 ∧ i≤j ≤i + 8\} \\ e_{2} : & h_{e_{2}}(i, j) = (i − 2, j + 2) & P_{e_{2}} = \{(i, j) : 2≤i≤5 ∧ i≤j ≤i + 5\} \\ e_{3} : & h_{e_{3}}(i) = (i, i) & P_{e_{3}} = \{i : 0≤i≤5\} \end{matrix}

In Abbildung 6 ist der Indexraum beider Zuweisungen jeweils mit den verschiedenen Abhängigkeiten dargestellt.

Wir stellen zunächst die Prototypen der Zuweisungen gemäß Gleichung (3) auf:

Θ_{S_{1}}(i, j) = μ_{1}^{0} + μ_{1}^{1}i + μ_{1}^{2}(5 − i) + μ_{1}^{3}(j − i) + μ_{1}^{4}(i − j + 9)

Θ_{S_{2}}(i) = μ_{2}^{0} + μ_{2}^{1}i + μ_{2}^{2}(5 − i)

Abbildung 6: Indexraum und Abhängigkeiten aus Beispiel 4.2.4.

Mit Anmerkung 4.2.3 und der Uniformität der Abhängigkeiten e_{1} und e_{2} vereinfacht sich die Konstruktion der Gleichungssysteme. Für diese beiden Abhängigkeiten erhält man direkt die Gleichungen

\begin{matrix}e_{1} : & μ_{1}^{1} − μ_{1}^{2} − μ_{1}^{3} + μ_{1}^{4} − 1 = λ_{1}^{0} \\ e_{2} : & 2μ_{1}^{1} − 2μ_{1}^{2} − 4μ_{1}^{3} + 4μ_{1}^{4} − 1 = λ_{2}^{0} \end{matrix}

Für die dritte Abhängigkeit gilt aufgrund von Gleichung (4):

μ_{2}^{0} + μ_{2}^{1}i + μ_{2}^{2}(5 − i)
− (μ_{1}^{0} + μ_{1}^{1}i + μ_{1}^{2}(5 − i) + μ_{1}^{3}(i − i) + μ_{1}^{4}(i − i + 9)) − 1
= λ_{3}^{0} + λ_{3}^{1}i + λ_{3}^{2}(5 − i)

Der Koeffizientenvergleich für i,j und den konstanten Anteil liefert folgendes Gleichungssystem für e_{3}:

μ_{2}^{1} − μ_{2}^{2} − μ_{1}^{1} + μ_{1}^{2} = λ_{3}^{1} − λ_{3}^{2}
0 = 0
μ_{2}^{0} + 5μ_{2}^{2} − μ_{1}^{0} − 5μ_{1}^{2} − 9μ_{1}^{4} − 1 = λ_{3}^{0} + 5λ_{3}^{2}

Die Lösungen der Vereinigung dieser Gleichungssysteme, die je Abhängigkeit erzeugt wurden, ergeben damit eine Menge von Werten für die Farkas-Multiplikatoren. Setzt man eine solche Lösung in die Prototypen ein, erhält man einen affinen, nicht negativen Schedule für das Eingabe-Programm. Die Menge von Schedules, die durch das Einsetzen aller Lösungen entsteht, bezeichnen wir mit G.

Wie bei der Hyperebenen-Methode sind wir auch hier an laufzeitoptimalen Schedules interessiert. Wir suchen also einen ”minimalen“ Schedule in G, wobei der Begriff der Minimalität noch festzulegen ist. Um diesen Prozeß zu vereinfachen, formen wir das erhaltene Gleichungssystem weiter um. Nach einer Heuristik [4] existieren in dem hergeleiteten Gleichungssystem etwa zweimal soviele Unbekannte wie Gleichungen. Somit bildet der Lösungsraum, falls er existiert, eine sehr große Menge, aus welcher der minimale Schedule zu wählen ist. Es ist jedoch möglich die Anzahl von Unbekannten zu reduzieren, ohne den Lösungsraum zu verändern. Die Grundlage dafür bildet die Nicht-Negativitäts-Bedingung für Farkas-Multiplikatoren. Wir führen dazu eine Abwandlung der Elimination nach Gauss-Jordan[7] durch. Wie im bekannten Verfahren löst man nach einer Unbekannten x, hier den Farkas-Multiplikatoren, auf und ersetzt diese durch den gewonnenen Ausdruck e in den übrigen Gleichungen. Da die Farkas-Multiplikatoren nicht negativ sind, kann man nun wegen der Gleichheit x = e eine Ungleichung der Form e ≥ 0 erzeugen. Im einzelnen:

1. Wähle eine Gleichung des Systems und daraus eine Unbekannte x mit Koeffizient ungleich null. Löse die Gleichung nach x auf. Es ergibt sich eine Gleichung der Form x = e.

2. Da x ein Farkas-Multiplikator ist, gilt x≥0. Erzeuge also eine Ungleichung e≥0 und streiche die ursprüngliche Gleichung aus dem System.

3. Ersetze x durch e in allen noch existierenden Gleichungen und in allen bereits erzeugten Ungleichungen.

4. Falls noch Gleichungen vorhanden sind, setze mit 1. fort.

Durch dieses Verfahren wird jede Gleichung durch eine Ungleichung ersetzt, wobei sich bei jeder Ersetzung die Zahl von Unbekannten um eins reduziert.

Anmerkung 4.2.5 Bei der Auswahl von Farkas-Multiplikatoren im ersten Schritt sollten bevorzugt die Farkas-Multiplikatoren λ_{e}^{i} der Verzögerung eliminiert werden. Zur späteren Konstruktion der Schedules aus den Prototypen sind lediglich die Werte der Farkas-Multiplikatoren μ_{S}^{i} von Interesse. Werden jedoch letztere eliminiert, so müssen diese später wieder aus den Werten der restlichen μ_{S}^{i}’s und λ_{e}^{i}’s berechnet werden. In einer Implementierung des Verfahrens erfordert dies also zusätzlichen Aufwand.

Es kann jedoch vorkommen, daß durch den Koeffizientenvergleich Gleichungen entstehen, die keine λ_{e}^{i}’s enthalten. In diesem Fall ergäbe sich für die Gleichung die Form e = 0, wobei e ein Ausdruck in den μ_{S}^{i}’s ist. Will man auch hier den Mehraufwand einsparen, so formt man aus dieser Gleichung die beiden Ungleichungen

e≥0 e≤0

und fügt diese ebenso in das Ungleichungssystem ein. Nachteilig ist dabei jedoch, daß damit die Anzahl von Ungleichungen steigt.

Das erzeugte Ungleichungssystem besitzt denselben Lösungsraum wie das vorherige Gleichungssystem, abgesehen von den eliminierten Unbekannten. Im folgenden sei es durch

G ( μ⃗ ; λ⃗ )≥h⃗ (5)

notiert. G bezeichne dabei die Koeffizientenmatrix, μ⃗ den Vektor der Farkas-Multiplikatoren aller Prototypen, λ⃗ den Vektor der Farkas-Multiplikatoren der Verzögerungen Δ_{e} aller Abhängigkeiten und schließlich h⃗ den Vektor der Konstanten in jeder Ungleichung.

Beispiel 4.2.6 Mit dem beschriebenen Verfahren zum Erzeugen der Ungleichungen erhält man für die Vereinigung der Gleichungssysteme in Beispiel 4.2.4 etwa (in Abhängigkeit der Wahl der zu eliminierenden Farkas-Multiplikatoren) folgendes Ungleichungssystem:

μ_{1}^{1} − μ_{1}^{2} − μ_{1}^{3} + μ_{1}^{4} − 1 ≥ 0
2μ_{1}^{1} − 2μ_{1}^{2} − 4μ_{1}^{3} + 4μ_{1}^{4} − 1 ≥ 0
μ_{1}^{1} − μ_{1}^{2} − μ_{2}^{1} + μ_{2}^{2} + λ_{3}^{1} ≥ 0
−μ_{1}^{0} − 5μ_{1}^{1} − 9μ_{1}^{4} + μ_{2}^{0} + 5μ_{2}^{1} − 5λ_{3}^{1} − 1 ≥ 0

Um die Matrix G und den Konstantenvektor h⃗ zu erhalten, stellen wir das hergeleitete Unlgeichungssystem in der folgenden Form dar:

\underbrace{( \begin{matrix}0 & 1 & −1 & −1 & 1 & 0 & 0 & 0 & 0 \\ 0 & 2 & −2 & −4 & 4 & 0 & 0 & 0 & 0 \\ 0 & 1 & −1 & 0 & 0 & 0 & −1 & 1 & 1 \\ −1 & −5 & 0 & 0 & −9 & 1 & 5 & 0 & −5 \end{matrix} )}_{G} ( \begin{matrix}μ_{1}^{0} \\ μ_{1}^{1} \\ μ_{1}^{2} \\ μ_{1}^{3} \\ μ_{1}^{4} \\ μ_{2}^{0} \\ μ_{2}^{1} \\ μ_{2}^{2} \\ λ_{3}^{1} \end{matrix} )≥\underbrace{( \begin{matrix}1 \\ 1 \\ 0 \\ 1 \end{matrix} )}_{h⃗}

Anmerkung 4.2.7 Ist das Ungleichungssystem (5) bei Einhaltung der Nicht-Negativitäts-Bedingung für die Farkas-Multiplikatoren nicht lösbar, so existiert kein affiner eindimensionaler Schedule für das gegebene Eingabe-Pro-gramm. Das Programm kann also nicht in linearer Zeit parallel ausgeführt werden. In diesem Fall weicht man auf mehrdimensionale Schedules aus. Eine Methode zu deren Berechnung ist in [2] beschrieben.

Im folgenden soll der Begriff der Minimalität festgelegt werden. Wir führen dazu eine Ordnung auf Schedules ein. Die naheliegende punktweise Ordnung wird wie folgt festgelegt:

Definition 4.2.8 (Ordnung auf G)

(∀ Θ_{1}, Θ_{2} : Θ_{1}, Θ_{2} ∈ G : Θ_{1} ≺_{G} Θ_{2} ⇔ (∀ S : S ∈ Ω_{L} :
(∀ i⃗ : i⃗ ∈ D_{S} : Θ_{1}(S, i⃗)≤Θ_{2}(S, i⃗))))

Es stellt sich jedoch als Problem heraus, daß diese Ordnung eine partielle Ordnung ist, was zur Folge hat, daß in ihr i. a. kein Minimum existiert. Einen Ausweg aus dieser Situation erreichen wir mit folgendem Satz, der eine Aussage über das punktweise Minimum zweier Schedules macht:

Satz 4.2.9 (Abgeschlossenheit der Scheduleeigenschaft) Sind Θ_{1} und Θ_{2} Schedules eines Schleifenprogramms L, so gilt dies auch für die Abbildung Θ_{3}, definiert durch:

(∀ S : S ∈ Ω_{L} : (∀ i⃗ : i⃗ ∈ D_{S} : Θ_{3}(S, i⃗) = min(Θ_{1}(S, i⃗), Θ_{2}(S, i⃗))))

Beweis:

(∀ S, T : S, T ∈ Ω_{L} : (∀ i⃗, j⃗ : i⃗ ∈ D_{S} ∧ j⃗ ∈ D_{T} ∧ (S, i⃗)δ(T, j⃗) :

Θ_{1}(S, i⃗) + 1≤Θ_{1}(T, j⃗) ∧ Θ_{2}(S, i⃗) + 1≤Θ_{2}(T, j⃗))) ⇒ { Anwendung des min-Operators und Monotonie des min-Operators }

(∀ S, T : S, T ∈ Ω_{L} : (∀ i⃗, j⃗ : i⃗ ∈ D_{S} ∧ j⃗ ∈ D_{T} ∧ (S, i⃗)δ(T, j⃗) :

min(Θ_{1}(S, i⃗) + 1, Θ_{2}(S, i⃗) + 1)≤min(Θ_{1}(T, j⃗), Θ_{2}(T, j⃗))) ⇒ { min(a + 1, b + 1) = min(a, b) + 1 }

(∀ S, T : S, T ∈ Ω_{L} : (∀ i⃗, j⃗ : i⃗ ∈ D_{S} ∧ j⃗ ∈ D_{T} ∧ (S, i⃗)δ(T, j⃗) : min(Θ_{1}(S, i⃗), Θ_{2}(S, i⃗)) + 1≤min(Θ_{1}(T, j⃗), Θ_{2}(T, j⃗)))

Dieser Satz stellt die mathematische Grundlage des Verfahrens dar. Es ist also möglich, den minimalen Schedule aus affinen Schedules in G zu konstruieren, wobei der dadurch entstehende Schedule nicht notwendig aus G ist, sondern eine stückweise affine Kombination derartiger Schedules ist. Die Auswahl der Kombination wird durch den Indexraum bestimmt, da durch den min-Operator jeweils der in bestimmten Bereichen des Indexraums zeitoptimale Schedule gewählt wird. Als Folgerung dieses Satzes ergibt sich der laufzeitoptimale Schedule Θ durch

(∀ S : S ∈ Ω_{L} : (∀ i⃗ : i⃗ ∈ D_{S} : Θ(S, i⃗) = (min t : t ∈ G : t(S, i⃗))))

Damit sind wir in der Lage das Optimierungsproblem für diese Methode zu formulieren:

Optimierungsproblem der Methode von Feautrier

(6)  \begin{aligned}Θ(S, i⃗) & = & min( μ_{S}^{0} + (∑ k : 1≤k ≤m_{S} : μ_{S}^{k}(a_{S}^{k} ( i⃗ ; n⃗ ) + b_{S}^{k})) : \\ & & G ( μ⃗ ; λ⃗ )≥h⃗ ∧ μ⃗≥0 ∧ λ⃗≥0 ) \end{aligned}

(Zur Erinnerung: m_{S} bezeichnet die Anzahl von Ungleichungen in der Beschreibung des Indexraums der Zuweisung S.) Die Lösung dieses Problems für eine feste Zuweisung liefert somit den entsprechenden zeitoptimalen Schedule der Zuweisung.

Beispiel 4.2.10 Für das Beispiel 4.2.4 ergeben sich folgende Optimierungsprobleme. Der Schedule der ersten Zuweisung ist:

Θ(S_{1}, (i, j)) = min( μ_{1}^{0} + μ_{1}^{1}i + μ_{1}^{2}(5 − i) + μ_{1}^{3}(j − i) + μ_{1}^{4}(i − j + 9) :

G ( μ⃗ ; λ⃗ )≥h⃗ ∧ μ⃗≥0 ∧ λ⃗≥0 )

Für die zweite Zuweisung ergibt sich:

Θ(S_{2}, i) = min( μ_{2}^{0} + μ_{2}^{1}i + μ_{2}^{2}(5 − i) : G ( μ⃗ ; λ⃗ )≥h⃗ ∧ μ⃗≥0 ∧ λ⃗≥0 )

wobei G ( μ⃗ ; λ⃗ )≥h⃗ folgendes Ungleichungssystem ist:

( \begin{matrix}0 & 1 & −1 & −1 & 1 & 0 & 0 & 0 & 0 \\ 0 & 2 & −2 & −4 & 4 & 0 & 0 & 0 & 0 \\ 0 & 1 & −1 & 0 & 0 & 0 & −1 & 1 & 1 \\ −1 & −5 & 0 & 0 & −9 & 1 & 5 & 0 & −5 \end{matrix} ) ( \begin{matrix}μ_{1}^{0} \\ μ_{1}^{1} \\ μ_{1}^{2} \\ μ_{1}^{3} \\ μ_{1}^{4} \\ μ_{2}^{0} \\ μ_{2}^{1} \\ μ_{2}^{2} \\ λ_{3}^{1} \end{matrix} )≥( \begin{matrix}1 \\ 1 \\ 0 \\ 1 \end{matrix} )

Eine implementierbare Lösung dieser Optimierungsprobleme mit Hilfe eines Programms zur Berechnung des lexikographischen Minimums wird in Abschnitt 4.2.3 vorgestellt.

Anmerkung 4.2.11 Der berechnete Schedule Θ ist i. a. eine rationale Funktion. Um den Schedule in eine Funktion Θ′ mit Bild in ℕ umzuwandeln, ist es ausreichend, rationale Bildpunkte auf die nächste ganze Zahl aufzurunden:

Satz 4.2.12 (Runden von Schedules) Ist Θ : (Ω_{L}, D_{S∈Ω_{L}}) → ℚ ein Schedule für ein Schleifenprogramm L so gilt dies auch für die Abbildung Θ′ : (Ω_{L}, D_{S∈Ω_{L}}) → ℕ, die durch Aufrunden rationaler Bildpunkte von Θ konstruiert wird:

Θ′(i⃗) = ⌈Θ(i⃗)⌉

Beweis: Wir zeigen dies durch einen Widerspruchsbeweis. Wir nehmen also an, Θ′ besitze nicht die Scheduleeigenschaft. Für zwei abhängige Instanzen (S, i⃗) und (T, j⃗) seien dazu

α = Θ′(S, i⃗) − Θ(S, i⃗)
β = Θ′(T, j⃗) − Θ(T, j⃗)

die rationalen Differenzen der beiden Abbildungen für diese Instanzen.

(∃ S, T : S, T ∈ Ω_{L} : (∃ i⃗, j⃗ : i⃗ ∈ D_{S} ∧ j⃗ ∈ D_{T} : Θ′(T, j⃗) − Θ′(S, i⃗) < 1)) ⇒ { Θ′ ist eine Abbildung mit Bild in ℕ }

(∃ S, T : : (∃ i⃗, j⃗ : : Θ′(T, j⃗) − Θ′(S, i⃗) = 0)) ⇔ { Festlegung von α und β }

(∃ S, T : : (∃ i⃗, j⃗ : : Θ(T, j⃗) + β − Θ(S, i⃗) − α = 0)) ⇔ { Arithmetische Umformung }

(∃ S, T : : (∃ i⃗, j⃗ : : Θ(T, j⃗) − Θ(S, i⃗) = α − β)) ⇒ { 0 ≤ α < 1 }

(∃ S, T : : (∃ i⃗, j⃗ : : Θ(T, j⃗) − Θ(S, i⃗) < 1 − β)) ⇒ { 0 ≤ β < 1 }

(∃ S, T : : (∃ i⃗, j⃗ : : Θ(T, j⃗) − Θ(S, i⃗) < 1))

Dies führt also zu einem Widerspruch zur Scheduleeigenschaft von Θ.

4.2.2 Bewertung der Methode

Die berechneten Schedules sind jeweils zeitoptimal in der Menge der stückweise affinen Abbildungen, die aus Schedules in G zusammengesetzt sind. Zur Bewertung dieser Schedules ist ein Vergleich mit dem sogenannten freien Schedule notwendig. Der freie Schedule ist der Schedule, der jeder Instanz unter Berücksichtigung der Abhängigkeiten den frühest möglichen Zeitpunkt der Ausführung zuordnet.

Wir stellen zunächst fest, daß der von diesem Verfahren berechnete Schedule nicht immer dem freien Schedule entspricht. Wir benötigen dazu einige Begriffe der linearen Algebra [7]:

Definition 4.2.13 (konvexe Menge) Sei V ein ℝ-Vektorraum. Für die Vektoren a⃗, b⃗ ∈ V heißt die Menge

S_{a⃗,b⃗} := { λa⃗ + (1 − λ)b⃗ : λ ∈ [0, 1] }

Strecke von a⃗ bis b⃗. Eine Menge A ⊆ V heißt konvex, wenn gilt

(∀ a⃗, b⃗ : a⃗, b⃗ ∈ A : S_{a⃗,b⃗} ⊆ A)

Da ein Indexraum eine konvexe Menge ist, bildet ein Schedule also Elemente einer konvexen Menge nach ℚ ab. Folgende Definition ermöglicht uns die Charakterisierung von Schedules, die mit der Methode von Feautrier berechnet wurden.

Definition 4.2.14 (konkave Abbildung) Sei A eine konvexe Teilmenge eines ℝ-Vektorraums. Eine Abbildung f : A → ℝ heißt konvex, wenn gilt:

(∀ a⃗, b⃗ : a⃗, b⃗ ∈ A : (∀ λ : λ ∈ [0, 1] : f(λa⃗+(1−λ)b⃗)≤λf(a⃗)+(1−λ)f(b⃗)))

f ist konkav, wenn −f konvex ist.

Geometrisch betrachtet besitzt eine konkave Abbildung die Eigenschaft, in keiner Komponente eine Linkskrümmung aufzuweisen, also nie in der Steigung zuzunehmen.

Schedules, die mit der Methode von Feautrier berechnet wurden, besitzen die Eigenschaft, konkave Funktionen zu sein. Das liegt daran, daß die Menge G affine Abbildungen enthält, die trivialerweise konkav und konvex sind. Das Bilden des punktweisen Minimums affiner Abbildungen ergibt eine konkave Funktion, wie der nächste Satz zeigt.

Satz 4.2.15 Sind f und g affine Abbildungen von einem ℝ-Vektorraum V nach ℝ, so ist die Abbildung

h(x⃗) = min(f(x⃗), g(x⃗))

für x⃗ ∈ V konkav.

Beweis: Für den Beweis benötigen wir folgende Eigenschaft des Minimums, die sich durch eine Fallunterscheidung leicht herleiten läßt :

(7) min(a + b, c + d) ≥ min(a, c) + min(b, d)

Seien also f und g affine Abbildungen:

(∀ a⃗, b⃗ : a⃗, b⃗ ∈ V : (∀ λ : λ ∈ [0, 1] : h(λa⃗ + (1 − λ)b⃗) =

min(f(λa⃗ + (1 − λ)b⃗), g(λa⃗ + (1 − λ)b⃗))))

⇔ { f und g sind affine Abbildungen }

(∀ a⃗, b⃗ : a⃗, b⃗ ∈ V : (∀ λ : λ ∈ [0, 1] : h(λa⃗ + (1 − λ)b⃗) =

min(λf(a⃗) + (1 − λ)f(b⃗), λg(a⃗) + (1 − λ)g(b⃗))))

⇒ { Ungleichung (7) }

(∀ a⃗, b⃗ : a⃗, b⃗ ∈ V : (∀ λ : λ ∈ [0, 1] : h(λa⃗ + (1 − λ)b⃗) ≥ min(λf(a⃗), λg(a⃗)) + min((1 − λ)f(b⃗), (1 − λ)g(b⃗))))

⇔ { Homogenität des min-Operators } (∀ a⃗, b⃗ : a⃗, b⃗ ∈ V : (∀ λ : λ ∈ [0, 1] : h(λa⃗ + (1 − λ)b⃗) ≥

λ min(f(a⃗), g(a⃗)) + (1 − λ) min(f(b⃗), g(b⃗))))

⇔ { Definition von h }

(∀ a⃗, b⃗ : a⃗, b⃗ ∈ V : (∀ λ : λ ∈ [0, 1] : h(λa⃗ + (1 − λ)b⃗) ≥

λh(a⃗) + (1 − λ)h(b⃗)))

⇔ { Multiplikation mit −1 }

(∀ a⃗, b⃗ : a⃗, b⃗ ∈ V : (∀ λ : λ ∈ [0, 1] : (−h)(λa⃗ + (1 − λ)b⃗) ≤

λ(−h)(a⃗) + (1 − λ)(−h)(b⃗)))

⇔ { Definition 4.2.14 } h ist konkav

Somit wird von diesem Verfahren der sogenannte beste konkave Schedule berechnet. Es gibt jedoch Programme, deren freier Schedule nicht konkav ist:

Beispiel 4.2.16 Gegeben sei folgendes Schleifenprogramm:

for i := 0 to 2n do a[i] := a[2n − i];

end

Der freie Schedule dafür ist:

t(i) = if i > n then 1 else 0

Diese Abbildung ist nicht konkav. Somit kann ein Schedule nach der Methode von Feautrier für dieses Programm nicht zeitoptimal sein. Der von dieser Methode gelieferte Schedule ist

Θ(i) = ⌈\frac{i}{2}⌉

Für den Fall, daß im Eingabe-Programm nur uniforme Abhängigkeiten enthalten sind, läßt sich jedoch zeigen, daß die Differenz zwischen der Laufzeit des freien Schedules und der des besten konkaven Schedules eine von den Strukturparametern meistens unabhängige Konstante ist. Nach einer Heuristik gilt diese Aussage^{3} für den besten affinen Schedule und den freien Schedule [9]. Da die Laufzeit des besten konkaven Schedule aufgrund seiner Konstruktion kleiner oder gleich der des besten affinen Schedules ist, folgt dies also auch für die Schedules dieses Verfahrens. Auf nicht uniforme Abhängigkeiten läßt sich diese Aussage allerdings nicht erweitern, wie Beispiel 4.2.16 zeigt, da der beste konkave Schedule hier die Laufzeit O(n) besitzt.

Eine weiteres Problem stellen Programme dar, für die kein affiner konkaver Schedule existiert. Für derartige Programme ist die Berechnung eines multidimensionalen Schedules erforderlich [2].

Abgesehen von diesen Einschränkungen gibt diese Methode jedoch eine schnelle und effiziente Möglichkeit zur Berechnung von Schedules.

Aus theoretischer Sicht ist damit die Berechnung eines Schedules für ein Eingabe-Programm durch die Lösung des genannten Optimierungsproblems möglich.

In der Praxis stellt sich jedoch das Problem, daß derzeit kein Verfahren in der Lage ist, dieses Optimierungsproblem zu lösen. Die Schwierigkeit liegt darin, daß in der Zielfunktion Parameter auftreten, nämlich die Vektoren i⃗ und n⃗. Somit scheidet beispielsweise der klassische Simplex-Algorithmus [14] aus.

PIP (Parametric Integer Programming) ist ein Programm, das unter bestimmten Voraussetzungen auch mit Parametern arbeiten kann [15]. Eine Einschränkung betrifft den Begriff der Affinität. Da ein Strukturparameter k im mathematischen Sinn eine Konstante darstellt, ist ein Ausdruck der Form k i mit i als Variable nach üblicher Definition affin. Dies gilt allerdings nicht bei PIP. Hier dürfen Strukturparameter lediglich additiv mit Variablen verknüpft werden. Im obigen Optimierungsproblem (6) sind die Farkas-Multiplikatoren die Variablen, i⃗ und n⃗ die Parameter der Zielfunktion. Diese sind jedoch multiplikativ verknüpft. Eine zweite Schwierigkeit besteht darin, daß PIP nicht mit Zielfunktionen arbeitet, sondern lediglich in der Lage ist, für ein gegebenes parametrisiertes Ungleichungssystem das lexikografische Minimum aller Lösungen zu berechnen. Somit besitzt das angegebene Optimierungsproblem (6) keine dieser beiden erforderlichen Eigenschaften.

3In [13] wird von einem Test mit 25 Algorithmen mit uniformen Abhängigkeiten berichtet. Dabei ergab sich für 23 dieser Algorithmen eine identische Laufzeit für freien Schedule und affinen Schedule. Bei zwei Algorithmen betrug die Laufzeit für einen affinen Schedule eine Zeiteinheit mehr als für den freien Schedule. Diese Größen blieben auch bei einer Änderung der Strukturparameter gleich.

Eine Notlösung für dieses Problem ist es, die angegebene Zielfunktion zu verwerfen, und nach dem lexikografischen Minimum der Lösungen von

G ( μ⃗ ; λ⃗ )≥h⃗

zu suchen:

Beispiel 4.2.17 Für das lexikografische Minimum des Gleichungssystems in Beispiel 4.2.10 ergibt sich für den Vektor μ⃗:

(μ_{1}^{0}, …, μ_{1}^{4}, μ_{2}^{0}, μ_{2}^{1}, μ_{2}^{2}, λ_{3}^{1}) = (0, 0, 0, 0, 1, 0, 2, 2, 0)

Dies entspricht nach dem Einsetzen der Farkas-Multiplikatoren in die Prototypen folgendem Schedule:

Θ_{S_{1}}(i, j) = i − j + 9
Θ_{S_{2}}(i) = 10

Unnötigerweise werden jedoch bei diesem Schedule alle Instanzen der Zuweisung S_{2} erst ausgeführt, nachdem alle Instanzen von Zuweisung S_{1} ausgeführt sind. Auch in Zuweisung S_{1} könnten einige Instanzen früher ausgeführt werden. Wir werden sehen, daß die Lösung des Optimierungsproblems zu einem anderen Resultat führt.

Problematisch ist bei dieser Notlösung, daß das Ergebnis von der bei der Konstruktion des Gleichungssystems gewählten Reihenfolge der Abhängigkeiten und Zuweisungen abhängt. Außerdem erhält man im Idealfall lediglich den besten affinen und nicht den besten stückweise affinen Schedule.

Eine weitere Möglichkeit besteht darin, für die Parameter i⃗ und n⃗ des Problems konkrete Werte einzusetzen, um somit ein nicht-parametrisiertes Problem zu erhalten, das mit einem Verfahren wie Simplex gelöst werden kann. Bei Zuweisungen mit großem Indexraum ist diese Lösung aber nicht tragbar. Außerdem liefert diese Möglichkeit lediglich die Werte des Schedules für einen isolierten Indexpunkt, jedoch keinen geschlossenen Ausdruck in Form einer mathematischen Abbildung.

Wir werden aber im nächsten Abschnitt sehen, daß durch einer Umwandlung des Optimierungsproblems diese Schwierigkeiten beseitigt werden können.

4.2.3 Anpassung an PIP

Die Grundlage für eine Umwandlung des Optimierungsproblems dieser Methode liefert folgendes aus der linearen Programmierung bekannte Theorem [7]:

Satz 4.2.18 (Dualitätstheorem) Für ein Minimierungsproblem mit nichtleerer Lösungsmenge gilt:

min( k⃗μ⃗ : μ⃗≥0 ∧ Gμ⃗≥h⃗) = max( ν⃗h⃗ : ν⃗ ≥0 ∧ ν⃗G≤k⃗)

Anmerkung 4.2.19 Die Voraussetzung zur Anwendung des Dualitätstheorem ist die Lösbarkeit des Ungleichungssystems. In einer Implementierung dieses Verfahrens muß also vor der Umformung des Problems das Ungleichungssystem Gμ⃗ ≥ h⃗ auf seine Lösbarkeit getestet werden. Da in diesem Ungleichungssystem keine Parameter auftreten, muß dazu nicht unbedingt PIP benutzt werden. Diese Aufgabe kann auch von einem Algorithmus (wie Simplex mit beliebiger Zielfunktion) gelöst werden. Wir gehen im folgenden, davon aus, daß das Ungleichungssystem Gμ⃗ ≥ h⃗ lösbar ist, und folglich ein Minimum dieser Lösungen existiert. Als Konsequenz des Satzes gilt dann auch die Existenz des Maximums.

Die Anwendung des Dualitätstheorems ergibt damit folgendes umgeformtes Optimierungsproblem:

Θ(S, i⃗) = max({ ν⃗h⃗ : ν⃗ ≥0 ∧

ν⃗G≤(\overbrace{0, …, 0}^{N}, 1, a_{S}^{1} ( i⃗ ; n⃗ ) + b_{S}^{1}, …, a_{S}^{m_{S}} ( i⃗ ; n⃗ ) + b_{S}^{m_{S}}, \overbrace{0, …, 0}^{N′}) \})

wobei für N und N′ gilt:

N = (∑ R : R ∈ Ω_{L} ∧ R <_{L} S : m_{R} + 1),
N′ = (∑ T : T ∈ Ω_{L} ∧ S <_{L} T : m_{T} + 1) + |λ⃗|

In der Zielfunktion treten nun keine Parameter mehr auf, da die Einträge des Vektors h⃗ stets Konstanten sind. Lediglich in den Bedingungen für ν⃗ existieren die Parameter i⃗ und n⃗, doch hier nur in additiver Form. Damit ist die Problemstellung im Sinne von PIP affin.

Anmerkung 4.2.20 Bei der Überführung des aus Gleichung (4) erzeugten Gleichungssystems in das Unlgeichungssystem (5) werden Farkas-Multiplika-toren eliminiert. Die obige Formulierung des Optimierungsproblems gibt aber für alle Farkas-Multiplikatoren Koeffizienten an, weil die Wahl der zu eliminierenden Multiplikatoren auf die Anzahl der insgesamt eliminierten Multiplikatoren Einfluß nimmt. Diese Koeffizienten sind in den entsprechenden Vektoren und in G mit null besetzt.

Wie zuvor erwähnt, berechnet PIP für ein parametrisiertes Ungleichungssystem das lexikografische Minimum aller Lösungen. Somit stehen wir vor der Aufgabe, das durch die Anwendung des Dualitätstheorems erhaltene Maximierungsproblem mit Hilfe des lexikografischen Minimums zu lösen.

Zunächst soll gezeigt werden, wie die Zielfunktion Θ in PIP dargestellt werden kann. Hier benutzen wir eine Eigenschaft des lexikografischen Maximums. Wir führen dazu eine neue Variable m in das Ungleichungssystem ein, die durch die beiden Ungleichungen

m≥ν⃗h⃗
m≤ν⃗h⃗

mit der Zielfunktion identisch ist. Ordnet man diese Variable, die in allen anderen Ungleichungen den Koeffizienten 0 hat, als erste Variable in das System ein, so wird bei der Berechnung des lexikografischen Maximums gerade das Maximum von m und damit von Θ berechnet.

Es bleibt zu klären, wie mit PIP das lexikografische Maximum berechnet wird. Dabei hilft uns folgende Eigenschaft der lexikografischen Ordnung. Ist die Existenz eines lexikografischen Maximums gegeben, so läßt es sich durch das lexikografische Minimum für einen Polyeder P = {i⃗ : Ai⃗ ≥ b⃗} aufgrund folgender Beziehungen darstellen:

lexmax( i⃗ : Ai⃗≥b⃗) = −lexmin( −i⃗ : Ai⃗≥b⃗)
= −lexmin( i⃗ : A(−i⃗)≥b⃗)

Geometrisch betrachtet führt man dadurch eine Spiegelung aller Punkte des Polyeders an den Koordinatenachsen durch. PIP ist aber auf Lösungsvektoren beschränkt, die in allen Komponenten positiv sind. Für einen Polyeder wie z. B. {(i, j) : 0 ≤ i ≤ n ∧ 0 ≤ j ≤ n} mit lexikographischen Maximum (n, n) würde PIP bei Anwendung der genannten Umformung also (0, 0) als Ergebnis liefern. Abhilfe schafft hier die Möglichkeit, daß in PIP einer der Parameter mit der Eigenschaft ausgestattet werden kann, beliebige positive Werte annehmen zu können. Wir führen dazu einen neuen Parameter p ein und führen zusätzlich zur Spiegelung eine Translation um p aus:

−lexmin( i⃗ : A(−i⃗)≥b⃗)
(8)  = lim_{p→∞}( p1 − lexmin( i⃗ : A(p1 − i⃗)≥b⃗))

Ist p groß genug, so ist damit der Lösungsvektor des lexikographischen Minimums, der nun von p abhängig ist, positiv.

Somit läßt sich mit PIP also auch das lexikographische Maximum für ein Ungleichungssystem berechnen, wobei auf das Resultat noch eine Rücktransformation, wie in (8) angegeben, angewendet werden muß.

Als Ergebnis der Anpassung an PIP ergibt sich insgesamt folgendes Ungleichungssystem als Eingabe:

(9) \begin{matrix}(p1 − ν⃗)G≤(\overbrace{0, …, 0}^{N}, 1, a_{S}^{1} ( i⃗ ; n⃗ ) + b_{S}^{1}, …, a_{S}^{m_{S}} ( i⃗ ; n⃗ ) + b_{S}^{m_{S}}, \overbrace{0, …, 0}^{N′}) \\ p − m≥(p1 − ν⃗)h⃗ \\ p − m≤(p1 − ν⃗)h⃗ \end{matrix}

mit der Variablenreihenfolge m, ν_{1}, …, ν_{k}, wobei k die Anzahl der Zeilen von G ist. Als Lösung ergibt sich somit p − m.

Beispiel 4.2.21 Die Transformation des Optimierungsproblems aus Beispiel 4.2.10 in eine für PIP zugängliche Form liefert mit (9) folgendes System von Ungleichungen als Eingabe für PIP:

Für Zuweisung S_{1}:

−ν_{1}^{4} + 1 + p ≥ 0
ν_{1}^{1} + 2ν_{1}^{2} + ν_{1}^{3} − 5ν_{1}^{4} + i + p ≥ 0
−ν_{1}^{1} − 2ν_{1}^{2} − ν_{1}^{3} − i + 5 + 4p ≥ 0
−ν_{1}^{1} − 4ν_{1}^{2} − i + j + 5p ≥ 0
ν_{1}^{1} + 4ν_{1}^{2} − 9ν_{1}^{4} + i − j + 9 + 4p ≥ 0
ν_{1}^{4} − p ≥ 0
−ν_{1}^{3} + 5ν_{1}^{4} − 4p ≥ 0
ν_{1}^{3} − p ≥ 0
ν_{1}^{3} − 5ν_{1}^{4} + 4p ≥ 0
−m_{1} + ν_{1}^{1} + ν_{1}^{2} + ν_{1}^{4} − 2p ≥ 0
m_{1} − ν_{1}^{1} − ν_{1}^{2} − ν_{1}^{4} + 2p ≥ 0

Für Zuweisung S_{2}:

−ν_{1}^{4} + 1 + p ≥ 0
ν_{1}^{1} + 2ν_{1}^{2} + ν_{1}^{3} − 5ν_{1}^{4} + p ≥ 0
−ν_{1}^{1} − 2ν_{1}^{2} − ν_{1}^{3} + 4p ≥ 0
−ν_{1}^{1} − 4ν_{1}^{2} + 5p ≥ 0
ν_{1}^{1} + 4ν_{1}^{2} − 9ν_{1}^{4} + 4p ≥ 0
ν_{1}^{4} + 1 − p ≥ 0
−ν_{1}^{3} + 5ν_{1}^{4} + i − 4p ≥ 0
ν_{1}^{3} + 5 − i − p ≥ 0
ν_{1}^{3} − 5ν_{1}^{4} + 4p ≥ 0
−m_{2} + ν_{1}^{1} + ν_{1}^{2} + ν_{1}^{4} − 2p ≥ 0
m_{2} − ν_{1}^{1} − ν_{1}^{2} − ν_{1}^{4} + 2p ≥ 0

[figure: 9; 0; Θ_{S_{2}}(i); i; j; 4; Θ_{S_{1}}(i, j); 5; 3; 2; 1; 5; 4; 5; 3; 2; 1; 6; 0]

Abbildung 7: Der Schedule aus Beispiel 4.2.21.

Als Lösung für p − m_{1} und p − m_{2} und damit für Θ_{S_{1}} und Θ_{S_{2}} ergibt sich schließlich:

Θ_{S_{1}}(i, j) = if j ≤9 then i else i − j + 9
Θ_{S_{2}}(i) = i + 1

In Abbildung 7 ist der Verlauf beider Abbildungen im Indexraum dargestellt. Der berechnete Schedule ist hier sogar der freie Schedule.

5 Berechnung von Allokationen

In diesem Abschnitt soll eine Methode zur Berechnung von Allokationen für Eingabe-Programme vorgestellt werden [3]. Im ersten Teil werden wir ein Verfahren kennenlernen, das in der Lage ist, für jedes solche Programm eine eindimensionale Allokation zu konstruieren, im zweiten Teil wird dieses Verfahren auf mehrdimensionale Allokationen erweitert.

5.1 Eindimensionale Allokation

Wie im ersten Abschnitt erwähnt, stellt eine Allokation eine Abbildung dar, die jeder Instanz einer Zuweisung einen virtuellen Prozessor zur Abarbeitung zuweist. Die Berechnung einer Allokation ist damit einerseits abhängig vom verwendeten Prozessorlayout, andererseits vom Speichertyp des Parallelrechners. Wir machen deshalb für das vorzustellende Verfahren folgende Annahmen:

• Jedem Prozessor steht ein eigener Speicher zur Verfügung; es gibt keinen gemeinsamen Speicher. Es gilt die owner computes rule, die besagt, daß die Ergebnisse einer Berechnung stets im Speicher des Prozessors abgelegt werden, der die Berechnung durchgeführt hat. Werden zur Berechnung Ergebnisse benötigt, die nicht im Speicher des ausführenden Prozessors liegen, so können diese über ein Kommunikationsnetz angefordert werden.

• Das Kommunikationsnetz stellt direkte Verbindungen zwischen allen Prozessoren her. Damit soll sichergestellt werden, daß die Kosten für die Kommunikation unabhängig von der Entfernung der beteiligten Prozessoren sind. Fordert ein Prozessor einen Variablenwert an, so entstehen keine Kosten, falls dieser Wert in seinem eigenen Speicher liegt. Ist dies nicht der Fall, fallen für die entsprechende Kommunikation sehr hohe, aber konstante Kosten an.

• Für die Eingabe-Programme dieser Methode gelten die im ersten Abschnitt genannten Einschränkungen.

Während bei der Berechnung eines Schedules meistens ein laufzeit-opti-males Programm angestrebt wird, steht bei der Berechnung einer Allokation die Reduktion der Anzahl benötigter Prozessoren oder der Kommunikation zwischen Prozessoren im Vordergrund. Das im folgenden beschriebene Verfahren hat letzteres Ziel.

Nach der owner computes rule wird eine Kommunikation immer dann erforderlich, wenn die auszuführende Instanz von einer anderen Instanz abhängig ist, die nicht auf demselben Prozessor ausgeführt wurde, und deren Resultat damit nicht in dessen Speicher vorhanden ist. Die Grundlage für die vorzustellende Methode besteht darin, die Allokation so zu konstruieren, daß abhängige Zuweisungen auf demselben Prozessor ausgeführt werden. Formal bedeutet das für eine Allokation Π und jede Abhängigkeit e:

(∀ j⃗ : j⃗ ∈ P_{e} : Π(δ(e), j⃗) = Π(σ(e), h_{e}(j⃗))) (10)

Die einfachste Lösung dieses Gleichungssystems ist die triviale, die alle Instanzen von Zuweisungen auf denselben Prozessor abbildet, aber damit keinerlei Parallelität ausnutzt. Es stellt sich jedoch heraus, daß es Programme gibt, die nur diese Lösung zulassen, wie z. B. folgendes Schleifenprogramm:

for i := 1 to 10 do for j := i − 9 to i

a[i, j] = a[i − 1, j] + a(i − 2, j + 2) end

end

Hier müsste die Allokation folgende Gleichung erfüllen:

Π(i, j) = Π(i − 1, j) = Π(i − 2, j + 2)

Eine Lösung ist offensichtlich nur für ein konstantes Π möglich. Deshalb versucht man eine Allokation zu finden, welche die Anzahl von Instanzen maximiert, für die Gleichung (10) erfüllt ist. Diese Aufgabe stellt sich jedoch als sehr schwierig heraus, da es allgemein aufwendig ist, die Größe eines Polyeders, also die Anzahl der enthaltenen Punkte, zu berechnen. Aus diesem Grund vereinfacht man die Problemstellung dadurch, daß man nach Allokationen sucht, welche die Anzahl von Abhängigkeiten maximieren, deren zugehörige Indexpunkte Gleichung (10) erfüllen und damit keine Kommunikation verursachen. Wir bezeichnen Abhängigkeiten mit dieser Eigenschaft bzgl. einer gegebenen Allokation im weiteren als Trennkanten.

Unser Ziel ist es nun, die Anzahl von Trennkanten zu maximieren. Dabei gilt die anfangs genannte Forderung der Nicht-Trivialität der Allokation bezüglich jeder Zuweisung. Jeder Abhängigkeit wird dabei ein Gewicht entsprechend des verursachten Kommunikationsaufwands zugeteilt. Abhängigkeiten mit großem Gewicht sollen bevorzugt eliminiert werden. Die geeignete Wahl einer Gewichtsfunktion wird am Ende dieses Abschnitts vorgestellt.

5.1.1 Vorstellung des Verfahrens

Wie bei der Methode von Feautrier zur Berechnung von Schedules bildet man auch hier zunächst Prototypen für die Allokation jeder Zuweisung, die jeweils affine Abbildungen mit unbekannten Koeffizienten darstellen:

Π(S, i⃗) = a⃗_{S}i⃗ + b⃗_{S}n⃗ + c_{S} (11)

Mit n⃗ bezeichnen wir wieder den Vektor der Strukturparameter.

Ausgangspunkt des Verfahrens ist die Tatsache, daß sich wegen Gleichung (10) jeder Abhängigkeit eine Distanzfunktion d_{e} : P_{e} → ℤ zuordnen läßt:

d_{e}(j⃗) = Π(δ(e), j⃗) − Π(σ(e), h_{e}(j⃗)) (12)

Anmerkung 5.1.1 In [3] werden die Strukturparameter nicht explizit in der Formulierung der Prototypen angegeben. Die Prototypen haben hier folgendes Format:

Π(S, i⃗) = a⃗_{S}i⃗ + b_{S}

Da aber eine Berücksichtigung von Strukturparametern erforderlich ist, etwa bei der h-Transformation h_{e}(i) = n − i, werden diese wohl in den konstanten Anteil b_{S} verlagert. Für eine Implementierung des Verfahrens bedeutet das unnötigen Aufwand zur Speicherung von Konstanten. Deshalb verwenden wir im folgenden die oben eingeführte Definition von Prototypen.

Eine Abhängigkeit ist offensichtlich eliminiert, wenn für alle Indexpunkte j⃗ ∈ P_{e} gilt: d_{e}(j⃗) = 0. Ist P_{e} volldimensional, so kann dies nur der Fall sein, wenn sämtliche Koeffizienten der Distanzfunktion d_{e}, die als Differenz affiner Abbildungen wieder affin ist, null sind. Dies gilt nicht, wenn P_{e} nicht volldimensional ist. Da wir in der weiteren Berechnung die Koeffizienten von d_{e} gleich null setzen wollen, stehen wir vor der Aufgabe entweder P_{e} in einen volldimensionalen Polyeder umzuformen, oder die Distanzfunktion so zu definieren, daß nur die Schleifenindizes abgebildet werden, die eine Dimension in P_{e} erzeugen. Dazu benötigen wir einige Begriffe aus der linearen Programmierung [7]:

Definition 5.1.2 (implizite Gleichung) Sei Ax⃗ ≤ b⃗ ein System von Ungleichungen. Dann heißt eine Ungleichung a⃗x⃗ ≤ β dieses Systems implizite Gleichung, wenn gilt:

(∀ x⃗ : Ax⃗≤b⃗ : a⃗x⃗ = β)

Mit A^{=}x⃗≤b⃗^{=} wird das System der impliziten Gleichungen in Ax⃗≤b⃗ beschrieben, A^{+}x⃗≤b⃗^{+} bezeichne das System der übrigen Ungleichungen.

Satz 5.1.3 (Dimension eines Polyeders) Sei Ax⃗≤b⃗ ein System von Ungleichungen in n Variablen, das einen Polyeder P beschreibt. Dann ist die Dimension von P gleich n minus dem Rang der Matrix, die die impliziten Gleichungen beschreibt:

dim(P) = n − rg(A^{=})

In [3] ist eine Methode angesprochen, wie P_{e} zu einem volldimensionalen Polyeder erweitert werden kann. Die Idee ist dabei, diesen Polyeder parametrisiert darzustellen, wobei die Parameter gerade die Schleifenindizes sind, die in P_{e} keine Dimension aufspannen. Wir werden hier die zweite Möglichkeit besprechen, welche die Distanzfunktion umformt.

Besitzt P_{e} die Dimension k und ist l die Anzahl umgebender Schleifenindizes i_{1}, …, i_{l} der Zuweisung δ(e), so existieren also l −k explizite Gleichungen der Form a⃗x⃗ − β = 0 in der Beschreibung von P_{e}, die jeweils einen Schleifenindex evtl. in Abhängigkeit von umgebenden Schleifenindizes festlegen. Löst man die impliziten Gleichungen nach diesen Schleifenindizes auf und ersetzt sie in der Distanzfunktion durch den erhaltenen Ausdruck, so erhält man für die Distanzfunktion eine Abbildung in ausschließlich den k Schleifenindizes, die eine Dimension erzeugen. Mit der neuen Festlegung der Distanzfunktion gilt dann wieder, daß diese Abbildung genau dann überall in P_{e} null ist, wenn ihre Koeffizienten null sind.

Für jede Abhängigkeit entsteht eine Abbildung d_{e}, deren Koeffizienten sich aus additiven Kombinationen von Koeffizienten der Prototypen zusammensetzen. Nach den eben genannten eventuellen Umformungen bzgl. d_{e}, kann d_{e} nur null sein, wenn die Koeffizienten von d_{e} null sind.

Führt man nun in d_{e} einen Koeffizientenvergleich bezüglich aller Variablen, Strukturparameter und des konstanten Anteils durch, so ergibt sich damit folgendes System von Gleichungen:

C_{e}α⃗_{S} = 0

wobei α⃗_{S} den Vektor aller unbekannten Koeffizienten a⃗_{S}, b⃗_{S} und c_{S} aus (11) darstellt. Die Vereinigung dieser Systeme, die jeweils durch eine Abhängigkeit erzeugt werden, notieren wir mit

Cα⃗ = 0

Nehmen wir im folgenden an, daß C n Zeilen und m Spalten besitzt.

Der Lösungsraum dieses Gleichungssystems für den Koeffizientenvektor α⃗besteht i. a. nur aus dem Nullvektor. Wie bereits erwähnt, ist aber unser Ziel die Berechnung von nicht-trivialen Allokationen. Eine Allokation wird als nicht-trivial betrachtet, wenn mindestens einer der Koeffizienten, die vor den Schleifenindizes stehen, nicht null ist. Zur Lösung des Gleichungssystems mit der Zusatzbedingung der nicht-trivialen Allokationen benutzen wir folgende Abwandlung der Gauss-Jordan-Elimination zur Berechnung der Koeffizienten [3]. Die zugrundeliegende Idee ist, nach jedem Ersetzungsschritt eines Koeffizienten durch eine additive Kombination anderer Koeffizienten zu testen, ob dadurch einer der Prototypen trivial wird. Ist dies der Fall, wird diese Ersetzung verworfen.

Algorithmus 5.1.4 (Bestimmung der Prototypen-Koeffizienten)

Eingabe: Ein Gleichungssystem Cα⃗ = 0 und eine Menge P von Prototypen.

Ausgabe: Substitution für α⃗.

Methode: Der Algorithmus benutzt folgende Strukturen und Funktionen:

– S und S′ bezeichnen Substitutionen von Variablen, hier den Koeffizienten α_{1}, …, α_{m}, durch lineare Ausdrücke. Eine Substitution wird mit [x ← expr] notiert, wodurch ausgedrückt wird, daß die Variable x durch den Ausdruck expr substituiert wird.

– Das Array eliminated vom Typ boolean speichert, welche Variablen eliminiert wurden.

– Die Funktion get equation(i, Cα⃗ = 0) liefert die i-te Gleichung des Systems.

– Die Funktion subst equation(S, eq) wendet eine Substitution S auf die Gleichung eq an und gibt als Resultat die umgeformte Gleichung zurück.

– zero equation bezeichnet die Gleichung 0 = 0.

– Die Funktion coeff (j, eq) hat den Koeffizienten der j-ten Variable in der Gleichung eq als Resultat.

– Die Funktion solve(j, eq) löst eine Gleichung eq nach der j-ten Variable auf und gibt die rechte Seite der umgeformten Gleichung zurück.

– Die Funktion prototypes trivial(S′, P) testet, ob sich die Prototypen in P durch die Anwendung einer Substitution S′ auf die Koeffizienten zu trivialen Funktionen auswertet, deren Koeffizienten alle null sind. Gilt dies für mindestens einen der Prototypen, so gibt diese Funktion true als Resultat zurück.

Initialisierung: S = leere Substitution (∀ i : 1≤i≤m : eliminated [i] = false)

for i := 1 to n do eq := get equation(i, Cα⃗ = 0); if subst equation(S, eq) ≠ zero equation then

j := 1;

while coeff (j, eq) = 0 or eliminated [j] do

j := j + 1; end expr := solve(j, eq); eliminated(j) := true; τ := [α_{j} ← expr]; S′ := τ ∘ S; if not prototypes trivial(S′, P) then

S := S′; else

eliminated(j) := false; fi

fi end

Alle nicht eliminierten Koeffizienten von α⃗ können nun auf Werte ungleich null gesetzt werden, so daß alle eliminierten Koeffizienten von α⃗ durch die berechnete Substitution ebenfalls ungleich null werden. Um dies zu sehen, stellen wir die Substition in allgemeiner Form dar. Wir führen dazu folgende Bezeichnungen ein:

– m bezeichne wieder die Anzahl der Koeffizienten.

– k ist die Anzahl nicht eliminierter Koeffizienten.

– l ist die Anzahl der Koeffizienten, die durch die Substitution explizit gleich null gesetzt werden.

– r ist die Anzahl aller übrigen Koeffizienten, deren Werte sich in Abhängigkeit der nicht eliminierten Koeffizienten ergeben.

– Der Vektor a′ bezeichne das Resultat einer Substitution des Vektors a⃗.

Nach entsprechender Sortierung hat dann die errechnete Substition die Form:

[( a′_{1} ; a′_{2} ; ⋮ ; a′_{n} ) ← ( \begin{matrix}S \\ N \\ \begin{matrix}0 & ⋯ & 0 \\ ⋮ & ... & ⋮ \\ 0 & ⋯ & 0 \end{matrix} & I_{k} \end{matrix} ) ( a_{1} ; a_{2} ; ⋮ ; a_{n} )] (13)

Die l × m-Matrix N hat ausschließlich Einträge gleich null. Sie repräsentiert diejenigen Substitutionen, durch die Koeffizienten gleich null gesetzt werden. Der untere Teil der Matrix stellt die nicht eliminierten Koeffizienten dar. Dies entspricht der Substitution eines Koeffizienten durch sich selbst. Die r × m-Matrix S gibt die Substitution eliminierter Koeffizienten wieder, die von nicht eliminierten Koeffizienten abhängig sind. Bezeichnet s⃗_{i} den i-ten Zeilenvektor von S, so kann die obige Aussage bzgl. der Existenz geeigneter Werte für nicht eliminierte Koeffizienten nun formal ausgedrückt werden:

(∃ a⃗ : a⃗ ∈ ℚ^{m} : (∀ i : 1≤i≤r : s⃗_{i}a⃗ ≠ 0))

Wir zeigen dies durch die Herleitung eines Widerspruch für die Negation:

(∀ a⃗ : a⃗ ∈ ℚ^{m} : (∃ i : 1≤i≤r : s⃗_{i}a⃗ = 0))

Dies führt sofort zum Widerspruch, weil dazu für eine Zeile der Matrix S eine Nullzeile benötigt wird und S nach Konstruktion keine Nullzeilen besitzt.

Anmerkung 5.1.5 Wie bei der Berechnung von Schedules im zweiten Verfahren kann es auch hier vorkommen, daß die berechnete Allokation rational ist. Aus theoretischer Sicht stellt dies kein Problem dar, da hier anders als für einen Schedule kein Mindestabstand für die Prozessornummern gefordert ist. In der Praxis ist es jedoch üblich Prozessoren ganzzahlig aufzuzählen. Eine Möglichkeit, dies zu erreichen, besteht darin, die Allokationen aller Zuweisungen mit dem gemeinsamen Nenner zu multiplizieren.

Anmerkung 5.1.6 Im vorgestellten Algorithmus wird durch das Akzeptieren einer Substitution erreicht, daß für eine Abhängigkeit die in Gleichung (12) angegebene Distanzfunktion in einer Komponente auf null gesetzt wird. Wird dagegen eine Substitution verworfen, hat das zur Folge, daß die Abhängigkeit, aus der diese Substitution mit Gleichung (12) erzeugt wurde, keine Trennkante sein kann. Alle anderen zu dieser Abhängigkeit gehörenden Substitutionen werden aber dennoch getestet und eventuell akzeptiert. Damit erreicht man, daß eine Nicht-Trennkante eine Distanzfunktion besitzt, die wenigstens in einigen Komponenten null ist. Diese Optimierung der Distanzfunktion ist aber unnötig, da wir bei den am Anfang des Abschnitts vorgestellten Annahmen von einer direkten Verbindung aller Prozessoren ausgehen. Deshalb ist folgende Änderung des Algorithmus naheliegend. Falls eine der zu einer Abhängigkeit gehörenden Substitutionen verworfen wird, so verwerfe alle bisherigen und folgenden Substitutionen dieser Abhängigkeit. Wir ermöglichen dadurch anderen Abhängigkeiten, zu Trennkanten zu werden.

Beispiel 5.1.7 Gegeben sei folgendes Schleifenprogramm:

for i := 0 to 20 do

for j := 0 to 20 do

for k := 0 to 20 do a[i, j] := b[i, j − 6] + d[i − 1, j + 3]; {S_{1}} b[i + 1, j − 1] := c[i + 2, j + 5]; {S_{2}} d[i, j − 1] := a[i, j − 1] {S_{3}}

end end

end

Nach Gleichung (11) ergibt sich für die Prototypen der Allokation jeder Zuweisung:

Π(S_{1}, (i, j, k)) = a_{1}^{1}i + a_{1}^{2}j + a_{1}^{3}k + c_{1}
Π(S_{2}, (i, j, k)) = a_{2}^{1}i + a_{2}^{2}j + a_{2}^{3}k + c_{2}
Π(S_{3}, (i, j, k)) = a_{3}^{1}i + a_{3}^{2}j + a_{3}^{3}k + c_{3}

Es existieren folgende drei Abhängigkeiten:

e_{1} = S_{1}δS_{3} : h_{e_{2}}(i, j, k) = (i, j − 1, k)
P_{e_{2}} = {(i, j, k) : 0≤i, k ≤20 ∧ 1≤j ≤20 }
e_{2} = S_{3}δS_{1} : h_{e_{3}}(i, j, k) = (i − 1, j + 4, k)
P_{e_{3}} = {(i, j, k) : 1≤i≤20 ∧ 0≤j ≤16 ∧ 1≤k ≤20 }
e_{3} = S_{2}δS_{1} : h_{e_{1}}(i, j, k) = (i − 1, j − 5, k)
P_{e_{1}} = {(i, j, k) : 1≤i≤20 ∧ 5≤j ≤20 ∧ 1≤k ≤20 }

Mit der Beschreibung der Abhängigkeiten stellen wir nun Gleichung (12) für e_{1}, e_{2} und e_{3} auf:

\begin{matrix}e_{1} : & a_{3}^{1}i + a_{3}^{2}j + a_{3}^{3}k + c_{3} − (a_{1}^{1}i + a_{1}^{2}(j − 1) + a_{1}^{3}k + c_{1}) = 0 \\ e_{2} : & a_{1}^{1}i + a_{1}^{2}j + a_{1}^{3}k + c_{1} − (a_{3}^{1}(i − 1) + a_{3}^{2}(j + 4) + a_{3}^{3}k + c_{3}) = 0 \\ e_{3} : & a_{1}^{1}i + a_{1}^{2}j + a_{1}^{3}k + c_{1} − (a_{2}^{1}(i − 1) + a_{2}^{2}(j − 5) + a_{2}^{3}k + c_{2}) = 0 \end{matrix}

Der Koeffizientenvergleich bzgl. aller Indizes und des konstanten Anteils liefert folgende Gleichungssysteme: e_{1}:

a_{3}^{1} − a_{1}^{1} = 0
a_{3}^{2} − a_{1}^{2} = 0
a_{3}^{3} − a_{1}^{3} = 0
c_{3} + a_{1}^{2} − c_{1} = 0

e_{2}:

a_{1}^{1} − a_{3}^{1} = 0
a_{1}^{2} − a_{3}^{2} = 0
a_{1}^{3} − a_{1}^{3} = 0
c_{1} + a_{3}^{1} − 4a_{3}^{2} − c_{3} = 0

e_{3}:

a_{1}^{1} − a_{2}^{1} = 0
a_{1}^{2} − a_{2}^{2} = 0
a_{1}^{3} − a_{2}^{3} = 0
c_{1} + a_{2}^{1} + 5a_{2}^{2} − c_{2} = 0

Faßt man die drei Gleichungssysteme nun zu einem zusammen und übergibt sie dem vorgestellten Algorithmus, so ergibt sich folgende Substitution:

[( \begin{matrix}a_{1}^{1} \\ a_{1}^{2} \\ a_{1}^{3} \\ c_{1} \\ a_{2}^{1} \\ a_{2}^{2} \\ a_{2}^{3} \\ c_{2} \\ a_{3}^{1} \\ a_{3}^{2} \\ a_{3}^{3} \\ c_{3} \end{matrix} ) ← ( \begin{matrix}0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 3 & 0 & 0 \\ 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 1 & 0 & 0 \\ 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 1 & 0 \\ 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 1 & 0 & 1 \\ 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 3 & 0 & 0 \\ 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 1 & 0 & 0 \\ 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 1 & 0 \\ 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 9 & 0 & 1 \\ 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 3 & 0 & 0 \\ 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 1 & 0 & 0 \\ 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 1 & 0 \\ 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 1 \end{matrix} ) ( \begin{matrix}a_{1}^{1} \\ a_{1}^{2} \\ a_{1}^{3} \\ c_{1} \\ a_{2}^{1} \\ a_{2}^{2} \\ a_{2}^{3} \\ c_{2} \\ a_{3}^{1} \\ a_{3}^{2} \\ a_{3}^{3} \\ c_{3} \end{matrix} )]

wobei a_{3}^{2}, a_{3}^{3} und c_{3} die nicht eliminierten Variablen sind. Es ist bereits ausreichend nur a_{3}^{2} gleich 1 zu setzen, um nicht-triviale Prototypen zu erhalten. Damit ergibt sich als Ergebnis:

Π_{S_{1}}(i, j, k) = 3i + j + 1
Π_{S_{2}}(i, j, k) = 3i + j + 9
Π_{S_{3}}(i, j, k) = 3i + j

5.1.2 Einführung einer Gewichtsfunktion

Im Algorithmus werden die Abhängigkeiten, deren Gleichungen als erste im Gleichungssystem auftreten, bevorzugt eliminiert. Dies kann man sich dadurch erklären, daß jede Substitution, die akzeptiert wird, zur Elimination einer Abhängigkeit beiträgt. Das Akzeptieren einer neuen Substitution ist davon abhängig, daß die neue Substitution in Verbindung mit den bereits bestehenden keine trivialen Prototypen erzeugt. Somit ist die Chance auf ein Akzeptieren von Substitutionen am Anfang des Algorithmus höher als im weiteren Verlauf.

Um die Kosten für die Kommunikation zu minimieren, sind wir daran interessiert, Abhängigkeiten zu eliminieren, die einen großen Kommunikationsaufwand erzeugen. Wir benötigen also ein Maß zur Angabe dieses Aufwands, das wir als Kommunikationsvolumen bezeichnen. Die Gleichungen der Abhängigkeiten werden dann nach diesem Maß absteigend geordnet dem vorgestellten Algorithmus übergeben.

Zunächst ist es naheliegended, die Anzahl von Instanzen, die Quellen für eine Abhängigkeit sind, als Maß zu wählen. Das entspricht wieder dem Problem, die Punkte eines Polyeders aufzuzählen, was i. a. sehr aufwendig ist. Deshalb beschränken wir uns auf die Dimension des Polyeders E_{e}, der die Menge dieser Instanzen beschreibt. Eine Voraussetzung für die Korrektheit dieses Ansatzes ist, daß der verwendete Parallelrechner über die Möglichkeit des broadcasting verfügt, daß also ein Wert, der von mehreren Prozessoren zur Ausführung benötigt wird, nur einmal an das Kommunikationsnetz gesendet werden muß.

Formal ist der Polyeder E_{e} isomorph zum Bild des Polyeders P_{e} bezüglich der h-Transformation h_{e}:

E_{e} = { h_{e}(j⃗) : j⃗ ∈ P_{e} }

Für die Berechnung der Dimension des Polyeders E_{e} benutzt man wieder Satz 5.1.3.

Anmerkung 5.1.8 Im Fall uniformer Abhängigkeiten mit gleichem Kommunikationsvolumen und identischen Indexräumen bzgl. ihrer Quellen und Ziele ist eine weitere Verfeinerung des genannten Masses dadurch möglich, daß man die Länge der Distanzvektoren betrachtet. Aus der Definition einer uniformen Abhängigkeit und den Festlegungen für die Polyeder P_{e} und E_{e} folgt unmittelbar, daß es anhand des Distanzvektors möglich ist, aus der Beschreibung des Indexraums sowohl den Polyeder P_{e} als auch den Polyeder E_{e} zu berechnen. Der Distanzvektor bestimmt dabei, wieweit sich die Grenzen des Indexraums bzgl. P_{e} und E_{e} verschieben. Je größer die Länge eines Distanzvektors ist, umso mehr wird P_{e} und E_{e} im Vergleich zum Indexraum eingeschränkt. Uniforme Abhängigkeiten mit kürzeren Distanzvektoren sind also solchen mit längeren Distanzvektoren in der Elimination vorzuziehen, sofern die Distanz keinen Einfluß auf die Kosten der Kommunikation hat.

Anmerkung 5.1.9 Der Nachteil bei der vorgestellten Festlegung der Gewichtsfunktion ist, daß die Auswahl der Ordnung von Abhängigkeiten mit gleichem Kommunikationsvolumen einen Einfluß auf die Elimination von Abhängigkeiten mit niedrigerem Kommunikationsvolumen haben kann. Diese Auswirkung ist aber erst nach der Ausführung des Algorithmus bekannt. Deshalb ist die Wahl von Gewichtsfunktionen noch immer Gegenstand der Forschung.

5.1.3 Bewertung der Methode

Obwohl diese Methode, wie aus der letzten Anmerkung zu sehen ist, noch nicht vollständig ausgereift ist, bietet sie dennoch eine schnelle und effiziente Möglichkeit zur Konstruktion von Allokationen.

Neben der Wahl der Gewichtsfunktion gibt es jedoch noch eine weitere Schwierigkeit. Das Ergebnis des Verfahrens ist stets eine eindimensionale Allokation. Dadurch entsteht der Nachteil, daß die in einem Programm enthaltene Parallelität nicht vollständig ausgenutzt werden kann. Eine Abhilfe für dieses Problem stellen mehrdimensionale Allokationen dar, die im nächsten Abschnitt behandelt werden.

Der vorgestellte Algorithmus ist zwar deterministisch, das Verfahren für sich jedoch nicht. Dies liegt an der Auswahl eines zu eliminierenden Koeffizienten aus einer Gleichung. Der Algorithmus wählt immer den erstmöglichen. Verwendet man beispielsweise im Algorithmus eine Auswahlregel, die immer den letztmöglichen Koeffizienten liefert, ergibt sich für die Allokation in Beispiel 5.1.7:

Π_{S_{1}}(i, j, k) = i + \frac{1}{3}j
Π_{S_{2}}(i, j, k) = i + \frac{1}{3}j + \frac{8}{3}
Π_{S_{3}}(i, j, k) = i + \frac{1}{3}j − \frac{1}{3}

In beiden Fällen werden hier alle Abhängigkeiten eliminiert, dies gilt aber nicht allgemein.

5.2 Mehrdimensionale Allokation

Moderne Rechnerarchitekturen, die eine große Anzahl von Prozessoren verwalten, ordnen diese häufig gemäß eines mehrdimensionalen Gitters an. Benutzt man für die Abbildung von Instanzen von Zuweisungen auf dieses Gitter eine eindimensionale Allokation, ist zusätzlich eine sogenannte Faltungsfunktion notwendig, die eine Transformation der virtuellen Prozessornummern in physikalische vornimmt. Zur Vermeidung dieses Zwischenschritts ist es deshalb naheliegend, eine Allokation mehrdimensional zu konstruieren. Die Komponenten einer mehrdimensionalen Allokation geben dann jeweils eine Koordinate des Prozessors in der jeweiligen Dimension an.

Neben der Anpassung der Allokation an die Rechnerarchitektur gibt es aber einen noch wesentlicheren Grund zur Berechnung mehrdimensionaler Allokationen. Bei der Verwendung von eindimensionalen Allokationen kann es vorkommen, daß die in einem Eingabe-Programm mögliche Parallelität nicht ausgeschöpft wird, oder daß die durch die Allokation festgelegte Anzahl von virtuellen Prozessoren kleiner ist, als die der physikalisch vorhandenen. Ergibt sich z. B. für eine Zuweisung mit n^{k}-dimensionalem Indexraum ein eindimensionaler Schedule, so können pro Zeiteinheit O(n^{k−1}) Instanzen der Zuweisung gleichzeitig ausgeführt werden. Existiert nun eine eindimensionale Allokation, die eine Menge von Instanzen stets auf O(n) virtuelle Prozessoren verteilt, so werden alle O(n^{k−1}) Instanzen auf diesen Prozessoren ausgeführt, während alle übrigen Prozessoren ohne Auftrag sind. Die vorhandene Rechenleistung wird also nicht ausgenutzt.

Wird für ein Schleifenprogramm L beispielsweise mit der Methode von Feautrier ein d_{Θ}-dimensionaler Schedule konstruiert, so besitzt das Bild der Allokation jeder Zuweisung S ∈ Ω_{L} im Idealfall die Dimension p_{S} −d_{Θ}, wobei p_{S} wieder die Dimension des Indexraums D_{S} von S bezeichnet. Diese Größen sind i. a. nicht für alle Zuweisungen gleich und damit auch nicht gleich der Dimension des Gitters zur Prozessoranordnung. Es können also Zuweisungen existieren, die nicht die notwendige Dimension im Indexraum aufweisen, um alle Prozessoren zu belegen.

Deshalb ist entweder ein System erforderlich, das es erlaubt, die Geometrie der Prozessoranordnung bzgl. jeder Zuweisung anzupassen oder man muß für alle Zuweisungen eine Allokation mit gleicher Dimension wählen und in Kauf nehmen, daß während der Abarbeitung des gesamten Programms manche Prozessoren vorübergehend keinen Abarbeitungsauftrag erhalten. Durch den ersten Ansatz geht keine Rechenleistung verloren, andererseits wird eine Anpassung der Geometrie von den meisten Rechnern nicht unterstützt.

Wir werden im folgenden deshalb den zweiten Ansatz verfolgen [3].

5.2.1 Vorstellung des Verfahrens

Unser Ziel ist es, eine Allokation zu berechnen, welche die Instanzen jeder Zuweisung nach ℚ^{d_{Π}} abbildet, wobei d_{Π} die maximal mögliche Dimension darstellt:

d_{Π} = max( p_{S} − d_{Θ} : S ∈ Ω_{L}) (14)

Um Allokationen mit dieser Eigenschaft zu berechnen ist lediglich eine kleine Erweiterung des Algorithmus zur Berechnung eindimensionaler Allokationen notwendig. Diese Erweiterung betrifft den Test auf die Trivialität der Prototypen. In der ursprünglichen Version wird getestet, ob mindestens einer der Koeffizienten eines Prototypen nicht eliminiert ist, und damit ungleich null gesetzt werden kann. Um nun d_{Π}-dimensionale Allokationen zu erhalten, muß für die Substitution von Koeffizienten eines Prototypen folgendes gelten. Es müssen so viele nicht-eliminierte Koeffizienten in dem Prototyp verbleiben, daß bei d_{Π}-maliger Besetzung aller nicht-eliminierten Variablen mit geeigneten Werten die dadurch entstehenden d_{Π} Komponenten der Allokation linear unabhängig sind. Ergibt sich beispielsweise nach der Substitution für einen Prototypen die Form

Π_{S}(i, j, k) = a_{S}^{1}i + 0j + a_{S}^{3}k + c_{S}

und ist die für die Allokation geforderte Dimension d_{Π} gleich 2, so kann beispielsweise durch eine Besetzung von a_{S}^{1},a_{S}^{3} und c_{S} mit (1 0 0) und (0 1 0) folgende zweidimensionale Allokation erzeugt werden:

Π_{S}(i, j, k) = ( \begin{matrix}i \\ k \end{matrix} )

Aufgrund der oben beschriebenen Eigenschaften ist dies nicht für alle Prototypen möglich, denn diejenigen Prototypen für deren Zuweisung p_{S} − d_{Θ} kleiner als d_{Π} ist, können nicht volldimensional konstruiert werden. Die Komponenten dieser Prototypen sind in diesem Fall linear abhängig.

Die Koeffizienten eines Prototypen sind nach einer Substitution, anders als im gerade gegebenen Beispiel, i. a. durch eine additive Kombination der nichteliminierten Koeffizienten gegeben. Um die Frage zu entscheiden, ob genügend nicht-elimierte Koeffizienten vorhanden sind, um eine d_{Π}-dimensionale Allokation zu konstruieren, betrachten wir folgende alternative Darstellung eines Prototypen:

Π_{S}(i⃗, α⃗_{S}) = M(α⃗_{S}i⃗)

α⃗_{S} bezeichne wieder den Vektor aller Koeffizienten von Π_{S}. M sei die Matrix, die aus den Zeilen der Substitutionsmatrix in (13) besteht, die Substitutionen für Koeffizienten dieses Prototypen darstellen. Mit einer Aussage der linearen Algebra gilt, daß die Dimension des Bildes von Π_{S} gleich dem Rang der Matrix M ist. Durch die Berechnung dieses Rangs läßt sich also die gestellte Frage beantworten.

Beispiel 5.2.1 In Beispiel 5.1.7 besitzen alle Zuweisungen denselben Indexraum mit Dimension 3. Nehmen wir an, daß für diese Zuweisungen jeweils ein Schedule mit Dimension 1 berechnet wurde, so ist d_{Π} = 2. Mit dem neuen Trivialitätstest ergibt sich dieselbe Substitution wie in Beispiel 5.1.7 mit a_{2}^{3}, a_{3}^{3} und c_{3} als nicht-eliminierten Variablen. Besetzt man diese Koeffizienten nun mit den Werten linear unabhängiger Vektoren, etwa (1 0 0) und (0 1 0) erhalten wir folgende zweidimensionale Allokation:

Π_{S_{1}}(i, j, k) = ( \begin{matrix}3i + j + 1 \\ k \end{matrix} )
Π_{S_{2}}(i, j, k) = ( \begin{matrix}3i + j + 9 \\ k \end{matrix} )
Π_{S_{3}}(i, j, k) = ( \begin{matrix}3i + j \\ k \end{matrix} )

Anmerkung 5.2.2 Die für eine Allokation maximal mögliche Dimension d_{Π} wird mit Gleichung (14) durch die Dimension des Schedules bestimmt. Steht ein Verfahren zur Verfügung, das es erlaubt neben ein- auch mehrdimensionale Schedules zu berechnen, so läßt sich durch die Wahl eines mehrdimensionalen Schedules die Dimension d_{Π} der Allokation verringern.

Ziel der Berechnung einer Allokation ist die Maximierung der Anzahl von Trennkanten und damit eine Verringerung der Kommunikationskosten. Bei empirischen Untersuchungen des vorgestellten Verfahrens zur Berechnung der Allokation fällt auf, daß eine Allokation mit geringerer Dimension mehr Trennkanten besitzt, als eine Allokation mit höherer Dimension. Das liegt daran, daß zur Einhaltung der geforderten Dimension im Algorithmus weniger Substitutionen akzeptiert werden können. Jede akzeptierte Substitution trägt aber zur Maximierung von Trennkanten bei.

Bei einer isolierten Betrachtung der Güte einer Allokation folgt also, daß für die Berechnung eine möglichst geringe Dimension zu wählen ist. Im Gegenzug verschlechtert sich aber die Laufzeit des Ziel-Programms, weil durch einen mehrdimensionalen Schedule die Parallelität im Eingabe-Programm i. a. nicht voll ausgenutzt wird. Hier muß also in Abwägung der Laufzeit und der Kosten für eine Kommunikation für ein konkretes Rechnersystem ein Mittelweg gefunden werden.

Anmerkung 5.2.3 In [3] wird keine Aussage über die Berechnung der Allokation für eine Zuweisung ohne Abhängigkeiten gemacht, wie sie z. B. in folgendem Schleifenprogramm mit Zuweisung S_{1} existiert: for i := 0 to 20 do

for j := 0 to 20 do a[i, j] := 0; {S_{1}} b[i, j] := b[i,j-1]; {S_{2}}

end end

Das Verfahren von Feautrier zur Berechnung eines Schedules weist solchen Zuweisungen die Nullabbildung als Schedule zu. Der Schedule ist also 0dimensional. Der Schedule Θ_{S_{2}}(i, j) = j für die Zuweisung S_{2} ist dagegen eindimensional. Für die Dimension der Allokation folgt also mit Gleichung (14), daß d_{Π} = 2 gilt.

5.2.2 Bewertung der Methode

Da sich die Berechnung mehrdimensionaler Allokationen direkt aus der für eindimensionale Allokation ableitet, gelten die gleichen Aussagen bezüglich der Wahl einer Gewichtsfunktion und des Nichtdeterminismus. Darüberhinaus entsteht aber neuer Nichtdeterminismus bei der Besetzung der nichteliminierten Koeffizienten. Bei mehrdimensionalen Allokationen wirkt sich dieser Nichtdeterminismus stark auf die Auslastung des Prozessorfeldes aus. Aus diesem Grund wird derzeit am Laboratoire MASI der Universität Versaille erforscht, wie die verwendete Rechnerarchitektur in das Berechnungsverfahren miteinzubeziehen ist.

Abgesehen von den beschriebenen Problemen stellt diese Methode aber eine sehr schnelle Möglichkeit zur Berechnung von Allokationen dar.

6 Erzeugung des parallelen Programms

In diesem Abschnitt wird eine allgemein anwendbare Methode zur Konstruktion des parallelen Programms gemäß Schedule und Allokation vorgestellt [3].

Die Idee des Verfahrens ist, sämtliche Kombinationen von Zeit- und Raumpunkten aufzustellen und zur Laufzeit für jede aufgezählte Instanz einer Zuweisung die Ausführung zu prüfen. Die Steuerung, die von Schedule und Allokation vorgegeben ist, wird damit explizit in den Zielcode integriert.

Dazu werden zunächst die Menge von Instanzen von Zuweisungen eines Schleifenprogramms L berechnet, die zu einem gegebenen Zeitpunkt ausgeführt werden sollen. Derartige Mengen werden als Front [3] bezeichnet. Für eine Front gilt also

F(t) = {(S, i⃗) : S ∈ Ω_{L} ∧ i⃗ ∈ D_{S} ∧ Θ(S, i⃗) = t}

Alle Instanzen einer Front können parallel ausgeführt werden und entsprechend der Allokation Π über das Prozessorfeld verteilt werden. Für die Ausführung von Instanzen einer Front auf einem Prozessor q ergibt sich also folgende Teilmenge einer Front F(t):

F(t, q) = {(S, i⃗) : S ∈ Ω_{L} ∧ i⃗ ∈ D_{S} ∧ Θ(S, i⃗) = t ∧ Π(S, i⃗) = q}

Ist l die Laufzeit des parallelen Programms bezüglich Θ und p die Anzahl von Prozessoren bezüglich Π, so kann daraus das parallele Programm wie in Abbildung 8 aufgebaut werden. Es kann durchaus vorkommen, daß auf einem Prozessor zum selben Zeitpunkt mehrere Instanzen von Zuweisungen abgearbeitet werden müssen, daß also die Mengen F(t, q) nicht einelementig sind. Das gilt insbesondere, wenn Schedule und Allokation, wie in den letzten beiden Abschnitten, unabhängig voneinander berechnet wurden. Deshalb wurde die innerste seqfor-Schleife eingefügt. Zu einem Konflikt führt dies allerdings nicht, da der Schedule die Unabhängigkeit garantiert. Die Reihenfolge der Abarbeitung ist dabei freigestellt. Lediglich die Prozessorauslastung des ganzen Systems verschlechtert sich in einem solchen Fall, da andere Prozessoren auf die Freischaltung durch barrier warten müssen. Die Zeitschleife kann dann allerdings nicht mehr als Systemzeit interpretiert werden. Sie stellt eine virtuelle, grobkörnige Zeit dar. Eine Implementierung dieses Verfahrens ist durch Prozesse möglich. Ein Prozeß besteht dann jeweils aus der Berechnung von Instanzen einer Menge F(t, q). Die Synchronisation erfolgt über parfor q := 0 to p do

seqfor t := 0 to l do {read remote data} seqfor

F(t, q) end {barrier}

end end

Abbildung 8: Explizite Steuerung durch F(t, q)

die Konstrukte read remote data und barrier jeweils vor und nach den zu berechnenden Instanzen.

Anmerkung 6.0.4 Die Berechnung der Werte p und l für eine Allokation Π und einen Schedule Θ kann mit Hilfe von PIP durchgeführt werden, indem man die in Abschnitt 4.2.3 vorgestellten Methoden zur Lösung eines Maximierungsproblems benutzt. Die Abbildungen Π und Θ stellen jeweils die Zielfunktion dar, während die Beschreibung des Indexraums das Ungleichungssystem Gμ⃗≥h⃗ bildet.

Die Frage, ob die Mengen F(t, q) einelementig sind, kann durch eine Betrachtung von Schedule und Allokation entschieden werden. Ist die Matrix

T = ( \begin{matrix}Θ(S, i⃗) \\ Π(S, i⃗) \end{matrix} )

für jede Zuweisung S invertierbar, so folgt, daß die Mengen F(t, q) einelementig sind. Die Matrix T ist eine Transformationsmatrix, wie sie in Abschnitt 4.1 eingeführt wurde. Für eine solche Matrix ist es möglich die Mengen F(t, q) durch die Anwendung der inversen Matrix T^{−1} bereits zur Compilezeit zu berechnen [6]. Gilt dies nicht, so wird in der konkreten Implementierung des Zielcodes die Berechnung der Mengen F(t, q) dadurch realisiert, daß zunächst der ursprüngliche Code anstelle von F(t, q) eingesetzt wird. Alle Schleifen des Eingabe-Programms werden jedoch durch parfor-Schleifen ersetzt. Durch ein umgebendes if-Konstrukt, das Schedule und Allokation beinhaltet, wird dann zur Laufzeit auf die Ausführung getestet.

Folgendes Beispiel zeigt die Programmkonstruktion nach diesem Verfahren:

Beispiel 6.0.5 Gegeben sei folgendes Schleifenprogramm:

for i := 1 to 10 do

for j := i to i + n + 5 do

b[i] := b[i] + 1 end

end

Mit den statementbasierten Methoden aus Abschnitt 2 und 3 berechnet man folgende Abbildungen als Schedule und Allokation:

Θ(i, j) = j − i
Π(i, j) = i

Damit entsteht folgender Zielcode:

parfor q := 0 to 10 do

seqfor t := 0 to n + 14 do {read remote data} parfor i := 1 to 10 do

parfor j := i to i + n + 5 do if j − i = t and i = q then

b[i] := b[i] + 1 fi end

end {barrier}

end end

7 Schlußbemerkungen

Mit den vorgestellten Methoden zur Berechnung eines Schedules und einer Allokation sind die Hauptbestandteile zur Parallelisierung eines sequentiellen Programms gegeben.

Wir haben für den Bereich Schedule zwei Verfahren kennengelernt. Dabei steht die Hyperebenen-Methode für Verfahren mit großen Einschränkungen an die Eingabe-Programme. Sie besitzt jedoch in der Form von Lamport eine sehr kurze Laufzeit. Die Hyperebenen-Methode wurde vielfach weiterentwickelt, etwa zum statement-orientierten Ansatz [12]. Ihre größte Einschränkung stellt jedoch die Beschränkung auf uniforme Abhängigkeiten und perfekte Verschachtelung dar.

Eine Methode ohne derartige Beschränkungen ist mit dem Verfahren von Feautrier gegeben. Es liefert unter den genannten Einschränkungen sogar zeitoptimale Schedules. Während die Methode für sich ebenfalls sehr schnell arbeitet, besteht der Großteil der benötigten Laufzeit aus der Lösung der Optimierungsprobleme mit PIP. Ein Problem beim Einsatz von PIP ist außerdem, daß es bei größeren Zahlen in der Eingabe nicht mehr stabil arbeitet. Dieses Programm ist aber bislang als einziges in der Lage, parametrisierte Probleme zu lösen.

Für den Bereich der Allokation wurde eine Methode vorgestellt, die ebenfalls geringe Einschränkungen für die Eingabe-Programme fordert. Diese Methode ist allerdings noch nicht vollständig ausgereift, da sie beispielsweise die verwendete Rechnerarchitektur nicht in die Berechnung miteinbezieht. Ebenfalls ist die Wahl einer geeigneten Gewichtsfunktion noch zu untersuchen.

Die vorliegende Arbeit sollte einerseits die notwendige Theorie der einzelnen Verfahren vermitteln, andererseits aber auch die Grundlage für eine Implementierung geben. Die vorgestellten Methoden wurden im Rahmen des Projekts LooPo (Loop Parallelization in the Polytope Model) implementiert. LooPo ist ein Tool zur Parallelisierung von sequentiellen Programmen und soll auch einen kleinen Teil zur Entwicklung des Gebiets der Parallelisierung beitragen.

Literatur

[1] Paul Feautrier. Some efficient solutions to the affin scheduling problem, Part I, One-dimensional Time. International Journal of Parallel Programming, Vol. 21, No. 5, pp. 313-347, 1992.

[2] Paul Feautrier. Some efficient solutions to the affin scheduling problem, Part II, Multidimensional Time. International Journal of Parallel Programming, Vol. 21, No. 6, pp. 389-420, 1992.

[3] Paul Feautrier. Toward Automatic Distribution. Parallel Processing Letters, Vol. 4, No. 3, pp. 233-244, 1994.

[4] Paul Feautrier. Dataflow Analysis of Array and Scalar References. International Journal of Parallel Programming, Vol. 20, pp. 23-53, 1991.

[5] Uptal Banerjee. Loop Transformations for Restructuring Compilers: The Foundations. Kluwer Academic Publishers, Boston/Dordrecht/London, 1993.

[6] Christian Lengauer. Loop Parallelization in the Polytope Model . Invited paper at CONCUR ’93, Lecture Notes in CS, Springer-Verlag, Rep. MIP-9303, Fakultät für Mathematik und Informatik, Universität Passau.

[7] Alexander Schrijver. Theory of Linear and Integer Programming. Series in Discrete Mathematics. John Wiley & Sons, 1986.

[8] Leslie Lamport. The Parallel Execution of DO Loops. Communications of the ACM, Vol. 17, No. 2, pp. 83-93, February 1974.

[9] Alain Darte, Leonid Khachiyan, Yves Robert. Linear Scheduling is nearly optimal. Parallel Processing Letters, Vol. 1, No. 2, pp. 73-81, 1991.

[10] Sabine Wetzel. Codegenerierung im Polyedermodell. Diplomarbeit, Fakultät für Mathematik und Informatik, Universität Passau, 1995, in Vorbereitung.

[11] E. W. Dijkstra, C. S. Scholten. Predicate Calculus and Program Semantics. Texts and Monographs in Computer Science. Springer-Verlag, 1990.

[12] Alain Darte, Yves Robert. Constructive Methods for Scheduling Uniform Loop Nests. IEEE Transactions on Parallel and Distributed Systems, Vol. 5,No. 8, August 1994.

[13] J. A. B. Fortes, F. Parisi-Presicce, Optimal linear schedules for the parallel execution of algorithms. Int. Conf. Parallel Processing, August 1984, pp. 322-328..

[14] G. B. Dantzig. Maximization of a linear function of variables subject to linear inequalities, in: Activity Analysis of Production and Allocation, Wiley, New York, 1951, pp. 339-347.

[15] Paul Feautrier. Résolution de systèmes d’inéquations linéaire; mode d’emploi du logiciel PIP. Laboratoire MASI, Paris, 1989.