(完整word版)矩阵乘法的OpenMP实现及性能分析
- 格式:doc
- 大小:210.50 KB
- 文档页数:10
一、实验目的1. 理解矩阵乘法的概念和运算规则。
2. 掌握矩阵乘法的编程实现方法。
3. 通过实验验证矩阵乘法的正确性。
二、实验环境1. 操作系统:Windows 102. 编程语言:Python3. 库:NumPy三、实验原理矩阵乘法是指两个矩阵相乘的运算。
设矩阵A为m×n的矩阵,矩阵B为n×p的矩阵,则它们的乘积C为一个m×p的矩阵。
矩阵乘法的运算规则如下:C[i][j] = Σ(A[i][k] B[k][j]),其中k为1到n的整数。
四、实验步骤1. 导入NumPy库。
```pythonimport numpy as np```2. 定义矩阵A和B。
```pythonA = np.array([[1, 2], [3, 4]])B = np.array([[5, 6], [7, 8]])```3. 计算矩阵A和B的乘积C。
```pythonC = np.dot(A, B)```4. 打印结果。
```pythonprint("矩阵A:")print(A)print("矩阵B:")print(B)print("矩阵C(A乘B):")print(C)```五、实验结果与分析1. 运行实验程序,得到以下结果:```矩阵A:[[1 2][3 4]]矩阵B:[[5 6][7 8]]矩阵C(A乘B):[[19 22][43 50]]```2. 分析结果:- 矩阵A为2×2的矩阵,矩阵B为2×2的矩阵,它们的乘积C为2×2的矩阵。
- 根据矩阵乘法的运算规则,我们可以计算出矩阵C的每个元素。
- 实验结果与理论计算相符,说明矩阵乘法的编程实现是正确的。
六、实验总结1. 本实验成功实现了矩阵乘法的编程,验证了矩阵乘法的正确性。
2. 通过实验,加深了对矩阵乘法概念和运算规则的理解。
3. NumPy库在矩阵运算方面具有强大的功能,为编程提供了便利。
DenseMulMatrixMPI.c#include <stdio.h>#include <stdlib.h>#include <mpi.h>#include <time.h>//#include <windows.h>#define Length 1000int *A,*B,*C;//A、B、C一维数组分别记录矩阵A、B、C,以行优先的方式存储int *buffer,*ans;int temp;//临时变量,存储乘法过程中的中间值double startTime,endTime,totalTime;//startTime endTime totalTime分别表示矩阵乘之前的时间、之后的时间,以及矩阵乘所花费的时间int procsID,procsNum;//进程标识、进程个数int i,j,k;int line;//每个进程可以分配的矩阵行数//初始化矩阵/*void InitMatrix(int *M,int len){//产生0到49之间的随机数//int count=0;srand((unsigned)time( NULL));for(i=0; i < len; i++){for(j = 0; j < len; j ++){M[i * len + j] = rand() % 50;//if (M[i * len + j] < 20 )// M[i * len + j] = 0;}}//统计矩阵中非元素的个数,占的个数达总个数的50%以上的话即可认为是稠密矩阵,因为矩阵元素是随机产生的所以只是大概估计//for(i=0; i < len*len; i++){// if(M[i]!=0)// count++;//}//printf("%d\n",count);}*/void InitMatrix(int *M,int len){//int count=0;srand((unsigned)time( NULL));for(i=0; i < len*len; i++){M[i] = rand() % 2;}//统计矩阵中非0元素的个数,占的个数达总个数的50%以上的话即可认为是稠密矩阵,因为矩阵元素是随机产生的所以只是大概估计//for(i=0; i < len*len; i++){// if(M[i]!=0)// count++;//}//printf("%d\n",count);}int main(int argc,char *argv[]){MPI_Status status;MPI_Init(&argc,&argv);MPI_Comm_rank(MPI_COMM_WORLD,&procsID);//获得当前进程号MPI_Comm_size(MPI_COMM_WORLD,&procsNum);//获得进程个数line = Length/procsNum;//将数据分为(进程数)个块,主进程也要处理数据A = (int*)malloc(sizeof(int)*Length*Length);B = (int*)malloc(sizeof(int)*Length*Length);C = (int*)malloc(sizeof(int)*Length*Length);//缓存大小大于等于要处理的数据大小,大于时只需关注实际数据那部分buffer = (int*)malloc(sizeof(int)*Length*line);//数据分组大小ans = (int*)malloc(sizeof(int)*Length*line);//保存数据块结算的结果//主进程对矩阵赋初值,并将矩阵B广播到各进程,将矩阵A分组广播到各进程if (procsID==0){//矩阵赋初值InitMatrix(A,Length);//Sleep(1000);//停顿1s,使得矩阵A和B产生的随机元素不一样(但是在命令窗口不识别它的头文件windows.h)InitMatrix(B,Length);startTime = MPI_Wtime();//将矩阵B发送给其他从进程for (i=1;i<procsNum;i++){MPI_Send(B,Length*Length,MPI_INT,i,0,MPI_COMM_WORLD);}//依次将A的各行发送给各从进程for (i=1;i<procsNum;i++){MPI_Send(A+(i-1)*line*Length,Length*line,MPI_INT,i,1,MPI_COMM_WORLD);}//接收从进程计算的结果for (k=1;k<procsNum;k++){MPI_Recv(ans,line*Length,MPI_INT,k,3,MPI_COMM_WORLD,&status);//将结果传递给数组Cfor (i=0;i<line;i++){for (j=0;j<Length;j++){C[((k-1)*line+i)*Length+j]=ans[i*Length+j];}}}//计算A剩下的数据for (i=(procsNum-1)*line;i<Length;i++){for (j=0;j<Length;j++){temp=0;for (k=0;k<Length;k++)temp += A[i*Length+k]*B[k*Length+j];C[i*Length+j]=temp;}}endTime = MPI_Wtime();totalTime=endTime-startTime;printf("并行稠密矩阵乘法过程总共花的时间:%.4fs\n",totalTime);}//其他进程接收数据,计算结果后,发送给主进程else{//接收广播的数据(矩阵N)MPI_Recv(B,Length*Length,MPI_INT,0,0,MPI_COMM_WORLD,&status);MPI_Recv(buffer,Length*line,MPI_INT,0,1,MPI_COMM_WORLD,&status);//计算乘积结果,并将结果发送给主进程for (i=0;i<line;i++){for (j=0;j<Length;j++){temp=0;for(k=0;k<Length;k++)temp += buffer[i*Length+k]*B[k*Length+j];ans[i*Length+j]=temp;}}//将计算结果传送给主进程MPI_Send(ans,line*Length,MPI_INT,0,3,MPI_COMM_WORLD);}MPI_Finalize();//结束return 0;}DenseMulMatrixSerial.c#include <stdio.h>#include <stdlib.h>#include <time.h>//#include <windows.h>#define Length 1000int *A,*B,*C;//A、B、C一维数组分别记录矩阵A、B、C,以行优先的方式存储int i,j,k;clock_t startTime, endTime;double totalTime;//初始化矩阵/*void InitMatrix(int *M,int len){//产生0到49之间的随机数//int count=0;srand((unsigned)time( NULL));for(i=0; i < len; i++){for(j = 0; j < len; j ++){M[i * len + j] = rand() % 50;//if (M[i * len + j] < 20 )// M[i * len + j] = 0;}}//统计矩阵中非元素的个数,占的个数达总个数的50%以上的话即可认为是稠密矩阵,因为矩阵元素是随机产生的所以只是大概估计//for(i=0; i < len*len; i++){// if(M[i]!=0)// count++;//}//printf("%d\n",count);}*/void InitMatrix(int *M,int len){//int count=0;srand((unsigned)time( NULL));for(i=0; i < len*len; i++){M[i] = rand() % 2;}//统计矩阵中非0元素的个数,占的个数达总个数的50%以上的话即可认为是稠密矩阵,因为矩阵元素是随机产生的所以只是大概估计//for(i=0; i < len*len; i++){// if(M[i]!=0)// count++;//}//printf("%d\n",count);}int main(){//矩阵A乘以矩阵B,结果为矩阵C,A、B、C均用一维的数组存储,有Length*Length个单元//动态为矩阵分配空间A = (int *)malloc(sizeof(int)*Length*Length);B = (int *)malloc(sizeof(int)*Length*Length);C = (int *)malloc(sizeof(int)*Length*Length);//初始化矩阵A BInitMatrix(A,Length);//Sleep(1000);//停顿1s,使得矩阵A和B中随机产生的元素不一样InitMatrix(B,Length);startTime = clock();//矩阵乘法过程for(i = 0; i < Length; i ++){for(j = 0; j < Length; j ++){C[i * Length + j] = 0;for (k = 0; k < Length; ++k){C[i * Length + j] += A[i * Length + k] * B[k * Length + j];}}}endTime = clock();//串行矩阵乘法总共用的时间totalTime = (double)(endTime - startTime) / CLOCKS_PER_SEC;printf("串行稠密矩阵乘法过程总共花的时间:%.4fs\n",totalTime);return 0;}SparseMulMatrixMPI.c#include <stdio.h>#include <stdlib.h>#include <mpi.h>#include <time.h>//#include <windows.h>#define Length 1000int *A,*B,*C;//A、B、C一维数组分别记录矩阵A、B、C,以行优先的方式存储int *buffer,*ans;int temp;//临时变量,存储乘法过程中的中间值double startTime,endTime,totalTime;//startTime endTime totalTime分别表示矩阵乘之前的时间、之后的时间,以及矩阵乘所花费的时间int procsID,procsNum;//进程标识、进程个数int i,j,k;int m,n;//统计矩阵中1时用int line;//每个进程可以分配的矩阵行数//初始化矩阵/*void InitMatrix(int *M,int len){//产生0 48 49 三个数//int count=0;srand((unsigned)time( NULL));for(i=0; i < len; i++){for(j = 0; j < len; j ++){M[i * len + j] = rand() % 50;if (M[i * len + j] < 48 )M[i * len + j] = 0;}}//统计矩阵中0元素的个数,占的个数小于总个数的5%的话即为稀疏矩阵,因为矩阵元素是随机产生的所以只是大概估计//for(i=0; i < len*len; i++){// if(M[i]==0)// count++;//}//printf("%d\n",count);}*/void InitMatrix(int *M,int len){//根据稀疏矩阵非0元素少于总元素个数的5%的原则设置下面的程序,随机产生元素位置,令其位置上的元素为1,其余为0//int count=0;for(i=0;i<len*len;i++)M[i]=0;srand((unsigned)time( NULL));for(m=0;m<224;m++){for(n=0;n<224;n++){i=rand()%1000;j=rand()%1000;M[i*len+j]=1;}}//统计矩阵中0元素的个数,占的个数小于总个数的5%的话即为稀疏矩阵,因为矩阵元素是随机产生的所以只是大概估计//for(i=0; i < len*len; i++){// if(M[i]==0)// count++;//}//printf("%d\n",count);}int main(int argc,char *argv[]){MPI_Status status;MPI_Init(&argc,&argv);MPI_Comm_rank(MPI_COMM_WORLD,&procsID);//获得当前进程号MPI_Comm_size(MPI_COMM_WORLD,&procsNum);//获得进程个数line = Length/procsNum;//将数据分为(进程数)个块,主进程也要处理数据A = (int*)malloc(sizeof(int)*Length*Length);B = (int*)malloc(sizeof(int)*Length*Length);C = (int*)malloc(sizeof(int)*Length*Length);//缓存大小大于等于要处理的数据大小,大于时只需关注实际数据那部分buffer = (int*)malloc(sizeof(int)*Length*line);//数据分组大小ans = (int*)malloc(sizeof(int)*Length*line);//保存数据块结算的结果//主进程对矩阵赋初值,并将矩阵B广播到各进程,将矩阵A分组广播到各进程if (procsID==0){//矩阵赋初值InitMatrix(A,Length);//Sleep(1000);//停顿1s,使得矩阵A和B产生的随机元素不一样(但是在命令窗口不识别它的头文件windows.h)InitMatrix(B,Length);startTime = MPI_Wtime();//将矩阵B发送给其他从进程for (i=1;i<procsNum;i++){MPI_Send(B,Length*Length,MPI_INT,i,0,MPI_COMM_WORLD);}//依次将A的各行发送给各从进程for (i=1;i<procsNum;i++){MPI_Send(A+(i-1)*line*Length,Length*line,MPI_INT,i,1,MPI_COMM_WORLD);}//接收从进程计算的结果for (k=1;k<procsNum;k++){MPI_Recv(ans,line*Length,MPI_INT,k,3,MPI_COMM_WORLD,&status);//将结果传递给数组Cfor (i=0;i<line;i++){for (j=0;j<Length;j++){C[((k-1)*line+i)*Length+j]=ans[i*Length+j];}}}//计算A剩下的数据for (i=(procsNum-1)*line;i<Length;i++){for (j=0;j<Length;j++){temp=0;for (k=0;k<Length;k++)temp += A[i*Length+k]*B[k*Length+j];C[i*Length+j]=temp;}}endTime = MPI_Wtime();totalTime=endTime-startTime;printf("并行稀疏矩阵乘法过程总共花的时间:%.4fs\n",totalTime);}//其他进程接收数据,计算结果后,发送给主进程else{//接收广播的数据(矩阵N)MPI_Recv(B,Length*Length,MPI_INT,0,0,MPI_COMM_WORLD,&status);MPI_Recv(buffer,Length*line,MPI_INT,0,1,MPI_COMM_WORLD,&status);//计算乘积结果,并将结果发送给主进程for (i=0;i<line;i++){for (j=0;j<Length;j++){temp=0;for(k=0;k<Length;k++)temp += buffer[i*Length+k]*B[k*Length+j];ans[i*Length+j]=temp;}}//将计算结果传送给主进程MPI_Send(ans,line*Length,MPI_INT,0,3,MPI_COMM_WORLD);}MPI_Finalize();//结束return 0;}SparseMulMatrixSerial.c#include <stdio.h>#include <stdlib.h>#include <time.h>//#include <windows.h>#define Length 1000int *A,*B,*C;//A、B、C一维数组分别记录矩阵A、B、C,以行优先的方式存储int i,j,k;int m,n;//统计矩阵中1时用clock_t startTime, endTime;double totalTime;//初始化矩阵/*void InitMatrix(int *M,int len){//产生0 48 49 三个数//int count=0;srand((unsigned)time( NULL));for(i=0; i < len; i++){for(j = 0; j < len; j ++){M[i * len + j] = rand() % 50;if (M[i * len + j] < 48 )M[i * len + j] = 0;}}//统计矩阵中0元素的个数,占的个数小于总个数的5%的话即为稀疏矩阵,因为矩阵元素是随机产生的所以只是大概估计//for(i=0; i < len*len; i++){// if(M[i]==0)// count++;//}//printf("%d\n",count);}*/void InitMatrix(int *M,int len){//根据稀疏矩阵非0元素少于总元素个数的5%的原则设置下面的程序,随机产生元素位置,令其位置上的元素为1,其余为0//int count=0;for(i=0;i<len*len;i++)M[i]=0;srand((unsigned)time( NULL));for(m=0;m<224;m++){for(n=0;n<224;n++){i=rand()%1000;j=rand()%1000;M[i*len+j]=1;}}//统计矩阵中0元素的个数,占的个数小于总个数的5%的话即为稀疏矩阵,因为矩阵元素是随机产生的所以只是大概估计//for(i=0; i < len*len; i++){// if(M[i]==0)// count++;//}//printf("%d\n",count);}int main(){//矩阵A乘以矩阵B,结果为矩阵C,A、B、C均用一维的数组存储,有Length*Length个单元//动态为矩阵分配空间A = (int *)malloc(sizeof(int)*Length*Length);B = (int *)malloc(sizeof(int)*Length*Length);C = (int *)malloc(sizeof(int)*Length*Length);//初始化矩阵A BInitMatrix(A,Length);//Sleep(1000);//停顿1s,使得矩阵A和B中随机产生的元素不一样InitMatrix(B,Length);startTime = clock();//矩阵乘法过程for(i = 0; i < Length; i ++){for(j = 0; j < Length; j ++){C[i * Length + j] = 0;for (k = 0; k < Length; ++k){C[i * Length + j] += A[i * Length + k] * B[k * Length + j];}}}endTime = clock();//串行矩阵乘法总共用的时间totalTime = (double)(endTime - startTime) / CLOCKS_PER_SEC;printf("串行稀疏矩阵乘法过程总共花的时间:%.4fs\n",totalTime);return 0;}附:另一种矩阵初始化的办法//稠密矩阵的生成方法void InitMatrix(int *M,int *N,int len){srand((unsigned)time( NULL));for(i=0; i < len*len; i++){M[i] = rand() % 2;}for(i=0;i<len;i++){for(j=0;j<len;j++){N[i*len+j]=M[j*len+i];}}}//稀疏矩阵的生成方法void InitMatrix(int *M, int *N, int len){for(i=0;i<len*len;i++)M[i]=0;srand((unsigned)time( NULL));for(m=0;m<224;m++){for(n=0;n<224;n++){i=rand()%len;j=rand()%len;M[i*len+j]=1;}}for(i=0;i<len;i++){for(j=0;j<len;j++){N[i*len+j]=M[j*len+i];}}}。
多核软件设计实验指导――OpenMP矩阵相乘开发者:开发时间:版本号:一. 问题描述矩阵相乘是线性代数中最常见的问题之一,它在数值计算中有广泛的应用,在计算机的世界里,矩阵相乘扮演着一个不可或缺的角色。
设A 和B 是2个 n n ⨯ 矩阵,,它们的乘积AB 同样是一个n n ⨯矩阵。
A 和B 乘积矩阵C 中元素C ]][[]][[1]][[j k B k i A nk j i ∑==二、串行算法描述若以上面的定义来计算A 和B 的乘积矩阵C ,则每计算C 的一个元素]][[j i C ,需要作n 次乘法运算和n -1次加法运算。
因此,算出矩阵C 的2n 个元素所需的计算时间为)(3n O 。
从程序运行的效率上来分析,如果n 比较大的时候,会耗费大量的空间和时间。
20世纪60年代末期,Strassen 采用了类似在大整数乘法中用过的分治技术,将计算2个n 阶矩阵乘积所需的计算时间改进到)()(81.27log n n O =O ,其基本思想还是使用分治法。
显然,这也没能够大幅度地提高运行速度和效率。
将矩阵a 与矩阵b 相乘得到矩阵c 的代码如下:for(i=0;i<dim;i++){ for(j=0;j<dim;j++){ c(i,j)=0;for(k=0;k<dim;k++){c(i ,j)+=a(i ,k)*b(k,j);}}}三、并行算法下面使用了OpenMP ,修改后的矩阵相乘程序使最外面的迭代在多个线程间静态分割,如图(a )所示:#pragma omp paralllel default(private) shared(a,b,c,dim) num_threads(2) #pragma omp for schedule(static)for(i=0;i<dim;i++){ for(j=0;j<dim;j++){ c(i,j)=0;for(k=0;k<dim;k++){c(i ,j)+=a(i ,k)*b(k,j);}}}Static调度类的一般形式为schedule(static[,chunk-size])。
mali opencl 矩阵乘法Mali OpenCL 矩阵乘法Mali OpenCL 是一种高性能计算框架,可用于在 Mali GPU 上执行并行计算。
利用 OpenCL,开发人员可以编写高效的代码,充分利用 GPU 的并行处理能力。
矩阵乘法是一种常见的线性代数运算,在图像处理、机器学习和科学计算等众多领域有广泛的应用。
本文将重点介绍如何使用 Mali OpenCL 实现矩阵乘法的并行化。
OpenCL 内核OpenCL 内核是并行执行的代码块,由工作组中的多个工作项组成。
对于矩阵乘法,我们可以定义一个内核,该内核执行以下操作:- 每个工作项负责计算乘法矩阵中单个元素的结果。
- 工作项获取其在矩阵中的行号和列号。
- 工作项迭代乘法矩阵的两个输入矩阵,累加部分和,并将其存储在输出矩阵中。
工作组和工作项工作组是一组同时执行的工作项。
Mali GPU 通常使用 32 个工作项组成的工作组。
工作组的数量由全局工作大小确定,而每个工作组中工作项的数量由局部工作大小确定。
数据结构在 OpenCL 中,矩阵表示为一维数组,其中数组索引对应于矩阵的行和列。
为了实现矩阵乘法,我们需要定义三个一维数组来存储输入矩阵 (A 和 B) 和输出矩阵 (C)。
代码示例以下是如何使用 Mali OpenCL 实现矩阵乘法的代码示例:```// 定义输入/输出矩阵的维度const int matrixSize = 1024;// 创建输入/输出矩阵的 OpenCL 缓冲区cl_mem aBuffer = clCreateBuffer(context,CL_MEM_READ_ONLY | CL_MEM_ALLOC_HOST_PTR, sizeof(float) matrixSize matrixSize, NULL, &err);cl_mem bBuffer = clCreateBuffer(context,CL_MEM_READ_ONLY | CL_MEM_ALLOC_HOST_PTR, sizeof(float) matrixSize matrixSize, NULL, &err);cl_mem cBuffer = clCreateBuffer(context,CL_MEM_WRITE_ONLY | CL_MEM_ALLOC_HOST_PTR, sizeof(float) matrixSize matrixSize, NULL, &err);// 将输入矩阵从主机复制到设备clEnqueueWriteBuffer(commandQueue, aBuffer, CL_FALSE, 0, sizeof(float) matrixSize matrixSize, a, 0, NULL, NULL);clEnqueueWriteBuffer(commandQueue, bBuffer, CL_FALSE, 0, sizeof(float) matrixSize matrixSize, b, 0, NULL, NULL);// 设置内核参数clSetKernelArg(kernel, 0, sizeof(cl_mem), &aBuffer);clSetKernelArg(kernel, 1, sizeof(cl_mem), &bBuffer);clSetKernelArg(kernel, 2, sizeof(cl_mem), &cBuffer);// 设置全局/局部工作大小const size_t globalWorkSize[] = { matrixSize, matrixSize };const size_t localWorkSize[] = { 32, 32 };// 提交内核执行clEnqueueNDRangeKernel(commandQueue, kernel, 2, NULL, globalWorkSize, localWorkSize, 0, NULL, NULL);// 将输出矩阵从设备复制到主机clEnqueueReadBuffer(commandQueue, cBuffer, CL_FALSE, 0, sizeof(float) matrixSize matrixSize, c, 0, NULL, NULL);```优化可以通过多种技术来优化 Mali OpenCL 矩阵乘法实现:- 使用本地内存:局部工作组中的工作项可以访问共享的本地内存。
1. 实验目的1.1掌握集群的使用方法。
1.2掌握以并行的方式分析问题、设计并行程序的方法。
1.3掌握如何对并行程序进行简单的性能分析2. 实验要求2.1使用MPI、OpenMp等并行程序设计方法设计矩阵乘法的并行程序。
2.2随机产生所需的矩阵元素,数据项不得少于1000*1000。
2.3尽量设计较高的加速比3. 实验环境3.1硬件环境:两个集群节点blade13、blade15。
3.2软件环境:Linux、gcc、Win7、VC++6.0。
3.3连接方式:Xmanager Enterprise4.0远程桌面连接211.69.198.203。
4. 实验程序4.1随机算法产生矩阵:srand((unsigned int)time(NULL));for (i=0; i<N; i++){for (j=0; j<N; j++){A[i][j] = rand() % 10;B[i][j] = rand() % 10;C[i][k] = 0;}}4.2串行程序设计time(&start);for (i=0; i<M; i++){for (k=0; k<M; k++){C[i][k] = 0;for (j=0; j<M; j++){C[i][k] += A[i][j]*B[j][k];}}}4.3并行程序设计MPI_Init(&argc,&argv) 和MPI_Finalize() MPI_Init用来初始化MPI执行环境,建立多个MPI 进程之间的联系,为后续通信做准备。
而MPI_Finalize则是结束MPI执行环境。
这两个函数就是定义MPI程序的并行区的,除了检测是否初始化的函数之外,不应该在这两个函数定义的区域外调用其它MPI函数。
这两个函数都返回整型值,标识函数是否调用成功。
intMPI_Comm_rank(MPI_Comm comm, int *rank) MPI_Comm_rank函数用来标识各个MPI进程,获取调用该函数进程的进程号,将自身与其他进程区分。
矩阵乘法是线性代数中一个基础且重要的运算,广泛应用于科学计算、数据分析和工程领域。
为了深入了解矩阵乘法的并行化实现,提高计算效率,本实验旨在通过MPI(Message Passing Interface)并行编程技术,实现矩阵乘法的并行计算,并分析其性能。
二、实验内容与方法1. 实验环境操作系统:Ubuntu Linux编译器:gcc并行计算平台:632核CPU、400GB内存的分布内存并行计算平台2. 实验方法(1)矩阵乘法算法本实验采用经典的矩阵乘法算法,即按行优先顺序进行计算。
具体步骤如下:① 将矩阵A、B、C划分为p个块,每个块包含m/p行和n/p列。
② 每个进程负责计算C的一个子块,即计算A的m/p行与B的n/p列的乘积。
③ 进程间通过MPI通信进行数据交换,实现并行计算。
(2)MPI编程本实验采用MPI编程实现矩阵乘法。
主要使用以下MPI语句:① MPI_Init:初始化MPI环境。
② MPI_Comm_size:获取进程总数。
③ MPI_Comm_rank:获取进程编号。
④ MPI_Send:发送数据。
⑤ MPI_Recv:接收数据。
⑥ MPI_Finalize:结束MPI环境。
1. 矩阵乘法结果验证通过比较串行计算和并行计算的结果,验证了程序的正确性。
2. 性能分析(1)执行时间在固定矩阵规模n=1000的情况下,分别测试进程数取1、2、4、8、16、32、64时的执行时间。
结果表明,随着进程数的增加,执行时间逐渐减少,且呈近似线性关系。
(2)加速比加速比是指并行计算时间与串行计算时间的比值。
本实验计算了不同进程数下的加速比,发现随着进程数的增加,加速比逐渐提高,且在进程数达到一定数量后,加速比趋于稳定。
(3)并行效率并行效率是指实际加速比与理论加速比之比。
本实验计算了不同进程数下的并行效率,发现随着进程数的增加,并行效率逐渐提高,但存在一个峰值,之后逐渐降低。
四、实验结论与展望1. 结论本实验通过MPI并行编程技术实现了矩阵乘法的并行计算,验证了程序的正确性。
一. 实验目的1) 用OpenMP 实现最基本的数值算法“矩阵乘法” 2) 掌握for 编译制导语句 3) 对并行程序进行简单的性能二. 实验环境1) 硬件环境:32核CPU 、32G 内存计算机;2) 软件环境:Linux 、Win2003、GCC 、MPICH 、VS2008;4) Windows 登录方式:通过远程桌面连接192.168.150.197,用户名和初始密码都是自己的学号。
三. 实验内容1. 用OpenMP 编写两个n 阶的方阵a 和b 的相乘程序,结果存放在方阵c 中,其中乘法用for 编译制导语句实现并行化操作,并调节for 编译制导中schedule 的参数,使得执行时间最短,写出代码。
方阵a 和b 的初始值如下:⎥⎥⎥⎥⎥⎥⎥⎥⎦⎤⎢⎢⎢⎢⎢⎢⎢⎢⎣⎡-++++=12,...,2,1,..2,...,5,4,31,...,4,3,2,...,3,2,1n n n n n n n a ⎥⎥⎥⎥⎥⎥⎥⎥⎦⎤⎢⎢⎢⎢⎢⎢⎢⎢⎣⎡=1,...,1,1,1..1,...,1,1,11,...,1,1,11,...,1,1,1b 输入:方阵的阶n 、并行域的线程数 输出:c 中所有元素之和、程序的执行时间 提示:a,b,c 的元素定义为int 型,c 中所有元素之各定义为long long 型。
Windows 计时:用<time.h>中的clock_t clock( void )函数得到当前程序执行的时间 Linux 计时: #include <sys/time.h> timeval start,end;gettimeofday(&start,NULL);gettimeofday(&end,NULL);cout<<"execution time:"<< (__sec)+(double)(__usec)/ 1000000<<"seconds" <<endl;答:在windows下使用Microsofe Visual Studio编程,源代码如下:#include<omp.h>#include<stdio.h>#include<time.h>#define NN 2000int a[NN][NN], b[NN][NN];long long c[NN][NN];void solve(int n, int num_thread){int i, j, t, k, time;clock_t startTime, endTime;long long sum;omp_set_num_threads(num_thread);for(i=0;i<n;i++)//对矩阵a和矩阵b进行初始化{t=i+1;for(j=0;j<n;j++){a[i][j]=t++;b[i][j]=1;}}startTime=clock();sum=0;#pragma omp parallel shared(a,b,c) private(i,j,k){#pragma omp for schedule(dynamic)for(i=0;i<n;i++){for(j=0;j<n;j++){c[i][j]=0;for(k=0;k<n;k++){c[i][j]+=a[i][k]*b[k][j];}}}}for(i=0;i<n;i++)for(j=0;j<n;j++) sum+=c[i][j];endTime=clock();time=endTime-startTime;printf("sum=%lld time=%dms\n",sum,time);}int main(){int n, num_thread;while(scanf("%d%d",&n,&num_thread)!=EOF){solve(n,num_thread);}return 0;}2.分析矩阵相乘程序的执行时间、加速比和效率:方阵阶固定为1000,节点数分别取1、2、4、8、16和32时,为减少误差,每项实验进行5次,取平均值作为实验结果。
答:串行执行时程序的执行时间为:T = 15.062s加速比=顺序执行时间/并行执行时间效率=加速比/节点数表1 不同节点数下程序的执行时间(秒)第1次 16.640 8.172 4.078 2.125 1.093 0.594 第2次 16.422 8.156 4.172 2.141 1.078 0.578 第3次 16.406 8.266 4.078 2.125 1.094 0.563 第4次 16.7818.172 4.079 2.109 1.094 0.563 第5次 16.4228.1714.078 2.125 1.093 0.578 平均值16.5342 8.18744.09702.12501.09040.5752图1 不同节点数下程序的执行时间图2 不同节点数下程序的加速比图3 不同节点数下程序的效率执行时间的分析:随着节点数的增加,程序的执行时间减少,大概可以从结果中得出,随着节点书的增加一倍,执行时间减少一半加速比的分析:随着节点数的增加,程序的加速比增加,大概可以从结果中得出,随着节点书的增加一倍,加速相应的增加接近一倍效率的分析:随着节点数的增加,程序的效率逐渐减少3.分析矩阵相乘程序的问题规模与效率的关系:固定节点数为4,让方阵阶从200到1600之间变化,每隔100取一个值。
(为了减少时间,每项实验可只执行1次)答:表2 相同节点数下不同问题规模程序的执行时间与效率方阵阶数并行执行时间串行执行时间效率200 0.015 0.047 0.783333 300 0.016 0.109 1.703125 400 0.063 0.297 1.178571 500 0.156 0.657 1.052885 600 0.406 1.64 1.009852 700 0.907 3.578 0.986218 800 1.609 6.36 0.988191 900 2.578 10.109 0.980314 1000 3.812 14.891 0.976587 1100 5.39 21.032 0.97551 1200 7.344 28.734 0.978145 1300 9.688 37.937 0.978969 1400 12.422 48.64 0.978908 1500 15.656 60.938 0.973077 1600 19.234 74.829 0.972614图3.1 不同问题规模下程序的效率问题规模与效率的关系分析:随着问题规模的增加,程序的效率趋于稳定,但是略微有点下降。
嵌套循环中,如果外层循环迭代次数较少时,如果将来CPU核数增加到一定程度时,创建的线程数将可能小于CPU核数。
另外如果内层循环存在负载平衡的情况下,很难调度外层循环使之达到负载平衡。
下面以矩阵乘法作为例子来讲述如何将嵌套循环并行化,以满足上述扩展性和负载平衡需求。
一个串行的矩阵乘法的函数代码如下:/**矩阵串行乘法函数@param int*a -指向要相乘的第个矩阵的指针@param int row_a -矩阵a的行数@param int col_a -矩阵a的列数@param int *b –指向要想成的第个矩阵的指针@param int row_b -矩阵b的行数@param int col_b -矩阵b的列数@param int *c -计算结果的矩阵的指针@param int c_size -矩阵c的空间大小(总元素个数)@return void –无*/void Martrix_Multiply(int *a, int row_a,int col_a,int*b,int row_b,int col_b,int*c,int c_size){If(col_a!=row_b||c_size<row_a*col_b){return;}int i,j,k;//#pragma omp for private(i,j,k)for(i = 0;i<row_a;i++){int row_i=i*col_a;int row_c=i*col_b;for(j=0;j<col_b;j++){c[row_c+j]=0;for(k=0;k<row_b;k++){c[row_c+j]+=a[row_i+k]*b[k*col_b+j];}}}}如果在外层循环前面加上OpenMP的for语句时,它就变成了一个并行的矩阵乘法函数,但是这样简单地将其并行化显然无法满足前面所述的扩展性需求。
其实可以采用一个简单地方法将最外层循环和第2层循环合并成一个循环,下面便是采用合并循环后的并行实现。
void Parallel_Matrix_Multiply(int *a,int row_a,int col_a,int *b,int row_b,int col_b,int *c,int c_size){If(col_a!=row_b){return;}int i,j,k;int index;int border=row_a*col_b;i=0;j=0;//#pragma omp parallel private(i,j,k) num_threads(dtn(border,1))for(index = 0;index<border;index++){i=index/col_b;j=index%col_b;int row_i=i*col_a;int row_c=i*col_b;c[row_c+j]=0;for(k=0;k<row_b;k++){c[row_c+j]+=a[row_i+k]*b[k*col_b+j];}}}从上面代码可以看出,合并后的循环便捷border=row_a*col_b;即等于原来的两个循环边界之积,然后再循环中计算出原来的外层循环和第2层循环的迭代变量i和j,采用除法和取余来求出i和j的值。
需要值得注意的是,上面求i和j的值必须要保证循环迭代的独立性,即不能有循环迭代间的依赖关系。
不能讲求i和j的值得过程优化成如下的形式if(j==col_b){j=0;i++;}//.......此处代表实际的矩阵乘法代码j++;上面这种优化,省去了除法,效率高,但是只能在串行代码中使用,因为它存在循环迭代间的依赖关系,无法将其正确地并行。