fdm restart for mpi

mpi communication established
gather4 only performs inverse transform
2014-08-22 14:19:00 +09:00 · 2014-08-22 13:59:56 +09:00 · 2014-08-22 13:43:31 +09:00 · 2014-08-22 13:41:11 +09:00 · 2014-08-22 12:01:12 +09:00 · 2014-08-22 11:51:42 +09:00
4 changed files with 398 additions and 215 deletions
--- a/4
+++ b/4
@ -8,8 +8,8 @@ ifeq ($(OMPI), 1)
        FFTW_HOME = ../libs/fftw-3.2.2/lib
        OPENMPI_HOME = /home/comb/ignis/openmpi/
        MPIF90 = $(OPENMPI_HOME)/bin/mpif90
-        FCFLAGS = -i8 -r8 -c -openmp $(MPI_COMPILE_FLAGS) -I$(FFTW_HOME)/include -I$(OPENMPI_HOME)/include
-        LDFLAGS = -i8 -r8 -openmp $(MPI_LD_FLAGS) -L$(FFTW_HOME)/lib -lfftw3 -L$(OPENMPI_HOME)/lib -lfftw3 -lm  -shared-intel
+        FCFLAGS = -i8 -r8 -c $(MPI_COMPILE_FLAGS) -I$(FFTW_HOME)/include -I$(OPENMPI_HOME)/include
+        LDFLAGS = -i8 -r8 $(MPI_LD_FLAGS) -L$(FFTW_HOME)/lib -lfftw3 -L$(OPENMPI_HOME)/lib -lfftw3 -lm  -shared-intel
        FCFLAGS_F77 = -Mextend

        MODULES += m_openmpi.o
--- a/gather_4.f90
+++ b/gather_4.f90
@ -45,18 +45,15 @@ subroutine gather_4
  ! velocities
  call xFFT3d(-1,1)
  fname = 'u.'//file_ext
-  tmp4(1:nx,1:ny,1:nz) = wrk(1:nx,1:ny,1:nz,1)
-  call gather_tmp4(u_)
+  u_(1:nx,1:ny,1:nz) = wrk(1:nx,1:ny,1:nz,1)

  call xFFT3d(-1,2)
  fname = 'v.'//file_ext
-  tmp4(1:nx,1:ny,1:nz) = wrk(1:nx,1:ny,1:nz,2)
-  call gather_tmp4(v_)
+  v_(1:nx,1:ny,1:nz) = wrk(1:nx,1:ny,1:nz,2)

  call xFFT3d(-1,3)
  fname = 'w.'//file_ext
-  tmp4(1:nx,1:ny,1:nz) = wrk(1:nx,1:ny,1:nz,3)
-  call gather_tmp4(w_)
+  w_(1:nx,1:ny,1:nz) = wrk(1:nx,1:ny,1:nz,3)

  ! scalars
  if (int_scalars) then
--- a/m_fdm_calc.f90
+++ b/m_fdm_calc.f90
@ -1,6 +1,7 @@
 module m_fdm_calc

  use m_parameters
+  use m_io
  use m_compact

      implicit none
@ -15,10 +16,17 @@ module m_fdm_calc
      integer :: fdmcyc,fdmsavecount
      real*8 :: fdmtime,fdmdt
      real*8 :: sumc,oldsumc,time_int,sum_wrate
+      real*8 :: sumc1,sum_wrate1
      real*8 :: visdt, convdt
      integer :: fdmstep
      real*8 :: coe,tt1,tt2,tt3,tt4
      real*8 :: umax,umin,vmax,vmin,wmax,wmin  ! J. Kwon
+      real*8, dimension(3) :: velmax, velmin, velmax1, velmin1
+
+   integer :: vtype1, vtype2
+   integer, dimension(128) :: scnt1, sdisp1, stype1, rcnt1, rdisp1, rtype1
+   integer, dimension(128) :: scnt2, sdisp2, stype2, rcnt2, rdisp2, rtype2
+

 !===========================================================================
 !===========================================================================
@ -26,134 +34,26 @@ module m_fdm_calc
  contains


-  subroutine fdm_restart_write
-      write(*,*) '======================================================='
-      write(*,*) 'Full results are being written',fullsavenum
-      OPEN (fullsavenum,form='unformatted',status='unknown')
-      write (fullsavenum) fdmtime,nx,ny,nz_all,oldsumc,time_int
-      write (fullsavenum) fdmcyc,DT,dummyu_
-      write (fullsavenum) dt_fdmsave,dt_fullsave
-      write (fullsavenum) t_fdmsave,t_fullsave
-      write (fullsavenum) in_yr,out_yr
-      write (fullsavenum) u_+dummyu_,v_,w_,y1
-      CLOSE (fullsavenum)
-      write(*,*) '======================================================='
-      fullsavenum=fullsavenum+1
-  end subroutine fdm_restart_write
+  subroutine fdm_restart_read
+    integer(kind=MPI_INTEGER_KIND) :: fh
+    integer(kind=MPI_OFFSET_KIND)  :: offset

-  subroutine prepare_fdm
+    write(fname,"('fort.',i4)") restartnum

-      implicit none
-   
-      integer :: i,j,ii,k
-      real*8 :: fl_location,wrate,yr,c
-      real*8 :: tar_sumc
-
-      allocate(u_(nx,ny,nz_all))
-      allocate(v_(nx,ny,nz_all))
-      allocate(w_(nx,ny,nz_all))
-
-      u_=0.0
-      v_=0.0
-      w_=0.0
-
-      if (fdm_sw .eq. 0) then
-         return
-      endif
-
-! DQ initializing
-      fdmsavecount=1        !FDM save count
-      sum_wrate=0.
-      sumc=0.       !for adjusting mean velocity
-      oldsumc=0.
-      time_int=0.
-      visdt=9999.
-      convdt=9999.
-      fdmcyc=0
-      fdmtime=0.
-      fullsavenum=1000      !full save file
-
-      allocate(y1(2,nx,ny,nz_all))
-      allocate(y2(2,nx,ny,nz_all))
-      allocate(yf(2,nx,ny,nz_all))
-
-      CALL ludcmp(nx,ny,nz_all,1,0,0)
-
-      refwr=pre*1.*exp(-ac/(1.+bc*c_ref))  ! Kwon
-
-
-
-!FDM normal start==================================
-      if (restartnum==0) then
-
-       OPEN(305,FILE='sfield.bin',form='unformatted',status='unknown')
-         DO i=1,nx
-           READ (305) y1(2,i,1,1)   ! Yr
-         ENDDO
-       CLOSE (305)
-
-       in_yr=y1(2,1,1,1)         ! inlet_Yr
-       out_yr=y1(2,nx,1,1)       ! outlet_Yr
- 
-       do i=1,ny 
-        do j=1,nz_all 
-         do ii=1,nx 
-           y1(1,ii,i,j)=1.            ! rho initializing
-           y1(2,ii,i,j)=y1(2,ii,1,1)  ! Yr  initializing
-         enddo 
-        enddo 
-       enddo 
- 
-
-! flame location setup
-         sum_wrate=0.; sumc=0.
-         DO ii=1,nz_all
-          DO j=1,ny
-           DO i=1,nx
-            yr=y1(2,i,j,ii)/y1(1,i,j,ii)
-            c=1.-yr
-            IF (c.lt.0.) c=0.
-            wrate=pre*yr*exp(-ac/(1.+bc*c))
-            ! cold boundary difficulty treatment
-            IF(c.le.c_ref) THEN
-               wrate=min_wr
-               IF(c.gt.c_cut) wrate=((refwr-min_wr)*exp(prof_wr*(c-c_ref))+ &
-                  min_wr-refwr*exp(prof_wr*(c_cut-c_ref)))/(1.-exp(prof_wr*(c_cut-c_ref)))
-            ENDIF
-            sum_wrate=sum_wrate+wrate*(hx*hy*hy)
-            sumc=sumc+(1.-yr)
-           ENDDO
-          ENDDO
-         ENDDO
-
-         sum_wrate=sum_wrate/(hy*hy*ny*nz_all)
-         write(*,633) sum_wrate
-
-         fl_location=(hx*(nx-1.))*(1.-(sumc/(nx*ny*nz_all)))
-         sumc=sumc*(hx*hy*hy)/(hy*(ny-1.)*hy*(ny-1.))
-         oldsumc=sumc
-         write(*,634) fl_location/(REAL(nx-1)*hx)*100.
-
-633     format (' ** Consumption Speed, Sc = ',f7.4)
-634     format (' ** Flame Location = ',f7.3, ' % point of x-domain.')
-
-! Restart setup ================================================================
-      else
+    if (myid.eq.master) then
         write (6,*) '********************************************************'
         write (6,*) ' FDM field is being initialized for restarting.'
         write (6,457) restartnum
 457     format('  Restart file number : ',i6)

-         OPEN (restartnum,form='unformatted',status='unknown')
+         OPEN (restartnum,file=fname,form='unformatted',status='unknown')
         read (restartnum) fdmtime,tt1,tt2,tt3,oldsumc,time_int
-         read (restartnum) fdmcyc,tt4,dummyu_
+         read (restartnum) fdmcyc,fdmdt,dummyu_
         read (restartnum) tt1,tt2    !dt_fdmsave,dt_fullsave
         read (restartnum) t_fdmsave,t_fullsave
         read (restartnum) in_yr,out_yr
-         read (restartnum) u_,v_,w_,y1
         CLOSE (restartnum)

-
         if(fdmtime.ge.t_fdmsave) t_fdmsave=t_fdmsave+dt_fdmsave
         if(fdmtime.ge.t_fullsave) t_fullsave=t_fullsave+dt_fullsave

@ -176,6 +76,273 @@ module m_fdm_calc
 455      format('  !! Spectral Time : ',f10.5,' / FDM Time : ',f10.5)
         endif

+    end if
+
+    count = 1
+
+    call MPI_BCAST(fdmtime, count, MPI_REAL8, 0, MPI_COMM_TASK, mpi_err)
+    call MPI_BCAST(tt1, count, MPI_INTEGER8, 0, MPI_COMM_TASK, mpi_err)
+    call MPI_BCAST(tt2, count, MPI_INTEGER8, 0, MPI_COMM_TASK, mpi_err)
+    call MPI_BCAST(tt3, count, MPI_INTEGER8, 0, MPI_COMM_TASK, mpi_err)
+    call MPI_BCAST(oldsumc, count, MPI_REAL8, 0, MPI_COMM_TASK, mpi_err)
+    call MPI_BCAST(time_int, count, MPI_REAL8, 0, MPI_COMM_TASK, mpi_err)
+
+    call MPI_BCAST(fdmcyc, count, MPI_INTEGER8, 0, MPI_COMM_TASK, mpi_err)
+    call MPI_BCAST(fdmdt, count, MPI_REAL8, 0, MPI_COMM_TASK, mpi_err)
+    call MPI_BCAST(dummyu_, count, MPI_REAL8, 0, MPI_COMM_TASK, mpi_err)
+
+    !call MPI_BCAST(dt_fdmsave, count, MPI_REAL8, 0, MPI_COMM_TASK, mpi_err)
+    !call MPI_BCAST(dt_fullsave, count, MPI_REAL8, 0, MPI_COMM_TASK, mpi_err)
+
+    !call MPI_BCAST(t_fdmsave, count, MPI_REAL8, 0, MPI_COMM_TASK, mpi_err)
+    !call MPI_BCAST(t_fullsave, count, MPI_REAL8, 0, MPI_COMM_TASK, mpi_err)
+
+    call MPI_BCAST(in_yr, count, MPI_REAL8, 0, MPI_COMM_TASK, mpi_err)
+    call MPI_BCAST(out_yr, count, MPI_REAL8, 0, MPI_COMM_TASK, mpi_err)
+
+
+  ! opening the file
+    call MPI_INFO_CREATE(mpi_info, mpi_err)
+    call MPI_FILE_OPEN(MPI_COMM_TASK,fname,MPI_MODE_RDONLY,mpi_info,fh,mpi_err)
+
+    count = nx * ny * nz
+    offset = 2*4+6*8 + 2*4+3*8 + 3*(2*4+2*8) + 4 + 0*nx*ny*nz_all*8 + myid*nx*ny*nz*8
+    call MPI_FILE_READ_AT_ALL(fh, offset, u_, count, MPI_REAL8, mpi_status, mpi_err)
+
+    offset = 2*4+6*8 + 2*4+3*8 + 3*(2*4+2*8) + 4 + 1*nx*ny*nz_all*8 + myid*nx*ny*nz*8
+    call MPI_FILE_READ_AT_ALL(fh, offset, v_, count, MPI_REAL8, mpi_status, mpi_err)
+
+    offset = 2*4+6*8 + 2*4+3*8 + 3*(2*4+2*8) + 4 + 2*nx*ny*nz_all*8 + myid*nx*ny*nz*8
+    call MPI_FILE_READ_AT_ALL(fh, offset, w_, count, MPI_REAL8, mpi_status, mpi_err)
+
+    count = 2 * nx * ny * nz
+    offset = 2*4+6*8 + 2*4+3*8 + 3*(2*4+2*8) + 4 + 3*nx*ny*nz_all*8 + myid*2*nx*ny*nz*8
+    call MPI_FILE_READ_AT_ALL(fh, offset, y1, count, MPI_REAL8, mpi_status, mpi_err)
+
+    call MPI_FILE_CLOSE(fh, mpi_err)
+    call MPI_INFO_FREE(mpi_info, mpi_err)
+
+  end subroutine fdm_restart_read
+
+  subroutine fdm_restart_write
+      integer(kind=MPI_INTEGER_KIND) :: fh
+      integer(kind=MPI_OFFSET_KIND)  :: offset
+  ! --------------- writing process ------------------
+
+      if (myid.eq.0) write(*,*) '======================================================='
+      if (myid.eq.0) write(*,*) 'Full results are being written',fullsavenum
+
+  !fname = 'fort.'//savenum
+      write(fname,"('fort.',i4)") fullsavenum
+
+  ! opening the file
+      call MPI_INFO_CREATE(mpi_info, mpi_err)
+      call MPI_FILE_OPEN(MPI_COMM_TASK,fname,MPI_MODE_WRONLY+MPI_MODE_CREATE,mpi_info,fh,mpi_err)
+
+  ! the master node writes the header with parameters
+      if (myid.eq.0) then
+         count = 1
+         call MPI_FILE_WRITE(fh,   8*6, count, MPI_INTEGER4, mpi_status, mpi_err)
+         call MPI_FILE_WRITE(fh,   fdmtime, count, MPI_REAL8, mpi_status, mpi_err)
+         call MPI_FILE_WRITE(fh,   nx, count, MPI_INTEGER8, mpi_status, mpi_err)
+         call MPI_FILE_WRITE(fh,   ny, count, MPI_INTEGER8, mpi_status, mpi_err)
+         call MPI_FILE_WRITE(fh,   nz_all, count, MPI_INTEGER8, mpi_status, mpi_err)
+         call MPI_FILE_WRITE(fh,   oldsumc, count, MPI_REAL8, mpi_status, mpi_err)
+         call MPI_FILE_WRITE(fh,   time_int, count, MPI_REAL8, mpi_status, mpi_err)
+         call MPI_FILE_WRITE(fh,   8*6, count, MPI_INTEGER4, mpi_status, mpi_err)
+
+         call MPI_FILE_WRITE(fh,   8*3, count, MPI_INTEGER4, mpi_status, mpi_err)
+         call MPI_FILE_WRITE(fh,   fdmcyc, count, MPI_INTEGER8, mpi_status, mpi_err)
+         call MPI_FILE_WRITE(fh,   DT, count, MPI_REAL8, mpi_status, mpi_err)
+         call MPI_FILE_WRITE(fh,   dummyu_, count, MPI_REAL8, mpi_status, mpi_err)
+         call MPI_FILE_WRITE(fh,   8*3, count, MPI_INTEGER4, mpi_status, mpi_err)
+
+         call MPI_FILE_WRITE(fh,   8*2, count, MPI_INTEGER4, mpi_status, mpi_err)
+         call MPI_FILE_WRITE(fh,   dt_fdmsave, count, MPI_REAL8, mpi_status, mpi_err)
+         call MPI_FILE_WRITE(fh,   dt_fullsave, count, MPI_REAL8, mpi_status, mpi_err)
+         call MPI_FILE_WRITE(fh,   8*2, count, MPI_INTEGER4, mpi_status, mpi_err)
+
+         call MPI_FILE_WRITE(fh,   8*2, count, MPI_INTEGER4, mpi_status, mpi_err)
+         call MPI_FILE_WRITE(fh,   t_fdmsave, count, MPI_REAL8, mpi_status, mpi_err)
+         call MPI_FILE_WRITE(fh,   t_fullsave, count, MPI_REAL8, mpi_status, mpi_err)
+         call MPI_FILE_WRITE(fh,   8*2, count, MPI_INTEGER4, mpi_status, mpi_err)
+
+         call MPI_FILE_WRITE(fh,   8*2, count, MPI_INTEGER4, mpi_status, mpi_err)
+         call MPI_FILE_WRITE(fh,   in_yr, count, MPI_REAL8, mpi_status, mpi_err)
+         call MPI_FILE_WRITE(fh,   out_yr, count, MPI_REAL8, mpi_status, mpi_err)
+         call MPI_FILE_WRITE(fh,   8*2, count, MPI_INTEGER4, mpi_status, mpi_err)
+
+         call MPI_FILE_WRITE(fh,   8*nx*ny*nz_all*5, count, MPI_INTEGER4, mpi_status, mpi_err)
+
+         offset = 2*4+6*8 + 2*4+3*8 + 3*(2*4+2*8) + 4 + 5*nx*ny*nz_all*8
+         call MPI_FILE_WRITE_AT(fh, offset, 8*nx*ny*nz_all*5, count, MPI_INTEGER4, mpi_status, mpi_err)
+
+      end if
+
+
+      count = nx * ny * nz
+      offset = 2*4+6*8 + 2*4+3*8 + 3*(2*4+2*8) + 4 + 0*nx*ny*nz_all*8 + myid*nx*ny*nz*8
+      call MPI_FILE_WRITE_AT(fh, offset, u_, count, MPI_REAL8, mpi_status, mpi_err)
+
+      offset = 2*4+6*8 + 2*4+3*8 + 3*(2*4+2*8) + 4 + 1*nx*ny*nz_all*8 + myid*nx*ny*nz*8
+      call MPI_FILE_WRITE_AT(fh, offset, v_, count, MPI_REAL8, mpi_status, mpi_err)
+
+      offset = 2*4+6*8 + 2*4+3*8 + 3*(2*4+2*8) + 4 + 2*nx*ny*nz_all*8 + myid*nx*ny*nz*8
+      call MPI_FILE_WRITE_AT(fh, offset, w_, count, MPI_REAL8, mpi_status, mpi_err)
+
+      count = 2 * nx * ny * nz
+      offset = 2*4+6*8 + 2*4+3*8 + 3*(2*4+2*8) + 4 + 3*nx*ny*nz_all*8 + myid*2*nx*ny*nz*8
+      call MPI_FILE_WRITE_AT(fh, offset, y1, count, MPI_REAL8, mpi_status, mpi_err)
+
+
+      call MPI_FILE_CLOSE(fh, mpi_err)
+      call MPI_INFO_FREE(mpi_info, mpi_err)
+
+      if (myid.eq.0) write(*,*) '======================================================='
+      fullsavenum=fullsavenum+1
+  end subroutine fdm_restart_write
+
+  subroutine prepare_fdm
+
+      implicit none
+   
+      integer :: i,j,ii,k
+      real*8 :: fl_location,wrate,yr,c
+      real*8 :: tar_sumc
+      real*8 :: yy(nx)
+
+      if (fdm_sw .eq. 0) then
+         return
+      endif
+
+      allocate(u_(nx,ny,nz))
+      allocate(v_(nx,ny,nz))
+      allocate(w_(nx,ny,nz))
+
+      u_=0.0
+      v_=0.0
+      w_=0.0
+
+! DQ initializing
+      fdmsavecount=1        !FDM save count
+      sum_wrate=0.
+      sum_wrate1=0.
+      sumc=0.       !for adjusting mean velocity
+      sumc1=0.       !for adjusting mean velocity
+      oldsumc=0.
+      time_int=0.
+      visdt=9999.
+      convdt=9999.
+      fdmcyc=0
+      fdmtime=0.
+      fullsavenum=1000      !full save file
+
+      allocate(y1(2,nx,ny,nz))
+      allocate(y2(2,nx,ny,nz))
+      allocate(yf(2,nx,ny,nz))
+
+      y1=0.0
+      y2=0.0
+      yf=0.0
+
+      CALL ludcmp(nx,ny,nz_all,1,0,0)
+
+      CALL MPI_TYPE_VECTOR (nz, nx*nz, nx*ny, MPI_REAL8, vtype1, mpi_err)
+      CALL MPI_TYPE_COMMIT (vtype1, mpi_err)
+
+      CALL MPI_TYPE_VECTOR (nz, 2*nx*nz, 2*nx*ny, MPI_REAL8, vtype2, mpi_err)
+      CALL MPI_TYPE_COMMIT (vtype2, mpi_err)
+
+
+      do i = 1, (ny/nz)
+         scnt1(i) = 1
+         rcnt1(i) = nx * nz * nz
+
+         sdisp1(i) = (i-1) * nx * nz  * 8
+         rdisp1(i) = (i-1) * nx * nz  * nz * 8
+
+         stype1(i) = vtype1
+         rtype1(i) = MPI_REAL8
+
+         scnt2(i) = 1
+         rcnt2(i) = 2 * nx * nz * nz
+
+         sdisp2(i) = (i-1) * 2 * nx * nz  * 8
+         rdisp2(i) = (i-1) * 2 * nx * nz  * nz * 8
+
+         stype2(i) = vtype2
+         rtype2(i) = MPI_REAL8
+
+      end do
+
+      refwr=pre*1.*exp(-ac/(1.+bc*c_ref))  ! Kwon
+
+
+
+!FDM normal start==================================
+      if (restartnum==0) then
+
+       if (myid.eq.master) then
+       OPEN(305,FILE='sfield.bin',form='unformatted',status='unknown')
+         DO i=1,nx
+           READ (305) yy(i)   ! Yr
+         ENDDO
+       CLOSE (305)
+       end if
+
+       call MPI_BCAST(yy, nx, MPI_REAL8, 0, MPI_COMM_TASK, mpi_err)
+
+       in_yr=yy(1)         ! inlet_Yr
+       out_yr=yy(nx)       ! outlet_Yr
+ 
+       do i=1,ny 
+        do j=1,nz
+         do ii=1,nx 
+           y1(1,ii,i,j)=1.            ! rho initializing
+           y1(2,ii,i,j)=yy(ii)  ! Yr  initializing
+         enddo 
+        enddo 
+       enddo 
+ 
+
+! flame location setup
+         sum_wrate=0.; sumc=0.
+         DO ii=1,nz
+          DO j=1,ny
+           DO i=1,nx
+            yr=y1(2,i,j,ii)/y1(1,i,j,ii)
+            c=1.-yr
+            IF (c.lt.0.) c=0.
+            wrate=pre*yr*exp(-ac/(1.+bc*c))
+            ! cold boundary difficulty treatment
+            IF(c.le.c_ref) THEN
+               wrate=min_wr
+               IF(c.gt.c_cut) wrate=((refwr-min_wr)*exp(prof_wr*(c-c_ref))+ &
+                  min_wr-refwr*exp(prof_wr*(c_cut-c_ref)))/(1.-exp(prof_wr*(c_cut-c_ref)))
+            ENDIF
+            sum_wrate=sum_wrate+wrate*(hx*hy*hy)
+            sumc=sumc+(1.-yr)
+           ENDDO
+          ENDDO
+         ENDDO
+
+         sum_wrate1=sum_wrate/(hy*hy*ny*nz_all)
+         CALL MPI_ALLREDUCE(sum_wrate1,sum_wrate,1,MPI_REAL8,MPI_SUM,MPI_COMM_TASK,mpi_err)
+         if (myid.eq.0) write(*,633) sum_wrate
+
+         fl_location=(hx*(nx-1.))*(1.-(sumc/(nx*ny*nz_all)))
+         sumc1=sumc*(hx*hy*hy)/(hy*(ny-1.)*hy*(ny-1.))
+         CALL MPI_ALLREDUCE(sumc1,sumc,1,MPI_REAL8,MPI_SUM,MPI_COMM_TASK,mpi_err)
+         oldsumc=sumc
+         if (myid.eq.0) write(*,634) fl_location/(REAL(nx-1)*hx)*100.
+
+633     format (' ** Consumption Speed, Sc = ',f7.4)
+634     format (' ** Flame Location = ',f7.3, ' % point of x-domain.')
+
+! Restart setup ================================================================
+      else
+
+        CALL fdm_restart_read
        fullsavenum=restartnum+1

      endif
@ -193,9 +360,12 @@ module m_fdm_calc
      real*8 :: fl_location,delu,wrate,yr
      real*8 :: c
      real*8 :: uvel,vvel,wvel  ! J. Kwon
+
+
+
      fdmcyc=itime
 ! Mean velocty setup  
-      do k=1,nz_all
+      do k=1,nz
       do i=1,nx
        do j=1,ny
          u_(i,j,k)=u_(i,j,k)+dummyu_
@ -203,32 +373,34 @@ module m_fdm_calc
       enddo
      enddo

-       write(*,*)
-       write(*,*) '=========================================================='
+       if (myid.eq.0) write(*,*)
+       if (myid.eq.0) write(*,*) '=========================================================='

-      WRITE(*,932) TIME,DT
+      if (myid.eq.0) WRITE(*,932) TIME,DT
 932   format(' Spectral results at time = ',f10.5,', dT = ',f7.5)
-      write(*,933) dummyu_
+      if (myid.eq.0) write(*,933) dummyu_
 933   format(' ** Mean U = ',f7.4)

 !----------------------------------------------------------------
 ! FDM DNS BGN
 !----------------------------------------------------------------
-          umax=-1.0d300;umin=1.0d300;vmax=-1.0d300;vmin=1.0d300;wmax=-1.0d300;wmin=1.0d300
+          velmax1(1)=maxval(u_)
+          velmin1(1)=minval(u_)
+          velmax1(2)=maxval(v_)
+          velmin1(2)=minval(v_)
+          velmax1(3)=maxval(w_)
+          velmin1(3)=minval(w_)

-!$omp parallel do private(i,j,k) reduction(max:umax,vmax,wmax) reduction(min:umin,vmin,wmin)  NUM_THREADS(6) schedule(static)
-          do k=1,nz_all
-          do j=1,ny
-          do i=1,nx
-            umax=max(umax,u_(i,j,k))
-            umin=min(umin,u_(i,j,k))
-            vmax=max(vmax,v_(i,j,k))
-            vmin=min(vmin,v_(i,j,k))
-            wmax=max(wmax,w_(i,j,k))
-            wmin=min(wmin,w_(i,j,k))
-          enddo
-          enddo
-          enddo
+          ! MPI_REDUCE [uvw](max|min)
+          CALL MPI_ALLREDUCE(velmax1,velmax,3,MPI_REAL8,MPI_MAX,MPI_COMM_TASK,mpi_err)
+          CALL MPI_ALLREDUCE(velmin1,velmin,3,MPI_REAL8,MPI_MIN,MPI_COMM_TASK,mpi_err)
+
+          umax=velmax(1)
+          umin=velmin(1)
+          vmax=velmax(2)
+          vmin=velmin(2)
+          wmax=velmax(3)
+          wmin=velmin(3)

          uvel = max(abs(umax), abs(umin))
          vvel = max(abs(vmax), abs(vmin))
@ -241,12 +413,14 @@ module m_fdm_calc

        fdmdt=min(DT,visdt,convdt)        
        DT=fdmdt
-        call solve(nx,ny,nz_all,u_,v_,w_,y1,y2,yf)
+        call solve(nx,ny,nz,u_,v_,w_,y1,y2,yf)
        fdmtime=time+fdmdt
-        write(*,'(a30,3x,4f12.7)')' ** DT, visdt, convdt, fdmdt =' , DT,visdt,convdt,fdmdt
+        if (myid.eq.0) write(*,'(a30,3x,4f12.7)')' ** DT, visdt, convdt, fdmdt =' , DT,visdt,convdt,fdmdt

-!$omp parallel do private(i,j,k,yr,wrate) reduction(+:sum_wrate,sumc) NUM_THREADS(6) schedule(static)
-        do k=1,nz_all
+        sumc=0.
+        sum_wrate=0.
+
+        do k=1,nz
        do j=1,ny
        do i=1,nx
           yr=y1(2,i,j,k)/y1(1,i,j,k)
@ -263,21 +437,24 @@ module m_fdm_calc
        enddo
        enddo

+
 ! Real time results for Sc and flame location.
-        sum_wrate=sum_wrate/(hy*hy*REAL(ny*nz_all))
-        write(*,633) sum_wrate
+        sum_wrate1=sum_wrate/(hy*hy*REAL(ny*nz_all))
+        CALL MPI_ALLREDUCE(sum_wrate1,sum_wrate,1,MPI_REAL8,MPI_SUM,MPI_COMM_TASK,mpi_err)
+        if (myid.eq.0) write(*,633) sum_wrate
 633     format (' ** Consumption Speed, Sc = ',f7.4)

        fl_location=(hx*REAL(nx-1))*(1.-(sumc/(REAL(nx*ny*nz_all))))
-        sumc=sumc*(hx*hy*hy)/(hy*(REAL(ny)-1.)*hy*(REAL(ny)-1.))
-        write(*,634) fl_location/(REAL(nx-1)*hx)*100.
+        sumc1=sumc*(hx*hy*hy)/(hy*(REAL(ny)-1.)*hy*(REAL(ny)-1.))
+        CALL MPI_ALLREDUCE(sumc1,sumc,1,MPI_REAL8,MPI_SUM,MPI_COMM_TASK,mpi_err)
+        if (myid.eq.0) write(*,634) fl_location/(REAL(nx-1)*hx)*100.
 634     format (' ** Flame Location = ',f7.3, ' % point of x-domain.')


-        write(St_data,'(f10.5,10e25.15)')fdmtime,dummyu_,sum_wrate,fl_location/(REAL(nx-1)*hx)*100.
+        if (myid.eq.0) write(St_data,'(f10.5,10e25.15)')fdmtime,dummyu_,sum_wrate,fl_location/(REAL(nx-1)*hx)*100.

 ! Control the inflow mean velocity, dummyu_
-        write(*,'(a20,2f12.8)') ' ** oldsumc, sumc = ',oldsumc,sumc
+        if (myid.eq.0) write(*,'(a20,2f12.8)') ' ** oldsumc, sumc = ',oldsumc,sumc
        if(swadtv.ne.0.and.mod((fdmcyc+1),swadtv).eq.0) then
          time_int=fdmtime-time_int
          delu=(sumc-oldsumc)/time_int
@ -291,11 +468,11 @@ module m_fdm_calc
 ! FDM DNS END
 !----------------------------------------------------------------------------
       
-      write(*,*) 'Number of whole calculation = ',fdmcyc
-      write(*,*) '=========================================================='
-      write(*,*)
+      if (myid.eq.0) write(*,*) 'Number of whole calculation = ',fdmcyc
+      if (myid.eq.0) write(*,*) '=========================================================='
+      if (myid.eq.0) write(*,*)

-       call write_vel_max ( umax,vmax,wmax,umin,vmin,wmin )
+      if (myid.eq.0) call write_vel_max ( umax,vmax,wmax,umin,vmin,wmin )

     return 

@ -315,8 +492,63 @@ module m_fdm_calc
      real*8 :: uz(4,yy),duz(4,yy),d2uz(yy)
      real*8 :: Ly(yy,zz),Dy

+      real*8 :: r1_t(2,xx,zz,yy),f_t(2,xx,zz,yy)
+      real*8 :: ww_t(xx,zz,yy)
+
+
+   !CALL MPI_ALLTOALLW (ww_, scnt1, sdisp1, stype1, ww_t, rcnt1, rdisp1, rtype1, MPI_COMM_TASK, mpi_err)
+
+   !CALL MPI_ALLTOALLW (r1_, scnt2, sdisp2, stype2, r1_t, rcnt2, rdisp2, rtype2, MPI_COMM_TASK, mpi_err)
+
+    do i = 1,numprocs
+        CALL MPI_ISEND(ww_(1,zz*(i-1)+1,1), scnt1(i), stype1(i), (i-1), myid, MPI_COMM_TASK, mpi_request, mpi_err)
+    end do
+    do i = 1,numprocs
+        CALL MPI_RECV(ww_t(1,1,zz*(i-1)+1), rcnt1(i), rtype1(i), (i-1), (i-1), MPI_COMM_TASK, mpi_status, mpi_err)
+    end do
+
+    do i = 1,numprocs
+        CALL MPI_ISEND(r1_(1,1,zz*(i-1)+1,1), scnt2(i), stype2(i), (i-1), myid, MPI_COMM_TASK, mpi_request, mpi_err)
+    end do
+    do i = 1,numprocs
+        CALL MPI_RECV(r1_t(1,1,1,zz*(i-1)+1), rcnt2(i), rtype2(i), (i-1), (i-1), MPI_COMM_TASK, mpi_status, mpi_err)
+    end do
+
+!! z-direction
+      DO j=1,zz
+       DO i=1,xx
+        DO k=1,yy
+         uz(1,k)=r1_t(1,i,j,k)                            ! 1:rho
+         uz(2,k)=r1_t(2,i,j,k)/r1_t(1,i,j,k)               ! 2:Y
+         uz(3,k)=uz(1,k)*ww_t(i,j,k)                       ! 3:rho*w
+         uz(4,k)=uz(3,k)*uz(2,k)                         ! 4:rho*w*Y
+        ENDDO
+
+       CALL dfp(yy,hy,uz(1:4,:),duz(1:4,:),4,3)
+       CALL d2fp(yy,hy,uz(2,:),d2uz(:),1,3)
+
+        DO k=1,yy
+! -( d(rho*w)/dz )
+         f_t(1,i,j,k)=-duz(3,k)         ! continuity
+
+! -( d(rho*w*Yr)/dz ) + d(rho*D* d(Yr)/dz)/dz
+! = -( d(rho*w*Yr)/dz )
+!    + D* (rho* d2(Yr)/dz2 + d(rho)/dz * d(Yr)/dz )
+         f_t(2,i,j,k)=-duz(4,k) + diff*(uz(1,k)*d2uz(k)+duz(1,k)*duz(2,k))   ! species conserv.
+        ENDDO
+       ENDDO
+      ENDDO
+
+   !CALL MPI_ALLTOALLW (f_t, rcnt2, rdisp2, rtype2, f_, scnt2, sdisp2, stype2, MPI_COMM_TASK, mpi_err)
+
+    do i = 1,numprocs
+        CALL MPI_ISEND(f_t(1,1,1,zz*(i-1)+1), rcnt2(i), rtype2(i), (i-1), myid, MPI_COMM_TASK, mpi_request, mpi_err)
+    end do
+    do i = 1,numprocs
+        CALL MPI_RECV(f_(1,1,zz*(i-1)+1,1), scnt2(i), stype2(i), (i-1), (i-1), MPI_COMM_TASK, mpi_status, mpi_err)
+    end do
+
 ! x-direction
-!$omp parallel do private(i,j,k,ux,dux,d2ux,wrate,uy,duy,d2uy) NUM_THREADS(6) schedule(static)
      DO k=1,zz
       DO j=1,yy
        DO i=1,xx
@ -339,12 +571,12 @@ module m_fdm_calc
 

 ! -( d(rho*u)/dx )
-         f_(1,i,j,k)=-dux(3,i)   ! continuity
+         f_(1,i,j,k) = f_(1,i,j,k) - dux(3,i)   ! continuity

 ! -( d(rho*u*Yr)/dx ) + d(rho*D* d(Yr)/dx)/dx
 ! = -( d(rho*u*Yr)/dx )
 !    + D* (rho* d2(Yr)/dx2 + d(rho)/dx * d(Yr)/dx )
-         f_(2,i,j,k)=-dux(4,i) + diff*(ux(1,i)*d2ux(i)+dux(1,i)*dux(2,i)) - wrate  ! species conservation
+         f_(2,i,j,k) = f_(2,i,j,k) - dux(4,i) + diff*(ux(1,i)*d2ux(i)+dux(1,i)*dux(2,i)) - wrate  ! species conservation

        ENDDO
       ENDDO
@ -373,34 +605,8 @@ module m_fdm_calc
       ENDDO
      ENDDO

-!! z-direction
-!$omp parallel do private(i,j,k,uz,duz,d2uz) NUM_THREADS(6) schedule(static)
-      DO j=1,yy
-       DO i=1,xx
-        DO k=1,zz
-         uz(1,k)=r1_(1,i,j,k)                            ! 1:rho
-         uz(2,k)=r1_(2,i,j,k)/r1_(1,i,j,k)               ! 2:Y
-         uz(3,k)=uz(1,k)*ww_(i,j,k)                       ! 3:rho*w
-         uz(4,k)=uz(3,k)*uz(2,k)                         ! 4:rho*w*Y
-        ENDDO
-
-       CALL dfp(zz,hy,uz(1:4,:),duz(1:4,:),4,3)
-       CALL d2fp(zz,hy,uz(2,:),d2uz(:),1,3)    
-
-        DO k=1,zz
-! -( d(rho*w)/dz ) 
-         f_(1,i,j,k)=f_(1,i,j,k)-duz(3,k)         ! continuity
-
-! -( d(rho*w*Yr)/dz ) + d(rho*D* d(Yr)/dz)/dz
-! = -( d(rho*w*Yr)/dz ) 
-!    + D* (rho* d2(Yr)/dz2 + d(rho)/dz * d(Yr)/dz )
-         f_(2,i,j,k)=f_(2,i,j,k)-duz(4,k) + diff*(uz(1,k)*d2uz(k)+duz(1,k)*duz(2,k))   ! species conserv.
-        ENDDO
-       ENDDO
-      ENDDO

 !    Boundary condition
-!$omp parallel do private(k,j,i) NUM_THREADS(6) schedule(static)
      DO k=1,zz
       DO j=1,yy
        DO i=1,yrsw
@ -473,7 +679,6 @@ module m_fdm_calc
      IF(istage<5) THEN
        at=a(istage)*fdmdt
        bt=(b(istage)-a(istage))*fdmdt
-!$omp parallel do private(k,j,i,nv) NUM_THREADS(6) schedule(static)
        DO k=1,zz
         DO j=1,yy
          DO i=1,xx
@ -486,9 +691,6 @@ module m_fdm_calc
        ENDDO
      ELSE
        bt=b(istage)*fdmdt
-        sumc=0.
-        sum_wrate=0.
-!$omp parallel do private(k,j,i,nv) NUM_THREADS(6) schedule(static)
        DO k=1,zz
         DO j=1,yy
          DO i=1,xx
--- a/main.f90
+++ b/main.f90
@ -67,10 +67,7 @@ program x_code
  ! Initializing the LES stuff
  if (les) call m_les_begin

-     if (task.eq.'hydro' .and. myid.eq.master) then
-         call prepare_fdm
-     end if
-     call MPI_BARRIER(MPI_COMM_TASK, mpi_err)
+     call prepare_fdm

     call stat_main

@ -134,20 +131,7 @@ program x_code

        if (fdm_sw.ne.0) then  ! Reacting
            call gather_4
-            if (myid.eq.master) then
               call fdm_exe
-               do id_to=1,numprocs-1
-                 call MPI_Isend (dummy, 1, MPI_REAL, id_to, 0, MPI_COMM_TASK, mpi_request, mpi_err)
-               end do
-            else
-               call MPI_Irecv (dummy, 1, MPI_REAL, master, 0, MPI_COMM_TASK, mpi_request, mpi_err)
-               do
-                 call MPI_Test(mpi_request, flag, mpi_status, mpi_err)
-                 if (flag) exit
-                 call SLEEPQQ (200)
-               end do
-            end if
-            call MPI_BARRIER(MPI_COMM_TASK, mpi_err)
        endif

        ! RHS for scalars
@ -202,18 +186,18 @@ program x_code
        ! write the restart file if it's the time
        if (mod(itime,IPRINT2).eq.0) then
           call restart_write_parallel
-           !remove hit restart file
-           call system ('rm -f '//run_name//'.64.'//ex_file_ext)
           if (fdm_sw.ne.0) then  ! Reacting
               call gather_4
-               if (myid.eq.master) then
-                  !write itime, savenum pair
-                  if(TIME.ge.bgn_save) call fdm_restart_write
-                  write (restart_pair, *) file_ext, fullsavenum-1
-                  call write_tp_field
+               if (TIME.ge.bgn_save) then
+                 call fdm_restart_write
+                 if (myid.eq.master) then
+                   write (restart_pair, *) file_ext, fullsavenum-1, '\n'
+                   !call write_tp_field
+                 end if
               end if
-               call MPI_BARRIER(MPI_COMM_TASK, mpi_err)
           endif
+           !remove HIT3D restart file except the latest
+           call system ('rm -f '//run_name//'.64.'//ex_file_ext)
           ex_file_ext = file_ext
        endif

@ -266,7 +250,7 @@ program x_code
           if (task_split) call fields_to_stats
           if (mod(itime,iprint1).eq.0) call stat_main
           if (fdm_sw.eq.0) then  ! non-reacting
-              if (mod(itime,iprint2).eq.0) call io_write_4
+              !if (mod(itime,iprint2).eq.0) call io_write_4
           endif
        end if
     end if stats
Author	SHA1	Message	Date
ignis	9854238fb3	fdm restart for mpi	2014-08-22 14:19:00 +09:00
ignis	c3345ba9c2	mpi communication established	2014-08-22 13:59:56 +09:00
ignis	344a5b085a	gather4 only performs inverse transform	2014-08-22 13:43:31 +09:00
ignis	5843ddccf0	every process call prepare_fdm and fdm_exe, do not work correctly	2014-08-22 13:41:11 +09:00
ignis	fb5a66c648	remove omp	2014-08-22 12:01:12 +09:00
ignis	cc56c9603b	remove old HIT3D save after writing fdm save	2014-08-22 11:51:42 +09:00
ignis	2e3bbc1163	disable io	2014-08-10 16:42:27 +09:00