[Fortran] 纯文本查看 复制代码
PROGRAM T5_2D
IMPLICIT REAL*8 (A-H,O-Z)
INCLUDE 'mpif.h'
PARAMETER (KK=20,NN=120,MM=160, KM=3,NN1=NN-1,MM1=MM-1)
PARAMETER (JP=2, IP=4, N=NN/JP, M=MM/IP, NP=JP*IP)
DIMENSION TT(KM,NN,MM)
DIMENSION U1(KK,0:N+1,0:M+1),V1(KK,0:N+1,0:M+1),PS1(0:N+1,0:M+1)
COMMON/BLK4/F1(KM,0:N+1,0:M+1),F2(KM,0:N+1,0:M+1),HXU(0:N+1,0:M+1),HXV(0:N+1,0:M+1),HMMX(0:N+1,0:M+1),HMMY(0:N+1,0:M+1)
COMMON/BLK5/VECINV(KK,KK),AM7(KK)
DIMENSION D7(0:N+1,0:M+1),D8(0:N+1,0:M+1),D00(KK,0:N+1,0:M+1)
PARAMETER (NDIM=2)
INTEGER ISTATUS(MPI_STATUS_SIZE), COMM2D,L_NBR, R_NBR, B_NBR, T_NBR, MY_CID, MY_COORD(NDIM),VECT_2D, VECT_3D
COMMON/BLKMPI/MYID, NPROC, ISTATUS, L_NBR, R_NBR, B_NBR, T_NBR,MY_CID, STRIDE2D, STRIDE3D,ISTART, ISTART2, IEND, IEND1, ISTARTM1, IENDP1,JSTART, JSTART2, JEND, JEND1, JSTARTM1, JENDP1, ISTARTG(0:NP),IENDG(0:NP),JSTARTG(0:NP),JENDG(0:NP)
CALL MPI_INIT (IERR)
CALL MPI_COMM_SIZE (MPI_COMM_WORLD, NPROC, IERR)
IF(NPROC.NE.NP) THEN
PRINT *,' NPROC NOT EQUAL TO ', NP,' PROGRAM WILL STOP'
STOP
ENDIF
CALL MPI_COMM_RANK (MPI_COMM_WORLD, MYID, IERR)
CALL NBR2D(COMM2D,MY_CID,MY_COORD,L_NBR,R_NBR,B_NBR,T_NBR,JP,IP)
CALL MPI_BARRIER (COMM2D,IERR)
CLOCK=MPI_WTIME()
CALL STARTEND (JP,1,NN,JSTARTG, JENDG, JCOUNTG)
CALL STARTEND (IP,1,MM,ISTARTG, IENDG, ICOUNTG)
ISTART=1
IEND=M
JSTART=1
JEND=N
! for DO I=x,MM1 (MM-1)
! for DO J=x,NN1 (NN-1)
IEND1=IEND
JEND1=JEND
IF( MY_COORD(2).EQ.IP-1) IEND1=IEND1-1
IF( MY_COORD(1).EQ.JP-1) JEND1=JEND1-1.
! for DO I=2,x
! for DO J=2,x
ISTART2=ISTART
JSTART2=JSTART
IF( MY_COORD(2).EQ.0) ISTART2=2
IF( MY_COORD(1).EQ.0) JSTART2=2
! for ghost point
ISTARTM1=ISTART-1
IENDP1=IEND+1
JSTARTM1=JSTART-1
JENDP1=JEND+1
! 原始数据的产生
!DO 10 I=1,MM1
DO 10 I=1,IEND1
II=I+ISTARTG(MY_COORD(2))-1
! DO 10 J=1,NN
DO 10 J=1,JEND
JJ=J+JSTARTG(MY_COORD(1))-1
DO 10 K=1,KK
U1(K,J,I)=1.D0/DFLOAT(II)+1.D0/DFLOAT(JJ)+1.D0/DFLOAT(K)
10 CONTINUE
! DO 20 I=1,MM
DO 20 I=1,IEND
II=I+ISTARTG(MY_COORD(2))-1
! DO 20 J=1,NN1
DO 20 J=1,JEND1
JJ=J+JSTARTG(MY_COORD(1))-1
DO 20 K=1,KK
V1(K,J,I)=2.D0/DFLOAT(II)+1.D0/DFLOAT(JJ)+1.D0/DFLOAT(K)
20 CONTINUE
! DO 30 I=1,MM
DO 30 I=1,IEND
II=I+ISTARTG(MY_COORD(2))-1
! DO 30 J=1,NN
DO 30 J=1,JEND
JJ=J+JSTARTG(MY_COORD(1))-1
PS1(J,I)=1.D0/DFLOAT(II)+1.D0/DFLOAT(JJ)
HXU(J,I)=2.D0/DFLOAT(II)+1.D0/DFLOAT(JJ)
HXV(J,I)=1.D0/DFLOAT(II)+2.D0/DFLOAT(JJ)
HMMX(J,I)=2.D0/DFLOAT(II)+1.D0/DFLOAT(JJ)
HMMY(J,I)=1.D0/DFLOAT(II)+2.D0/DFLOAT(JJ)
30 CONTINUE
DO 40 K=1,KK
AM7(K)=1.D0/DFLOAT(K)
DO 40 KA=1,KK
VECINV(KA,K)=1.D0/DFLOAT(KA)+1.D0/DFLOAT(K)
40 CONTINUE
! 开始计算
N2=N+2
N2KK=N2*KK
CALL MPI_TYPE_VECTOR (M, KK, N2KK, MPI_REAL8, IVECT_3D, IERR)
CALL MPI_TYPE_COMMIT (IVECT_3D, MPI_ERR)
CALL MPI_TYPE_VECTOR (M, 1, N2, MPI_REAL8, IVECT_2D, IERR)
CALL MPI_TYPE_COMMIT (IVECT_2D, IERR)
CALL MPI_BARRIER (COMM2D, IERR)
ITAG=10
CALL MPI_SENDRECV (U1(1,1,IEND ), N2KK, MPI_REAL8, T_NBR, ITAG,U1(1,1,ISTARTM1), N2KK, MPI_REAL8, B_NBR, ITAG,COMM2D, ISTATUS, IERR)
ITAG=20
CALL MPI_SENDRECV(V1(1,JEND,1 ), 1, IVECT_3D, R_NBR, ITAG,V1(1,JSTARTM1,1), 1, IVECT_3D, L_NBR, ITAG,COMM2D, ISTATUS, IERR)
ITAG=30
CALL MPI_SENDRECV (PS1(1,ISTART), N, MPI_REAL8, B_NBR, ITAG,PS1(1,IENDP1), N, MPI_REAL8, T_NBR, ITAG,COMM2D, ISTATUS, IERR)
ITAG=40
CALL MPI_SENDRECV (PS1(JSTART,1), 1, IVECT_2D, L_NBR, ITAG,PS1(JENDP1,1), 1, IVECT_2D, R_NBR, ITAG,COMM2D, ISTATUS, IERR)
! DO 210 I=1,MM
! DO 210 J=1,NN
DO 210 I=ISTART,IEND
DO 210 J=JSTART,JEND
DO 210 K=1,KM
F1(K,J,I)=0.0D0
F2(K,J,I)=0.0D0
210 CONTINUE
! DO 220 I=1,MM1
! DO 220 J=2,NN1
DO 220 I=ISTART,IEND1
DO 220 J=JSTART2,JEND1
D7(J,I)=(PS1(J,I+1)+PS1(J,I))*0.5D0*HXU(J,I)
220 CONTINUE
! DO 230 I=2,MM1
! DO 230 J=1,NN1
DO 230 I=ISTART2,IEND1
DO 230 J=JSTART,JEND1
D8(J,I)=(PS1(J+1,I)+PS1(J,I))*0.5D0*HXV(J,I)
230 CONTINUE
CALL MPI_BARRIER(COMM2D, MPIERROR)
ITAG=50
CALL MPI_SENDRECV (D7(1,IEND), N, MPI_REAL8, T_NBR, ITAG,D7(1,ISTARTM1), N, MPI_REAL8, B_NBR, ITAG,COMM2D, ISTATUS, IERR)
ITAG=60
CALL MPI_SENDRECV (D8(JEND,1), 1, IVECT_2D, R_NBR, ITAG,D8(JSTARTM1,1), 1, IVECT_2D, L_NBR, ITAG,COMM2D, ISTATUS, IERR)
! DO 240 I=2,MM1
! DO 240 J=2,NN1
DO 240 I=ISTART2,IEND1
DO 240 J=JSTART2,JEND1
DO 240 K=1,KK
D00(K,J,I)=(D7(J,I)*U1(K,J,I)-D7(J,I-1)*U1(K,J,I-1))*HMMX(J,I)+(D8(J,I)*V1(K,J,I)-D8(J-1,I)*V1(K,J-1,I))*HMMY(J,I)
240 CONTINUE
! DO 260 I=2,MM1
DO 260 I=ISTART2,IEND1
DO 260 KA=1,KK
! DO 260 J=2,NN1
DO 260 J=JSTART2,JEND1
DO 260 K=1,KM
F1(K,J,I)=F1(K,J,I)-VECINV(K,KA)*D00(KA,J,I)
260 CONTINUE
SUMF1=0.D0
SUMF2=0.D0
! DO 270 I=2,MM1
DO 270 I=ISTART2,IEND1
! DO 270 J=2,NN1
DO 270 J=JSTART2,JEND1
DO 270 K=1,KM
F2(K,J,I)=-AM7(K)*PS1(J,I)
SUMF1=SUMF1+F1(K,J,I)
SUMF2=SUMF2+F2(K,J,I)
270 CONTINUE
! 输出数据用来验证
CALL MPI_BARRIER (COMM2D,IERR)
IROOT=0
CALL MPI_REDUCE (SUMF1, GSUMF1, 1, MPI_REAL8, MPI_SUM, IROOT, COMM2D, IERR)
CALL MPI_REDUCE (SUMF2, GSUMF2, 1, MPI_REAL8, MPI_SUM, IROOT, COMM2D, IERR)
KOUNT=KM*(N+2)*(M+2)
ITAG=70
IF (MY_CID.NE.0) THEN
CALL MPI_SEND (F2, KOUNT, MPI_REAL8, IROOT, ITAG, COMM2D, IERR)
ELSE
CALL COPY1(MY_CID, F2, TT,ISTARTG,JSTARTG)
DO ISRC=1,NPROC-1
CALL MPI_RECV (F2, KOUNT, MPI_REAL8, ISRC, ITAG,COMM2D, ISTATUS, IERR)
CALL COPY1 (ISRC, F2, TT,ISTARTG,JSTARTG)
ENDDO
ENDIF
301 FORMAT(8E10.3)
IF(MY_CID.EQ.0) THEN
PRINT *,'SUMF1,SUMF2=', GSUMF1,GSUMF2
PRINT *,' F2(2,2,I),I=1, 160,5'
PRINT 301,(TT(2,2,I),I=1, 160,5)
ENDIF
CLOCK=MPI_WTIME() - CLOCK
PRINT *,' MY_CID, CLOCK TIME=', MY_CID,CLOCK
CALL MPI_FINALIZE (IERR)
STOP
END
SUBROUTINE NBR2D(COMM2D, MY_CID, MY_COORD, L_NBR, R_NBR,
& B_NBR, T_NBR, JP, IP)
INCLUDE 'mpif.h'
PARAMETER (NDIM=2)
INTEGER COMM2D, MY_CID, MY_COORD(NDIM), L_NBR, R_NBR,
& B_NBR, T_NBR, JP, IP
INTEGER IPART(2), SIDEWAYS, UPDOWN, RIGHT, UP
LOGICAL PERIODS(2), REORDER
!
IPART(1)=JP
IPART(2)=IP
PERIODS(1)=.FALSE.
PERIODS(2)=.FALSE.
REORDER=.TRUE.
SIDEWAYS=0
UPDOWN=1
RIGHT=1
UP=1
CALL MPI_CART_CREATE ( MPI_COMM_WORLD, NDIM,I PART, PERIODS,REORDER, COMM2d, MPI_ERR)
CALL MPI_COMM_RANK ( COMM2D, MY_CID, MPI_ERR)
CALL MPI_CART_COORDS( COMM2D, MY_CID, NDIM, MY_COORD, MPI_ERR)
CALL MPI_CART_SHIFT( COMM2D,SIDEWAYS,RIGHT,L_NBR,R_NBR,MPI_ERR)
CALL MPI_CART_SHIFT( COMM2D, UPDOWN, UP, B_NBR, T_NBR, MPI_ERR)
PRINT *,'MY_CID=',MY_CID,', COORD=(',MY_COORD(1),',',MY_COORD(2),')',', L_,R_,T_,B_NBR=', L_NBR,R_NBR,T_NBR,B_NBR
RETURN
END
SUBROUTINE COPY1(ISRC,FF,TT,ISTARTG,JSTARTG)
INCLUDE 'mpif.h'
! Copy partitioned array FF to global array TT
PARAMETER (KK=20,NN=120,MM=160, KM=3,MM1=MM-1,NN1=NN-1)
PARAMETER (JP=2, IP=4, N=NN/JP, M=MM/IP, NP=IP*JP)
REAL*8 FF(KM,0:N+1,0:M+1),TT(KM,NN,MM),ISTARTG(0:NP),JSTARTG(0:NP)
IF(ISRC.LT.IP) THEN
JJ=0
II=ISRC
ELSE
JJ=ISRC/IP
II=ISRC-JJ*IP
ENDIF
DO I=1,M
IG=ISTARTG(II)+I-1
DO J=1,N
JG=JSTARTG(JJ)+J-1
DO K=1,KM
TT(K,JG,IG)=FF(K,J,I)
ENDDO
ENDDO
ENDDO
RETURN
END
SUBROUTINE STARTEND(NPROC,IS1,IS2,GSTART,GEND,GCOUNT)
INTEGER ID,NPROC,IS1,IS2,NP,GSTART(0:31),GEND(0:31),GCOUNT(0:31)
LENG=IS2-IS1+1
IBLOCK=LENG/NPROC
IR=LENG-IBLOCK*NPROC
DO I=0,NPROC-1
IF(I.LT.IR) THEN
GSTART(I)=IS1+I*(IBLOCK+1)
GEND(I)=GSTART(I)+IBLOCK
ELSE
GSTART(I)=IS1+I*IBLOCK+IR
GEND(I)=GSTART(I)+IBLOCK-1
ENDIF
IF(LENG.LT.1) THEN
GSTART(I)=1
GEND(I)=0
ENDIF
GCOUNT(I)=GEND(I)-GSTART(I)+1
ENDDO
END