高性能计算和网格技术实验报告实验题目OpenMP和MPI编程姓名学号专业计算机系统结构指导教师助教所在学院计算机科学与工程学院论文提交日期一、实验目的本实验的目的是通过练习掌握OpenMP和MPI并行编程的知识和技巧。1、熟悉OpenMP和MPI编程环境和工具的使用;2、掌握并行程序编写的基本步骤;3、了解并行程序调试和调优的技巧。二、实验要求1、独立完成实验内容;2、了解并行算法的设计基础;3、熟悉OpenMP和MPI的编程环境以及运行环境;4、理解不同线程数,进程数对于加速比的影响。三、实验内容3.1、矩阵LU分解算法的设计:参考文档sy6.doc所使用的并行算法:在LU分解的过程中,主要的计算是利用主行i对其余各行j,(ji)作初等行变换,各行计算之间没有数据相关关系,因此可以对矩阵A按行划分来实现并行计算。考虑到在计算过程中处理器之间的负载均衡,对A采用行交叉划分:设处理器个数为p,矩阵A的阶数为n,pnm/,对矩阵A行交叉划分后,编号为i(i=0,1,…,p-1)的处理器存有A的第i,i+p,…,i+(m-1)p行。然后依次以第0,1,…,n-1行作为主行,将其广播给所有处理器,各处理器利用主行对其部分行向量做行变换,这实际上是各处理器轮流选出主行并广播。若以编号为my_rank的处理器的第i行元素作为主行,并将它广播给所有处理器,则编号大于等于my_rank的处理器利用主行元素对其第i+1,…,m-1行数据做行变换,其它处理器利用主行元素对其第i,…,m-1行数据做行变换。根据上述算法原理用代码表示如下(关键代码):for(k=0;kN;k++){for(i=0;iTHREADS_NUM;i++){thread_data_arrray[i].thread_id=i;thread_data_arrray[i].K_number=k;thread_data_arrray[i].chushu=a[k][k];//创建线程rc=pthread_create(&pid[i],NULL,work,(void*)&thread_data_arrray[i]);…}for(i=0;iTHREADS_NUM;i++){//等待线程同步rc=pthread_join(pid[i],&ret);…}}void*work(void*arg){structthread_data*my_data;my_data=(structthread_data*)arg;intmyid=my_data-thread_id;//线程IDintmyk=my_data-K_number;//外层循环计数Kfloatmychushu=my_data-chushu;//对角线的值ints,e;inti,j;s=(N-myk-1)*myid/THREADS_NUM;//确定起始循环的行数的相对位置e=(N-myk-1)*(myid+1)/THREADS_NUM;//确定终止循环的行数的相对位置for(i=s+myk+1;ie+myk+1;i++)//由于矩阵规模在缩小,找到偏移位置{a[i][myk]=a[i][myk]/mychushu;for(j=myk+1;jN;j++)a[i][j]=a[i][j]-a[i][myk]*a[myk][j];}//printMatrix(a);returnNULL;}第一部分为入口函数,其创建指定的线程数,并根据不同的线程id按行划分矩阵,将矩阵的不同部分作为参数传递给线程,在多处理器电脑上,不同的线程并行执行,实现并行计算LU分解。在LU分解的过程中,主要的计算是利用主行i对其余各行j,(j)i)做初等行变换,由于各行计算之间没有数据相关关系,因此可以对矩阵按行划分来实现并行算法。考虑到计算过程中处理器负载的均衡,对矩阵采用行交叉划分;假设处理器个数为p,矩阵的阶数为n,则每个处理器处理的行数为pnm/。由于在OpenMP和MPI中并行算法的实现不太一样,所以接下来的两小节中我将分别针对两个编程环境设计LU分解的并行实现。3.2、OpenMP编程因为OpenMP是基于线程的编程模型,所以设计了一个基于多线程的OpenMP的LU分解算法,关键代码如下:for(k=0;kN;k++){omp_set_num_threads(THREADS_NUM);#pragmaompparallelprivate(tid){tid=omp_get_thread_num();//当前线程IDintmyid=tid;printf(helloworldfromOMPthread%d\n,tid);intmyk=k;floatmychushu=A[k][k];ints,e;inti,j;s=(N-myk-1)*myid/THREADS_NUM;//确定起始循环的行数的相对位置e=(N-myk-1)*(myid+1)/THREADS_NUM;//确定终止循环的行数的相对位置for(i=s+myk+1;ie+myk+1;i++)//由于矩阵规模在缩小,找到偏移位置{A[i][myk]=A[i][myk]/mychushu;for(j=myk+1;jN;j++)A[i][j]=A[i][j]-A[i][myk]*A[myk][j];//对行进行初等行变换}}}其主要思想为:外层设置一个列循环,在每次循环中开设THREAD_NUMS个线程,每个线程处理的矩阵A的行为上述的m,一次循环过后则完成对应列的变换,这样在N此循环过后便可完成矩阵A的LU分解。即L为A[k][j]中kj的元素,其对角线上元素为1.0,其它为0,U为A[k][j]中k=j的元素,其余为0。这里如果我们使用的是一般的多线程编程,则在开启THREAD_NUMS个线程后,在下次循环开始之前,需要手动配置等待线程同步,不然可能出现错误。但由于OpenMP使用Fork-Join并行执行模型,其会在线程队执行完以后才转到主线程执行,所以不需要等待线程同步。详细的代码请参看附带源程序。3.3、MPI编程设处理器个数为p,矩阵A的阶数为n,pnm/,对矩阵A行交叉划分后,编号为i(i=0,1,…,p-1)的处理器存有A的第i,i+p,…,i+(m-1)p行。然后依次以第0,1,…,n-1行作为主行,将其广播给所有处理器,各处理器利用主行对其部分行向量做行变换,这实际上是各处理器轮流选出主行并广播。若以编号为my_rank的处理器的第i行元素作为主行,并将它广播给所有处理器,则编号大于等于my_rank的处理器利用主行元素对其第i+1,…,m-1行数据做行变换,其它处理器利用主行元素对其第i,…,m-1行数据做行变换,计算完成后,编号为0的处理器收集各处理器中的计算结果,并从经过初等行变换的矩阵A中分离出下三角矩阵L和上三角矩阵U。关键代码如下:/*0号进程采用行交叉划分将矩阵A划分为大小m*M的p块子矩阵,依次发送给1至p-1号进程*/if(my_rank==0){for(i=0;im;i++)for(j=0;jM;j++)a(i,j)=A((i*p),j);for(i=0;iM;i++)if((i%p)!=0){i1=i%p;i2=i/p+1;MPI_Send(&A(i,0),M,MPI_FLOAT,i1,i2,MPI_COMM_WORLD);}}else{for(i=0;im;i++)MPI_Recv(&a(i,0),M,MPI_FLOAT,0,i+1,MPI_COMM_WORLD,&status);}for(i=0;im;i++)for(j=0;jp;j++){/*j号进程负责广播主行元素*/if(my_rank==j){v=i*p+j;for(k=v;kM;k++)f[k]=a(i,k);MPI_Bcast(f,M,MPI_FLOAT,my_rank,MPI_COMM_WORLD);}else{v=i*p+j;MPI_Bcast(f,M,MPI_FLOAT,j,MPI_COMM_WORLD);}/*编号小于my_rank的进程(包括my_rank本身)利用主行对其第i+1,…,m-1行数据做行变换*/if(my_rank=j){for(k=i+1;km;k++){a(k,v)=a(k,v)/f[v];for(w=v+1;wM;w++)a(k,w)=a(k,w)-f[w]*a(k,v);}}/*编号大于my_rank的进程利用主行对其第i,…,m-1行数据做行变换*/if(my_rankj){for(k=i;km;k++){a(k,v)=a(k,v)/f[v];for(w=v+1;wM;w++)a(k,w)=a(k,w)-f[w]*a(k,v);}}}/*0号进程从其余各进程中接收子矩阵a,得到经过变换的矩阵A*/if(my_rank==0){for(i=0;im;i++)for(j=0;jM;j++)A(i*p,j)=a(i,j);}if(my_rank!=0){for(i=0;im;i++)MPI_Send(&a(i,0),M,MPI_FLOAT,0,i,MPI_COMM_WORLD);}else{for(i=1;ip;i++)for(j=0;jm;j++){MPI_Recv(&a(j,0),M,MPI_FLOAT,i,j,MPI_COMM_WORLD,&status);for(k=0;kM;k++)A((j*p+i),k)=a(j,k);}}3.4、程序调优:OpenMP和MPI混合编程我们知道OpenMP是基于线程的并行编程模型,一个共享存储的进程由多个线程组成,OpenMP就是基于已有线程的共享编程模型;而MPI属于消息传递的并行编程模型,这个从前两小节中可以看到,因为在LU的MPI实现中,我们对矩阵采用交叉划分,根据p(处理器号)划分行,因此可以对每个划分出来的矩阵采用多线程并行算法,也即可以采用OpenMP计算。在MPI的编号大于等于my_rank的处理器利用主行元素对其第i+1,…,m-1行数据做行变换,其它处理器利用主行元素对其第i,…,m-1行数据做行变换部分采用OpenMP计算,所以混合编程的核心代码如下:/*编号小于my_rank的进程(包括my_rank本身)利用主行对其第i+1,…,m-1行数据做行变换*/if(my_rank=j){inttid;omp_set_num_threads(THREADS_NUM);#pragmaompparallelprivate(tid){tid=omp_get_thread_num();intmyid=tid;intmyk=i+1;floatmychushu=f[v];ints,e;intc,d;s=(m-myk-1)*myid/THREADS_NUM;e=(m-myk-1)*(myid+1)/THREADS_NUM;for(c=s+myk+1;ce+myk+1;c++){a(c,v)=a(c,v)/mychushu;for(d=v+1;dM;d++){a(c,d)=a(c,d)-f[d]*a(c,v);}}}}四、程序运行效果及分析这里将通过测试在确定的线程下,对于随机生成的不同大小的矩阵,串行及OpenMP&MPI的运行时间来分析程序性能。并通过改变线程的数量,在不同情况下多次测量,测试编译参数对程序性能的影响。4.1、固定线程数时,不同矩阵大小下性能测试a.当线程数为5时:测试得到的不同矩阵大小下的加速比数据如图1所示:图1根据图1中所示测得的数据,生成折线图如下:b.当线程数为100时,测得的实验数据如图2所示:图2根据图2中所示测得实验数据生成的折线图如下:结果分析:从以上的结果可以看出在不同的线程下,OpenMP的性能都会随着矩阵规模的增大而变好,也就是在矩阵规模变大时,OpenMP的加速比变大,虽然在实验中有出现下降,但总体来说加速比在增大。同时,可以看到在矩阵很小时,加速比几乎为零,也即此时的OpenMP运行时间