按行分配
思路和MPI基本类似，不过OpenMP是共享内存的，不必做分发和聚集，申请的矩阵空间就不必是完全连续的。
1 #include<stdio.h> 2 #include<p.h> 3 #include<stdlib.h> 4 5 #define N //规模(方针的阶数) 6 int i,j;//通用游标 7 double **mat=NULL;//矩阵对象 8 int *vec=NULL;//列向量对象 9 double *result=NULL;//结果向量对象 10 11 int main(int argc,char* argv[]) 12 36 37 //并行for:为结果向量赋初始值,这样能避免calloc时间或多余的if判断 38 # pragma p parallel for num_threads(thrdCnt) 39 for(i=0;i<N;i++) 40 result[i]=mat[i][0]*vec[0]; 41 42 //并行for:计算结果 43 # pragma p parallel for num_threads(thrdCnt) private(j) 44 for(i=0;i<N;i++) 45 for(j=1;j<N;j++) 46 result[i]+=mat[i][j]*vec[j]; 47 48 //采样输出结果看一下 49 for(i=0;i<N;i+=N/11) 50 printf("%f\n",result[i]); 51 52 return 0; 53 }输出
1 [lzh@hostlzh OpenMP]$ gcc fopenmp o p_roo p_roc 2 [lzh@hostlzh OpenMP]$ ./p_roo 7 3 79800.000000 4 79800.000000 5 79800.000000 6 79800.000000 7 79800.000000 8 79800.000000 9 79800.000000 10 79800.000000 11 79800.000000 12 79800.000000 13 79800.000000 14 79800.000000 15 [lzh@hostlzh OpenMP]$按列分配
按列分配有很多细节要注意，如果不对result数组保护则可能会发生冲突，如果用critical或者atic会导致计算过程实际是串行的，虽然atic的加解锁过程是critical的7倍。
1 #include<stdio.h> 2 #include<p.h> 3 #include<stdlib.h> 4 5 #define N //规模(方针的阶数) 6 int i,j;//通用游标 7 double **mat=NULL;//矩阵对象 8 int *vec=NULL;//列向量对象 9 double *result=NULL;//结果向量对象 10 11 int main(int argc,char* argv[]) 12 36 37 //并行for:为结果向量赋初始值,这样能避免calloc时间或多余的if判断 38 # pragma p parallel for num_threads(thrdCnt) 39 for(i=0;i<N;i++) 40 result[i]=mat[i][0]*vec[0]; 41 42 //创建存放每个线程临时结果的数组,初始化为0 43 //这里calloc还能像前面那样改进,但为了程序可读性暂时不做改进 44 double* sum; 45 46 //OpenMP默认块划分给每个进程(除了最后一个进程)的循环次数 47 //是总的循环次数处以进程数向上取整,需要先计算出这个数 48 int needPlus;//记录是否向上+1 49 needPlus=((N1)%thrdCnt==0)?0:1; 50 //每个进程(除了最后一个进程)的循环次数 51 int numThrdFor=(N1)/thrdCnt+needPlus; 52 53 //并行for:计算结果,按列分配时这个大的外层for用j,且从1开始 54 # pragma p parallel for num_threads(thrdCnt) \ 55 private(i) firstprivate(sum)//对i只要每个线程私有,对sum数组还需初值 56 for(j=1;j<N;j++) 57 65 66 //对于这其中的每一短行 67 //(这个i被private保护,为每个线程单独创建了私有的i) 68 for(i=0;i<N;i++) 69 77 78 //仅当本线程即将结束前,把算好的sum数组加上来 79 //判断条件:从1开始计数下能整除次数,或当前是最后一次循环 80 if(j%numThrdFor==0 || j==N1) 81 86 } 87 88 //采样输出结果看一下 89 for(i=0;i<N;i+=N/11) 90 printf("%f\n",result[i]); 91 92 return 0; 93 }输出
1 [lzh@hostlzh OpenMP]$ gcc fopenmp o p_col.o p_col 2 [lzh@hostlzh OpenMP]$ ./p_col.o 12 3 79800.000000 4 79800.000000 5 79800.000000 6 79800.000000 7 79800.000000 8 79800.000000 9 79800.000000 10 79800.000000 11 79800.000000 12 79800.000000 13 79800.000000 14 79800.000000 15 [lzh@hostlzh OpenMP]$按块分配
每次循环都会重新开个线程，虽然操作了共享变量但是测试时没出现问题(为啥)。
1 #include<stdio.h> 2 #include<p.h> 3 #include<stdlib.h> 4 5 #define N //规模(方针的阶数) 6 int i,j;//通用游标 7 double **mat=NULL;//矩阵对象 8 int *vec=NULL;//列向量对象 9 double *result=NULL;//结果向量对象 10 int blkSqrt=1;//块数的开算数平方 11 12 //自己的求平方根的函数,因为math里的sqrt报错 13 int mysqrt(int k) 14 22 23 int main(int argc,char* argv[]) 24 33 //矩阵一级挂载空间 34 mat=(double**)malloc(N*sizeof(double *)); 35 //存向量的空间 36 vec=(int *)malloc(N*sizeof(int)); 37 //存结果的空间 38 result=(double *)malloc(N*sizeof(double)); 39 40 //并行for:申请存矩阵的空间 41 //因为OpenMP不需要分发聚集等,所以不必做连续空间申请 42 # pragma p parallel for num_threads(thrdCnt) 43 for(i=0;i<N;i++) 44 mat[i]=(double*)malloc(N*sizeof(double)); 45 46 //并行for:为矩阵和列向量赋值(实际做时是从文件读入) 47 # pragma p parallel for num_threads(thrdCnt) private(j) 48 for(i=0;i<N;i++) 49 54 55 //并行for:为结果向量赋初始值,这样能避免calloc时间或多余的if判断 56 # pragma p parallel for num_threads(thrdCnt) 57 for(i=0;i<N;i++) 58 result[i]=mat[i][0]*vec[0]; 59 60 //并行for:按块分配计算结果,即行列分别分开 61 # pragma p parallel for num_threads(blkSqrt) 62 for(i=0;i<N;i++) # pragma p parallel for num_threads(blkSqrt) 64 for(j=1;j<N;j++) 65 result[i]+=mat[i][j]*vec[j]; 66 67 //采样输出结果看一下 68 for(i=0;i<N;i+=N/11) 69 printf("%f\n",result[i]); 70 71 return 0; 72 }输出
1 [lzh@hostlzh OpenMP]$ gcc fopenmp o p_blk.o p_blk 2 [lzh@hostlzh OpenMP]$ ./p_blk.o 9 3 79800.000000 4 79800.000000 5 79800.000000 6 79800.000000 7 79800.000000 8 79800.000000 9 79800.000000 10 79800.000000 11 79800.000000 12 79800.000000 13 79800.000000 14 79800.000000 15 [lzh@hostlzh OpenMP]$上一篇：Qt 5.14.2配置QVTKOpenGLNativeWidget控件并在QT中显示VTK窗口
下一篇：共聚焦显微镜针孔效应
OpenMP
