Cannon乘法的MPI实现解析

《并行算法实践》要求学生在KD60实验平台上设计并行算法并实现。实验平台由一块处理板、一块监控板和一块背板等组成。处理板逻辑结构如图1所示。处理板承载4个处理单元，每个处理单元包括一个龙芯3号四核CPU、2GBDDR2内存、RTL8110千兆以太网卡芯片、BIOSFlash、串口收发芯片以及电源变换电路等。四个龙芯3号处理器通过HT总线实现互连。监控电路检测4个处理单元的状态，并实现对其控制。图1处理板逻辑结构实验平台的系统软件以开源软件为主（见图2），具有兼容性强、易维护、易升级、易使用等特点。处理单元操作系统为DebianGNU/Linux无盘系统，采用稳定高效的2.6.27内核。图2软件系统结构要求选修《并行算法实践》同学在下面表1中选一个题目，（1）阐述基本原理（包括对算法改进和优化方法）；(2)根据KD60实验平台设计实验方法和步骤（包括主要调试过程要求拷屏）。(3)数据及结果分析：根据不同的实验内容，记录具体的实验数据或程序运行结果（要求拷屏）。实验数据量较大时，最好制成表格形式。附上程序功能、模块说明、完整源代码，源代码中尽量多带注释；(4)分析和总结：对程序与算法性能改进结论，总结和体会。表1《并行算法实践》题目序号题目名称基本方法和内容要求1LU分解的OpenMP实现编写LU分解的OpenMP程序2KMP算法的OpenMP实现编写KMP算法的OpenMP程序3高斯消元法解线性方程组的OpenMP实现编写高斯消元法解线性方程组的OpenMP程序4高斯消元法解线性方程组的MPI实现编写高斯消元法解线性方程组的MPI程序5高斯-塞德尔迭代解线性方程组的MPI实现编写高斯-塞德尔迭代解线性方程组的MPI程序6Cannon乘法的MPI实现编写Cannon乘法的MPI程序7LU分解的MPI实现编写LU分解的MPI程序8随机串匹配算法的MPI实现编写随机串匹配算法的MPI程序9单源最短路径Dijkstra算法的MPI实现编写单源最短路径Dijkstra算法的MPI程序10快速排序算法的MPI实现编写快速排序算法的MPI程序11KMP串匹配的MPI实现编写KMP串匹配算法的MPI程序图2软件系统结构Cannon乘法的MPI实现及性能分析摘要：cannon算法是矩阵的并行乘法，属于数值并行算法MPI编程实现一篇，其中关于数值并行算法MPI编程由于要处理的数据量巨大，程序循环次数多，对于串行而言，处理时间将非常长，将其并行化非常必要。本文将矩阵数据进行棋盘划分成多个子矩阵，再分别指派给多个处理器，使个处理器并行运算。关键字：cannon乘法并行计算数据划分一、Cannon乘法的MPI实现基本原理Cannon乘法属于数值并行算法MPI编程实现一篇，其中关于数值并行算法MPI编程由于要处理的数据量巨大，程序循环次数多，对于串行而言，处理时间将非常长，使其并行化的一般方法有：1）数据相关分析2）数据划分和处理器指派3）循环重构对原有程序并行化，首先要分析计算程序中所有语句间的依赖关系，这称之为相关分析。本项目Cannon乘法的mpi实现，是矩阵运算，阶往往都很高，而且行列之间数据依赖关系也不强，所以就对矩阵进行划分，然后指派给不同的处理器进行处理。最常用的矩阵划分有带状划分和块状划分。1．带状划分方法带状划分又叫行列划分，就是将矩阵整行或整列地分成若干组，各组指派给一个处理器。也可以将若干行或列指派给一个处理器，而且这些行和列可以是连续的，也可以是等间距的，前者称为块带状的，后者称为循环带状的。2.块状划分方法块状划分又叫棋盘划分，就是将矩阵划分成若干个子矩阵，每个子矩阵指派给一个处理器，此时任意处理器均不含整行或整列。和带状划分类似，棋盘划分也可分为块棋盘划分和循环棋盘划分。棋盘划分比带状划分可开发更高的并行度，Cannon乘法的mpi实现也正是基于棋盘划分的并行实现。循环重构是指在数据分解之后，相应地将串行程序循环部分进行重构，以实现这种划分所确定的并行计算，主要方法有1）循环交换2）拉伸法3）分裂法4）轮转法5）并列法在三种程序并行化的方法中，数据相关分析和循环重构目的都是挖掘语句间的并行性，而数据划分和处理器指派则重在策略，宏观上挖掘并行性。Cannon算法是一种存储有效的算法，设矩阵Ann和Bnn相乘。为了使两矩阵下标满足相乘的要求，和带状的并行分块乘法不同，不是仅仅让B矩阵的各列块循环移动，而是有目的地让A的各行块以及B的各列块皆施行循环移位，从而实现对C的子块的计算。将矩阵A和B分成p个方块Aij和Bij，)1,0(pji，每块大小为pn/pn/，并将它们分配给pp个处理器),...,,(110100ppPPP。开始时处理器Pij存放块Aij和Bij，并负责计算块Cij，然后算法开始执行：⑴将块Aij向左循环移动i步；将块Bij向上循环移动j步；⑵Pij执行乘加运算后将块Aij向左循环移动1步，块Bij向上循环移动1步；⑶重复第⑵步，总共执行p次乘加运算和p次块Aij和Bij的循环单步移位。二、Cannon乘法的MPI实现内容和步骤实验涉及内容主要有：1）数据划分和指派处理器最常用的矩阵数据划分有带状划分和块状划分。棋盘划分比带状划分可开发更高的并行度，Cannon乘法的mpi实现也正是基于棋盘划分的并行实现。设有P个处理器，将矩阵A和B分成p个方块Aij和Bij，)1,0(pji，每块大小为pn/pn/，并将它们分配给pp个处理器),...,,(110100ppPPP。2）子矩阵的循环移动处理器Pij存放块Aij和Bij，并负责计算块Cij，在使A矩阵的左右循环移动和B矩阵的上下循环移动时，为了避免在通信过程中发生死锁，奇数号及偶数号处理器的收发顺序被错开，使偶数号处理器先发送后接收；而奇数号处理器先将子矩阵块存于缓冲区Buffer中，然后接收编号在其后面的处理器所发送的子矩阵块，最后再将缓冲区中子矩阵块发送给编号在其前面的处理器。基本算法如下：Begin(1)if(j=0)then/*最左端的子块*/(1.1)将所存的A的子块发送到同行最右端子块所在的处理器中(1.2)接收其右邻处理器中发来的A的子块endif(2)if((j=sqrt(p)-1)and(jmod2=0))then/*最右端子块处理器且块列号为偶数*/(2.1)将所存的A的子块发送到其左邻处理器中(2.2)接收其同行最左端子块所在的处理器发来的A的子块endif(3)if((j=sqrt(p)-1)and(jmod2≠0))then/*最右端子块处理器且块列号为奇数*/(3.1)将所存的A的子块在缓冲区buffer中做备份(3.2)接收其同行最左端子块所在的处理器发来的A的子块(3.3)将在缓冲区buffer中所存的A的子块发送到其左邻处理器中endif(4)if((j≠sqrt(p)-1)and(jmod2=0)and(j≠0))then/*其余的偶数号处理器*/(4.1)将所存的A的子块发送到其左邻处理器中(4.2)接收其右邻处理器中发来的A的子块endif(5)if((j≠sqrt(p)-1)and(jmod2=1)and(j≠0))then/*其余的奇数号处理器*/(5.1)将所存的A的子块在缓冲区buffer中做备份(5.2)接收其右邻处理器中发来的A的子块(5.3)将在缓冲区buffer中所存的A的子块发送到其左邻处理器中endifEnd实验步骤1）登陆KD-60图2.1KD-60登陆界面2）转至node80节点，上传程序输入命令：sshloongson@node80和密码进入图界面图2.2转到节点80的界面再命令vim，进入vim编辑器加入程序,保存为cannon.c3）编译程序输入命令：mpicccannon.c–ocannon–lm在目录中查看，已成功。如下图图2.3将程序保存并编译后界面4）运行程序输入：mpirun–np4cannon4,其中第一个4是指定的处理器个数，第二个4是产生随机矩阵的维数，这两个参数在实验过程中可以调整，但要求第一个参数即处理器的个数必须是一个数的平方数。输出：图2.4cannon乘法运行结果图2.4并行程序运行界面两个参数都是4，分别输出两个随机矩阵和矩阵的乘积三、数据及结果1.下面列出了两组数据，分别是用一个处理器进行串行运算和四个处理器进行并行运算矩阵维数为200的计算时间比较。四个处理器处理阶数为200的矩阵相乘时，所花时间为：1.705844秒。单个处理器处理阶数为200的矩阵相乘时，所花时间为：3.727210秒。如图3.1和图3.2所示。图3.1四个处理器并行执行结果图图3.2单个处理器串行执行结果图附：1.程序模块伪代码：输入：An×n，Bn×n输出：Cn×nBegin对所有处理器my_rank(my_rank=0,…,p-1)同时执行如下的算法:(1)计算子块的行号i=my_rank/sqrt(p)计算子块的列号j=my_rankmodsqrt(p)(2)fork=0top-1doif(ik)thenLeftmoveonestep(a)endif/*a循环左移至同行相邻处理器中*/if(jk)thenUpmoveonestep(b)endif/*b循环上移至同列相邻处理器中*/endfor(3)fori=0tom-1doforj=0tom-1doc[i,j]=0endforendfor(4)fork=0top-1dofori=0tom-1doforj=0tom-1dofork1=0tom-1doc[i,j]=c[i,j]+a[i,k1]*b[k1,j]endforendforendforLeftmoveonestep(a)/*子块a循环左移至同行相邻的处理器中*/Upmoveonestep(b)/*子块b循环上移至同列相邻的处理器中*/endforEndLeftmoveonestep(a)见实验内容处附2.程序源码#includestdlib.h#includestring.h#includempi.h#includetime.h#includestdio.h#includemath.h/*全局变量声明*/float**A,**B,**C;/*总矩阵,C=A*B*/float*a,*b,*c,*tmp_a,*tmp_b;/*a、b、c表分块，tmp_a、tmp_b表缓冲区*/intdg,dl,dl2,p,sp;/*dg:总矩阵维数;dl:矩阵块维数;dl2=dl*dl;p:处理器个数;sp＝sqrt(p)*/intmy_rank,my_row,my_col;/*my_rank:处理器ID;(my_row,my_col):处理器逻辑阵列坐标*/MPI_Statusstatus;floatstarttime;floattime1;/**函数名:get_index*功能：处理器逻辑阵列坐标至rank号的转换*输入：坐标、逻辑阵列维数*输出：rank号*/intget_index(introw,intcol,intsp){return((row+sp)%sp)*sp+(col+sp)%sp;}/**函数名：random_A_B*功能：随机生成矩阵A和B*/voidrandom_A_B(){inti,j;srand((unsignedint)time(NULL));/*设随机数种子*//*随机生成A,B,并初始化C*/for(i=0;idg;i++)for(j=0;jdg;j++){A[i][j]=rand();B[i][j]=rand();C[i][j]=0.0;}}/*函数名：scatter_A_B*功能：rank为0的处理器向其他处理器发送A、B矩阵的相关块*/voidscatter_A_B(){inti,j,k,l;intp_imin,p_imax,p_jmin,p_jmax;for(k=0;kp;k++){/*计算相应处理器所分得的矩阵块在总矩阵中的坐标范围*/p_jmin=(k%sp)