OpenMP

[Ben04]

Inhaltsverzeichnis

• 1 OpenMP - Was ist das?
• 2 Compilerunterstützung
• 3 Syntax
• 3.1 Parallel
• 3.2 Thread ID, single und master
• 3.3 Sections
• 3.4 Barrier
• 3.5 For
• 3.5.1 Verschachtelte Schleifen mit und ohne Collapse
• 3.5.2 Partielle sequentielle Abarbeitung mit ordered
• 3.6 Kurzschreibweise
• 4 Gemeinsamer Speicher
• 4.1 Critcal
• 4.2 Flush
• 4.3 Reduction
• 4.4 Locks/Mutex
• 5 Zeitmessung
• 6 Suchen in einem ungeordneten Array
• 7 Akkumulation
• 7.1 STL
• 8 Gibt es noch mehr?
• 9 Quellen

1 OpenMP - Was ist das?

OpenMP ist eine standardisierte Schnittstelle für Fortran, C und C++, welche die Entwickelung von Programmen, welche mehrere Prozessor-Kerne nutzen, vereinfachen soll. Ich werde in diesem Artikel aber nur auf die C++- Schnittstelle eingehen.

2 Compiler Unterstützung

Es gibt mehrere Versionen von OpenMP. Dieser Artikel beschäftigt sich hauptsächlich mit OpenMP 2.5. Ich werde allerdings auch auf ein paar Neuerungen in der Version 3.0 eingehen. OpenMP wird von einer Reihe von Compilern unterstützt, jedoch meistens nicht in den Basisausstattungen dieser. Mit dem GCC gibt es aber mindestens ein Compiler, der jedem Leser zugänglich ist.

Folgende Tabelle liefert einen kleinen Überblick.

• GCC >= 4.2 (mit dem -fopenmp Parameter)
• MinGW Technology-Preview-Version (mit dem -fopenmp Parameter)
• VC 2005 Professional und Team System (mit dem /openmp Parameter, Standard- und Express-Versionen gehen nicht)

Auf der OpenMP-Webseite findet man eine detailliertere Tabelle.

3 Syntax

OpenMP besteht zum Großteil aus Compiler-Erweiterungen. Es gibt eine Reihe von Pragma-Anweisungen, welche die OpenMP-Funktionalität bereitstellen. Jede von ihnen hat folgende Syntax:

#pragma omp Anweisung
C++ Code, auf den sich die Anweisung bezieht

Neben diesen pragma-Konstrukten gibt es noch eine Reihe von normalen Funktionen. Diese werden im Header omp.h definiert. Alle haben das Präfix omp_.

Jeder Compiler, der sich an den C++-Standard hält, muss unbekannte Pragma-Anweisungen ignorieren und höchstens mit einer Warnung versehenen. Fehler sind nicht standardkonform. Durch das Benutzen eines Dummy-Headers ist es deswegen in vielen Fällen möglich, eine OpenMP-Anwendung sequentiell auf nicht OpenMP-fähigen Compilern zu übersetzen. Dies ist auch empfehlenswert, wenn der Compiler zwar OpenMP unterstützt, allerdings der Zielprozessor keine echte Parallelisierung unterstützt.

In diesem Artikel werde ich versuchen, die Anweisungen, welche kein sequentielles Fallback-Konstrukt haben, zu ignorieren. Sie bieten auch eher selten Funktionalität an, welche man nicht auch anders erhalten kann. Daher ist dies an sich kein großer Verlust.

3.1 Parallel

Die grundlegende Anweisung, ohne die bei OpenMP nichts auskommt, ist die parallel-Anweisung. Sie erzeugt eine Reihe von Threads und führt den nachfolgenden C++ Block mit jedem Thread aus. Dies bedeutet, dass der gleiche Code mehrmals gleichzeitig ausgeführt werden kann, und zwar von jedem Thread einmal. Wie viele Threads erstellt werden, ist der Implementierung von OpenMP freigestellt. Meistens werden aber so viele Threads erzeugt wie es im Rechner Kerne gibt.

Es ist also in der Regel nicht möglich, vorherzusagen, wie oft ein parallel-Block ausgeführt wird.

int main(){
	printf("Hier ist noch alles einfach und sequentiell.\n");
	#pragma omp parallel
	{
		printf("Hallo ");
		printf("Welt!\n");
	}
	printf("Alle Threads haben Hallo gesagt.");
}

Im parallel-Block ist ganz normaler C++-Code erlaubt. Man kann C++-Objekte anlegen und auch Funktionen aufrufen. Sogar Ausnahmen sind erlaubt, solange sie nicht aus dem parallel-Block rausfliegen.

Solange sich die Threads nicht in die Quere kommen, ist das auch schon alles, was man machen muss, um sämtliche Kerne seines Rechners zu benutzen. Allerdings ist dies leider nur sehr selten der Fall. Bereits in diesem einfachen Beispiel ist es möglich, dass Threads sich streiten. Wenn zwei Threads gleichzeitig etwas ausgeben, ist es nämlich nicht klar, welcher zuerst dran kommt. Alles von Funktionieren-wie-erwartet bis hin zu Abstürzen ist eine mögliche Konsequenz.

Was das Problem noch komplizierter macht, ist, dass das Verhalten nicht vollständig deterministisch ist. Wie viele Threads erzeugt werden, kann von Lauf zu Lauf unterschiedlich sein, und welcher Thread wann drankommt, kann man noch schlechter vorhersagen. Dies führt dazu, dass manche Bugs nur mit einer gewissen Wahrscheinlichkeit eintreten. Das Debuggen solcher Fehler ist sehr schwer.

An sich müsste ich hier sicherstellen, dass zu jedem Zeitpunkt immer höchstens nur ein Thread etwas auf die Konsole schreibt. Bei printf kann man aber normalerweise davon ausgehen, dass lediglich die Reihenfolge der Ausgaben etwas durcheinander geraten kann. Abstürze gehören eher in den Bereich der Theorie. printf leidet normalerweise weniger unter diesem Problem als std::cout, deswegen werde ich in diesem Artikel ausnahmsweise printf bevorzugen. Im Folgenden werde ich davon ausgehen, dass Ausgaben mit printf einfach funktionieren. Alles andere würde die Beispiele zu kompliziert machen.

Da man nicht wirklich weiß, wie oft ein parallel-Block ausgeführt wird, ist eine Parallelisierung recht umständlich ohne weitere OpenMP-Konstrukte. Er eignet sich aber in Regel vorzüglich, um threadspezifische Variablen zu initialisieren.

Mit Hilfe der num_threads(n)-Anweisung kann man festlegen, wie viele Threads erstellt werden. n muss eine Ganzzahl sein und braucht nicht zum Zeitpunkt der Übersetzung bekannt zu sein.

int main(){
	printf("Wie viele Threads sollen es sein?\n");
	int n;
	scanf("%d", &n);
	#pragma omp parallel num_threads(n)
	{
		printf("Hallo ");
		printf("Welt!\n");
	}
	printf("Alle Threads haben Hallo gesagt.");
}

Wenn man für n eine große Zahl eingibt, kann man mit großer Wahrscheinlichkeit beobachten, wie sich die Threads in die Quere kommen. Es folgt eine Beispielausgabe, so wie ich sie beobachten konnte.

Hallo Hallo Welt!
Hallo Welt!
Welt!
Hallo Welt!

Während der eine Thread noch am Schreiben ist, funkt ein anderer ihm einfach dazwischen. Versuchen Sie, nachzuvollziehen, wie diese Ausgabe zustande kam. Es ist sehr wichtig, dieses Problem zu verstehen, denn sonst werden Sie mit OpenMP oder auch Multithread-Programmierung im Allgemeinen nur Ärger haben.

omp_get_max_threads gibt, an wie viele Threads im nächsten parallel-Block erzeugt werden.

Bei den meisten OpenMP-Implementierungen werden die Threads nur einmal erstellt. Dies gilt auch, wenn es mehrere parallel-Blocks gibt. Zwischen den Blocks gehen die zusätzlichen Threads in einen Wartezustand. Sie werden erst am Ende des Programms zerstört. Die in der Regel recht hohen Threaderstellungskosten muss man also nur genau einmal zahlen.

parallel-Blöcke dürfen verschachtelt werden. In der Standardeinstellung erzeugt der verschachtelte Block aber keine neuen Threads.

3.2 Thread ID, single und master

Jeder Thread besitzt eine eindeutige Identifikationsnummer. Es handelt sich dabei um fortlaufende Nummern und nicht etwa um gecastete Zeiger oder ähnliches. Man kann die ID also als Index in einem Array benutzen. Den Thread mit der ID 0 nennt man auch noch Master.

Mit Hilfe von omp_get_thread_num kann man die ID eines Threads bestimmen und durch einen Aufruf von omp_get_num_threads kann man die Gesamtanzahl der Threads, welche vom parallel-Block erzeugt wurden, herausfinden.

int main(){
	#pragma omp parallel
	{
		printf("Thread %d sagt Hallo \n", omp_get_thread_num());	
		if(omp_get_thread_num() == 0)
			printf("Es gibt %d Threads.\n", omp_get_num_threads());		
	}

}

Wie man an diesem Beispiel sieht, ist es möglich, auch in einem parallel-Block Instruktionen nur einmal auszuführen. Allerdings gibt es dafür bessere OpenMP-Konstrukte.

int main(){
	#pragma omp parallel
	{
		printf("Thread %d sagt Hallo \n", omp_get_thread_num());	

		#pragma omp master
		{
			printf("Es gibt %d Threads.\n", omp_get_num_threads());	
		}

		#pragma omp single
		{
			printf("Thread %d war hier.\n", omp_get_thread_num());	
		}

		#pragma omp single
		{
			printf("Danach war Thread %d da.\n", omp_get_thread_num());	
		}
	}
}

Ein master-Block wird nur vom Master-Thread ausgeführt. Der Code ist äquivalent zur if-Abfrage aus dem vorherigen Beispiel. Ein single-Block wird genau einmal ausgeführt, allerdings ist es egal, von welchem Thread.

Sämtliche Threads warten am Ende eines single-Blocks, bis auch der letzte es bis dahin geschafft hat. Dadurch ist sichergestellt, dass der zweite single Block nicht ausgeführt wird, ehe der erste vollständig abgearbeitet wurde. Durch eine nowait-Anweisung kann man das Warten unterdrücken.

int main(){
	#pragma omp parallel
	{
		#pragma omp single nowait
		{
			printf("Manchmal bin ich als erstes da.\n");	
		}

		#pragma omp single nowait
		{
			printf("Hin und wieder gewinne aber auch ich das Rennen.\n");	
		}
	}
}

Wenn man nowait verwendet, dann wird nur noch garantiert, dass die Ausführung des ersten Blocks vor der des zweiten beginnt. Es ist aber gut möglich, dass der zweite vor dem ersten beendet wird. Es ist auch möglich, dass die Ausführung des ersten Threads noch vor der ersten Anweisung im Block unterbrochen wird und es dadurch erscheint, als wenn der zweite Block vor dem ersten gestartet werden würde. In der Praxis kann man also nicht darauf bauen, dass der erste vor dem zweiten gestartet wird, auch wenn dies mit einer hohen Wahrscheinlichkeit so sein wird.

Am Ende eines Master-Blocks wird nicht gewartet.

Bei verschachtelten parallel-Blöcken werden die IDs durch den innersten Block neu verteilt. In der Standardeinstellung erhalten alle Threads die ID 0 und die Anzahl der Threads wird als 1 gemeldet. Die Idee dahinter ist, dass sich alles auf den innersten Block bezieht. Dieser hat aber keinen Thread erzeugt und folglich gibt es nur den, der den Block angestoßen hat. Eine weitere Konsequenz ist, dass single- und master-Blöcke mehrfach ausgeführt werden. Ein verschachtelter parallel-Block wird demnach also in der Standardeinstellung sequentiell abgearbeitet. Mit num_threads(n) kann man das Erzeugen weiterer Threads erzwingen.

Die Problematik wird durch folgendes Beispiel gut illustriert.

int main(){
	#pragma omp parallel
	{
		#pragma omp parallel
		{
			printf("omp_get_thread_num() = %d,"
				" omp_get_num_threads() = %d\n", 
				omp_get_thread_num(), 
				omp_get_num_threads());	
		}
	}
}

3.3 Sections

Alternativ zu den single-Blöcken kann man auch sections verwenden.

int main(){
	#pragma omp parallel
	{
		#pragma omp sections
		{
			{
				printf("Thread %d war hier.\n", omp_get_thread_num());			
			}
			#pragma omp section
			{
				printf("Aber hier war Thread %d.\n", omp_get_thread_num());		
			}
			#pragma omp section
			{
				printf("und eine weitere section");		
			}
		}
		printf("Thread %d sagt Tschüss \n", omp_get_thread_num());
	}
}

Jede der drei Sektionen wird genau einmal ausgeführt. Am Ende des sections-Block warten die Threads, ehe sie fortfahren. Auch hier kann man das Warten durch nowait unterdrücken.

Man könnte obiges Beispiel auch mit single ausdrücken.

int main(){
	#pragma omp parallel
	{
		#pragma omp single nowait
		{
			printf("Thread %d war hier.\n", omp_get_thread_num());			
		}
		#pragma omp single nowait
		{
			printf("Aber hier war Thread %d.\n", omp_get_thread_num());		
		}
		#pragma omp single nowait
		{
			printf("und eine weitere section");		
		}
		printf("Thread %d sagt Tschüss \n", omp_get_thread_num());
	}
}

Bei der single-Variante ist der erste Block mit hoher Wahrscheinlichkeit der erste, der gestartet wird. Bei der sections-Variante sind alle Blöcke gleich wahrscheinlich.

3.4 Barrier

Wenn an einer beliebigen Stelle gewartet werden soll, kann man eine barrier-Anweisung benutzen.

int main(){
	#pragma omp parallel
	{
		printf("Thread %d sagt Hallo \n", omp_get_thread_num());
		#pragma omp barrier
		printf("Thread %d sagt Tschüss \n", omp_get_thread_num());
	}
}

In diesem Beispiel ist sichergestellt, dass alle Threads "Hallo" sagen, ehe sich der erste verabschiedet.

3.5 For

Mit Hilfe von OpenMP kann man auch einfache Schleifen parallelisieren.

int main(){
	#pragma omp parallel
	{
		#pragma omp for
		for(int i=0; i<10; ++i)
			printf("%d^2 = %d\n", i, i*i);
	}
}

Anders als bei einer normalen Schleife, werden die Durchläufe hier nicht unbedingt nacheinander abgearbeitet, sondern möglicherweise gleichzeitig und auch in zufälliger Reihenfolge. Dies führt dazu, dass die Ausgabe des obigen Codes nicht immer die gleiche ist. Die Reihenfolge der Zeilen kann beliebig permutiert sein.

Um dies zu bewerkstelligen, werden die Durchläufe möglichst gleichmäßig zwischen den Threads verteilt. Bei zwei Threads könnte z.B. der erste die Durchläufe mit 0<=i<=4 abarbeiten und der andere die restlichen mit 5<=i<=9. Die genaue Verteilung ist dem Compiler freigestellt. Sie hängt jedoch nur von der Anzahl der Threads und der Anzahl der Schleifendurchläufe ab, das heißt, wenn aus irgendeinem Grund ein Thread massive Verspätung hat, dann wird nicht kurzerhand die Arbeit umgeschichtet.

int main(){
	#pragma omp parallel
	{
		// Sorgt für massive Verspätung eines Threads
		#pragma omp single nowait
		{
			for(int i=0; i<INT_MAX; ++i)
				for(int j=0; j<INT_MAX; ++j)
					{}
		}

		// Tja hier warten dann die anderen...
		#pragma omp for
		for(int i=0; i<10; ++i)
			printf("%d^2 = %d\n", i, i*i);
	}
}

Obwohl der eine Thread sehr lange mit dem single-Block beschäftigt ist, kommt keiner der wartenden Threads auf die Idee, bei der for-Schleife Überstunden zu leisten. Dies ist nur mit einer dynamischen Arbeitsverteilung möglich.

int main(){
	#pragma omp parallel
	{
		#pragma omp single nowait
		{
			for(int i=0; i<INT_MAX; ++i)
				for(int j=0; j<INT_MAX; ++j)
					{}
		}

		#pragma omp for schedule(dynamic)
		for(int i=0; i<10; ++i)
			printf("%d^2 = %d\n", i, i*i);

		// Hier wird erst gewartet wenn die Schleife vollständig abgearbeitet wurde.
	}
}

Wenn ein Thread fertig ist, schaut er, ob es noch Arbeit gibt, und übernimmt diese dann auch gleich. Ein ähnliches Problem ergibt sich, wenn die Durchläufe unterschiedlich lang dauern.

Am Ende eines for-Blocks wird gewartet. Dies kann man wie üblich mit nowait unterdrücken.

Da schedule(dynamic) größere Verwaltungskosten als die normale Version hat, kann es sein, dass die normale Version in manchen Fällen schneller ist, obwohl nicht alle Threads immer beschäftigt sind. Die normale Version kommt nämlich ganz ohne Kontakt zwischen den einzelnen Threads zurecht und ist daher sehr schnell. Ich werde weiter unten auf konkrete Implementierungsmöglichkeiten eingehen.

Damit eine Schleife überhaupt parallelisiert werden kann muss sie eine Reihe von sehr strikten Kriterien erfüllen.

• Die Schleifenvariable muss eine Ganzzahl sein. In OpenMP 3.0 sind auch Random-Access-Iteratoren und Zeiger erlaubt.
• Die Abbruchbedingung muss ein Ordnungsvergleich (<,>, <= oder sein und auf der linken Seite muss sich die Schleifenvariable und nur die befinden. Die rechte Seite darf nicht während der Ausführung der Schleife verändert werden.
• Die Inkrementierung muss in immer gleich großen Schritten stattfinden. Dies bedeutet, dass man nur --, ++, += oder -= verwenden darf, um die Schleifenvariable zu verändern. Die Größe der Schritte darf nicht während der Ausführung der Schleife verändert werden.
• break und return sind verboten. continue ist erlaubt. goto darf nur verwendet werden, wenn man nicht aus der Schleife raus springt.

Folgende Schleifen lassen sich nicht mit OpenMP parallelisieren.

int n = foo();

// Auf der einen Seite der Bedingung muss ein einfaches i stehen.
for(int i=0; i*i<n; ++i)
	printf("%d", i);

// Die Inkrementierungsschritte sind nicht gleichgroß
for(int i=0; i<n; i*=3)
	printf("%d", i);

// Dies ist erlaubt, allerdings ist unspezifiziert, wie oft bar und foobar aufgerufen werden.
for(int i=-10; i<n*n*bar(); i+=foobar(n))
	printf("%d", i);

// Man darf nicht aus der Schleife herausspringen.
int a[10] = { ... }
for(int j=0; j<10; ++j)
	if(a[j] == n)
		return i;

// n darf nicht in der Scheife verändert werden.
for(int i=0; i<n*10; i+=n)
	++n;

3.5.1 Verschachtelte Schleifen mit und ohne Collapse

OpenMP-for-Schleifen dürfen nicht verschachtelt werden. Man kann allerdings mit einem zusätzlichen parallel-Block Abhilfe schaffen. Dies ist jedoch selten sinnvoll, da die Schleife zwar eine OpenMP-Schleife ist, aber dennoch sequentiell abgearbeitet wird.

#pragma omp parallel
{
	#pragma omp for
	for(int i=0; i<10; ++i)
	{
		#pragma omp parallel
		{
			#pragma omp for
			for(int j=0; j<10; ++j)
			{
				printf("%d*%d = %d\n", i, j, i*j);
			}
		}
	}
}

Es wird aber nur die äußere Schleife parallelisiert und die innere sorgt nur für Probleme und Kosten. Besser ist:

#pragma omp parallel
{
	#pragma omp for
	for(int i=0; i<10; ++i)
	{
		// Dieses for ist ein ganz normales C/C++ for
		for(int j=0; j<10; ++j)
		{
			printf("%d*%d = %d\n", i, j, i*j);
		}

	}
}

In OpenMP 3.0 gibt es eine spezielle collapse-Schleife, um über ein mehrdimensionales Array zu iterieren. Dies entspricht einem der häufigsten Einatzgebiete von mehrfach verschachtelten Schleifen. Die Schleifen müssen allerdings sehr einfach gehalten sein, damit eine Parallelisierung überhaupt erfolgen kann. Alle Schleifenvariablen müssen völlig unabhängig voneinander sein. Es darf nur die innerste Schleife beliebigen Code enthalten. Alle anderen dürfen nur verschachtelte Schleifen enthalten. Das Argument der collapse-Anweisung gibt an, wie tief die Schleifen verschachtelt sind.

Dies sind sehr starke Einschränkungen und deswegen deckt das folgende Beispiel an sich auch quasi jede erlaubte Variante einer collapse-Schleife ab.

int main(){
	#pragma omp parallel
	{
		int foo[10][3][7];
		#pragma omp for collapse(3)
		for(int l=0; l<10; ++l)
			for(int j=0; j<3; ++j)
				for(int k=0; k<7; ++k){
					// Nur hier darf beliebiger Code stehen!
					foo[l][j][k] = 0;
				}

}

3.5.2 Partielle sequentielle Abarbeitung mit ordered

Man kann in eine OpenMP-for Schleife einen ordered-Block packen. Dieser sorgt dafür, dass ein Teil der Schleife sequentiell abgearbeitet wird.

int main(){
	#pragma omp parallel
	{
		#pragma omp for ordered
		for(int i=0; i<=10; ++i){
			double v = i;
			for(int j=1; j<10; ++j)
				v *= i;
			#pragma omp ordered
			{
				printf("%d^10 = %.0f\n", i, v);
			}
		}
	}
}

Die Berechnung der Werte kann in beliebiger Reihenfolge und parallel erfolgen. Die Ausgabe ist jedoch sequentiell, das heißt, die Ausgabe ist immer die gleiche.

0^10 = 0
1^10 = 1
2^10 = 1024
3^10 = 59049
4^10 = 1048576
5^10 = 9765625
6^10 = 60466176
7^10 = 282475249
8^10 = 1073741824
9^10 = 3486784401
10^10 = 10000000000

Ein ordered-Block darf nur in einer speziellen ordered-for-Schleife auftreten. Jeder Durchlauf darf nur maximal einem ordered-Block begegnen. Folgendes Beispiel ist demnach falsch.

int main(){
	#pragma omp parallel
	{
		#pragma omp for ordered
		for(int i=0; i<10; ++i){
			#pragma omp ordered
			printf("%d^10 =", i);

			double v = i;
			for(int j=1; j<10; ++j)
				v *= i;

			// Hier wird etwas schiefgehen, da 
			// dies der zweite Block ist.
			#pragma omp ordered
			printf("%.0f\n", v);
		}
	}
}

Dies bedeutet jedoch nicht, dass man nur einen ordered-Block pro Schleife haben darf.

int main(){
	#pragma omp parallel
	{
		#pragma omp for ordered
		for(int i=0; i<10; ++i){
			if(i%3 == 0){
				#pragma omp ordered
				printf("%d ist ein Vielfaches von 3\n", i);
			}else{
				#pragma omp ordered
				printf("%d ist kein ein Vielfaches von 3\n", i);
			}
		}
	}
}

In dem Beispiel ist sichergestellt, dass kein Durchlauf zwei verschiedene Blöcke zu sehen bekommt.

3.6 Kurzschreibweise

Wenn ein parallel-Block nur eine for-Schleife enthält, dann gibt es eine Kurzschreibweise.

#pragma omp parallel for
for(int i=0; i<10; ++i)
{
	for(int j=0; j<10; ++j)
	{
		printf("%d*%d = %d\n", i, j, i*j);
	}

}

Analog gibt es auch einen parallel sections-Block.

Folgendes Beispiel berechnet alle Primzahlen kleiner als eine gegebene Zahl.

bool is_prime(unsigned n){
	for(unsigned i=2; i<n; ++i)
		if(n%i == 0)
			return false;
	return true;
}

int main(){
	unsigned n;

	printf("Geben sie eine Ganzzahl ein: ");
	scanf("%u", &n);

	// schedule(dynamic) ist notwendig da is_prime keine
	// konstante Laufzeit besitzt. 
	#pragma omp parallel for schedule(dynamic)
	for(int i=2; i<=n; ++i)
		if(is_prime(i))
			printf("%u ist prim\n", i);
}

4 Gemeinsamer Speicher

Die Zugriffsregeln für Variablen sind die, die man erwarten würde. Variablen leben so lange wie der Block, in dem sie definiert wurden. Sie können verdeckt werden und man kann auf Variablen einer höheren Ebene zugreifen.

int a;

int main(){
	int b;
	#pragma omp parallel
	{
		int c;
		// a und b werden geteilt, jeder Thread hat sein eigenens c
	}
}

Wenn ein Thread schreibend auf a oder b zugreift, dann ergeben sich sämtliche vom Multithreading her bekannte Probleme. OpenMP bietet aber leider nur wenig Revolutionäres, um dieses altbekannte Problem zu lösen.

4.1 Critical

Ein critical-Block ist ein Block, in dem sich immer nur genau ein Thread befindet. Wenn bereits ein Thread im Block ist, wartet jeder Thread, der hinein will, bis der Block wieder frei wird.

Betrachten Sie folgendes Beispiel:

int main(){
	#pragma omp parallel for
	for(int i=0; i<10; ++i){
		int n;
		#pragma omp critical
		{
			printf("Geben sie eine Zahl ein: ");
			scanf("%d", &n);
		}
		printf("%d * %d = %d\n", i, n, i*n);
	}
}

Hier wird ein critical-Block verwendet, um sicherzustellen, dass der Benutzer nie aufgefordert wird, zwei oder mehr Zahlen gleichzeitig einzugeben. Nimmt man das critical weg, dann kann es zu folgender Ausgabe kommen.

Geben sie eine Zahl ein: Geben sie eine Zahl ein:

Jeder critical-Block besitzt einen Namen. Wenn zwei Blöcke denselben Namen besitzen, kann nur genau ein Thread in einem der beiden sein. Zwei Blöcke mit verschiedenen Namen beeinflussen sich nicht. Das heißt, es ist durchaus möglich, dass zwei Threads in critical-Blöcken sind, sofern sie verschiedene Namen haben.

Wenn bei einem critical-Block kein Name angegeben wird, so wie im ersten Beispiel, besteht der Name aus dem leeren Wort. Alle leeren Wörter werden als gleich angesehen. Im Klartext bedeutet dies, dass keine zwei Threads gleichzeitig in irgendeinem namenlosen Block sein können.

Leider sind die Namen global und liegen in keinem Namenraum. Dadurch sind Namenskonflikte möglich. Es ist allerdings sehr unwahrscheinlich, dass dies andere Auswirkungen hat als einfach nur ein langsameres Programm, da lediglich ein Thread warten muss, obwohl dies nicht nötig gewesen wäre. Es ist sichergestellt, dass sich normale Symbolen wie Funktionsnamen nicht in die Quere kommen.

int main(){
	FILE*file = fopen("out.txt", "w");
	#pragma omp parallel for
	for(int i=0; i<10; ++i){
		int n;

		#pragma omp critical(console)
		{
			printf("Geben sie eine Zahl ein : ");
			scanf("%d", &n);
		}

		#pragma omp critical(console)
		{
			printf("%d * %d = %d\n", n, n, n*n);
		}

		#pragma omp critical(file)
		{
			fprintf(file, "%d * %d = %d\n", n, n, n*n);
			fprintf(file, "%d + %d = %d\n", n, n, n+n);
		}

		#pragma omp critical(console)
		{
			printf("%d + %d = %d\n", n, n, n+n);
		}
	}
	fclose(file);
}

In diesem Beispiel ist sichergestellt, dass nichts auf die Konsole geschrieben wird, falls der Benutzer aufgefordert wurde, etwas einzugeben. In der Datei müssen beide zueinander passenden Rechnungen hintereinander ausgegeben werden. Auf der Konsole ist dies nicht garantiert.

4.2 Flush

Viele Rechner besitzen Register und mehrere Level von Speicher-Caches, um den Zugriff auf Variablen zu beschleunigen. Anstatt direkt mit der Variablen im Hauptspeicher zu arbeiten, wird diese zuerst in einen schnelleren Cache-Speicher kopiert. Anschließend wird mit der Kopie gerechnet und nach einiger Zeit wird die Variable in den Hauptspeicher zurück kopiert. Dies führt dazu, dass es mehrere Versionen von der gleichen Variable gibt, welche nicht immer gleich sind.

In einem sequentiellen Programm stellt dies kein Problem dar, da zu jedem Zeitpunkt bekannt ist, wo sich die aktuelle Version der Variable befindet. Jedoch ist dies bei mehreren Threads nicht mehr so einfach. Der eine weiß nicht, wann der andere fertig ist und wann er die Werte synchronisieren muss.

Um dieses Problem zu lösen, bietet OpenMP die flush-Anweisung an. Sie synchronisiert die lokale Version des Threads mit der des Hauptspeichers. Dies bedeutet, dass jeder Thread, nachdem er fertig ist, eine geteilte Variable zu verändern, diese auch flushen muss. Ein Thread, der lesend auf eine Variable zugreift will, muss die Variable auch flushen, damit er sicher ist, die neuste Version zu Gesicht zu bekommen. Wenn ein Thread auch mit einer veralteten Version arbeiten kann, braucht er die Variable nicht zu flushen. Wenn sicher ist, dass eine Variable in keinem parallel-Block verändert wird, dann kann man sich das flush auch sparen. Es ist nämlich egal, wann ein Thread seine lokale Version aus dem Hauptspeicher lädt, da er immer die gleichen Daten lesen würde.

Folgendes Programm berechnet die Fakultät einer Zahl.

int main(){
	unsigned n, fact = 1;
	printf("Geben sie eine Zahl ein: ");
	scanf("%u", &n);

	#pragma omp parallel
	{
		unsigned local_fact = 1;

		// Auf n kann man ohne flush zugreifen,
		// da es nicht verändert wird.
		#pragma omp for
		for(int i=1; i<=n; ++i){
			// local_fact braucht nicht geflusht 
			// zu werden, da es nicht geteilt wird.
			local_fact *= i;
		}

		#pragma omp critical
		{
			// Zuerst holt der Thread sich die neuste Version
			// aus dem Hauptspeicher.
			#pragma omp flush(fact)

			// Danach wird die Variable verändert.
			fact *= local_fact;

			// Anschließend wird die Variable wieder in den
			// Hauptspeicher geschrieben, damit andere Threads
			// sich die neuste Version von dort laden können.
			#pragma omp flush(fact)
		}	
	}

	printf("%u! = %u\n", n, fact);
}

Ein fehlendes flush bedeutet nicht, dass die Variablen nicht synchronisiert werden. Es steht der Implementierung frei, dies zu tun, wann immer sie will. Dies führt leider oft zu nur sehr schwer reproduzierbaren Bugs. Zum Beispiel ist es sehr wahrscheinlich, dass das Betriebssystem die Rechnerzeit anders auf die Threads verteilt, wenn die Hardware unter voller Last läuft, als wenn es nur mit kleinen Texteingaben konfrontiert wird. Es ist daher sinnvoll, im Zweifelsfall lieber ein flush zuviel zu schreiben.

4.3 Reduction

Mit reduction lassen sich Berechnungen vom Typ s += f(1)+f(2)+...+f(n) parallelisieren. Jeder Thread rechnet ein paar Terme aus und summiert die auf, welche er berechnet hat. Nachdem alle Threads fertig sind, zählt jeder Thread sein Zwischenergebnis zur Variable s hinzu.

Diese Zwischenergebnis wird in einer Variable mit dem gleichen Namen wie s gespeichert. Diese Hilfsvariable wird im reduction-Block angelegt und mit dem neutralen Element (0 im Fall der Addition) initialisiert. Man sollte sich aber nicht auf dieses Verhalten verlassen. Wenn ein Programm sequentiell übersetzt wird, das heißt, ohne OpenMP-Unterstützung, dann gibt es keine Zwischenergebnisse und deswegen hat s beim ersten Durchlauf auch nicht unbedingt den Wert des neutralen Elements.

Neben + lassen sich auch -, *, &, |, ^, && und || auf diese Weise parallelisieren. Es muss sich bei den Operatoren jedoch um ihre Buildin-Version handeln.

• +, -, | und ^ besitzen das neutrale Element 0.
• * besitzt das neutrale Element 1.
• & besitzt das neutrale Element ~0. (Bei ~0 sind alle Bits gesetzt.)
• && besitzt das neutrale Element true.
• || besitzt das neutrale Element false.

Bei Gleitkommazahlen hängt die Genauigkeit von der Auswertungsreihenfolge ab. Die Reihenfolge ist nicht spezifisiert, also auch nicht die Genauigkeit des Ergebnisses. Es ist nicht einmal garantiert, dass es immer das gleiche Ergebnis sein wird. Gleitkommazahlen in Kombination mit einer Reduktion haben einen gewissen nichtdeterministischen Teil im Ergebnis.

In der Fortran-Version von OpenMP gibt es zusätzlich noch Min- und Max-Reduktionen. Die wurden aus mir unerklärlichen Gründen in der C/C++ Version vergessen.

Folgendes Beispiel bestimmt, ob eine Zahl prim ist.

int main(){
	int n;

	printf("Geben sie eine Ganzzahl ein: ");
	scanf("%d", &n);

	// Der Initialwert von is_prime ist wichtig, da dieser auch
	// in die Berechnung einbezogen wird.
	bool is_prime = true;

	#pragma omp parallel reduction(&&: is_prime)
	{
		// jeder bekommt sein eigenes is_prime,
		// welches mit true, dem neutralen Element von &&,
		// initialisiert wurde.

		#pragma omp for
		for(int i=2; i<n; ++i)
			if(n%i == 0)
				is_prime = false;

		// Alle is_prime (auch das aus der main Funktion) werden 
		// mittels && verknüpft. Im Klartext bedeutet dies, dass, 
		// falls eines false ist, das Resultat false ist.
	}

	if(is_prime)
		printf("%d ist prim\n", n);
	else
		printf("%d ist nicht prim\n", n);

}

Am Ende des parallel Blocks wird folgende Rechnung ausgeführt:

is_prime = is_prime && thread[0].is_prime && thread[1].is_prime && ...  && thread[thread_num-1].is_prime

Bei der thread[j]-Notation handelt es sich legendlich um Pseudocode.

Auch hier gibt es eine Kurzschreibweise. Folgendes Beispiel berechnet die Fakultät einer Zahl.

int main(){
	unsigned n;

	printf("Geben sie eine Ganzzahl ein: ");
	scanf("%u", &n);

	unsigned fact = 1;

	#pragma omp parallel for reduction(*: fact)
	for(int i=2; i<n; ++i)
		fact *= i;

	printf("%u! = %u\n", n, fact);
}

reduction erlaubt es dem Compiler, eine Reihe von Optimierungen durchzuführen, welche bei der äquivalenten Formulierungen mit critical und flush nur sehr schwer möglich sind.

4.4 Locks/Mutex

OpenMP bietet auch Mutexe an. Hier weicht OpenMP an sich nicht von allgemein bekannten Konzepten ab. Es gibt einen Lock-Type namens omp_lock_t. Jeder Mutex ist entweder offen oder geschlossen. Folgende Funktionen werden bereitgestellt, um damit zu arbeiten:

• void omp_init_lock(omp_lock_t*)
  Initialisiert den Mutex. Der Mutex ist danach offen.
• void omp_destroy_lock(omp_lock_t*)
  Zerstört den Mutex.
• void omp_set_lock(omp_lock_t*)
  Versucht, den Mutex zu schließen. Wurde der Mutex bereits von einem anderen Thread geschlossen, so wird gewartet, bis dieser ihn wieder geöffnet hat.
• void omp_unset_lock(omp_lock_t*)
  Öffnet den Mutex wieder
• int omp_test_lock(omp_lock_t*)
  Versucht, den Mutex zu schließen. Wurde der Mutex bereits von einem anderen Thread geschlossen, so wird 0 zurückgegeben, ansonsten 1.

omp_lock_t ist ein nicht rekursives Mutex. Das heißt, ein Thread darf einen Mutex nicht mehrmals schließen. Er würde in dem Fall nämlich auf sich selbst warten. Es gibt aber auch eine rekursive Variante und die heißt omp_nest_lock_t. Die Funktionen, um damit umzugehen, heißen ähnlich und funktionieren analog.

• void omp_init_nest_lock(omp_next_lock_t*)
• void omp_destroy_nest_lock(omp_next_lock_t*)
• void omp_set_nest_lock(omp_next_lock_t*)
• void omp_unset_nest_lock(omp_next_lock_t*)
• int omp_test_nest_lock(omp_next_lock_t*)

Da diese Schnittstelle alles anderes als modernes C++ ist, bietet es sich an, diese Funktionen zu wrappen.

class Mutex
{
public:
	Mutex()				{omp_init_lock(&lock);}
	~Mutex()			{omp_destroy_lock(&lock);}
	void lock()			{omp_set_lock(&lock);}
	void unlock()			{omp_unset_lock(&lock);}
	bool try_to_lock()		{return omp_test_lock(&lock);}
private:
	Mutex(const Mutex&);
	Mutex&operator=(const Mutex&);
	mp_lock_t lock;
};

class ScopedLock
{
public:
	ScopedLock(Mutex&lock)		:lock(lock){lock.lock();}
	~ScopedLock()			{lock.unlock();}
private:
	ScopedLock(const ScopedLock&);
	ScopedLock&operator=(const ScopedLock&);
	Mutex&lock;
};

5 Zeitmessung

Neben den Funktionen, welche direkt mit der Multithread-Programmierung zusammenhängen, bietet OpenMP noch zwei nützliche Funktionen zur Zeitmessung an.

• double omp_get_wtime();
  Gibt die momentane Zeit in Sekunden zurück. Es ist nicht definiert, wann der Zeitpunkt 0 ist. Es ist legendlich bekannt, dass es nicht während dem Lauf eines Programms verändert wird und in der Vergangenheit liegt. Dies ist aber nicht weiter schlimm, wenn man die Funktion nur für Zeitmessungen einsetzen will. Die Genauigkeit dieser Funktion ist nicht spezifisiert und darf leicht von Thread zu Thread abweichen. Man kann die Genauigkeit aber mit omp_get_wtick in Erfahrung bringen.
• double omp_get_wtick();
  Gibt den kleinstmöglichen Abstand zwischen zwei verschiedenen von omp_get_wtime zurückgegebenen Zeitpunkten an. Der Rückgabewert von omp_get_wtick darf von Thread zu Thread variieren.

6 Suchen in einem ungeordneten Array

Das Suchen in einem unsortierten Array lässt sich sehr leicht mit OpenMP beschleunigen.

// Ohne OpenMP
int find_first(int arr[], int size, int value){
	for(int j=0; j<size; ++j)
		if(arr[j] == value)
			return j;
	return -1;
}

// Mit OpenMP
int find_any(int arr[], int size, int value){
	int i = -1;
	#pragma omp parallel for
	for(int j=0; j<size; ++j)
		if(arr[j] == value){
			i = j;
			#pragma omp flush(i)		
		}
	return i;
}

Das aktuelle i braucht nicht zuerst per flush-Anweisung geladen zu werden, da der Wert ohnehin nur überschrieben werden würde.

Wie man bereits am Namen erkennt, ist bei der parallelen Version nicht definiert, welcher Wert gefunden wird, wenn es mehrere gleiche Werte im Array gibt. Wenn dies jedoch nötig ist, kann man dieses Problem leicht mit ordered beheben.

int find_first(int arr[], int size, int value){
	int i = -1;
	#pragma omp parallel for ordered
	for(int j=size-1; j>=0; --j)
		if(arr[j] == value){
			#pragma omp ordered
			i = j;
			#pragma omp flush(i)		
		}
	return i;
}

Eines der Nachteile der parallelen Version ist, dass immer alle Elemente angeschaut werden. Dies ist aber nicht immer notwendig, denn man kann die Schleife abbrechen, sobald ein Element gefunden wurde. OpenMP bietet allerdings keine Möglichkeit, dies zu tun. Man kann versuchen, die Schleife nachzubauen.

int find_any(int arr[], int size, int value){
	int i = -1;
	#pragma omp parallel
	{
		const int begin = size*omp_get_thread_num() / omp_get_num_threads();
		const int end = size*(omp_get_thread_num()+1) / omp_get_num_threads();

		for(int j=begin; j<end; ++j){
			if(arr[j] == value)
				i = j;
			#pragma omp flush(i)
			if(i != -1)
				break;
		}
	}
	return i;
}

Dieser Code dürfte recht nahe an dem sein, was der Compiler aus einer OpenMP-Schleife mit statischer Verteilung macht. In der Regel sollte man aber die OpenMP-for-Schleife bevorzugen. Diese ist weit weniger fehleranfällig. So kann es beim Beispiel oben zum Beispiel zu Überläufen bei der Berechnung von begin und end kommen.

Ein weiteres Probleme ist, dass jeder Durchlauf den Wert von i laden muss. Dies kann je nach Architektur sehr langsam sein und somit den potentiellen Geschwindigkeitsgewinn wieder zunichtemachen. Um diesm zweifelsfrei festzustellen, muss man die Zeit messen.

7 Akkumulation

Die Akkumulation ist die Standardanwendung der reduction-Anweisung. Es geht darum, alle Elemente einer Liste auf irgendeine Weise zu kombinieren. Das Addieren aller Zahlen einer Liste ist ein Beispiel, welches mit einer solchen Anweisung gelöst werden kann.

int accumulate(int arr[], int size){
	int sum = 0;
	#pragma omp parallel for reduction(+:sum)
	for(int i=0; i<size; ++i)
		sum += arr;
	return sum;
}

Leider handelt es sich jedoch nicht immer um einen Buildinoperator. Sei [i]Matrix* eine Klasse, welche eine 2-dimensionale Matrix beschreibt und für welche der *= Operator entsprechend überladen wurde. Der Defaultkonstruktor initialisiert die Matrix auf die Einheitsmatrix.

class Matrix{
private:
	int 
		a,b,
		c,d
	;
public:
	Matrix():a(1), b(0), c(0), d(1){}
	Matrix&operator*=(const Matrix&right){
		Matrix left = *this;
		a = left.a*right.a + left.b*right.c;
		b = left.a*right.b + left.b*right.d;
		c = left.c*right.a + left.d*right.c;
		d = left.c*right.b + left.d*right.d;
		return *this;
	}
};

Matrix multiply(Matrix arr[], int size){
	const unsigned thread_count = omp_get_max_threads();
	Matrix*sub_prod = new Matrix[thread_count]; 
	#pragma omp parallel num_threads(thread_count)
	{
		const unsigned id = omp_get_thread_num();
		const int begin = size*id/thread_count;
		const int end = size*(id+1)/thread_count;

		for(int j=begin; j<end; ++j)
			sub_prod[id] *= arr[j];
	}
	Matrix prod;
	for(int j=0; j<thread_count; ++j)
		prod *= sub_prod[j];
	delete[]sub_prod;
	return prod;
}

Die Idee ist, die Assoziativität der Multiplikation auszunutzen, um das Produkt in mehrere Teilprodukte zu zerlegen und diese zuerst von jedem Thread einzeln auswerten zu lassen.

Zuerst wird bestimmt, wie viele Threads überhaupt zur Verfügung stehen. Anschießend wird ein Array angelegt, welches einen Eintrag für jeden Thread besitzt. In diesem Array werden die Zwischenresultate abgelegt. Die Matrizenmultiplikation ist nicht kommutativ, daher verwende ich zur Berechnung der Teilprodukte keine OpenMP-for-Schleife. Diese garantiert nämliche keine bestimmte Reihenfolge. Eine ordered OpenMP-for-Schleife ist zu steif, da nicht alle Multiplikationen sequentiell ablaufen müssen. Nachdem alle Teilprodukte berechnet wurden, werden sie multipliziert und das Ergebnis wird zurückgegeben.

Ich habe keinen STL-vector verwendet, da dieser in der Regel nicht dafür ausgelegt ist, mit mehreren Threads klar zu kommen. In diesem speziellen Fall ist es jedoch höchst unwahrscheinlich, dass dies ein Problem darstellen würde. Man könnte natürlich einen Smartpointer verwenden, um die Speicherfreigabe zu garantieren. Allerdings ist das im Fall einer Ausnahme das kleinere Problem.

Die einzigen Operationen, die eventuell eine Ausnahme werfen könnten, sind die Speicherallokierung selbst und der *= Operator. Wenn die Speicherallokierung fehlschlägt, geht kein Speicher verloren. Der *= Operator unserer Matrixklasse wirft natürlich keine Ausnahmen, allerdings ist dies bei einer allgemeinen Matrixklasse nicht mehr selbstverständlich. In diesem Fall muss man allerdings sicherstellen, dass die Ausnahme sicher aus dem parallel-Block hinaus propagiert wird. Dies ist ein wesentlich größeres Problem als das kleine Speicherleck.

7.1 STL

Im Forum kommt es immer wieder zu Diskussionen zum Thema Parallelisierung und wie man die einsetzt, ohne Programme von Grund auf neu zu bauen. Oft kommt dann die Idee auf, die Algorithmen der STL zu parallelisieren und den Benutzer-Code unverändert zu lassen. Dies ist allerdings schwer umzusetzen. Im Folgenden werde ich zeigen, wie man std::accumulate parallelisieren kann und welche Probleme es dabei gibt.

Es sei folgende recht einfache Funktion gegeben:

template<class Iter, class T>
T accumulate(Iter begin, Iter end, T sum){
	for(; begin!=end; ++begin)
		sum += *begin;
	return sum;
}

Unser Ziel ist es, diese Funktion durch Benutzung von OpenMP zu beschleunigen. Auf den ersten Blick erscheint dies wie eine sehr einfache Aufgabe, dies ist sie aber nicht!

Zuerst müssen wir sicherstellen, dass die zu Grunde liegende +-Operation parallel ausgeführt werden kann. Speichert der Operator zum Beispiel Zwischenwerte in einem statischen Array, dann können wir eine Parallelisierung, ohne die Operation selbst zu verändern, bereits vergessen. Es ist auch möglich, dass Werte vorberechnet werden müssen und der Operator auf diese zugreift, das heißt, es muss etwas initialisiert werden. Zum Beispiel sind eine Primzahlliste oder das Pascalsche Dreieck zur Bestimmung großer Binomialkoeffizienten gute Kandidaten. Dies muss vor dem ersten Aufruf passiert sein. Es darf nicht erst bis zur ersten Verwendung gewartet werden, ehe man Initialisierungsarbeit erledigt.

Die +-Operation darf keine Ausnahme werfen, da C++ leider keine Möglichkeit bietet, diese über Threadgrenzen hinaus zu propagieren.

Des Weiteren benötigt man die Assoziativität der +-Operation. Die Assoziativität ist notwendig, damit die Summe in Teilsummen aufgespalten werden kann. Die Kommutativität ist notwendig, da nicht klar ist, in welcher Reihe die Threads fertig werden. Dies erscheint offensichtlich, allerdings ist das eine zusätzliche Einschränkung gegenüber der sequentiellen Version. Oft ist sie auch nicht gegeben. Zum Beispiel ist die +-Operation von Gleitkommazahlen nicht assoziativ da die Genauigkeit von der Auswertungsreihenfolge abhängt. Normalerweise sagt man bei Gleitkommazahlen jedoch, dass die Genauigkeit nicht wichtig ist und nimmt einfach an, dass die Assoziativität gegeben ist.

Ohne die Kommutativität kann man auskommen, wenn es sein muss, wie das vorherige Beispiel zeigt. Wenn sie allerdings zur Verfügung steht wird die Implementierung einfacher. Ich gehe im Folgenden davon aus, dass sie gegeben ist.

C++-Templates bieten leider praktisch keine Möglichkeiten, diese Bedingungen sicherzustellen. Die Dokumentation der Funktion ist also die einzige Stelle, wo man sie angeben kann. Es bleibt nur zu hoffen, dass der Benutzer der Funktion die Dokumentation auch gründlich liest und versteht...

Ist sichergestellt, dass die +-Operation parallel ausgeführt werden kann, müssen wir zusehen, dass wir wahlfreien Zugriff auf die Liste von Elementen haben. Dies ist nötig, da es einem OpenMP-for freisteht, in welcher Reihenfolge er die Schleifenduchläufe ausführen will. Ist dies nicht gegeben, kann man entweder die sequentielle Version benutzen oder die Elemente zu erst kopieren. Diese Kopie muss sequentiell erfolgen, da es keinen wahlfreien Zugriff gibt. Was besser ist, hängt wieder von der +-Operation ab.

Bei Ganzzahlen ist die sequentielle Version deutlich überlegen, da die Addition wahrscheinlich so schnell ist wie das Kopieren einer Zahl. Das einfache Ausrechnen ist also in etwa genauso schnell wie das Vorbereiten der parallelen Berechnung. Bei der Matrizenmultiplikation ist es allerdings besser, zu kopieren. Werden Matrizen unterschiedlicher Dimensionen multipliziert, so lohnt es sich unter Umständen viel mehr, zuerst dynamisch die optimale Auswertungsreihenfolge zu bestimmen als einfach nur zu parallelisieren. Ich werde im Folgenden davon ausgehen, dass es sich lohnt, zu kopieren.

Der Code sieht bis hierhin in etwa folgendermaßen aus:

template<class Iter, class T>
T accumulate(Iter begin, Iter end, T init){
	return accumulate_check_category(begin, end, init, std::iterator_traits<Iter>::iterator_category());
}

template<class Iter, class T>
T accumulate_check_category(Iter begin, Iter end, T init, std::random_access_iterator_tag){
	return accumulate_parallel(begin, end-begin);
}

template<class Iter, class T>
T accumulate_check_category(Iter begin, Iter end, T init, std::input_iterator_tag){
	typename std::iterator_traits<Iter>::value_type T;
	std::vector<T>buffer(begin, end);
	return accumulate_parallel(buffer.begin(), buffer.size());
}

template<class Iter, class T>
T accumulate_parallel(Iter begin, int len, T init){
	// ... 
}

Die eigentliche Parallelisierung ist ein typischer Reduktionsfall. Die entsprechende OpenMP-Anweisung kann aber nur im Fall eines buildin-Typen verwendet werden. Da durch die Parallelisierung Kosten entstehen, lohnt es sich nicht, unter einer Gewissen Grenze zu parallelisieren. Wie viele Element es geben muss, damit sich das Ganze lohnt, hängt von der +-Operation, vom Prozessor und dem verwendeten Compiler ab.

Die OpenMP-Anweisung dürfte bei einigen Compilern schneller sein als ihr Nachbau. Beim Nachbau muss man zwischen Nachbau mit Kommutativität und ohne Kommutativität unterscheiden. Bei der OpenMP-Anweisung braucht man dies nicht, da es keinen nicht kommutativen buildin +-Operator gibt. Es gibt also mindestens drei Konstanten, welche man praktisch nur empirisch bestimmen kann, die man in den Code einbauen muss. Da ich die Kommutativität angenommen habe, sind es bei mir nur zwei. Ich nehme einfach mal 50000 für beide an.

Ein weiteres Problem, das sich stellt, ist die Frage, ob die Laufzeit des Operator+ von seinen Operanden abhängt oder nicht. Ob die for-Anweisung ein schedule(dynamic) benötigt oder nicht, hängt davon ab. Will man alle Fälle abdecken, so muss man jetzt schon sechs Konstanten bestimmen. Ich gehe mal davon aus, dass es in etwa immer gleich schnell ist. Ein Gegenbeispiel wäre die Multiplikation von Matrizen unterschiedlicher Dimensionen.

template<class Iter, class T>
T accumulate_parallel(Iter begin, int len, T init){
	accumulate_parallel_reduction(begin, len, init, boost::is_arithmetic<typename std::iterator_traits<T>::value_type>::type());
}

template<class Iter, class T>
T accumulate_parallel_reduction(Iter begin, int len, T sum, boost::true_type){
	if(len < 50000){
		for(int i=0; i<len; ++i)
			sum += *(begin+i);
	}else{
		#pragma omp parallel reduction(+:ret) for
		for(int i=0; i<len; ++i)
			sum += *(begin+i);
	}
	return sum;
}

template<class Iter, class T>
T accumulate_parallel_reduction(Iter begin, int len, T sum, boost::false_type){
	if(len < 50000){
		for(int i=0; i<len; ++i)
			sum += *(begin+i);
	}else{
		#pragma omp parallel 
		{
			T sub_sum = T();
			for(int i=0; i<len; ++i)
				sub_sum += *(begin+i);
			#pragma omp critical
			{
				#pragma omp flush(sum)
				sum += sub_sum;
				#pragma omp flush(sum)
			}
		}
	}
	return sum;
}

Ich habe auf meinem 4-Kern-Prozessor und dem GCC-4.2 ein paar Tests gemacht, um zu bestimmen, wie groß bei mir diese Konstanten liegen. Unter etwa 50000 int-Elementen ist die sequentielle Version bei mir im Vorteil. Bei mehr Elementen ist die parallelisierte Version besser. Erstaunlich ist, dass die Reduktions-Answeisung in etwa genau so schnell ist wie ihr Nachbau. Ich gehe mal davon aus, dass die GCC-Entwickler den OpenMP-Code noch nicht wirklich optimiert haben und in dieser Version vor allem auf korrektes Verhalten Wert gelegt haben.

8 Gibt es noch mehr?

In diesem Artikel habe ich nicht alle OpenMP-Anweisungen vorgestellt. Ich habe mich legendlich auf die mir am wichtigsten erscheinenden beschränkt. Wer OpenMP wirklich produktiv einsetzen will, dem empfehle ich, die Spezifikation durchzulesen. Für eine Spezifikation ist diese erstaunlich leicht verständlich und es sind auch viele Beispielanwendungen darin enthalten.

9 Quellen

• http://bisqwit.iki.fi/story/howto/openmp/
• http://www.openmp.org/
• http://algo2.iti.uni-karlsruhe.de/singler/mcstl/

[paracuter]

hey vielen vielen dank dafür ist echt sehr interessant!

[Kenner des Autor]

und wieder mal ein Spitzenartikel! Weiter so!

[SAn]

Excellent!

Can I use this material to teach my students (in russian?)
Can I translate this into russian and put a reference to the original article?

[Ben04]

SAn schrieb:

Excellent!

Can I use this material to teach my students (in russian?)
Can I translate this into russian and put a reference to the original article?

Yes, as long as you state the original sources. (Link to http://magazin.c-plusplus.net/artikel/OpenMP)

[Nanyuki]

Ausgezeichneter Artikel. Da g++ jetzt auch OpenMP unterstützt, habe ich bereits nach etwas gesucht, was mir einen leichten Einstieg verschaffen könnte.

[DeLoghi]

Geniale Aufarbeitung des Themas. Weiter so... Hier noch der direktlink zum Download: http://openmp.org/wp/openmp-specifications/

[Petrechen_88]

Ohne den Artikel ganz gelesen zu haben, muss wenn ich verschachtelte for Schleifen parallelisieren will die Anzahl der Iterationen schon zur kompilierzeit feststehen oder reicht zur Laufzeit?

[Walli]

Sie muss vor Eintritt in die Schleife feststehen.

[Peterchen_88]

Hallo,

ich will ueber ein 2D Array itertieren und die Werte addieren, aber das Ergebnis st falsch, lässt sich das überhaupt mit OpenMP parallelisieren?

void foo(int dim1, int dim2)
{
    double sum = 0;
    #pragma omp parallel reduction(+: sum) 
    { 
		#pragma omp parallel 
		{	
			#pragma omp for 
			for(int i=0;i<dim1;i++)
			{
				for(int j=0;j<dim2;j++)
				{
					sum += a[i][j];
				}
			}
		}
	}

[Ben04]

Lass das innere parallel weg

[freaksdd]

hi,
ich versuche mich gerade in openMP rein zu lesen aber zu meinem Problem komme ich auf keine wirkich gute Lösung.

ich habe eine gemeinsame Liste die ich für 2 Prozesse nutzen möchte. Nach Programmstart kommen Funktionsaufrufe. Die Aufrufe sollen unabhängig von einander ausgeführt werden und während dessen auf die gemeinsame Liste zugreifen.

ich habe es schon mit "section/-s" und "single nowait" probiert und mir immer die Prozess nummer ausgeben lassen. es kam nur selten vor das die Aufteilung von zwei verschiedenen Prozessen funktioniert hat. ansonsten hat immer der Master-Prozess (0) die Funktionen ausgeführt.

Wie bekomme ich eine Trennung auf sichere Art und Weise hin?

Mfg Peter