摘要：在阵列信号处理中需要大量的矩阵运算，而其中最基本的就是矩阵相乘运算。本文就矩阵相乘的行划分并行实现进行了改进，将A矩阵的一行和整个B矩阵传输到每个工作进程，其中第一个工作进程指定在某台机器，其余工作进程由PVM选择在合适的机器上产生。该并行实现基于PVM环境，采用主从编程模式，然后给出了改进的行划分矩阵相乘在DSP上实现的方案。
一、引言
并行算法可以根据算法的特点进行算法分解,首先需要分析数据的依存关系和依赖关系,寻找任务的并行性,将多个操作合并成一个任务,然后将整个程序分解为单个任务,同时还可以分析目标的结构通信连接和时序关系,进行并行处理。而工作站机群能充分利用现有的工作站资源和局域网，在工作站机群上，对传统的并行数值计算算法进行优化和改进，可以大幅度提高计算的效率，缩短计算的时间。
数字信号处理对信号在时域及变换域的特性进行分析、处理，能使我们对信号的特性和本质有更清楚的认识和理解，在阵列信号处理中需要大量的矩阵运算，而其中最基本的就是矩阵乘法运算。本文就矩阵相乘在PVM环境下进行了并行化改进实现，并在数字信号处理器中得以实现。
二、问题的提出
在文献[1]和[2]中介绍了一种矩阵相乘的算法。为了提高矩阵相乘的效率，我们将其进行并行化改进，其中矩阵A采用按行划分，矩阵B整个的传给每个工作进程。而PVM程序采用Master/Slave编程模式。首先，主进程产生n个工作进程，其中第一个工作进程指定在某台机器，其余工作进程由PVM选择合适的机器产生。每个工作进程完成矩阵乘积后，将计算结果返回Master进程。
假设矩阵A、B、C均是n阶方阵，其中A被划分为n个n阶的子矩阵，A、B是将进行乘法运算的初始矩阵，C存放运算结果，C在运算前为零矩阵。现在假设有n个处理机Pi,n(1
而在矩阵相乘的直接实现中，将A矩阵的一行和B矩阵的一列传送到一台处理机，这样就需要n2台处理机，实现框图和过程如下所示：
矩阵的串行实现代码如下[2,3]：
for(i←1) to n do
for(j←1) to n do
t←0
for(k←1) to j do
t←t+ai,k*bk,j
repeat
ci,j←t
在单个DSP中实现需要三个循环,使得算法复杂度为N (O 3) 。
矩阵相乘的直接实现图
三、矩阵相乘任务分配及并行PVM程序
3.1 子任务分配策略
现在考虑矩阵的分配问题。假设有n个处理机，且在工作站机群上，采用PVM并行编程环境实现矩阵相乘的并行算法，必须将n个处理机对应于n个并行子任务，每个子任务可以独立地进行局部矩阵运算并相互传送数据。这n个子任务将被分配到工作站机群的有限个工作站上进行，此时，一个工作站上可能同时运行多个子任务。一般情况下，工作站机群中所包含的工作站个数不会恰好等于n。但是，若工作站个数大于n,则并行子任务个数小于工作站个数，机群的全部资源不能得到充分利用，这时应该增加对矩阵的划分数，既增加n的值，使工作站个数小于或等于n。此时，就存在如何将n个子任务分配到小于n个工作站上去的问题，具体分配的策略必须满足下面两个规则：
(1) 各工作站的计算任务负载尽量均衡，减少工作站的空闲等待时间。
(2) 因为机群中各工作站的网络带宽有限，应尽量减少工作站间的数据传输量，以减少子任务间的数据传输时间开销。
3.2 PVM并行程序描述
PVM即并行虚拟机，它是一个自包含的自由软件，源代码公开，可在同构/异构型网络环境下模拟实现一个通用的基于分布存储的并行计算机系统，提供基于消息传递的并行程序设计接口，使网上的用户能够集中地使用众多的机器资源来求解大规模计算问题。而且支持多用户、多任务，多个用户可将系统配置成相互重叠的虚拟机，每个用户可执行多个应用程序和子任务。
如前所述，该矩阵相乘的并行算法矩阵A采用了按行划分，矩阵B整个的传给每个工作进程。而PVM程序采用Master/Slave编程模式。首先，主进程产生n个工作进程，其中第一个工作进程指定在某台机器，其余工作进程由PVM选择合适得机器产生。每个工作进程完成矩阵乘积后，将计算结果返回Master进程。
Master程序主要代码描述：
// 在特定的机器上生成子任务，并将子任务号保留在wtids数组里，但是在此并不返回出错的信息
rCODe=pvm_spawn("mm1.w",NULL,PvmTaskDefault,"",1,&wtids[i]);
// 初始化矩阵
a[i][j]=i+j;
b[i][j]=i*j;
// 数据打包并向其余工作站发送数据
rcode=pvm_iNItsend(PvmDataDefault);
rcode=pvm_pkint(&offset,1,1);
rcode=pvm_pkint(&rows,1,1);
rcode=pvm_pkdouble(&a[offset][0],rows*n,1);
rcode=pvm_pkdouble(&b[0][0],n*n,1);
rcode=pvm_send(wtids[i],mtype);
// 从各个工作站接收数据
rcode=pvm_recv(-1,mtype);
rcode=pvm_upkint(&offset,1,1);
rcode=pvm_upkint(&rows,1,1);
rcode=pvm_upkdouble(&c[offset][0],rows*n,1)
Slave主要代码描述：
// 从主机接收数据
rcode=pvm_recv(mtid,mtype);
rcode=pvm_upkint(&offset,1,1);
rcode=pvm_upkint(&rows,1,1);
rcode=pvm_upkdouble(&a[0][0],rows*n,1);
rcode=pvm_upkdouble(&b[0][0],n*n,1);
// 进行矩阵相乘
c[i][k]=c[i][k]+a[i][j]*b[j][k];
// 数据打包并向主节点发送结果
rcode=pvm_initsend(PvmDataDefault);
rcode=pvm_pkint(&offset,1,1);
rcode=pvm_pkint(&rows,1,1);
rcode=pvm_pkdouble(&c[0][0],rows*n,1);
rcode=pvm_send(mtid,mtype);
四、性能分析
以上所采用的并行算法，它的性能分析包括通讯时间tcomm 和计算时间tcomp 两部分的计算。第一是主从程序间的通讯时间( tcomm) 。将数据分发到wtids[i]台从处理器上,则每台处理器收到矩阵A的一行和整个B矩阵的元素(即n3个元素) 。依赖于具体的广播算法,首先根据矩阵A的行数来分配进程数, 再将矩阵B整体利用PVM中的多播方式传递给每个Slave 进程,然后每个Slave进程得到矩阵A的一行数据,和矩阵B的所有数据，产生矩阵C的一行数据,计算结束,将结果返回Master程序显示。修改后的算法通讯时间为:
tcomm = ( tstartup + ( k ×wtids[i] tdata) + n*wtids[i]*( tstartup + ktdata)
+ n*( tstartup + n*tdata)
其中，tstartup为启动时间, tdata表示发送一个数据字所需的传送时间。然后是每个处理器的计算时间tcomp。计算仅发生在从处理器上,每台从处理器并行执行k*k个乘法操作和k*k个加法操作,即tcomp = 2 k*k 。总的处理时间为tcomm + tcomp 。
五、DSP阵列实现
多DSP的并行执行,是一种结构上的并行。对于阵列信号的处理常用脉动矩阵结构。在脉动矩阵中,处理器计算出的结果送往其他处理器,只有最终处理器的运算结果送回到存储器。它具有一维,二维阵列结构,阵列配置的多个处理器有简单和相同的功能，处理器间数据流和
控制流是规则局域的，通过流水线处理和并行处理相结合能成正比于处理器数目提高处理速
度。首先要对变量进行流水化,分解后得到:
for ( i = 1 ; i < = n ; i + + )
a[ i ] = a[ i - 1 ] ;
for ( j = 1 ; j < = n ; j + + )
b[ j ] = b[ j - 1 ] ;
for ( j = 1 ; j < = N ; j + + )
c[ j ] = c[ j - 1 ] + a[ j - 1 ] * b[ j - 1 ] ;
为了使用脉动矩阵进行计算,需要把运算进行变换,即对运算进行一种影射。其中定义a’ = a , b’ = b, c’ = c + a*b。矩阵数据从阵列的上边和左边进入，计算结果从右边进入。这样，处理器的个数就减少到n。在DSP中实现框图如下所示：
阵列处理框图如下图所示：
在整个系统中总线分为XY两个总线, Y总线用于加载程序,信号阵列输入采样数据,即传输矩阵A的数据。X总线用于连接标准总线,负责传输矩阵B的数据,同时也用于扩展外部存储器。总线控制器包括了XY两个总线控制器,负责控制X总线和Y总线。还包括总线切换控制器,控制XY总线的数据交换。程序加载控制器控制程序存储器,给每一行DSP进行程序加载。在程序加载时,每行由第一个DSP控制使用链路口进行后续程序的加载,每行程序加载完成便进入等待,同时最后一个节点通过Flag标志位传输Flag标志给程序加载控制器,表明此行加载完毕。以此类推直到最后一行的最后一个节点加载完成加载后,再由程序加载控制器发出启动命令,DSP 阵列才启动。系统启动后每个节点按照节点处理方程进行处理。系统启动后每个节点按照节点处理方程a’ = a , b’ = b, c’ = c + a*b 进行处理。
在系统中还需要一些辅助电路,其中时钟驱动器驱动模块用于整个DSP 阵列的时钟同步驱动。电源管理系统负责整阵列系统的功率分配和控制,及系统的散热冷却管理[12]。我们按照不同的深度的颜色把阵列中的各节点处理器运行负荷表示出来,可以发现,在阵列边沿特别是XY总线相交处的处理器运行负荷较大。
六、结束语
DSP 的运算速度越来越高,一些复杂的算法可以更快的实现。但是随着通信系统的性能的提高,对计算能力的要求也更高。并行处理是提高系统运算能力最有效的方法。因此并
行处理在复杂的阵列信号处理中将有着广泛的应用[12]。
本文作者创新点为：在阵列信号处理中需要大量的矩阵运算，而其中最基本的就是矩阵相乘运算。本文就矩阵相乘的行划分并行实现进行了改进，该并行实现基于PVM环境，采用主从编程模式，然后给出了改进的行划分矩阵相乘在DSP上实现的方案。
七、参考文献
[1] 李晓梅,莫则尧, 胡庆丰等著.可扩展并行算法的设计与分析.北京:国防工业出版社,2000 李晓梅,
[2] 陈国良编著. 并行计算-结构.算法.编程. 北京: 高等教育出版社, 2003
[3] 余祥宣,崔国华,邹海明. 计算机算法基础[M] . 武汉:华中科技大学出版社,2000.
[4] 李小洲, 李庆华.矩阵相乘Cannon并行算法在工作站机群上的实现.计算机与数字工程,Vol 29(5), p5-8
[5] 张锦雄.矩阵相乘并行算法的MPI实现. 广西科学院学报. 2004,11,Vol20 (4)
[6] 唐俊奇. 多处理机中矩阵乘法的算法研究. 中国西部科技. 2007,2,p4-
[7] Geist A, Beguelin A, Dongarra J.PVM:tual machine USA:MIT Press 1994
[8] Wilkinson B ,Allen M. 并行程序设计[ M] . 陆鑫达等译. 北京:机械工业出版社,2002.
[9] 朱福喜,何炎祥, 并行分布计算中的调度算法理论与设计[M].武汉:武汉大学出版社,2003
[10] 张学波,李晓梅.基于对角划分的矩阵乘并行算法[J].计算机工程,2004,(3):42-43
[11] 温馀洪,沈美明.分布并行系统的并行程序设计环境.小型微型计算机系统.1995,16(2)
[12] 雷 晶,金心宇,王 锐.矩阵相乘的并行计算及其DSP实现.传感技术学报.19(3),2006,6
[13] 杨兴国,郭勇,马厚雪. 基于DSP的取样数字式平均器的设计与实现.微计算机信息.2007,23:2-2
