并行程序设计王刚矩阵乘法的并行化的设计与实现.docx
《并行程序设计王刚矩阵乘法的并行化的设计与实现.docx》由会员分享，可在线阅读，更多相关《并行程序设计王刚矩阵乘法的并行化的设计与实现.docx（19页珍藏版）》请在冰豆网上搜索。
并行程序设计王刚矩阵乘法的并行化的设计与实现
南开大学现代远程教育学院考试卷
2020年度春季学期期末（2020.9）《并行程序设计》
主讲教师：
王刚
学习中心：
____________________________专业：
_______________________
姓名：
_________________学号：
_______________成绩：
___________
一、请同学们在下列题目中任选一题，写成期末论文。
（一）并行算法研究类
对某一问题，研究其并行算法的设计、实现，分析其性能，进行实验验证，撰写研究论文。
例如：
1、对矩阵相乘问题，设计pthread多线程结合SSE/AVX的两层并行算法，实现并行程序。
讨论算法层面不同策略对性能的影响，例如多个线程间不同的任务分配方式、不同的线程同步策略等，讨论不同并行编程方法对性能的影响，例如SSE/AVX的对齐和不对齐内存访问等等。
对不同的矩阵规模、不同的线程数测试程序性能，撰写研究论文。
2、对高斯消去法问题（其串行算法伪代码示意如下面算法1所示），设计pthread多线程结合SSE/AVX的两层并行算法，实现并行程序。
讨论算法层面不同策略对性能的影响，例如多个线程间不同的任务分配方式、不同的线程同步策略等，讨论不同并行编程方法对性能的影响，例如SSE/AVX的对齐和不对齐内存访问等等。
对不同的矩阵规模、不同的线程数测试程序性能，撰写研究论文。
3、其他类似难度的问题。
（二）并行编程工具调研类
对某种并行编程工具进行调研，选取某个问题（例如矩阵相乘问题），用这种编程工具编写并行程序求解这个问题，进行实验验证，撰写研究论文介绍这种并行编程工具的特色、基本编程（使用）方法、如何用它解决实际问题（以你选定的问题为例）。
例如：
1、C++、Java等语言本身对并行编程提供的支持。
2、HadoopMapReduce编程工具。
3、其它并行编程工具。
二、论文写作要求
（一）并行算法研究类
1、论文应详细描述清楚所研究的问题，并行算法的设计。
2、鼓励大家选择课堂教学之外的问题，通过文献调研，研究其并行求解方法，甚至有自己提出新的方法。
3、最好能有求解一个问题的多种并行算法之间的对比分析。
（二）并行编程工具调研类
1、应调研较新的工具，避免调研太“古老”的工具。
2、不能只是工具相关资料的调研和文字的汇总、整理，重点仍是并行编程——用调研的工具编程解决一个具体问题。
3、鼓励大家进行不同并行编程工具间的对比，例如调研的工具与课堂讲授的工具之间的对比。
三、论文写作格式要求：
论文题目要求为宋体三号字，加粗居中；
正文部分要求为宋体小四号字，标题加粗，行间距为1.5倍行距；
应符合科技论文写作规范，题目、摘要、关键字、章节、参考文献等等完整、正确。
这方面可参考附件范文。
四、论文提交注意事项：
1、论文一律以此文件为封面，写明学习中心、专业、姓名、学号等信息。
论文保存为word文件，以“课程名+学号+姓名”命名。
2、论文一律采用线上提交方式，在学院规定时间内上传到教学教务平台，逾期平台关闭，将不接受补交。
3、不接受纸质论文。
4、与论文一同打包提交源程序，注意，是提交.cpp、.h等源程序，不要将工程文件、编译后的目标文件等打包提交。
5、如有抄袭雷同现象，将按学院规定严肃处理。
矩阵乘法的并行化的设计与实现
摘要：
矩阵乘法是最基本的矩阵运算之一,由于其计算密集的特点,适合于在FPGA上实现。
本文给定两个n
阶矩阵A
与B
，矩阵乘法是指计算C=A×B
，现在对两个矩阵乘法进行串行和并行的实验和分析。
关键词：
矩阵乘法；并行算法；实验；
一、算法原理:
1、串行算法
通常的O（n3）矩阵乘矩阵的串行计算过程如算法1所示，此外为计算矩阵相乘，还可以有对3层循环采用其他嵌套形式的串行算法。
算法1:
稠密矩阵相乘的i,j,k形式串行算法
2、并行算法
两个矩阵相乘的行列划分并行算法假设一共有P个进程，将矩阵A按行分成P个块，将矩阵B按列分成P个块：
每块包含连续若干个行．为使得负载平衡，应使得每块中的行数尽量相等．将Ak与Bk分别存储在进程Pk的A’与B’中．将C分为P×P块，且将Ci,j存储在
的p’i中，如算法2
算法2稠密矩阵乘C=A×B的行列划分并行算法
实验由MPICH2在VS2010上进行并行环境的配置来完成，单机情况下用进程数的个数模拟多处理器。
在实验中算法由以下几个函数实现：
voidreadData（）;此函数被rankID为0的进程调用，负责从dataIn.txt文件中A[M,K],B[P,N]两个相乘矩阵的数据，并为结果矩阵C[M,N]分配空间。
其中C[N,N]=A[M,K]*B[P,N]。
intgcd（intM,intN,intgroup_size）此函数用来返回两个整数的不大于group_size的最大公因子，即算法所用到的处理器个数，为了保证行划分和列划分可以平均的划分，通过求M，N不大于group_size的最大公因子来确定实际用到的处理器p。
voidprintResult（）;此函数被rankID为0的进程调用，用来将A,B,C矩阵打印输出给用户，并输出用于分发数据和并行计算的时间。
intmain（intargc,char**argv）;程序的主函数。
算法分析（可扩展性分析）：
在LogP模型上，算法2并行执行时间为：
由此可知，并行效率为：
因此，等效率函数为：
3、算法的MPI程序:
//matrix.cpp:
定义控制台应用程序的入口点。
//
#include"stdafx.h"
#include"stdio.h"
#include"stdlib.h"
#include"mpi.h"
#include
#defineintsizesizeof（int）
#definefloatsizesizeof（float）
#definecharsizesizeof（char）
#defineA（x,y）A[x*K+y]
#defineB（x,y）B[x*N+y]
#defineC（x,y）C[x*N+y]
#definea（x,y）a[x*K+y]
#defineb（x,y）b[x*n+y]
#definebuffer（x,y）buffer[x*n+y]/*此宏用来简化对标号为奇数的处理器内的缓冲空间的访问*/
#definec（l,x,y）c[x*N+y+l*n]
float*a,*b,*c,*buffer;
ints;
float*A,*B,*C;/*A[M,K],B[P,N].正确的情况下K应该等于P,否则无法进行矩阵相乘*/
intM,N,K,P;
intm,n;
intmyid;
intp;/*保存工作站集群中处理器数目，也即通信子大小*/
FILE*dataFile;/*用于读取输入文件内容和将计算结果输出到结果文件的临时文件指针*/
MPI_Statusstatus;
doubletime1;
doublestarttime,endtime;
/*
*函数名:
readData
*功能：
此函数被rankID为0的进程调用，负责从dataIn.txt文件中读入
*A[M,K],B[P,N]两个相乘矩阵的数据，并为结果矩阵C[M,N]分配空间。
*其中C[N,N]=A[M,K]*B[P,N]
*输入：
无
*返回值：
无
*/
voidreadData（）
{
inti,j;
starttime=MPI_Wtime（）;
dataFile=fopen（"dataIn.txt","r"）;
fscanf（dataFile,"%d%d",&M,&K）;/*读取矩阵A的行，列数M,K*/
A=（float*）malloc（floatsize*M*K）;/*为矩阵A分配空间*/
for（i=0;i{
for（j=0;j{
fscanf（dataFile,"%f",A+i*K+j）;
}
}
fscanf（dataFile,"%d%d",&P,&N）;/*读取矩阵B的行，列数P,N*/
if（K!
=P）/*K应该等于P,否则矩阵无法相乘*/
{
printf（"theinputiswrong\n"）;
exit
（1）;
}
B=（float*）malloc（floatsize*K*N）;/*为矩阵B分配空间*/
for（i=0;i{
for（j=0;j{
fscanf（dataFile,"%f",B+i*N+j）;
}
}
fclose（dataFile）;
printf（"Inputoffile\"dataIn.txt\"\n"）;
printf（"%d\t%d\n",M,K）;/*输出A矩阵的维数*/
for（i=0;i{
for（j=0;jprintf（"\n"）;
}
printf（"%d\t%d\n",K,N）;/*输出B矩阵的维数*/
for（i=0;i{
for（j=0;jprintf（"\n"）;
}
C=（float*）malloc（floatsize*M*N）;/*为结果矩阵C[M,N]分配空间*/
}
/*
*函数名:
gcd
*功能：
此函数用来返回两个整数的不大于group_size的最大公因子
*输入：
M,N:
要求最大公因数的两个整数
*group_size所求公因子必须小于此参数，此参数代表用户指定的通信子大小
*返回值：
M和N的不大于group_size的最大公因子
*/
intgcd（intM,intN,intgroup_size）
{
inti;
for（i=M;i>0;i--）
{
if（（M%i==0）&&（N%i==0）&&（i<=group_size））
returni;
}
return1;
}
/*
*函数名:
printResult
*功能：
此函数被rankID为0的进程调用，用来将A,B,C矩阵打印输出给用户，
*并输出用于分发数据和并行计算的时间
*输入：
无
*返回值：
无
*/
voidprintResult（）
{
inti,j;
printf（"\nOutputofMatrixC=AB\n"）;
for（i=0;i{
for（j=0;jprintf（"\n"）;
}
endtime=MPI_Wtime（）;
printf（"\n"）;
printf（"Wholerunningtime=%fseconds\n",endtime-starttime）;
printf（"Distributedatatime=%fseconds\n",time1-starttime）;
printf（"Parallelcomputetime=%fseconds\n",endtime-time1）;
}
/*
*函数名:
main
*功能：
程序的主函数
*输入：
argc为命令行参数个数；
*argv为每个命令行参数组成的字符串数组。
*输出：
返回0代表程序正常结束；其它值表明程序出错。
*/
intmain（intargc,char**argv）
{
inti,j,k,l,group_size,mp1,mm1;
MPI_Init（&argc,&argv）;
MPI_Comm_size（MPI_COMM_WORLD,&group_size）;
MPI_Comm_rank（MPI_COMM_WORLD,&myid）;
p=group_size;
//下面一段程序负责从dataIn.txt文件中读入A[M,K],B[P,N]两个相乘矩阵的数据，
//并为结果矩阵C[M,N]分配空间。
C[N,N]=A[M,K]*B[P,N]
//注意这段程序只有编号为0的处理器才执行此步操作
if（myid==0）
{
readData（）;
}
if（myid==0）/*由编号为0的进程将A,B两矩阵的行列维数M,K,N发送给所有其他进程*/
for（i=1;i{
MPI_Send（&M,1,MPI_INT,i,i,MPI_COMM_WORLD）;
MPI_Send（&K,1,MPI_INT,i,i,MPI_COMM_WORLD）;
MPI_Send（&N,1,MPI_INT,i,i,MPI_COMM_WORLD）;
}
else/*编号非0的进程负责接收A,B两矩阵的行列维数M,K,N*/
{
MPI_Recv（&M,1,MPI_INT,0,myid,MPI_COMM_WORLD,&status）;
MPI_Recv（&K,1,MPI_INT,0,myid,MPI_COMM_WORLD,&status）;
MPI_Recv（&N,1,MPI_INT,0,myid,MPI_COMM_WORLD,&status）;
}
p=gcd（M,N,group_size）;
m=M/p;/*m代表将矩阵按行分块后每块的行数*/
n=N/p;/*m代表将矩阵按列分块后每块的列数*/
if（myid{
a=（float*）malloc（floatsize*m*K）;/*a[m,K]用来存储本处理器拥有的矩阵A的行块*/
b=（float*）malloc（floatsize*K*n）;/*b[K,n]用来存储此时处理器拥有的矩阵B的列块*/
c=（float*）malloc（floatsize*m*N）;/*c[m,N]用来存储本处理器计算p-1次得到所有结果*/
if（myid%2!
=0）/*为标号为奇数的处理器分配发送缓冲空间*/
buffer=（float*）malloc（K*n*floatsize）;
if（a==NULL||b==NULL||c==NULL）/*如果分配空间出错，则打印出错信息*/
printf（"Allocatespacefora,borcfail!
"）;
if（myid==0）/*标号为0的处理器将应该它拥有的矩阵A,B的元素读入自己的a,b中*/
{
for（i=0;ifor（j=0;ja（i,j）=A（i,j）;
for（i=0;ifor（j=0;jb（i,j）=B（i,j）;
}
if（myid==0）/*标号为0的处理器将其他处理器的初始数据分别发给各处理器*/
{
for（i=1;i{
MPI_Send（&A（m*i,0）,K*m,MPI_FLOAT,i,i,MPI_COMM_WORLD）;
for（j=0;jMPI_Send（&B（j,n*i）,n,MPI_FLOAT,i,i,MPI_COMM_WORLD）;
}
free（A）;
free（B）;/*至此，A,B两矩阵的数据已经完全被分散到各处理器。
释放A,B所占空间*/
}
else/*标号非0的处理器从0处理器接受各自的初始矩阵数据*/
{
MPI_Recv（a,K*m,MPI_FLOAT,0,myid,MPI_COMM_WORLD,&status）;
for（j=0;jMPI_Recv（&b（j,0）,n,MPI_FLOAT,0,myid,MPI_COMM_WORLD,&status）;
}
if（myid==0）
time1=MPI_Wtime（）;/*标号为0的处理器记录开始矩阵相乘计算的时间*/
for（i=0;i{
l=（i+myid）%p;
for（k=0;kfor（j=0;jfor（c（l,k,j）=0,s=0;sc（l,k,j）+=a（k,s）*b（s,j）;
mm1=（p+myid-1）%p;/*计算本进程的前一个进程的标号*/
mp1=（myid+1）%p;/*计算本进程的后一个进程的标号*/
if（i!
=p-1）
{
if（myid%2==0）/*偶数号处理器先发送后接收*/
{
MPI_Send（b,K*n,MPI_FLOAT,mm1,mm1,MPI_COMM_WORLD）;
MPI_Recv（b,K*n,MPI_FLOAT,mp1,myid,MPI_COMM_WORLD,&status）;
}
else/*奇数号处理器先将B的列块存于缓冲区buffer中，然后接收编号在其后面的
处理器所发送的B的列块，最后再将缓冲区中原矩阵B的列块发送给编号
在其前面的处理器*/
{
for（k=0;kfor（j=0;jbuffer（k,j）=b（k,j）;
MPI_Recv（b,K*n,MPI_FLOAT,mp1,myid,MPI_COMM_WORLD,&status）;
MPI_Send（buffer,K*n,MPI_FLOAT,mm1,mm1,MPI_COMM_WORLD）;
}
}
}
if（myid==0）/*标号为0的进程直接将计算结果保存到结果矩阵C中*/
for（i=0;ifor（j=0;jC（i,j）=*（c+i*N+j）;
if（myid!
=0）/*标号非0的进程则要把计算结果发送到标号为0的处理器中去*/
MPI_Send（c,m*N,MPI_FLOAT,0,myid,MPI_COMM_WORLD）;
else/*标号为0的进程负责接收其他进程的计算结果并保存到结果矩阵C中*/
{
for（k=1;k{
MPI_Recv（c,m*N,MPI_FLOAT,k,k,MPI_COMM_WORLD,&status）;
for（i=0;ifor（j=0;jC（（k*m+i）,j）=*（c+i*N+j）;
}
}
if（myid==0）/*0号处理器负责将A,B,C矩阵打印输出给用户，并输出用于分发数据和并行计算的时间*/
printResult（）;
}
MPI_Finalize（）;
if（myid{
free（a）;
free（b）;
free（c）;
if（myid==0）/*只有0号进程要释放C*/
free（C）;
if（myid%2!
=0）/*只有奇数号进程要释放buffer*/
free（buffer）;
}
printf（"%d",p）;
getchar（）;
return（0）;
}
三、算例结果及分析
上面的dataIn.txt是两个10乘以10的矩阵，总的时间是0.197588秒，可以看到主要的时间用来通讯和分发数据，并行计算的时间其实较短。
换一个算例，即10乘以20矩阵乘以20乘以10的矩阵，下面是运行时间：
四、总结
数据量增大，但是并行时间并没有大量增加，而是在通讯时间上面消耗更多时间，跟算法可扩展的等效率分析一致。
参考文献
[1]都志辉，《高性能计算之并行编程技术——MPI并行程序设计》
[2]何红旗,邵仪.矩阵乘法的FPGA并行设计与实现[C]//计算机研究新进展（2010）——河南省计算机学会2010年学术年会论文集.2010.
升级会员 