Kniffliger algorithmus zu Parallelisierung

[Minimalist (Gast)]

Hallo,
ich sitze hier grad vor einem Problem und bräuchte etwas zusätzlichen 
Hirnschmalz.

Kurze Einführung:
Im meinem Projekt (geht um numerische Simulationen) wird eine 
pseudo-Diagonalisierung einer Matrix verwendet, da die vollständige 
Diagonalisierung zu lange dauert. Das ist ausgiebeig getestet und 
bestätigt. Wir verwenden eine Pseudodiagonalisierung, die im Grunde auch 
parallelisierbar wäre, wenn da nicht ein Problem wäre. Und zwar führt 
der Algorithmus immer eine Paarweise 2x2 Rotation von Eigenvektoren 
durch, dadurch parallelisiert der code nicht mehr.

Pseudocode:
für alle EigenvektorenA:
   für alle EigenvektorenB:
     rotiere A-B


dabei wird auf A und B geschrieben, d.h. wenn ein Tread bei A1-B4 und 
ein anderer bei A2-B4 ist, dann knallts. Es geht hier um >5000 
Eigenvektoren und um bis zu 28 Threads. Die Menge der Eigenvektoren A 
und B ist nicht identisch.
OMP ORDERED DO mit geschützten Zugriffen ist zu langsam und zerstört die 
Vektorisierung der innersten Loop, wo die Rotation dann durchgehführt 
wird.

Daher suche einen Algorithmus, mit dem ich die Paarungen "batchen" kann, 
so dass innerhalb einer Batch ein Index nie zweimal vorkommt.

Im Grunde eine Art Jeder-Gegen-Jeden Mannschafts Problem mit ungleich 
großen Mannschaften.
Wie ordne ich die Begegnungen an, wenn jeder in Mannschaft A mit jedem 
aus Mannschaft B gekämpft haben soll, und möglichst viele Kämpfe 
parallel laufen sollen. Das heißt ab Ende soll etwas rauskommen wie:
Batch 1:
1-1, 2-2, 3-3, 4-4
Batch 2:
1-2, 2-1, 3-4, 4-3
Batch 3:
1-3, 2-4, 3-1 ,4-2
Batch 4:
1-4, 2-3, 3-2, 4-1
Nur mit unterschiedlich (und beliebig) großen Mannschaften.

Hoffe, ich hab mich halbwegs verständlich ausgedrückt.
Ideal wäre natürlich wenn dieses Preescreening ebensfalls 
parallelisierbar bzw. Vektorisierbar ist.

Vielen Dank euch,
Minimalist

[minimal (Gast)]

Eine effizientes Ablaufschema für dieses Beispiel wird als eines der 
einfachsten Parallelisierungskonzepte in jedem Numeriklehrbuch 
abgehandelt. Es wird auf erläutert, dass diese Lösung nur theoretischen 
Wert hat, da sie praktisch dem Compiler überlassen wird. Einfach mal 
eins der Standardbücher lesen.

Wenn du allerdings den Ablauf strikt erzwingst widerspricht dies der 
Portierbarkeit. Ich weiß ja nicht, was du veranstaltest, aber ein 
optimierter Compiler löst solche Abläufe grundsätzlich nahezu perfekt 
auf, wenn man nicht bewusst Parallelisierungs-unfreundlich programmiert.

Schließlich, wenn du meinst, auch noch die Ablauferstellung 
parallelisieren zu müssen, hast du den grundlegenden Sinn der 
Parallelisierung noch nicht verstanden.

[Klaus P. (Gast)]

Minimalist schrieb:
> Wie ordne ich die Begegnungen an, wenn jeder in Mannschaft A mit jedem
> aus Mannschaft B gekämpft haben soll, und möglichst viele Kämpfe
> parallel laufen sollen.

Unter der Annahme, dass die "Kämpfe" annähernd gleich lange dauern, kann 
die kleinere "Mannschaft" ja bis zur Maximalanzahl der Threads parallel 
arbeiten. Die größere "Mannschaft" wird dabei nach jedem Lauf im Simme 
eines Ringspeichers rotiert.

Also im Prinzip so:

KleinereAB = EigenvektorenA.Count < EigenvektorenB.Count ? 
EigenvektorenA : EigenvektorenB;
GrößereAB = EigenvektorenA.Count < EigenvektorenB.Count ? EigenvektorenB 
: EigenvektorenB;

für (i = 0 ... (KleinereAB.Count / Anzahl Threads + 1))
   KleineAbTeilmenge = KleinereAB[i * Anzahl Threads]
      für alle GrößereAB
         rotiere parallel KleineAbTeilmenge - GrößereAB
         Schiebe GrößereAB im Ringspeicher

[Minimalist (Gast)]

minimal schrieb:
> Ich weiß ja nicht, was du veranstaltest, aber ein
> optimierter Compiler löst solche Abläufe grundsätzlich nahezu perfekt
> auf, wenn man nicht bewusst Parallelisierungs-unfreundlich programmiert.

Die Routine ist selbst eher überschaubar, und wurde von Menschen 
geschrieben, die wirklich wissen, was sie tun. Namentlich Pulay, der 
Pulay der auch das DIIS-Verfahren entwickelt hat. Beschrieben ist das 
verfahren hier: STEWART. J.J.P., CSASZAR, P., PULAY, P., J. COMP. 
CHEM.,3, 227, (1982)

Der Compiler löst solche Probleme von selbst gar nicht.
Mal etwas detaillierter

für alle EigenvektorenA:
   für alle EigenvektorenB:
     für alle elemente M:
               a=EigenvektorA(m)
               b=EigenvektorB(m)
               EigenvektorA(m)=alpha*a+beta*b
               EigenvektorB(m)=alpha*b-beta*a
Die innerste Schleife ist sehr gut Vektorisierbar. Aber in jedem 
Schleifendurchlauf wird gelesen und geschrieben.

Daher parallelisiert der Comiler erstmal nichts. OpenMP paralelisierung 
geht zwar prinzipiell wie oben beschrieben, muss aber LOCKS setzen, um 
die Speicherzugriffe zu schützen. Das zerstört die innere Vektorisierung 
und bremst den Code aus. Du hast schon Recht, das Compiler i.d.R sehr 
effizient optimieren, mir geht es um die Frage, ob es möglich (und 
Sinnvoll) ist hier Vektorisierung UND Parallelisierung effektiv zu 
nutzen. Außerdem wäre ich dir für konkrete Literaturhinweise zu diesem 
Thema sehr dankbar, denn die meisten Beispiele die ich gefunden habe 
lösen das klassische jeder-gegen-jeden Problem, aber nicht das oben 
beschriebene spezielle.

Trotzdem schon mal vielen Dank für die Antwort

[Minimalist (Gast)]

@Klaus P. (Gast)
Danke, da hatte ich einen Knoten im Hirn. Ich hatte mich darin verrannt 
die kleinere Menge zu rotieren und den Rest dann Rekrsiv zu verarbeiten 
bis kein Rest mehr bleibt.

Manchmal brauchts nur einen, der mal von Draußen drauf schaut...