Fortran结合MPI并行
该Bolg整理在Fortran用MPI编写并行程序，用在科学计算中速度提升肉眼可见，
{:.info}
代码解析
先给一个完整的代码，是用来计算自旋极化率的，实际上在Fortran结合MPI并行计算自旋极化率中已经出现过了，这里就是解读一下其中的并行部分是如何用MPI写的。
1 | module param |
MPI并行解读
首先是要在程序中调用MPI的参数库1
use mpi
我在代码中的并行思想就是将一个$[-kn,kn-1)$，长度为$2*kn$的循环进行分拆。首先是MPI环境的初始化1
2
3
4
5
6call MPI_INIT(ierr) ! 初始化进程
call MPI_COMM_RANK(MPI_COMM_WORLD, indcore, ierr) ! 得到本进程在通信空间中的rank值,即在组中的逻辑编号(该 rank值为0到p-1间的整数,相当于进程的ID。)
call MPI_COMM_SIZE(MPI_COMM_WORLD, numcore, ierr) !获得进程个数 size, 这里用变量p保存
call MPI_Barrier(MPI_COMM_WORLD,ierr)
nki = floor(indcore * (2.0 * kn)/numcore) - kn
nkf = floor((indcore + 1) * (2.0 * kn)/numcore) - kn - 1
{:.success}
并行会同时调用多个CPU进行计算，在调用MPI的时候，每个CPU都会有自己单独的编号，也就是上面的indcore这个参数来进行标识，而numcore就是在程序执行的时候调用CPU的总数量。此时可以看到对于每个CPU而言，nki和nkf都是不同的，这样就可以在不同的CPU上面执行$[-kn,kn-1)$这个区间的不同子区间。
这里为了保证并行的效率，需要做到不同的区间之间是没有交叠的，因为后续我们要对不同CPU中计算结果的收集，用的是一个求和的方式，所以如果此时计算区间是有重叠的，那么计算的结果中就会有一部分是重复叠加的，所以一定要保证并行区间是没有重叠，且它们的并集等于$[-kn,kn-1)$.
{:.warning}
下面就是分拆之后的循环在不同的CPU上计算不同的区间1
2
3
4
5
6
7
8
9
10
11
12
13
14do iky = nki,nkf
qy = pi * iky/kn
do ikx = -kn,kn - 1
qx = pi * ikx/kn
call chi0cal(qx,qy,chi0(iky,ikx,:,:)) ! 得到裸极化率
call inv(ones - matmul(chi0(iky,ikx,:,:),Umat),temp2) ! 矩阵求逆
chi(iky,ikx,:,:) = matmul(temp2,chi0(iky,ikx,:,:))
temp3 = sum(chi(iky,ikx,:,:))
local_rechi(iky,ikx,1) = qx
local_rechi(iky,ikx,2) = qy
local_rechi(iky,ikx,3) = real(temp3)
local_rechi(iky,ikx,4) = aimag(temp3)
end do
end do
每个CPU上都会有自己独立的local_rechi变量，下面就是要等所有的CPU计算结束，将结果收集起来1
2call MPI_Barrier(MPI_COMM_WORLD,ierr)
call MPI_Reduce(local_rechi, rechi, (2 * kn)**2 * 4, MPI_REAL, MPI_SUM, 0, MPI_COMM_WORLD,ierr)
这里调用call MPI_Barrier(MPI_COMM_WORLD,ierr)的目的就是等所有的CPU都计算结束，如果有一些CPU已经计算完循环了，而有一些并没有计算结束，这个时候就让程序在这里等候所有的CPU计算都接受，然后调用call MPI_Reduce(local_rechi, rechi, (2 * kn)**2 * 4, MPI_REAL, MPI_SUM, 0, MPI_COMM_WORLD,ierr)来收集数据。其中其中的rechi是一个维度和local_rechi都完全相同的变量，第三个参数(2 * kn)**2 * 4表示要收集的数据量大小。举个例子，你有一个矩阵$A(10,10)$，它一共有100个元素，那么此时需要收集的数据量就是100。第四个参数MPI_REAL则用来声明收集数据的类型，还有整型MPI_INTEGER和复数型MPI_COMPLEX。第五个参数MPI_SUM就是数据收集的方式。我这里就需要用求和即可，因为程序在设计的时候不同的CPU计算的结果实际上是存储在矩阵不同位置，没有存储的位置自然就是0，所以用求和就能达到收集数据的目的。第六个参数则是申明在哪个CPU核心上来执行数据收集这个操作，我这里就选择了root核心。剩下的就是一些公共参数了，没什么好解释的了。
计算结束之后总是要进行数据存储的，这个操作需要在一个特定的核心上进行1
2
3
4
5
6
7
8
9
10
11
12
13
14if (indcore .eq. 0) then
char1 = "fortran-chi-"
write(char2,"(I3.3)")2 * kn
filename = trim(char1)//trim(char2)
char1 = ".dat"
filename = trim(filename)//trim(char1)
open(12,file = filename)
do iky = -kn,kn - 1
do ikx = -kn,kn - 1
write(12,"(4F8.3)")rechi(iky,ikx,1),rechi(iky,ikx,2),rechi(iky,ikx,3),rechi(iky,ikx,4)
end do
end do
close(12)
end if
程序执行结束之后还申明结束MPI并行1
call MPI_Finalize(ierr)
上面实际上只完成了对一个循环的并行，实际情况可以更加丰富。比如我需要用到前面并行得到的结果进行后续计算，而且后续的计算也需要并行，此时就需要对计算得到的变量进行数据广播，然后再进行并行计算，这个方法在后续会继续更新，目前的代码还没有涉及到这个问题。
公众号
相关内容均会在公众号进行同步，若对该Blog感兴趣，欢迎关注微信公众号。
{:.info}
|
yxliphy@gmail.com |
wechat
alipay
