MPI程序例子 test_8_1_2.c -- 对等模式的MPI程序,Jacobi迭代 (MPI_Send、MPI_Recv)

NOTE:

         这里首先需要弄明白 Jacobi迭代是做什么的，怎么操作。

 网上找到的一篇讲解使用 MPI解决Jacobi迭代并行化的文章，这个与都志辉 并行程序一书的例子有相似之处。链接http://www.doc88.com/p-99830356011.html。下边程序算法主要步骤代码段为，这个步骤为矩阵每个数取相邻四个数的中值。

           数据分块：

                               

for(j=begin_col;j<=end_col;j++)
                 for(i=1;i<totalsize-1;i++)
                         b[i][j]=0.25*(a[i][j+1]+a[i][j-1]+a[i+1][j]+a[i-1][j]);
         

         for(j=begin_col;j<=end_col;j++)
                  for(i=1;i<totalsize-1;i++)
                         a[i][j]=b[i][j];

/*
        对等模式的MPI程序：test_8_1_2.c
        Jacobi迭代：迭代数据按列进行分割，
        并假设一共有4个进程同时并行计算。
*/
#include "mpi.h"
#include <stdio.h>

#define totalsize 16
#define mysize totalsize/4    //分成四块，每块大小
#define steps 10

void main(int argc,char* argv[])
{
	int myid,numprocs,n,i,j,rc;
	float a[totalsize][mysize+2],b[totalsize][mysize+2];  //除分块大小外，还包括左右两边各一列
	float temp[totalsize];  /* 临时数组 */
	int begin_col,end_col,ierr;
	MPI_Status status;
	MPI_Init(&argc,&argv);

	/* 得到当前进程标识和总的进程个数 */
	MPI_Comm_rank(MPI_COMM_WORLD,&myid);
	MPI_Comm_size(MPI_COMM_WORLD,&numprocs);
	fprintf(stderr,"Process %d of %d is alive.\n",myid,numprocs);
   
	/* 数组初始化 */
	//NOTE: 整个数组置值为0, 行数为 totalsize, 列数为 mysize+2
	for(j=0;j<mysize+2;j++)
		for(i=0;i<totalsize;i++)
			a[i][j]=0.0;
	if(myid==0)
		for(i=0;i<totalsize;i++)
			a[i][1]=8.0;

	if(myid==3)
		for(i=0;i<totalsize;i++)
			a[i][mysize]=8.0;
	for(i=1;i<mysize+1;i++){   
		a[0][i]=8.0;
		a[totalsize-1][i]=8.0;
	}    
	/* Jacobi迭代部分 */
	for(n=1;n<=steps;n++)
	{
		//这里从两边分别"获取数据"和"发送数据",每次获取或发送数据时先判断是否会越界
		/* 从右侧的邻居得到数据 */
		if(myid<3)
		{
			MPI_Recv(&temp[0],totalsize,MPI_FLOAT,myid+1,10,MPI_COMM_WORLD,&status);
			for(i=0;i<totalsize;i++)
					a[i][mysize+1]=temp[i];
			//a[i][mysize+1]为块最右一列，表示相邻右侧一块最左边一列
		}
		/* 向左侧的邻居发送数据 */
		if(myid>0)
		{
			for(i=0;i<totalsize;i++)
					temp[i]=a[i][1];
			MPI_Send(&temp[0],totalsize,MPI_FLOAT,myid-1,10,MPI_COMM_WORLD);
		}
		/* 向右侧的邻居发送数据 */
		if(myid<3)
		{
				for(i=0;i<totalsize;i++)
						temp[i]=a[i][mysize];
				MPI_Send(&temp[0],totalsize,MPI_FLOAT,myid+1,10,MPI_COMM_WORLD);
		}
		/* 从左侧的邻居得到数据 */
		if(myid<3)
		{
			MPI_Recv(&temp[0],totalsize,MPI_FLOAT,myid-1,10,MPI_COMM_WORLD,&status);
			for(i=0;i<totalsize;i++)
					a[i][0]=temp[i];
			//a[i][0] 为块最左一列，表示相邻左侧一块最右一列
		}
			 
		begin_col=1;
		end_col=mysize;
		if(myid==0)
				begin_col=2;
		if(myid==3)
				end_col=mysize-1;
		for(j=begin_col;j<=end_col;j++)
				for(i=1;i<totalsize-1;i++)
						b[i][j]=0.25*(a[i][j+1]+a[i][j-1]+a[i+1][j]+a[i-1][j]);
		for(j=begin_col;j<=end_col;j++)
				for(i=1;i<totalsize-1;i++)
						a[i][j]=b[i][j];
		//      MPI_Barrier(MPI_COMM_WORLD);
	}/*     迭代结束 */
	MPI_Barrier(MPI_COMM_WORLD);
	/* 输出结果 */
	fprintf(stderr,"\nProcess %d :\n",myid);
	begin_col=1;
	end_col=mysize;
	for(i=0;i<totalsize;i++)
	{
		for(j=begin_col;j<=end_col;j++)
			fprintf(stderr,"%.2fP%d\t",a[i][j],myid);
		fprintf(stderr,"\n");
	}
	fprintf(stderr,"\n");
	MPI_Finalize();
}