95MIMD异步通信模型的并行算法

19.5MIMD异步通信模型的并行算法（一）快速排序并行算法快速排序(QuickSort)是一种最基本的排序算法，它的基本思想是，在当前无序序列R[1,n]中取一个记录作为比较的基准，用此基准将当前的无序序列R[1,n]划分成两个无序序列R[1,i-1]和R[i,n](1≤i≤n)，且R[1,i-1]中记录的所有关键字均小于等于基准的关键字，R[i,n]中记录的所有关键字均大于等于基准的关键字；当R[1,i-1]和R[i,n]非空时，分别对他们重复上述的划分过程。对每次划分过后所得到的两个序列分别使用两个处理器完成递归排序。例如对一个长为n的序列，首先划分得到两个长为n/2的序列，将其交给两个处理器分别处理；而后进一步划分得到4个长为n/4的序列，在分别交给4个处理器处理；如此递归下去最终得到排序号的序列[1001]。该并行算法的描述如下：1234567891011121314151617算法9.5.1快速排序并行算法输入：无序数组data[1,n]，使用的处理器个数2m输出：有序数组data[1,n]SORT{para_quicksort(data,1,n,m,0);}para_quicksort(data,i,j,m,id){if((j-i)≤k||m==0){Pid:quicksort(data,i,j);}else{Pid:r=partition(data,i,j);Pidsenddata[r+q,m-1]toPid+2m-1;para_quicksort(data,I,r-1,m-1,id);para_quicksort(data,r+1,j,m-1,id+2m-1);Pid+2m-1senddata[r+1,m-1]toPid;}}在最优的情况下该并行算法形成一个高度为logn的排序树，其计算复杂度为O(n)，通信复杂度为O(n)；在最坏情况下其计算复杂度降为O(n2)；正常情况下该算法的计算复杂度为O(n)。（二）二维网孔上的矩阵转置并行算法网孔上的矩阵转置并行算法思路是，实现矩阵转置时，若处理器个数为p，且它们的编号依次是0,1,…,p-1，则将n阶矩阵A分成p个大小为m×m的子块，m=n/p。p个字块组成一个p×p的子块阵列。记其中第i行第j列的子块为Aij，它含有A的第(i-1)m+1至第im行中的第(j-1)m+1至第jm列的所有元素。对每一处理器按行主方式赋以二维下标，记编号为i的处理器的而为下标为(v,u)，2其中v=pi/，u=imodp，将A的子块存入下标为(v,u)表示的对应处理器中。这样，转置过程分两步进行：第一步，子块转置，具体过程如图9_22所示；第二步，处理器内部局部转置。图9_22子块转置为了避免对应子块交换数据时处理器发生死锁，可令下三角子块先向与之对应的上三角子块发送数据，然后从上三角子块接收数据；上三角子块先将数据存放在缓冲区buffer中，然后从与之对应的下三角子块接收数据；最后再将缓冲区中的数据发送给下三角子块[1001]。该并行算法的描述如下：1234567891011121314151617181920212223算法9.5.2矩阵转置并行算法输入：矩阵An×n输出：矩阵An×n的转置矩阵ATn×nTRANSPOSEDMATRIX{/*对所有处理器my_rank(my_rank=0,…,p-1)*/vmy_rank/sqrt(p);/*计算子块的行号*/umy_rankmodsqrt(p);/*计算子块的列号*/if(uv)/*对存放下三角块的处理器*/{send所存的子块to其对角块所在的处理器;receive其对角块所在的处理器中发来的子块;}else/*对存放上三角块的处理器*/{将所存的子块在缓冲区buffer中做备份;receive其对角块所在的处理器中发来的子块;sendbuffer中所存子块to其对角块所在的处理器;}for(i=1;i=m;i++)/*处理器内部局部转置*/{for(j=1;j=i;j++){交换a[i,j]和a[j,i];}}}3若记ts为发送启动时间，tw为单位数据传输时间，th为处理器间的延迟时间，则第一步由于每个子块有n2/p个元素，又由于通信过程中为了避免死锁，错开下三角子块与上三角子块的发送顺序，因此子块的交换时间为)/(22ptpntthws；第二步，假定一对数据的交换时间为一个单位时间，则局部转置时间为n2/2p。因此所需的并行计算时间ptpntptpnThwsp22222。（三）矩阵并行分块乘法算法矩阵并行分块乘法的基本思想：将矩阵A按行划分为p块(p为处理器个数)，设u=m/p，每块含有连续的u行向量，这些行块依次记为A0,A1,…,Ap-1，分别存放在标号为0,1,…,p-1的处理器中。对矩阵B按列划分为p块，记v=k/p，每块含有连续的v列向量，这些列块依次记为B0,B1,…,Bp-1，分别存放在标号为0,1,…,p-1的处理器中。将结果矩阵C也相应的同时进行行、列划分，得到p×p个大小为u×v的子矩阵，记第i行第j列的子矩阵为Cij，则有Cij=Ai×Bj，其中Ai大小为u×n，Bj大小为n×v。开始，各处理器的存储内容如图9_23(a)所示。此时各处理器并行计算Cij=Ai×Bi，其中i,j=0,1,…,p-1，此后第i号处理器将所存储的B的列块送至第i-1号处理器(第0号处理器将B的列块送至第p-1号处理器中，形成循环传送)，各处理器中存储内容如图9_23(b)所示。它们再次并行计算Cij=Ai×Bj，这里j=(i+1)modp。B的列块在处理器中以这样的方式循环传送p-1次，并做p次子矩阵相乘运算，就生成了矩阵C的所有子矩阵。编号为i的处理器的内部存储器存有子矩阵Ci0,Ci1,…,Ci(p-1)。为了避免在通信过程中发生死锁，奇数号及偶数号处理器的收发顺序被错开，使偶数号处理器先发送后接收；而奇数号处理器先将B的列块存于缓冲区buffer中，然后接受编号在其后面的处理器所发送的B的列块，最后再将缓冲区的原矩阵B的列块发送给编号在其前面的处理器[1001]。该并行算法的描述如下：处理器编号存储器内容0A0B01A1B12A2B2………p-1Ap-1Bp-1处理器编号存储器内容0A0B11A1B22A2B3………p-1Ap-1B0(a)(b)图9_23矩阵相乘并行算法中的数据交换12算法9.5.3矩阵并行分块乘法算法输入：矩阵Am×n,Bn×k输出：矩阵Cm×kBLOCK-MATRIXMUTIPLICATION{/*对所有处理器my_rank(my_rank=0,…,p-1)*/4345678910111213141516171819202122232425262728l(i+my_rank)%p;/*计算C的子块号*/for(z=0;z=u-1;z++){c[l,z,j]0；for(j=0;j=v-1;j++){c[l,z,j]c[l,z,j]+a[z,s]*b[s,j];}}mm1(p+my_rank-1)%p;/*计算左邻处理器的标号*/mm1(my_rank+1)%p;/*计算右邻处理器的标号*/if(i≠p-1){if(my_rank%2==0)/*编号为偶数的处理器*/{将所存的B的子块发送到其左邻处理器中接收其右邻处理器中发来的B的子块}if(my_rank%2≠0)/*编号为奇数的处理器*/{将所存的B子块在缓冲区buffer中做备份;receive其右邻处理器中发来的B的子块;sendbuffer中所存的B的子块to其左邻处理器;}}}矩阵并行分块乘法中，设一次乘法和加法运算时间为一个单位时间，由于每个处理器计算p个u×n与n×v阶的子矩阵相乘，因此计算时间为u×n×v×p；所有处理器交换数据p-1次，每次的通信量为v×n，通信过程中为了避免死锁，错开奇数号及偶数号处理器的收发顺序，通信时间为2(p-1)(ts+nvtw)=O(nk)，所以并行计算时间Tp=uvnp+2(p-1)(ts+nvtw)=mnk/p+2(p-1)(ts+nvtw)。（四）分布式矩阵求逆的并行算法矩阵求逆的并行算法的基本思想：在矩阵求逆的过程中，一次利用主行i(i=0,1,…,n-1)对其余各行j(j≠i)作初等变换，由于各行计算之间没有数据相关关系，因此对矩阵A按行划分来实现并行计算。考虑到在计算过程中处理器之间的负载均衡，对A采用行交叉划分：设处理器个数为p，矩阵A的阶数为n,m=n/p，对矩阵A行交叉划分后，编号为i(i=0,1,…,p-1)的处理器存有A的第i,i+p,…,i+(m-1)p行。在计算中，依次将第0,1,...,n-1行作为主行，将其广播给所有处理器，这实际上是各处理器轮流选出主行并广播。发送主行数据的处理器利用主行对其主行之外的m-1行行向量做行变换，其余处理器则利用主行对其m行行向量做行变换[1001]。该并行算法的描述如下：算法9.5.4矩阵求逆并行算法51234567891011121314151617181920212223242526272829303132333435363738394041输入：矩阵An×n输出：矩阵A-1n×nINVERSEMATRIX{/*对所有处理器my_rank(my_rank=0,…,p-1)*/for(i=0;i=m-1;i++){for(j=0;j=p-1;j++){if(my_rank=j)/*主元素在本处理器中*/{vi*p+j;/*v为主元素的行号*/a[i,v]1/a[i,v];for(k=0;k=n-1;k++){if(k≠v){a[i,k]a[i,k]*a[i,v];}}for(k=0;k=n-1;k++){f[k]a[i,k];}broadcast变换后的主行元素(存于数组f中)to所有处理器;else/*主元素不在本处理器中*/{vi*p+j;/*v为主元素行号*/reveive主行元素存于数组中;}if(my_rank≠j)/*主行元素不在本处理器中*/{for(k=0;k=m-1;k++)/*处理非主行、非主列元素*/{for(w=0;w=n-1;w++){if(w≠v){a[k,w]a[k,w]-f[w]*a[k,v];}}}for(k=0;k=m-1;k++)/*处理主行元素*/{64243444546474849505152535455565758596061626364656667686970a[k,v]-f[v]*a[k,v];}}else/*处理主行所在处理器中的其它元素*/{for(k=0;k=m-1;k++){if(k≠i){for(w=0;w=n-1;w++){if(w≠v){a[k,w]a[k,w]-f[w]*a[k,v];}}}}for(k=0;k=m-1;k++){if(k≠i){a[k,v]-f[v]*a[k,v];}}}}}}矩阵求逆的并行算法中，若一次乘法和加法运算或一次除法运算时间为一个单位时间，则所需的计算时间为mn2；又由于共有n行数据依次作为主行被广播，其通信时间为n(ts+ntw)logp，所以该算法并行计算时间为Tp=mn2+n(ts+ntw)logp。（五）分布式高斯消去并行算法高斯消去法是利用主行i对其余各行j(ji)，作初等变换，各行计算之间没有数据相关关系，因此可以对矩阵A按行划分。考虑到在计算过程中处理器之间的负载均衡，对A采用行交叉划分。设处理器个数为p，矩阵A的阶数为n,m=n/p，对矩阵A行交叉划分后，编号为i(i=0,1,…,p-1)的处理器含有A的第i,i+p,…,i+(m-1)p行和向量B的第i,i+p,…,i+(m-1)p一共m个元素。消去过程的并行是依次以第0,1,…n-1行为主行进行消去计算，由于对行的交叉划分与分布，这实际上是由各处理器轮流选出主行。在每次消去计算前，各处理器并行求其局部存储器中右下角子阵的最大元。若以编号为my_rank的处理器