并行程序设计之Openmp

2022-10-23 19:41:04

Openmp复习

下载安装Openmp
第一个omp程序

Hello World
编译指令
注意事项

Gauss Elimination

下载安装Openmp

操作系统; ubuntu 16.04
处理器: Core I3
命令: 不需安装, gcc 支持omp

第一个omp程序

Hello World

#include <stdio.h>
#include <stdlib.h>
#include <string.h>

#include <math.h>
#include "omp.h"

const int MAXN = 1000;

int main(int argc, char* argv[]) {
	int thread_count = strtol(argv[1], NULL, 10);

	int i, j, k;
	int a[MAXN][MAXN];
	int n = MAXN;
	double start, finish, duration;
	start = omp_get_wtime();
#	pragma omp parallel num_threads(thread_count) \
	shared(a) \
	private(i, j, k)
	for (i = 0; i < n; i++) {
		for (j = i; j < n; j++) {
#			pragma omp for
			for (k = i; k < n; k++) {
				a[i][j] = a[k][k] * a[i][i];
			}
		}
	}
	finish = omp_get_wtime();
	duration = finish - start;
	printf("the total time: %f\n", duration);
	return 0;
}

编译指令

openmp是POSIX标准, 许多库都支持它, 所以编译时我们只需要链接库就可以

gcc -fopenmp demo.c -o demo

运行时与常规c语言程序同样运行即可.

注意事项

openmp使用很方便, 但是功能单一, 只能对已有代码进行改动, 而且只能是for循环的代码, 在需要改动的地方添加parallel for预编译指令即可, 另一方面应当注意对于多重for循环, 不能在循环里面直接加预编译指令, 而是在最外层循环的外面加线程分配指令, 这样可以避免线程分配和销毁的额外开销. 最后, for循环内部最好不要用分支指令, 因为分支指令涉及到硬件层面的分支预测, 可能导致cache未命中等一些特殊情况, 导致openmp无法加速.

Gauss Elimination

高斯消元法在数学上的实现非常简单, 首先将参数矩阵上三角化, 再将上三角矩阵对角化, 就可以得到最后的结果, 但是我们应当注意几个问题, 在最后一步对角化时挖掘算法的并行性时有些类似于cannon算法, 是整体上的并行化, 在最后一步得到结果, 不能像正常的串行化算法一样, 一步得到结果. 另外, if条件句不要有, 条件的测试可以放到for循环的外面去做.

最后我们来分析以下高斯消元法的复杂度, 时间复杂度就是O(nnn/p), p是并行线程数, 对于我的电脑来说, I3的处理器p可以认为是2, 因为线程间通信采用的是细粒度线程的方式, 所以通信开销相对mpi来说要小一些, 可以认为加速比可以达到p. 另一方面, 空间复杂度就是O(n*n), 因为每个线程仅分配大矩阵的一部分.

代码如下

#include <stdlib.h>
#include <string.h>
#include <stdio.h>

#include <math.h>
#include "omp.h"

const int MAXN = 1000;

int main(int argc, char* argv[]) {
	int thread_count = strtol(argv[1], NULL, 10);

	int i, j, k; 
	int tmp;
	int n = MAXN;
	double a[MAXN][MAXN], b[MAXN];
	double start, finish, duration;
	start = omp_get_wtime();
#	pragma omp parallel num_threads(thread_count) \
	shared(a, b) \
	private(i, j, k, tmp)
	for (k = 0; k < n; k++) {
		for (tmp = k; tmp < n; tmp++) {
			a[k][tmp] /= a[k][k];
		} 
		b[k] /= a[k][k];
		a[k][k] = 1.0;
#		pragma omp for
		for (i = k + 1; i < n; i++) {
			for (j = k + 1; j < n; j++) {
				a[i][j] -= a[k][j] * a[i][k];
			}
			a[i][k] = 0.0;
		}
	}

#	pragma omp parallel num_threads(thread_count) \
	shared(a, b) \
	private(k, j) 
	for (k = n - 1; k >= 0; k--) {
#		pragma omp for
		for (j = k + 1; j < n; j++) {
			b[k] -= a[k][j] * b[j];
		}
	}
	finish = omp_get_wtime();
	duration = finish - start;
	printf("the total time: %f\n", duration);
	return 0; 
}

挖掘并行性主要在第一部分矩阵由常规的矩阵变为上三角矩阵, 另一部分在上三角矩阵变为对角矩阵, 最终试验结果如图

可以看出当我使用2个线程时, 基本达到了理想加速比, 这是大小为1000维的矩阵, 当我使用4个线程时, 运行速度反而降低了, 这是因为电脑是I3处理器, 最多只支持两个节点的并行.

码农公寓