基于MPICH并行实验报告

1、MPI并行程序实验报告王琨114712017MPI并行程序实验报告 王琨 114712017一、实验名称MPICH并行程序实验二、实验环境 Windows XP +MPICH+ CodeBlocks 10.5 实验室里的若干台PC机，分别为SOFT-F3，SOFT-F3，SOFT-F3，SOFT-F3，SOFT-F3，SOFT-F3，SOFT-F3，SOFT-F3，SOFT-F3。下面着重说明如何搭建MPICH并行程序设计环境：1. 系统要求安装 MPICH for Microsoft Windows对系统有如下要求： （1）Windows NT4/2000/XP 的Professional

2、或Server版（不支持 Windows 95/98） （2）所有主机必须能够建立TCP/IP 连接 MPICH支持的编译器有：MS VC+ 6.x，MS VC+.NET， Compaq Visual Fortran 6.x， Intel Fortran，gcc，以及g77。2. 创建账户安装MPICH必须以管理员的身份登录每台主机，在所有主机上建立一个同样的账户（当然也可以每个机器使用不同的用户名和账户，然后建立一个配置文件，使用命令行的方式运行程序），我们组统一创建管理员账户：215，如1-1图：1-1创建账户215点击下一步：1-2再点击创建账户，即可成功创建新账户：215最后，为新账户

3、创建密码：2151-33. 安装MPICH运行下载的安装文件mpich.nt.1.2.5.exe，将MPICH 安装到每台主机上。 打开“任务管理器”中的“进程”选项卡，查看是否有一个 mpd.exe的进程。如图1-4所示，如果有的话说明安装成功。以后每次启动系统，该进程将自动运行。1-4 mpd.exe进程4. 注册与配置安装好 MPICH 之后还必须对每台计算机进行注册和配置才能使用。其中注册必须每台计算机都要进行，配置只要在主控的计算机执行就行了。注册的目的是，将先前在每台计算机上申请的账号与密码注册到 MPICH 中去，这样MPICH 才能在网络环境中访问每台主机。配置方法：运行“mp

4、ichmpdbinMPIRegister.exe”首先会提示输入用户账号，然后会提示输入两边密码，之后会问你是否保持上面的设定。如果选择是，则上面的信息将写入硬盘，否则保存在内存中，再重新启动之后就不存在了。如图1-5所示：1-5注册为了让程序在许多主机上执行，而不需建立配置文件来给出相应的各个主机的信息，主控机必须直到当前可用的主机的信息。这时就需要运行MPICH的配置程序来进行配置了。我们选择july为主控计算机，配置过程如下：从“开始->程序->MPICH->mpd-> MPICH Configuration tool”启动。启动之后的界面如1-6图所示：1-6配

5、置整个界面分为三栏，在第一栏中点击Select（1号按钮），然后在跳出的对话框中选择安装了MPICH的主机名（SOFT-F8，SOFT-F7，SOFT-F33，SOFT-F9，SOFT-F10，SOFT-F3，SOFT-F4，SOFT-F2，SOFT-F6，）。之后在第一栏的编辑框中会显示出所有选择的主机。 检查无误后点击第二栏的Apply（2号按钮），这时下方的进度条会显示对各主机核查的情况，如果没问题整个进度条会变为蓝色。最后点击OK（3号按钮）。整个配置就完成了。5. MPICH与编译环境的整合MPICH提供了C语言和Fortran语言的接口。要编译一个MPI+C或MPI+Fortran

6、的程序必须对编译器进行设置。由于本人习惯使用Code Block软件进行编程，所以这里介绍在Code Block中如何进行编译环境的整合。在Code Block中编译一个MPI+C的程序的步骤如下：(1) 打开Code Block 10.5。 （2）设置settings->complier and debugger setting。 将MPICH2中lib目录下的文件都加到link libraries中。（3）将MPICH2中include目录分别添加到Complier、Linker、Resource complier选项卡中。到此，Code Block的MPI环境整合完毕，下面就可以编

7、译MPI的源码了。可以看出，配置Code Block环境要比VC+6.0方便的多，我们也来看看如何配置VC+6.0的。(1) 打开VC+ 6.0。(2) 新建一个工程（工程名为：compute如下图1-7），通常为Win32 Console Application。 1-7 新建工程(3) 在新的工程的编辑界面下，按Alt+F7打开工程设置对话框。(4) 切换到C/C+选项卡。 （如图1-8）首先选择“Win32 Debug”，再选择 “Code Generation”，再选择“Debug Multithreaded”。这时在“工程选项”的文本框中显示“/MT”表示设置成功。1-8 Visua

8、l C+设置然后选择“Win32Release”重复上述步骤，如下图1-91-9 Visual C+设置(5) 在C/C+选项卡中，选择“所有配置”。选择“预处理器”，在“附加包含路径”的文本框中输入MPICH所附带的头文件的目录。如下图1-10：1-10 Visual C+设置(6) 在“连接”选项卡中，选择“所有配置”。然后再选择“输入” ，在“附加库路径”的文本框中输入MPICH所附带的库文件的目录。如下图1-11：1-11 Visual C+设置(7) 在“连接”选项卡中，选择“所有配置” 。然后再选择“常规”，然后在 “对象/库模块” 的文本框中添加 “ws2_32.lib” 。点击

9、“确定” 。这时在“Common 选项”中会出现“ws2_32.lib” 。如图1-12：1-12 Visual C+设置(8) 在“连接”选项卡中，选择“Win32 Debug” 。然后再选择“常规” ，然后在“对象/库模块”的文本框中添加“mpichd.lib” 。点击“确定” 。这时在“工程选项”中会出现“mpichd.lib” 。如图1-13：1-13 Visual C+设置(9) 在“连接”选项卡中，选择“Win32 Release” 。然后再选择“常规” ，然后在“对象/库模块”的文本框中添加“mpich.lib” 。点击“确定” 。这时在“工程选项”中会出现“mpich.lib”

10、 。如图1-14：1-14 Visual C+设置(10) 关闭工程设置对话框。至此，MPICH与编译环境Visual C+6.0整合结束。同配置Code Block相比，VC+6.0要复杂的多，因此推荐用Code Block软件进行编译程序以节省时间。6.设置环境变量将MPICH的bin所在文件目录C:Program FilesMPICHmpdbin加入到系统变量path中，如下图所示：7. 关闭Windows防火墙 为了实验的顺利进行，实验之前先关闭Windows防火墙。我们都以bingxing账户的身份登陆主机进行实验。三、实验内容本人进行了两个并行实验：分别为PI的并行计算和矩阵相乘的

11、并行计算。实验1：MPICH软件中自带的计算PI的程序。算法分析：设计求值并行算法的关键是构造一个合适的函数f ( x) ,使得它计算起来既简便,误差又小。即使 这里取a = 0 , b = 1 , N 为子区间分割数,于是h = 1 / N。由于：故我们选取的函数是，可以直接进行积分，即比较方便。算式中的在计算各个子区间的小矩形面积时,按左边界计算取，按右边界计算为，按中间值计算为，按梯形计算时则要同时用到左、右两个边界，如下图所示，则最后的计算式子为：并行算法实现时，分配每一个进程或CPU 计算其中的若干个子区间,然后通过归约,得到最终结果 。计算过程中，精度越高，要求N 越大。对于足够大

12、的N ，要求在p 个进程或处理机上尽量均衡分配计算任务，通过归约完成总的计算。下面是实现的代码：/* This is an interactive version of cpi */#include "mpi.h"#include <stdio.h>#include <math.h>double f(double);double f(double a) return (4.0 / (1.0 + a*a);int main(int argc,char *argv) int done = 0, n, myid, numprocs, i; double P

13、I25DT = 3.141592653589793238462643; double mypi, pi, h, sum, x; double startwtime = 0.0, endwtime; int namelen; char processor_nameMPI_MAX_PROCESSOR_NAME; MPI_Init(&argc,&argv); MPI_Comm_size(MPI_COMM_WORLD,&numprocs); MPI_Comm_rank(MPI_COMM_WORLD,&myid); MPI_Get_processor_name(proce

14、ssor_name,&namelen); /* fprintf(stdout,"Process %d of %d is on %sn", myid, numprocs, processor_name); fflush(stdout); */ while (!done) if (myid = 0) fprintf(stdout, "Enter the number of intervals: (0 quits) "); fflush(stdout); if (scanf("%d",&n) != 1) fprintf( s

15、tdout, "No number entered; quittingn" );n = 0; startwtime = MPI_Wtime(); MPI_Bcast(&n, 1, MPI_INT, 0, MPI_COMM_WORLD); if (n = 0) done = 1; else h = 1.0 / (double) n; sum = 0.0; for (i = myid + 1; i <= n; i += numprocs) x = h * (double)i - 0.5); sum += f(x); mypi = h * sum; MPI_Redu

16、ce(&mypi, &pi, 1, MPI_DOUBLE, MPI_SUM, 0, MPI_COMM_WORLD); if (myid = 0) printf("pi is approximately %.16f, Error is %.16fn", pi, fabs(pi - PI25DT);endwtime = MPI_Wtime();printf("wall clock time = %fn", endwtime-startwtime); fflush( stdout ); MPI_Finalize(); return 0;实验步骤

17、如下：1. 在主控机（SOFT-F6）上，将计算PI的源代码复制到已建工程中，进行编译，组建，生成可执行文件cpi.exe。2. 分别在三台主机C盘根目录下创建文件夹test，并将cpi.exe拷贝到该文件夹下。3. 在主控机上运行：从“开始->程序->MPICH->mpd->MPIRun”启动图形方式的MPI环境。首先，选择待运行的程序，点击右上角的“.”按钮，选择相应的程序(C:testcpi.exe)；然后，选择需要用于计算的主机个数(9)，指定参加运算的主机名(SOFT-F3，SOFT-F3，SOFT-F3，SOFT-F3，SOFT-F3，SOFT-F3，SOF

18、T-F3，SOFT-F3，SOFT-F3)；最后，点击Run，运行程序。4. 我们也可以直接在Windows的命令行下使用mpirun np 进程数 程序名，即可启动该mpi程序。程序运行过程中要求输入的intervals是细分数目，数字越大，计算的PI的值越精确，当然计算量也相应的增大。由于采用9台主机并行运算，因此可以输入一个较大的数量级，否则并行运算的效果无法体现，经过多次试验，我们采取比该程序源码所限定的double双精度的数量级低一个数量级的数字来运算，即10的十五次方。如果到达10的十六次方，程序将会溢出出现错误。另外，我们分别采用不同数量的主机和不同线程数来运算进行比较。运行的结

19、果会显示在命令行窗口下，不同主机的输出信息会以不同的颜色标出。运行结果如下：实验一PI的并行计算试验共进行了多次，由于篇幅原因，上面只给出3张不同数量主机运行结果的截图，以下为全部试验的数据结果图表：CPUIntel® Pentium E5800 ( 2M Cache 3.0GHz ) 线程数主机数502510542117.3227.3177.3127.3257.32814.61851.4791.4861.4983.24190.8960.9511.467从上图可以分析得到以下几个结果：1 在主机数为1的实验类别中，可以观察到，当为1个线程单核心时，计算时间为14.6秒而多线程即双核同

20、时运作时，计算时间大致为1个线程单核心时的二分之一，这基本符合了之前的预想，根据Amdahl定律，能达到的最大加速比为2，试验结果与之相吻合。但是我们可以看到，当线程数继续增加时，计算的时间基本保持不变，经过分析，我认为当系统全部资源用来进行本程序计算时，线程的数量增加并不能增加系统加速比。由于本机只有2个CPU核心，在同一时间各自只能有两个激活线程真正的同时计算。从这个角度出发，系统线程数的增加应该不会降低运行时间，反而由于线程的切换消耗会略有降低。2 在主机数为5的实验类别中，计算的时间大致为主机数为1的五分之一，基本符合之前的预想，同理，当线程减到5时，各台主机的线程只有一个，此时为单核

21、心运作，因此计算时间增大到为之前的两倍。3 当使用更多的主机到达9台时，计算时间大大的减小了，让我们充分认识到并行计算的重要性。实验2-矩阵相乘算法分析：两个矩阵相乘问题的表达式为：C=A×B，计算定义如下：其中A和B分别是m×k和k×n矩阵，C是m×n矩阵。就上面矩阵看见，串行算法需要进行m×n×k次加法和乘法。随着维数的增加，其复杂度可想而知，为了降低矩阵相乘的运行时间，采用并行算法是非常必要的。对于矩阵相乘的并行算法，可以有三种：对矩阵按行划分、按列划分和棋盘式分块划分。和按行或列划分相比，棋盘式划分可以开发出更高的并行度。对于

22、一个n×n的方阵，棋盘划分最多可以使用n2个处理器进行并行计算，但使用按行或列分解最多可以使用n个。对矩阵相乘采用棋盘式划分的算法通常称作Cannon算法。1. 矩阵的划分 在实际工程数值计算中，矩阵一般阶都很高，元素数量都以万计，为了实现并行计算，必须得对矩阵进行加以划分，然后指派给不同的处理器进行处理。常见的矩阵划分有行列划分和棋盘划分。行列划分，又叫带状划分（Striped Partitioning），就是将矩阵整行或者整列分成若干个组，每个组指派给一个处理器。下图所例为4个CPU，8×8矩阵的带状划分。在带状划分情况下，每个CPU将会均匀分配到2行(列)数据。8&#

23、215;8矩阵变成了一个1×4或4×1的分块矩阵，每个CPU所属的分块矩阵大小为8×2或2×8。图1 矩阵按列划分图2 矩阵按行划分棋盘划分（Checker Board Partitioning）就是将矩阵分成若干个子矩阵，每个子矩阵指派给一个处理器，此时任一处理器均不包含整行或者整列。下图所示即为4个处理器情况下8×8矩阵的棋盘划分，其中处理器阵列为2×2，每个处理器分配到的子矩阵大小为4×4。2. 矩阵Cannon算法描述对于两矩阵相乘C=A×B，假设A大小为m×k，B大小为k×n，C大小为

24、m×n。使用p个处理器把矩阵A和B分解成sqrt(p)行sqrt(p)列，p为平方数，p个处理器被映射成二维的网孔，分别与A和B的分块矩阵对应。Ci,j的计算公式如下：如下图所示，进程Pi,j存储分块矩阵这一部分。块矩阵乘法要计算所有匹配的Ai,k和Bk,j，然而只有在主对角线的进程才是匹配的。因此需要循环移动分块矩阵的方法来使每个进程Pi,j都有一对可以直接相乘的匹配的块，具体方法如下：（1）将排第i行的Ai,j块循环左移i-1个位置；将第j列Bi,j块循环上移j-1个位置；（2）Pi,j进程执行乘-加运算；然后，将移动得到的Ai,j块循环左移1个位置；将移动得到的Bi,j块循环上

25、移1个位置（3）重复第2步(sqrt(p)-1) 次，每次移动后Pi,j进程执行乘-加运算；下图所示(假设处理器有p=16个)，可以看出，经过第1步循环移动后的矩阵元素变化情况。具体实现的代码如下：#include <stdlib.h>#include <string.h>#include <mpi.h>#include <time.h>#include <stdio.h>#include <math.h>/* 全局变量声明 */double starttime;double time1,time2;float *A, *

26、B, *C; /* 总矩阵,C = A * B */float *a, *b, *c, *tmp_a, *tmp_b; /* a、b、c表分块，tmp_a、tmp_b表缓冲区 */int dg, dl, dl2,p, sp; /* dg:总矩阵维数;dl:矩阵块维数;dl2=dl*dl;p:处理器个数;spsqrt(p) */int my_rank, my_row, my_col; /* my_rank:处理器ID;(my_row,my_col):处理器逻辑阵列坐标 */MPI_Status status;/* *函数名: get_index *功能：处理器逻辑阵列坐标至rank号的转换 *输

27、入：坐标、逻辑阵列维数 *输出：rank号 */int get_index(int row, int col, int sp) return (row+sp)%sp)*sp + (col+sp)%sp;/* *函数名：random_A_B *功能：随机生成矩阵A和B */void random_A_B() int i,j; srand(unsigned int)time(NULL); /*设随机数种子*/*随机生成A,B,并初始化C*/ for(i=0; i<dg ; i+) for(j=0; j<dg ; j+) Aij = 2+rand(); Bij = 2+rand(); C

28、ij = 0.0; /* 函数名：scatter_A_B * 功能：rank为0的处理器向其他处理器发送A、B矩阵的相关块 */void scatter_A_B() int i,j,k,l; int p_imin,p_imax,p_jmin,p_jmax; for(k=0; k<p; k+) /*计算相应处理器所分得的矩阵块在总矩阵中的坐标范围*/ p_jmin = (k % sp ) * dl; p_jmax = (k % sp + 1) * dl-1; p_imin = (k - (k % sp)/sp * dl; p_imax = (k - (k % sp)/sp +1) *dl

29、-1; l = 0; /*rank=0的处理器将A,B中的相应块拷至tmp_a,tmp_b，准备向其他处理器发送*/ for(i=p_imin; i<=p_imax; i+) for(j=p_jmin; j<=p_jmax; j+) tmp_al = Aij; tmp_bl = Bij; l+; /*rank=0的处理器直接将自己对应的矩阵块从tmp_a,tmp_b拷至a,b*/ if(k=0) memcpy(a, tmp_a, dl2 * sizeof(float); memcpy(b, tmp_b, dl2 * sizeof(float); else /*rank=0的处理器向

30、其他处理器发送tmp_a,tmp_b中相关的矩阵块*/ MPI_Send(tmp_a, dl2, MPI_FLOAT, k, 1, MPI_COMM_WORLD); MPI_Send(tmp_b, dl2, MPI_FLOAT, k, 2, MPI_COMM_WORLD); /* *函数名:init_alignment *功能:矩阵A和B初始对准 */void init_alignment() /*将A中坐标为(i,j)的分块A(i,j)向左循环移动i步*/ MPI_Sendrecv(a, dl2, MPI_FLOAT, get_index(my_row,my_col-my_row,sp),

31、1, tmp_a, dl2, MPI_FLOAT, get_index(my_row,my_col+my_row,sp), 1, MPI_COMM_WORLD, &status); memcpy(a, tmp_a, dl2 * sizeof(float) ); /*将B中坐标为(i,j)的分块B(i,j)向上循环移动j步*/ MPI_Sendrecv(b, dl2, MPI_FLOAT, get_index(my_row-my_col,my_col,sp), 1, tmp_b, dl2, MPI_FLOAT, get_index(my_row+my_col,my_col,sp), 1,

32、 MPI_COMM_WORLD, &status); memcpy(b, tmp_b, dl2 * sizeof(float) );/* *函数名：main_shift *功能：分块矩阵左移和上移，并计算分块c */void main_shift() int i,j,k,l; for(l=0; l<sp; l+) /*矩阵块相乘，c+=a*b */ for(i=0; i<dl; i+) for(j=0; j<dl; j+) for(k=0; k<dl; k+) ci*dl+j += ai*dl+k*bk*dl+j; /* 将分块a左移1位 */ MPI_Send

33、(a , dl2, MPI_FLOAT, get_index(my_row, my_col-1, sp), 1, MPI_COMM_WORLD); MPI_Recv(a , dl2, MPI_FLOAT, get_index(my_row, my_col+1, sp), 1, MPI_COMM_WORLD, &status); /* 将分块b上移1位 */ MPI_Send(b , dl2, MPI_FLOAT, get_index(my_row-1, my_col, sp), 1, MPI_COMM_WORLD); MPI_Recv(b , dl2, MPI_FLOAT, get_i

34、ndex(my_row+1, my_col, sp), 1, MPI_COMM_WORLD, &status); /* *函数名：collect_c *功能：rank为0的处理器从其余处理器收集分块矩阵c */void collect_C() int i,j,i2,j2,k; int p_imin,p_imax,p_jmin,p_jmax; /* 分块矩阵在总矩阵中顶点边界值 */ /* 将rank为0的处理器中分块矩阵c结果赋给总矩阵C对应位置 */ for (i=0;i<dl;i+) for(j=0;j<dl;j+) Cij=ci*dl+j; for (k=1;k<

35、;p;k+) /*将rank为0的处理器从其他处理器接收相应的分块c*/ MPI_Recv(c, dl2, MPI_FLOAT, k, 1, MPI_COMM_WORLD, &status); p_jmin = (k % sp ) *dl; p_jmax = (k % sp + 1) *dl-1; p_imin = (k - (k % sp)/sp *dl; p_imax = (k - (k % sp)/sp +1) *dl -1; i2=0; /*将接收到的c拷至C中的相应位置,从而构造出C*/ for(i=p_imin; i<=p_imax; i+) j2=0; for(j=

36、p_jmin; j<=p_jmax; j+) Cij=ci2*dl+j2; j2+; i2+; /*函数名：print *功能：打印矩阵 *输入：指向矩阵指针的指针，字符串 */void print(float *m,char *str) int i,j; printf("%s",str); /*打印矩阵m*/ for(i=0;i<dg;i+) for(j=0;j<dg;j+) printf("%15.0f ",mij); printf("n"); printf("n");/* *函数名：main

37、 *功能：主过程，Cannon算法，矩阵相乘 *输入：argc为命令行参数个数，argv为每个命令行参数组成的字符串数组 */int main(int argc, char *argv) int i; MPI_Init(&argc, &argv); /* 启动MPI计算 */ MPI_Comm_size(MPI_COMM_WORLD, &p); /* 确定处理器个数 */ MPI_Comm_rank(MPI_COMM_WORLD, &my_rank); /* 确定各自的处理器标识符 */ sp = sqrt(p); /* 确保处理器个数是完全平方数，否则打印错误

38、信息，程序退出 */ if (sp*sp != p) if (my_rank = 0) printf("Number of processors is not a quadratic number!n"); MPI_Finalize(); exit(1); if (argc != 2) if (my_rank = 0) printf("usage: mpirun -np ProcNum cannon MatrixDimensionn"); MPI_Finalize(); exit(1); starttime = MPI_Wtime(); dg = ato

39、i(argv1); /* 总矩阵维数 */ dl = dg / sp; /* 计算分块矩阵维数 */ dl2 = dl * dl; /* 计算处理器在逻辑阵列中的坐标 */ my_col = my_rank % sp ; my_row = (my_rank-my_col) / sp ; /* 为a、b、c分配空间 */ a = (float *)malloc( dl2 * sizeof(float) ); b = (float *)malloc( dl2 * sizeof(float) ); c = (float *)malloc( dl2 * sizeof(float) ); /* 初始化c

40、 */ for(i=0; i<dl2 ; i+) ci = 0.0; /* 为tmp_a、tmp_b分配空间 */ tmp_a = (float *)malloc( dl2 * sizeof(float) ); tmp_b = (float *)malloc( dl2 * sizeof(float) ); if (my_rank = 0) /* rank为0的处理器为A、B、C分配空间 */ A = (float *)malloc( dg * sizeof(float*) ); B = (float *)malloc( dg * sizeof(float*) ); C = (float

41、*)malloc( dg * sizeof(float*) ); for(i=0; i<dg; i+) Ai = (float *)malloc( dg * sizeof(float) ); Bi = (float *)malloc( dg * sizeof(float) ); Ci = (float *)malloc( dg * sizeof(float) ); random_A_B(); /* rank为0的处理器随机化生成A、B矩阵 */ scatter_A_B(); /* rank为0的处理器向其他处理器发送A、B矩阵的相关块 */ else /* rank不为0的处理器接收来自

42、rank为0的处理器的相应矩阵分块 */ MPI_Recv(a, dl2, MPI_FLOAT, 0 , 1, MPI_COMM_WORLD, &status); MPI_Recv(b, dl2, MPI_FLOAT, 0 , 2, MPI_COMM_WORLD, &status); time1=MPI_Wtime(); init_alignment(); /* A、B矩阵的初始对准 */ main_shift(); /* 分块矩阵左移、上移, cannon算法的主过程 */ if(my_rank = 0) collect_C(); /* rank为0的处理器从其余处理器收集分

43、块矩阵c */ print(A,"random matrix A : n"); /* 打印矩阵A */ print(B,"random matrix B : n"); /* 打印矩阵B */ print(C,"Matrix C = A * B : n"); /* 打印矩阵C */ else MPI_Send(c,dl2,MPI_FLOAT,0,1,MPI_COMM_WORLD); /* rank不为0的处理器向rank为0的处理器发送矩阵块c */ time2=MPI_Wtime(); if (my_rank=0) printf("n"); printf("Whole running time = %f secondn",time2-starttime);