我可以并行化我的程序吗？

Question

这是我的程序：

program test
implicit none
integer n,m,k,i,j,Errorflag
real :: Yabs(39,39),angle(39,39)
real ,dimension(67,1) :: deltaA,A
real :: V(1,39),d(1,39),v1(29,1),d1(38,1),Ps(1,38),Qs(1,39),Jac(67,67),invJac(67,67)
real :: B1(1,38),B2(1,29),MF(1,67),trnsMF(67,1),P0(1,39),Q0(1,39)
real, dimension(38,38) :: dia1,offdia1,J1
real, dimension(29,29) :: dia2,dia3,dia4,offdia4,J4
real,dimension(38,29) ::offdia2,J2
real,dimension(29,38) ::offdia3,J3
real p,p1,q,q1

n=39;m=9
MF(1,1)=10

open(unit=3,file="ybus.dat",status="old")
open(unit=4,file="angle.dat",status="old")
 do i=1,39
        read(3,*) Yabs(i,1:39)
        read(4,*)angle(i,1:39)
    end do
close(3)
close(4)

open(unit=5,file="activepower.dat",status="old")
open(unit=8,file="reactivepower.dat",status="old")


read(5,*)Ps(1,1:38)
read(8,*)Qs(1,1:29)


close(5)
close(8)

do i=1,67
deltaA(i,1)=0
end do
v1(1:29,1)=1
d1(1:38,1)=0 
A(1:38,1)=d1(1:38,1)
A(39:67,1)=v1(1:29,1)
!call cpu_time(t1)

do while(maxval(abs(MF))>0.0001)

V(1,1)=0.982
V(1,2:30)=v1(1:29,1)
V(1,31)=1.03
V(1,32)=0.9831           
V(1,33)=1.0123
V(1,34)=0.9972
V(1,35)=1.0493
V(1,36)=1.0635 
V(1,37)=1.0278
V(1,38)=1.0265
V(1,39)=1.0475
d(1,1)=0
d(1,2:39)=d1(1:38,1)
 ! % % % %------Active Power Calculation-----%
p1=0;p=0
do i=2,n
    do j=1,n
        p1=(V(i)*V(j)*Yabs(i,j)*cos(angle(i,j)-d(i)+d(j)))
        p=p1+p
     end do
     P0(i-1)=p
     p=0
end do
! % % % %------Reactive Power Calculation-----%  
p=0;p1=0
do i=2,(n-m)
      do j=1,n
p1=-(V(i)*V(j)*Yabs(i,j)*sin(angle(i,j)-d(i)+d(j)))
p=p1+p
      end do
      Q0(i-1)=p
      p=0
end do
!!!!!!!!!!!mismatch factor
do i=1,(n-1)
   B1(i)=Ps(i)-P0(i)
  end do
do i=1,(n-m-1)
   B2(i)=Qs(i)-Q0(i)
end do

MF(1,1:38)=B1(1,1:38)
MF(1,39:67)=B2(1,1:29)
!!!!!!!!jacobian calculation for preddictor step
!!!!!!!!!!!!!!!!!!!!!!dia of j1
p=0;p1=0
do i=2,n
    do j=1,n
       if(j .ne. i)then
 p1=V(i)*V(j)*Yabs(i,j)*sin(angle(i,j)-d(i)+d(j))
 !print*,p1
 p=p1+p
       end if

    end do
    i=i-1
    dia1(i,i)=p

    p=0
    i=i+1
end do

!!!!!!!!!!!!!!off dia. of j1
q=0;q1=0;
do k=2,n
 i=k
    do j=2,n
        if(j .ne. i)then
        q1=V(i)*V(j)*Yabs(i,j)*sin(angle(i,j)-d(i)+d(j))
        end if
        i=i-1;j=j-1
offdia1(i,j)=-q1
q1=0
i=i+1;j=j+1
    end do
end do
do i=1,38
do j=1,38
J1(i,j)=offdia1(i,j)+dia1(i,j)
end do
end do
!!!!!!!!!!!!!!!!!!!dia. of j2
p=0;p1=0
do i=2,(n-m)
    do j=1,n
       if(j .ne. i)then
         p1=V(j)*Yabs(i,j)*cos(angle(i,j)-d(i)+d(j))
 p=p1+p
       end if

    end do

    dia2(i-1,i-1)=p+(2*V(i)*Yabs(i,i)*cos(angle(i,i)))
    p=0;

end do
!!!!!!!!!!!!!!!!!!off dia. of j2
p1=0;
do k=2,n
 i=k
    do j=2,(n-m)
        if(j .ne. i)then
       p1=V(i)*Yabs(i,j)*cos(angle(i,j)-d(i)+d(j));
        end if
        i=i-1;j=j-1
offdia2(i,j)=p1
    p1=0;
    i=i+1;j=j+1
    end do
end do
do i=1,(n-m-1)

   offdia2(i,i)=dia2(i,i)

end do
J2=offdia2
!!!!!!!!!!!!!!!!!!!!dia. of j3
p=0;p1=0
do i=2,(n-m)
    do j=1,n
       if(j .ne. i)then
          p1=V(i)*V(j)*Yabs(i,j)*cos(angle(i,j)-d(i)+d(j))
 p=p1+p;
       end if
    end do
    i=i-1;
    dia3(i,i)=p
    p=0;
    i=i+1;
 end do
 !!!!!!!!!!!!!!off dia of j3
p=0;p1=0
do k=2,(n-m)
 i=k;
    do j=2,n
        if(j .ne. i)then
    p1=V(i)*V(j)*Yabs(i,j)*cos(angle(i,j)-d(i)+d(j))
        end if
     i=i-1;j=j-1
    offdia3(i,j)=-p1;
    p1=0;
    i=i+1;j=j+1
    end do
end do
do i=1,(n-m-1)
      offdia3(i,i)=dia3(i,i)

end do
J3=offdia3
!!!!!!!!!!dia of j4
p=0;p1=0
do i=2,(n-m)
    do j=1,n
       if(j .ne. i)then
           p1=V(j)*Yabs(i,j)*sin(angle(i,j)-d(i)+d(j))
 p=p1+p
       end if

    end do

    dia4(i-1,i-1)=-(2*V(i)*Yabs(i,i)*sin(angle(i,i)))-p
    p=0;p1=0
end do
 !!!!!!!!!!!!!!!off dia of j4
p1=0;p=0
do k=2,(n-m)
 i=k;
    do j=2,(n-m)
        if(j .ne. i)then
       p1=V(i)*Yabs(i,j)*sin(angle(i,j)-d(i)+d(j))
        end if
    i=i-1;j=j-1
    offdia4(i,j)=-p1
    p1=0;
    i=i+1;j=j+1
    end do
 end do
 do i=1,(n-m-1)
offdia4(i,i)=dia4(i,i);
 end do
J4=offdia4
!!!!!!!
!!!!!!!!!!!!!!!!!!!formation of final jacobian!!!!!!!!!!
Jac( 1:38, 1:38) = J1 (1:38,1:38)
Jac( 1:38,39:67) = J2 (1:38,1:29)
Jac(39:67, 1:38) = J3 (1:29,1:38)
Jac(39:67,39:67) = J4 (1:29,1:29)
!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
!print*,Jac(23,21)
CALL FindInv(Jac,invJac ,67, ErrorFlag)
trnsMF=transpose(MF)

deltaA=matmul( invJac, trnsMF)
do i=1,67
A(i)=A(i)+deltaA(i)
end do

!!!!!!!!!!!!updating values
do i=1,(n-1)
    d1(i)=A(i)
end do

k=0
do i=n,(2*n-2-m)
    k=1+k
    v1(k)=A(i)
end do
end do

end program test

数组“Ps”包含一些值。现在，如果我将 Ps(15) 的值增加 Ps(15)+1，那么对于这两个值，我可以并行化此代码以快速获得答案。

我正在为 CUDA FORTRAN 使用 PGI 编译器。

Answer 1

您的代码相当简单，包含许多独立的并行循环。这些并行循环似乎被包裹在一个外部收敛do while 循环中，因此只要您将数据保存在设备上以进行收敛循环的所有迭代，就不会受到传输的瓶颈。

我建议从该代码的编译器指令开始，而不是深入研究 CUDA Fortran。编译器指令适用于像这样的简单独立循环——它们是您放置在代码 cmets 中的简单提示，告诉编译器要并行化哪些循环、要复制哪些数据等。

您可以先尝试OpenMP 加速到多个 CPU 内核。然后，您可以使用 GPU 指令，例如 OpenACC，这将很快在 PGI、Cray 和 CAPS 的编译器中可用。为了抢占先机，您可以下载PGI compiler 的免费试用版并使用他们的"Accelerator" 指令。 Accelerator 在语法上与 OpenACC 非常相似。

Answer 2

是的，您可以使用PGI compiler 编写 CUDA 内核并进行 CUDA API 调用。

PGI Fortran CUDA Homepage

我想，你的意思是“应该我并行化这段代码吗？”

我的回答是，是的，您一眼就能看到并行化的一些轻微好处。

例如，像这样的片段：

 do i=2,n
    do j=1,n
       if(j .ne. i)then
 p1=V(i)*V(j)*Yabs(i,j)*sin(angle(i,j)-d(i)+d(j))
 !print*,p1
 p=p1+p
       end if

    end do
    i=i-1
    dia1(i,i)=p

    p=0
    i=i+1
end do

是 N^2 组独立计算（在这种情况下，您设置 n=39，但我认为它可能会改变）。因此，您至少要处理几百个计算。虽然理想情况下您希望在并行化方面进行更多计算，但至少您的许多循环似乎在每一步都在做相同的独立工作 - 非常适合线程应用程序。

因此，您可以看到编写 CUDA 内核来替换数据后处理算法中的循环代码段的一些轻微好处。请注意，PCI 总线在内存传输方面的延迟确实会抵消一些性能提升，尤其是对于小型系统。

因此我会说，是的，如果你是游戏玩家，无论如何你都可以并且应该尝试这个，但不要指望它会快 100 倍……如果你编码得好，可能会快 2-10 倍，具体取决于您的循环边界大小和特定循环内的发散程度。

最坏的情况是你看不到任何收益，甚至看到放缓，但至少你学到了一些东西！！