!********************************************************************** subroutine plotf(t,fout,myid,numprocs,comm,twoslice,subslice, > id1d,id2d, > u,v,w,temp,uf,vf,wf,tempf) !********************************************************************** ! !-Subroutine that dumps out 3-D fiels of (u,v,w,rho) in !-SINGLE PRECISION for future postprocessing. Uses MPI-2 !-in exactly the same fashion as subroutine *output* does. !-Created by PD on 8/8/05. !-NOTE: As in serial version of code data is outputted in !-Physical Space. Though this may entail some additional FFTs !-it makes the data more "shareable" with the community. include 'dim.h' include 'mpif.h' ! include 'comwave.h' include 'comleng.h' include 'comgridptxyz.h' include 'comsubd.h' include 'comsubdpar.h' !-3-D Arrays (GLOBAL) !-a) Phys. Space real, dimension(nxpp,nypl,nzpl) :: u,v,w,temp !-b) Fourier Space real, dimension(nxpl,ny,nzpl) :: uf,vf,wf,tempf !-Even in the 2-D case, the global array is considered 3-D !-with 3rd dimension equal to 1. integer, parameter :: ndims_io=3 character (len=60) :: fout !------------- !MPI Variables !------------- integer :: myid, numprocs, ierr, comm integer :: twoslice, subslice integer, dimension(nproch*nprocv,2) :: id2d integer, dimension(nproch,nprocv) :: id1d integer :: fh, filetype, local_array_size !-Note the special definition of the displacement variable integer (kind=MPI_OFFSET_KIND) disp integer status_array(MPI_STATUS_SIZE) integer sizeofreal, sizeofinteger, sizeoffloat !-Variables needed to set up subarray for MPI I/O integer, dimension(3) :: gsize ! Size of each dimension of global array integer, dimension(3) :: lsize ! Size of each dimension of local array integer, dimension(3) :: start_index ! Beginning index of each local array ! element in global array !-Some auxilliary variables integer, dimension(8) :: auxint real, dimension(4) :: auxreal !-Temporary output array !-NOTE: Kind = 4 is compiler dependent !! real(kind=4), dimension(:,:,:), allocatable :: tempout integer AllocateStatus !-Allocate local arrays allocate ( tempout(nxpp,nypl,nzpl), > stat = AllocateStatus) if (AllocateStatus /= 0) then stop "**Not Enough Memory - OUTPUT**" end if !-Before doing anything else convert flow field back to !-Physical Space call horfft(u,uf,1, > myid,numprocs,comm, > twoslice,subslice,id1d,id2d) call horfft(v,vf,1, > myid,numprocs,comm, > twoslice,subslice,id1d,id2d) call horfft(w,wf,1, > myid,numprocs,comm, > twoslice,subslice,id1d,id2d) call horfft(temp,tempf,1, > myid,numprocs,comm, > twoslice,subslice,id1d,id2d) !-Find the size in bytes of MPI_REAL8, MPI_REAL and MPI_INTEGER call MPI_TYPE_EXTENT(MPI_REAL,sizeoffloat,ierr) call MPI_TYPE_EXTENT(MPI_REAL8,sizeofreal,ierr) call MPI_TYPE_EXTENT(MPI_INTEGER,sizeofinteger,ierr) !-Size of each dimension of global array gsize(1) = nxpp gsize(2) = ny gsize(3) = nz !-Size of each dimension of local array lsize(1) = nxpp lsize(2) = nypl lsize(3) = nzpl !-Starting indices of the first element of the local array of !-this specific processor in the global array !-Need to first identify 2-D coordinates of processor hproc = id2d(myid+1,1) vproc = id2d(myid+1,2) !-Remember C-indexing is used here ! start_index(1) = 0 start_index(2) = (hproc-1)*nypl start_index(3) = (vproc-1)*nzpl !-Now create the subarray !-Data are outputted in REAL type for restart ! !-(REMEMBER, WE'RE USING SINGLE PRECISION !) !-Then commit it call MPI_TYPE_CREATE_SUBARRAY( ndims_io,gsize,lsize,start_index, > MPI_ORDER_FORTRAN,MPI_REAL,filetype, > ierr ) call MPI_TYPE_COMMIT(filetype,ierr) !-Open file call MPI_FILE_OPEN(comm,fout, > MPI_MODE_CREATE + MPI_MODE_WRONLY, > MPI_INFO_NULL,fh,ierr) disp = 0 !-Create view for individual processor !-to dump out 3-D arrays call MPI_FILE_SET_VIEW(fh, disp, MPI_REAL, filetype, > 'native', > MPI_INFO_NULL,ierr) !-Now dump out arrays (These are all in FOURIER space) local_array_size = lsize(1)*lsize(2)*lsize(3) !-The temporary copies are not done through some subroutine !-to avoid screw-ups with memory addressing (PD: Improve this asap !) !-Dump u-velocity do k=1,nzpl do j=1,nypl do i=1,nxpp tempout(i,j,k) = real(u(i,j,k)) enddo enddo enddo call MPI_FILE_WRITE_ALL(fh, tempout, local_array_size, > MPI_REAL, > MPI_STATUS_IGNORE,ierr) !-Dump v-velocity do k=1,nzpl do j=1,nypl do i=1,nxpp tempout(i,j,k) = real(v(i,j,k)) enddo enddo enddo call MPI_FILE_WRITE_ALL(fh, tempout, local_array_size, > MPI_REAL, > MPI_STATUS_IGNORE,ierr) !-Dump w-velocity do k=1,nzpl do j=1,nypl do i=1,nxpp tempout(i,j,k) = real(w(i,j,k)) enddo enddo enddo call MPI_FILE_WRITE_ALL(fh, tempout, local_array_size, > MPI_REAL, > MPI_STATUS_IGNORE,ierr) !-Dump Active Scalar do k=1,nzpl do j=1,nypl do i=1,nxpp tempout(i,j,k) = real(temp(i,j,k)) enddo enddo enddo call MPI_FILE_WRITE_ALL(fh, tempout, local_array_size, > MPI_REAL, > MPI_STATUS_IGNORE,ierr) !-Now dump out information on grid resolution/number of processors !----------------------------------------------------------------- disp = gsize(1)*gsize(2)*gsize(3)*sizeoffloat*4 call MPI_FILE_SET_VIEW(fh, disp, MPI_INTEGER, filetype, > 'native', > MPI_INFO_NULL,ierr) auxint(1) = nx auxint(2) = ny auxint(3) = nz auxint(4) = nzloc auxint(5) = nsubd auxint(6) = nsdtype auxint(7) = nproch auxint(8) = nprocv if (myid == 0) call MPI_FILE_WRITE(fh,auxint,8,MPI_INTEGER, >status_array,ierr) disp = disp + sizeofinteger*8 100 continue !-Dump out some physical constants etc. !-------------------------------------- call MPI_FILE_SET_VIEW(fh, disp, MPI_REAL8, filetype, 'native', > MPI_INFO_NULL,ierr) auxreal(1) = t auxreal(2) = xlen auxreal(3) = ylen auxreal(4) = zlen if (myid == 0) call MPI_FILE_WRITE(fh,auxreal,4,MPI_REAL8, > status_array,ierr) disp = disp + sizeofreal*4 200 continue !-Close output file call MPI_FILE_CLOSE(fh,ierr) !-Free Subarray datatype call MPI_TYPE_FREE(filetype,ierr) !-Convert back to Fourier space call horfft(u,uf,-1, > myid,numprocs,comm, > twoslice,subslice,id1d,id2d) call horfft(v,vf,-1, > myid,numprocs,comm, > twoslice,subslice,id1d,id2d) call horfft(w,wf,-1, > myid,numprocs,comm, > twoslice,subslice,id1d,id2d) call horfft(temp,tempf,-1, > myid,numprocs,comm, > twoslice,subslice,id1d,id2d) !-De-Allocate local array deallocate (tempout, > stat=AllocateStatus) if (AllocateStatus /= 0) then stop "**Error Deallocating - OUTPUT**" end if end subroutine plotf