rendered paste body USE MPI
include 'param-3'
dimension XVR3x (NVRM), XVR3y (NVRM), XVR3z (NVRM),&
VR3x (NVRM), VR3y (NVRM), VR3z (NVRM),&
UVVpx (NVRM), UVVpy (NVRM), UVVpz (NVRM),&
GVVpx (NVRM), GVVpy (NVRM), GVVpz (NVRM)
dimension XO(3),VRO(3),UVR(3),GVR(3)
real*4 TARRAY(2)
integer p,source,dest,ierr,my_rank,tag,NVR
integer status (MPI_Status_Size)
call MPI_Init(ierr)
call MPI_Comm_Rank(MPI_COMM_WORLD,my_rank,ierr)
call MPI_Comm_size(MPI_COMM_WORLD,p,ierr)
if (my_rank==0) then
open (15,file='VORTP450ok')
rewind(15)
nv=0
do i=1,50000!NVRM
nv=nv+1
read(15,*,end=30) XVR3x(nv), XVR3y(nv), XVR3z(nv),&
VR3x(nv) , VR3y(nv) , VR3z(nv)
UVVpx(i)=0.
UVVpy(i)=0.
UVVpz(i)=0.
GVVpx(i)=0.
GVVpy(i)=0.
GVVpz(i)=0.
enddo
write (*,*)XVR3y(1)
30 close(15)
NVR=nv-1
write (*,*) 'NVR=',NVR
do dest=1,p-1
tag=10
call MPI_SEND(NVR,1,MPI_INTEGER,dest,tag,MPI_COMM_WORLD,ierr)
tag=20
call MPI_SEND(XVR3x,NVR,MPI_REAL ,dest,tag,MPI_COMM_WORLD,ierr)
tag=30
call MPI_SEND(XVR3y,NVR,MPI_REAL ,dest,tag,MPI_COMM_WORLD,ierr)
tag=40
call MPI_SEND(XVR3z,NVR,MPI_REAL ,dest,tag,MPI_COMM_WORLD,ierr)
tag=50
call MPI_SEND(VR3x,NVR ,MPI_REAL ,dest,tag,MPI_COMM_WORLD,ierr)
tag=60
call MPI_SEND(VR3y,NVR ,MPI_REAL ,dest,tag,MPI_COMM_WORLD,ierr)
tag=70
call MPI_SEND(VR3z,NVR ,MPI_REAL ,dest,tag,MPI_COMM_WORLD,ierr)
enddo
! call MPI_BCAST(NVR ,1,MPI_INTEGER,0,MPI_COMM_WORLD,ierr)
! call MPI_BCAST(XVR3x(1),NVR,MPI_REAL,0,MPI_COMM_WORLD,ierr)
! call MPI_BCAST(XVR3y(1),NVR,MPI_REAL,0,MPI_COMM_WORLD,ierr)
! call MPI_BCAST(XVR3z(1),NVR,MPI_REAL,0,MPI_COMM_WORLD,ierr)
! call MPI_BCAST(VR3x (1),NVR,MPI_REAL,0,MPI_COMM_WORLD,ierr)
! call MPI_BCAST(VR3y (1),NVR,MPI_REAL,0,MPI_COMM_WORLD,ierr)
! call MPI_BCAST(VR3z (1),NVR,MPI_REAL,0,MPI_COMM_WORLD,ierr)
else
tag=10
call MPI_RECV(NVR,1,MPI_INTEGER,0,tag,MPI_COMM_WORLD,status,ierr)
tag=20
call MPI_RECV(XVR3x,NVR,MPI_REAL,0 ,tag,MPI_COMM_WORLD,status,ierr)
tag=30
call MPI_RECV(XVR3y,NVR,MPI_REAL,0 ,tag,MPI_COMM_WORLD,status,ierr)
tag=40
call MPI_RECV(XVR3z,NVR,MPI_REAL,0 ,tag,MPI_COMM_WORLD,status,ierr)
tag=50
call MPI_RECV( VR3x,NVR,MPI_REAL,0 ,tag,MPI_COMM_WORLD,status,ierr)
tag=60
call MPI_RECV( VR3y,NVR,MPI_REAL,0 ,tag,MPI_COMM_WORLD,status,ierr)
tag=70
call MPI_RECV( VR3z,NVR,MPI_REAL,0 ,tag,MPI_COMM_WORLD,status,ierr)
endif
I0=1+NVR/p * my_rank
I1=NVR/p * (my_rank+1)+(NVR -(NVR/p *p))*int(my_rank/max(1,(p-1)))
! TARRAY(1)=0.
! TARRAY(2)=0.
! dt1=dtime(TARRAY)
open(3,file='a')
write(3,*) "aaa",I0,I1,I1-I0+1
close (3)
do i = I0,I1
XO(1) = XVR3x(i)
XO(2) = XVR3y(i)
XO(3) = XVR3z(i)
VRO(1) = VR3x(i)
VRO(2) = VR3y(i)
VRO(3) = VR3z(i)
call DVORT (XVR3x,XVR3y,XVR3z,VR3x,VR3y,VR3z,&
UVVpx,UVVpy,UVVpz,GVVpx,GVVpy,GVVpz,XO,&
VRO, UVR , GVR ,NVR)
UVVpx(i)=UVR(1)
UVVpy(i)=UVR(2)
UVVpz(i)=UVR(3)
GVVpx(i)=GVR(1)
GVVpy(i)=GVR(2)
GVVpz(i)=GVR(3)
! write(*,*) UVVpx(i),GVVpx(i)
enddo
if (my_rank/=0) then
tag=1!*my_rank**5
dest=0
call MPI_SEND(UVVpx(I0),I1-I0+1,MPI_REAL,dest,tag,MPI_COMM_WORLD,ierr)
tag=2!*my_rank**5
call MPI_SEND(UVVpy(I0),I1-I0+1,MPI_REAL,dest,tag,MPI_COMM_WORLD,ierr)
tag=3!*my_rank**5
call MPI_SEND(UVVpz(I0),I1-I0+1,MPI_REAL,dest,tag,MPI_COMM_WORLD,ierr)
tag=4!*my_rank**5
call MPI_SEND(GVVpx(I0),I1-I0+1,MPI_REAL,dest,tag,MPI_COMM_WORLD,ierr)
tag=5!*my_rank**5
call MPI_SEND(GVVpy(I0),I1-I0+1,MPI_REAL,dest,tag,MPI_COMM_WORLD,ierr)
tag=6!*my_rank**5
call MPI_SEND(GVVpz(I0),I1-I0+1,MPI_REAL,dest,tag,MPI_COMM_WORLD,ierr)
endif
! dt1=dtime(TARRAY)
! write(*,*) 'TIME for BS:',dt1,TARRAY(1),TARRAY(2)
if(my_rank==0) then
do source=1,p-1
I0=1+NVR/p* source
I1=NVR/p * (source+1)+(NVR -(NVR/p *p))*int(source/(max(1,p-1)))
write (*,*) I0,I1
tag=1!*source**5
call MPI_RECV(UVVpx(I0),I1-I0+1,MPI_REAL,source,tag,MPI_COMM_WORLD,status,ierr)
tag=2!*source**5
call MPI_RECV(UVVpy(I0),I1-I0+1,MPI_REAL,source,tag,MPI_COMM_WORLD,status,ierr)
tag=3!*source**5
call MPI_RECV(UVVpz(I0),I1-I0+1,MPI_REAL,source,tag,MPI_COMM_WORLD,status,ierr)
tag=4!*source**5
call MPI_RECV(GVVpx(I0),I1-I0+1,MPI_REAL,source,tag,MPI_COMM_WORLD,status,ierr)
tag=5!*source**5
call MPI_RECV(GVVpy(I0),I1-I0+1,MPI_REAL,source,tag,MPI_COMM_WORLD,status,ierr)
tag=6!*source**5
call MPI_RECV(GVVpz(I0),I1-I0+1,MPI_REAL,source,tag,MPI_COMM_WORLD,status,ierr)
enddo
open(3,file='res.dat')
do i=1,nvr
write (3,*) UVVPx(i)
enddo
close(3)
endif
call MPI_Finalize(ierr)
end
!-----------------------------------------------------------------------
SUBROUTINE DVORT (XVR3x,XVR3y,XVR3z,VR3x,VR3y,VR3z,&
UVVpx,UVVpy,UVVpz,GVVpx,GVVpy,GVVpz,XO, VRO,&
UVV, GVV,NVR)
!-----------------------------------------------------------------------
include 'param-3.f'
dimension XVR3x (NVRM), XVR3y (NVRM), XVR3z (NVRM),&
VR3x (NVRM), VR3y (NVRM), VR3z (NVRM),&
UVVpx (NVRM), UVVpy (NVRM), UVVpz (NVRM),&
GVVpx (NVRM), GVVpy (NVRM), GVVpz (NVRM)
dimension XO(3),VRO(3),UVV(3),GVV(3),XS(3),VS(3)
do 100 K = 1,3
UVV (K) = 0.
GVV (K) = 0.
100 continue
if (NVR.eq.0) return
do I = 1,NVR
!jim do K = 1, 3
! XS (K) = XVR (K,I)
! VS (K) = VR (K,I)
! enddo
XS (1) = XVR3x (I)
XS (2) = XVR3y (I)
XS (3) = XVR3z (I)
VS (1) = VR3x (I)
VS (2) = VR3y (I)
VS (3) = VR3z (I)
EPSI = 0.05
REFLEN = 1000.
call DELEV (XO,VRO,XS,VS,GVV,UVV,EPSI,REFLEN)
enddo
! write(*,*) UVV(1)
return
end
!-----------------------------------------------------------------------
SUBROUTINE DELEV ( XO, VO, XJ, VJ, G, U, EPS, REFLEN )
!-----------------------------------------------------------------------
! transpose
parameter (CU = 0.079577471)
dimension XO(3),XJ(3),VO(3),VJ(3),G(3),R(3),U(3),GG(3),UU(3)
! CU = 1./(PI4)
R(1) = XO(1)-XJ(1)
R(2) = XO(2)-XJ(2)
R(3) = XO(3)-XJ(3)
DR = R(1)*R(1)+R(2)*R(2)+R(3)*R(3)
SDR = sqrt (DR)
if (SDR.gt.REFLEN) return
SDR3 = DR*SDR
if (SDR.eq.0.) return
R(1) = R(1)/SDR
R(2) = R(2)/SDR
R(3) = R(3)/SDR
UU(1) = -R(2)*VJ(3)+R(3)*VJ(2)
UU(2) = -R(3)*VJ(1)+R(1)*VJ(3)
UU(3) = -R(1)*VJ(2)+R(2)*VJ(1)
C2 = UU(1)*VO(1)+UU(2)*VO(2)+UU(3)*VO(3)
GG(1) = VO(2)*VJ(3)-VJ(2)*VO(3)
GG(2) = VO(3)*VJ(1)-VJ(3)*VO(1)
GG(3) = VO(1)*VJ(2)-VJ(1)*VO(2)
UU(1) = UU(1)*CU
UU(2) = UU(2)*CU
UU(3) = UU(3)*CU
if (SDR.lt.(3.*EPS)) go to 5
FUNG = 1./SDR3
FUNU = 1./DR
go to 6
5 if (SDR.lt.0.00001) return
FUN = -SDR3/(EPS*EPS*EPS)
FUNE = EXP(FUN)
FUNU = (1.-FUNE)/DR
FUNG = (1.-FUNE+FUN*FUNE)/SDR3
6 U(1) = UU(1)*FUNU+U(1)
U(2) = UU(2)*FUNU+U(2)
U(3) = UU(3)*FUNU+U(3)
G(1) = G(1)+CU*(FUNU*GG(1)/SDR-3.*FUNG*R(1)*C2)
G(2) = G(2)+CU*(FUNU*GG(2)/SDR-3.*FUNG*R(2)*C2)
G(3) = G(3)+CU*(FUNU*GG(3)/SDR-3.*FUNG*R(3)*C2)
return
END