OpenMP與MPI聯合程式設計

研究一下如何一起使用mpi和openmp

先上程式：

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

#ifdef _OPENMP
#include <omp.h>
#endif

#ifdef MPICH_CC
  #include <mpi.h>
#endif

#define SIZE 512
#define SIZEM SIZE
#define SIZEN SIZE
#define SIZEK SIZE

#define DISPLAYT 1

double mytic() 
{
	double result = 0.0;
	struct timeval tv; 
	gettimeofday( &tv, NULL );
	result = tv.tv_sec*1e6 + tv.tv_usec;
	return result;
}

int main(int argc, char** argv)
{
#ifdef MPICH_CC
	int mpirank;
	int mpisize;
	MPI_Status mpistatus;
#endif

	float *ma;
	float *mb;
	float *mc;
	int m, n, k;
	int r, s, t;
	float tsum;
	double tstart, tend;

#ifdef MPICH_CC
	MPI_Init(&argc, &argv);
	MPI_Comm_rank(MPI_COMM_WORLD, &mpirank);
#endif

#ifdef _OPENMP
	omp_set_num_threads(2);
#endif

	/*Matrx A, B, C allocation.*/
	m = SIZEM; n = SIZEN; k = SIZEK;

	ma = (float*)malloc(m*n*sizeof(float));
	if (!ma){return -1;}
	mb = (float*)malloc(n*k*sizeof(float));
	if (!mb){free(ma);return -1;}
	mc = (float*)malloc(m*k*sizeof(float));
	if (!mc){free(ma);free(mb);return -1;}

	/*Matrx A, B, C initialization.*/
	srand((unsigned int)time(0));
	tstart = mytic();
	for (r = 0; r < m; r++)
	{
		for (s = 0; s < n; s++)
		{
			ma[r*n+s] = (1.0*rand())/RAND_MAX;
		}
	}
#if DISPLAYT 
	printf("[T] ma initialization time : %.3fms\n", 1e-3*(mytic()-tstart));
#endif
	tstart = mytic();
	for (r = 0; r < n; r++)
	{
		for (s = 0; s < k; s++)
		{
			mb[r*k+s] = (1.0*rand())/RAND_MAX;
		}
	}
#if DISPLAYT 
	printf("[T] mb initialization time : %.3fms\n", 1e-3*(mytic()-tstart));
#endif
	/*Compute matrix C*/
	tstart = mytic();
#pragma omp parallel for private(r,t,s,tsum) shared(m,n,k,ma,mb,mc)
	for (r = 0; r < m; r++)
	{
		for (t = 0; t < k; t++)
		{
			tsum = 0.0f;
			for (s = 0; s < n; s++)
			{
				tsum += ma[r*n+s]*mb[s*k+t];
			}
			mc[r*m+t] = tsum;
		}
	}
#if DISPLAYT 
	printf("[T] mc computation time : %.3fms\n", 1e-3*(mytic()-tstart));
#endif

#ifdef MPICH_CC
	printf("[R] rank = %d with ", mpirank);
#pragma omp parallel 
	{
		if (0==omp_get_thread_num())
		{
			printf("omp_size = %d\n", omp_get_num_threads());
		}
	}
#endif
#ifdef MPICH_CC
#endif

#ifdef MPICH_CC
	printf("MPI -- parallel \n");
	MPI_Finalize();
#else
	printf("       serial   \n");
#endif
	free(ma);free(mb);free(mc);
	return 0;
}
           

使用如下指令編譯這個代碼

gcc -DMPICH_CC -O0 -fopenmp test.c `pkg-config --cflags --libs mpich` -o mpiomp.exe
           

執行代碼（我用的是mpich）

mpiexec -n 2  ./mpiomp.exe
           

上面啟動了兩個程序來執行代碼。

程式中，每個程序做的是一樣的事情，都是計算兩個矩陣的乘積，計算矩陣乘積的時候，我用openmp将矩陣乘積的計算配置設定給幾個不同的線程來計算，以增加并行度。

如果你想在執行的時候設定openmp的線程數，可以把代碼中第50行注釋掉，然後用下面的指令執行。我是在ubuntu12.04上測試的，用的gcc-4.6.3，如果先不注釋執行的話，再注釋後按下面的指令執行的話，可能會得到跟之前一樣的結果。這個原因我不是很清楚。有知道的，煩請告知一下呀。

mpiexec -n 2 -env OMP_NUM_THREADS 4 ./mpiomp.exe
           

注意幾點：

因為有很多新手，可能以上來執行不出來。是以給出一下幾個可能出問提的點。

* 是否安裝MPICH？ 沒有安裝的話，麻煩先安裝。源碼安裝的看這裡 點選打開連結

* pkg-config出問題。 檢查下 mpich.pc是否在pkg-config的搜尋路徑，或者man pkg-config研究下pkg-config