From e4c0c57f35f67322fc5660efd69d6309a04895e2 Mon Sep 17 00:00:00 2001 From: chowland Date: Fri, 7 Oct 2022 10:01:25 +0200 Subject: [PATCH] New HDF5 interface for movie writing --- Makefile | 7 +- src/flow_solver/MakeMovieXCut.F90 | 201 ---------- src/flow_solver/MakeMovieYCut.F90 | 281 ------------- src/flow_solver/MakeMovieZCut.F90 | 370 ------------------ .../HdfReadContinua.F90 | 0 src/{flow_solver => h5tools}/HdfRoutines.F90 | 0 src/h5tools/MakeMovieXCut.F90 | 70 ++++ src/h5tools/MakeMovieYCut.F90 | 96 +++++ src/h5tools/MakeMovieZCut.F90 | 96 +++++ src/h5tools/h5_tools.F90 | 129 ++++++ src/h5tools/mean_zplane.F90 | 228 +++++++++++ src/h5tools/means.F90 | 157 ++++++++ src/main.F90 | 10 +- 13 files changed, 791 insertions(+), 854 deletions(-) delete mode 100644 src/flow_solver/MakeMovieXCut.F90 delete mode 100644 src/flow_solver/MakeMovieYCut.F90 delete mode 100644 src/flow_solver/MakeMovieZCut.F90 rename src/{flow_solver => h5tools}/HdfReadContinua.F90 (100%) rename src/{flow_solver => h5tools}/HdfRoutines.F90 (100%) create mode 100644 src/h5tools/MakeMovieXCut.F90 create mode 100644 src/h5tools/MakeMovieYCut.F90 create mode 100644 src/h5tools/MakeMovieZCut.F90 create mode 100644 src/h5tools/h5_tools.F90 create mode 100644 src/h5tools/mean_zplane.F90 create mode 100644 src/h5tools/means.F90 diff --git a/Makefile b/Makefile index 1d118389..80df15ca 100644 --- a/Makefile +++ b/Makefile @@ -107,9 +107,12 @@ OBJS += obj/AddLatentHeat.o obj/DeallocatePFVariables.o obj/ExplicitTermsPhi.o \ obj/InterpTempMgrd.o obj/SolveImpEqnUpdate_Phi.o obj/CreateICPF.o \ obj/ImmersedBoundary.o +# Object files for plane writing +OBJS += obj/mean_zplane.o + # Module object files MOBJS = obj/param.o obj/decomp_2d.o obj/AuxiliaryRoutines.o obj/decomp_2d_fft.o \ - obj/HermiteInterpolations.o obj/GridModule.o + obj/HermiteInterpolations.o obj/GridModule.o obj/h5_tools.o obj/means.o #======================================================================= # Files that create modules: @@ -142,6 +145,8 @@ $(OBJDIR)/%.o: src/%.F90 $(MOBJS) $(FC) -c -o $@ $< $(LDFLAGS) $(OBJDIR)/%.o: src/flow_solver/%.F90 $(MOBJS) $(FC) -c -o $@ $< $(LDFLAGS) +$(OBJDIR)/%.o: src/h5tools/%.F90 $(MOBJS) + $(FC) -c -o $@ $< $(LDFLAGS) $(OBJDIR)/%.o: src/multires/%.F90 $(MOBJS) $(FC) -c -o $@ $< $(LDFLAGS) $(OBJDIR)/%.o: src/multires/IC_interpolation/%.F90 $(MOBJS) diff --git a/src/flow_solver/MakeMovieXCut.F90 b/src/flow_solver/MakeMovieXCut.F90 deleted file mode 100644 index 39006d8c..00000000 --- a/src/flow_solver/MakeMovieXCut.F90 +++ /dev/null @@ -1,201 +0,0 @@ -!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! -! ! -! FILE: MakeMovieXCut.F90 ! -! CONTAINS: subroutine Mkmov_xcut ! -! ! -! PURPOSE: Write down xcut snapshot for flow vars ! -! ! -!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! - -subroutine Mkmov_xcut - use param - use mpih - use hdf5 - use decomp_2d, only: xstart,xend,xstartr,xendr,DECOMP_2D_COMM_CART_X - use local_arrays, only: vz,vy,vx,temp - use mgrd_arrays, only: sal,phi - implicit none - character(70) :: filename - character(30) :: dsetname - character( 5) :: frame - - integer(HID_T) :: file_id, group_id, plist_id - integer(HID_T) :: filespace, slabspace, memspace - integer(HID_T) :: dset_vx, dset_vy, dset_vz, dset_temp, dset_sal - integer(HSIZE_T), dimension(2) :: dims, data_count, data_offset - integer :: ic, hdf_error, ndims, comm, info, ierror - - logical :: fexist, dsetexists - - ! Select plane - plane next to lower wall - ic = 1 - - ! Record frame number and filename as strings - write(frame,"(i5.5)")nint(time/tframe) - filename = trim("outputdir/flowmov/movie_xcut.h5") - - ! MPI - call MPI_CART_SUB(DECOMP_2D_COMM_CART_X, (/.true.,.true./), comm, ierror) - info = MPI_INFO_NULL - - ! Begin HDF5 - call h5pcreate_f(H5P_FILE_ACCESS_F, plist_id, hdf_error) - call h5pset_fapl_mpio_f(plist_id, comm, info, hdf_error) - - ! Check if movie file already exists, and create it if not - inquire(file=filename,exist=fexist) - if (fexist) then - call h5fopen_f(filename, H5F_ACC_RDWR_F, file_id, hdf_error, access_prp=plist_id) - else - call h5fcreate_f(filename, H5F_ACC_TRUNC_F, file_id, hdf_error,access_prp=plist_id) - call h5gcreate_f(file_id, "vx", group_id, hdf_error) - call h5gclose_f(group_id, hdf_error) - call h5gcreate_f(file_id, "vy", group_id, hdf_error) - call h5gclose_f(group_id, hdf_error) - call h5gcreate_f(file_id, "vz", group_id, hdf_error) - call h5gclose_f(group_id, hdf_error) - call h5gcreate_f(file_id, "temp", group_id, hdf_error) - call h5gclose_f(group_id, hdf_error) - if (salinity) then - call h5gcreate_f(file_id, "sal", group_id, hdf_error) - call h5gclose_f(group_id, hdf_error) - end if - if (phasefield) then - call h5gcreate_f(file_id, "phi", group_id, hdf_error) - call h5gclose_f(group_id, hdf_error) - end if - end if - - call h5pclose_f(plist_id, hdf_error) - - ! Define dimension - ndims = 2 - dims(1) = nym - dims(2) = nzm - - data_count(1) = xend(2) - xstart(2) + 1 - data_count(2) = xend(3) - xstart(3) + 1 - - data_offset(1) = xstart(2) - 1 - data_offset(2) = xstart(3) - 1 - - ! Create dataspace (file & memory) - call h5screate_simple_f(ndims, dims, filespace, hdf_error) - call h5screate_simple_f(ndims, data_count, memspace, hdf_error) - - !! Create datasets and write data to file - ! vx - dsetname = trim("vx/"//frame) - call h5lexists_f(file_id, dsetname, dsetexists, hdf_error) - if (dsetexists) call h5ldelete_f(file_id, dsetname, hdf_error) - call h5dcreate_f(file_id, dsetname, H5T_NATIVE_DOUBLE, filespace, dset_vx, hdf_error) - call h5dget_space_f(dset_vx, filespace, hdf_error) - call h5sselect_hyperslab_f (filespace, H5S_SELECT_SET_F,data_offset, data_count, hdf_error) - call h5pcreate_f(H5P_DATASET_XFER_F, plist_id, hdf_error) - call h5pset_dxpl_mpio_f(plist_id, H5FD_MPIO_COLLECTIVE_F, hdf_error) - call h5dwrite_f(& - dset_vx, H5T_NATIVE_DOUBLE,& - vx(ic,xstart(2):xend(2),xstart(3):xend(3)),& - data_count, hdf_error, file_space_id = filespace,& - mem_space_id = memspace, xfer_prp = plist_id) - call h5pclose_f(plist_id, hdf_error) - call h5dclose_f(dset_vx, hdf_error) - - ! vy - dsetname = trim("vy/"//frame) - call h5lexists_f(file_id, dsetname, dsetexists, hdf_error) - if (dsetexists) call h5ldelete_f(file_id, dsetname, hdf_error) - call h5dcreate_f(file_id, dsetname, H5T_NATIVE_DOUBLE, filespace, dset_vy, hdf_error) - call h5dget_space_f(dset_vy, filespace, hdf_error) - call h5sselect_hyperslab_f (filespace, H5S_SELECT_SET_F,data_offset, data_count, hdf_error) - call h5pcreate_f(H5P_DATASET_XFER_F, plist_id, hdf_error) - call h5pset_dxpl_mpio_f(plist_id, H5FD_MPIO_COLLECTIVE_F, hdf_error) - call h5dwrite_f(& - dset_vy, H5T_NATIVE_DOUBLE,& - vy(ic,xstart(2):xend(2),xstart(3):xend(3)),& - data_count, hdf_error, file_space_id = filespace,& - mem_space_id = memspace, xfer_prp = plist_id) - call h5pclose_f(plist_id, hdf_error) - call h5dclose_f(dset_vy, hdf_error) - - ! vz - dsetname = trim("vz/"//frame) - call h5lexists_f(file_id, dsetname, dsetexists, hdf_error) - if (dsetexists) call h5ldelete_f(file_id, dsetname, hdf_error) - call h5dcreate_f(file_id, dsetname, H5T_NATIVE_DOUBLE, filespace, dset_vz, hdf_error) - call h5dget_space_f(dset_vz, filespace, hdf_error) - call h5sselect_hyperslab_f (filespace, H5S_SELECT_SET_F,data_offset, data_count, hdf_error) - call h5pcreate_f(H5P_DATASET_XFER_F, plist_id, hdf_error) - call h5pset_dxpl_mpio_f(plist_id, H5FD_MPIO_COLLECTIVE_F, hdf_error) - call h5dwrite_f(& - dset_vz, H5T_NATIVE_DOUBLE,& - vz(ic,xstart(2):xend(2),xstart(3):xend(3)),& - data_count, hdf_error, file_space_id = filespace,& - mem_space_id = memspace, xfer_prp = plist_id) - call h5pclose_f(plist_id, hdf_error) - call h5dclose_f(dset_vz, hdf_error) - - ! temp - dsetname = trim("temp/"//frame) - call h5lexists_f(file_id, dsetname, dsetexists, hdf_error) - if (dsetexists) call h5ldelete_f(file_id, dsetname, hdf_error) - call h5dcreate_f(file_id, dsetname, H5T_NATIVE_DOUBLE, filespace, dset_temp, hdf_error) - call h5dget_space_f(dset_temp, filespace, hdf_error) - call h5sselect_hyperslab_f (filespace, H5S_SELECT_SET_F,data_offset, data_count, hdf_error) - call h5pcreate_f(H5P_DATASET_XFER_F, plist_id, hdf_error) - call h5pset_dxpl_mpio_f(plist_id, H5FD_MPIO_COLLECTIVE_F, hdf_error) - call h5dwrite_f(& - dset_temp, H5T_NATIVE_DOUBLE,& - temp(ic,xstart(2):xend(2),xstart(3):xend(3)),& - data_count, hdf_error, file_space_id = filespace,& - mem_space_id = memspace, xfer_prp = plist_id) - call h5pclose_f(plist_id, hdf_error) - call h5dclose_f(dset_temp, hdf_error) - - call h5sclose_f(filespace, hdf_error) - call h5sclose_f(memspace, hdf_error) - - if (salinity) then - !! Repeat on refined grid to save salinity - ! Select wall plane index for refined grid - ic = 1 - - ! Define dimension - ndims = 2 - dims(1) = nymr - dims(2) = nzmr - - data_count(1) = xendr(2) - xstartr(2) + 1 - data_count(2) = xendr(3) - xstartr(3) + 1 - - data_offset(1) = xstartr(2) - 1 - data_offset(2) = xstartr(3) - 1 - - ! Create dataspace (file & memory) - call h5screate_simple_f(ndims, dims, filespace, hdf_error) - call h5screate_simple_f(ndims, data_count, memspace, hdf_error) - - ! Create and write salinity data - dsetname = trim("sal/"//frame) - call h5lexists_f(file_id, dsetname, dsetexists, hdf_error) - if (dsetexists) call h5ldelete_f(file_id, dsetname, hdf_error) - call h5dcreate_f(file_id, dsetname, H5T_NATIVE_DOUBLE, filespace, dset_sal, hdf_error) - call h5dget_space_f(dset_sal, filespace, hdf_error) - call h5sselect_hyperslab_f (filespace, H5S_SELECT_SET_F,data_offset, data_count, hdf_error) - call h5pcreate_f(H5P_DATASET_XFER_F, plist_id, hdf_error) - call h5pset_dxpl_mpio_f(plist_id, H5FD_MPIO_COLLECTIVE_F, hdf_error) - call h5dwrite_f(& - dset_sal, H5T_NATIVE_DOUBLE,& - sal(ic,xstartr(2):xendr(2),xstartr(3):xendr(3)),& - data_count, hdf_error, file_space_id = filespace,& - mem_space_id = memspace, xfer_prp = plist_id) - call h5pclose_f(plist_id, hdf_error) - call h5dclose_f(dset_sal, hdf_error) - - call h5sclose_f(memspace, hdf_error) - call h5sclose_f(filespace, hdf_error) - end if - - call h5fclose_f(file_id, hdf_error) - -end subroutine Mkmov_xcut \ No newline at end of file diff --git a/src/flow_solver/MakeMovieYCut.F90 b/src/flow_solver/MakeMovieYCut.F90 deleted file mode 100644 index 8808dc1e..00000000 --- a/src/flow_solver/MakeMovieYCut.F90 +++ /dev/null @@ -1,281 +0,0 @@ -!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! -! ! -! FILE: MakeMovieYCut.F90 ! -! CONTAINS: subroutine Mkmov_ycut ! -! ! -! PURPOSE: Write down ycut snapshot for flow vars ! -! ! -!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! - -subroutine Mkmov_ycut - use param - use mpih - use hdf5 - use decomp_2d, only: xstart,xend,xstartr,xendr,DECOMP_2D_COMM_CART_X - use local_arrays, only: vz,vy,vx,temp - use mgrd_arrays, only: sal,phi - implicit none - character(70) :: filename - character(30) :: dsetname - character( 5) :: frame - - integer(HID_T) :: file_id, group_id, plist_id - integer(HID_T) :: filespace, slabspace, memspace - integer(HID_T) :: dset_vx, dset_vy, dset_vz, dset_temp, dset_sal - integer(HSIZE_T), dimension(2) :: dims, data_count, data_offset - integer :: ic, hdf_error, ndims, comm, info, ierror - - logical :: fexist, dsetexists - - ! Select mid-plane - ic = nym/2 + 1 - - ! Record frame number and filename as strings - write(frame,"(i5.5)")nint(time/tframe) - filename = trim("outputdir/flowmov/movie_ycut.h5") - - ! MPI - call MPI_CART_SUB(DECOMP_2D_COMM_CART_X, (/.false.,.true./), comm, ierror) - info = MPI_INFO_NULL - - if (ic.le.xend(2) .and. ic.ge.xstart(2)) then - - ! Begin HDF5 - call h5pcreate_f(H5P_FILE_ACCESS_F, plist_id, hdf_error) - call h5pset_fapl_mpio_f(plist_id, comm, info, hdf_error) - - ! Check if movie file already exists, and create it if not - inquire(file=filename,exist=fexist) - if (fexist) then - call h5fopen_f(filename, H5F_ACC_RDWR_F, file_id, hdf_error, access_prp=plist_id) - else - call h5fcreate_f(filename, H5F_ACC_TRUNC_F, file_id, hdf_error,access_prp=plist_id) - call h5gcreate_f(file_id, "vx", group_id, hdf_error) - call h5gclose_f(group_id, hdf_error) - call h5gcreate_f(file_id, "vy", group_id, hdf_error) - call h5gclose_f(group_id, hdf_error) - call h5gcreate_f(file_id, "vz", group_id, hdf_error) - call h5gclose_f(group_id, hdf_error) - call h5gcreate_f(file_id, "temp", group_id, hdf_error) - call h5gclose_f(group_id, hdf_error) - if (salinity) then - call h5gcreate_f(file_id, "sal", group_id, hdf_error) - call h5gclose_f(group_id, hdf_error) - end if - if (phasefield) then - call h5gcreate_f(file_id, "phi", group_id, hdf_error) - call h5gclose_f(group_id, hdf_error) - end if - end if - - call h5pclose_f(plist_id, hdf_error) - - ! Define dimension - ndims = 2 - dims(1) = nx - dims(2) = nzm - - data_count(1) = nx - data_count(2) = xend(3) - xstart(3) + 1 - - data_offset(1) = 0 - data_offset(2) = xstart(3) - 1 - - ! Create dataspace (file & memory) - call h5screate_simple_f(ndims, dims, filespace, hdf_error) - call h5screate_simple_f(ndims, data_count, memspace, hdf_error) - - !! Create datasets and write data to file - ! vx - dsetname = trim("vx/"//frame) - call h5lexists_f(file_id, dsetname, dsetexists, hdf_error) - if (dsetexists) call h5ldelete_f(file_id, dsetname, hdf_error) - call h5dcreate_f(file_id, dsetname, H5T_NATIVE_DOUBLE, filespace, dset_vx, hdf_error) - call h5dget_space_f(dset_vx, filespace, hdf_error) - call h5sselect_hyperslab_f (filespace, H5S_SELECT_SET_F,data_offset, data_count, hdf_error) - call h5pcreate_f(H5P_DATASET_XFER_F, plist_id, hdf_error) - call h5pset_dxpl_mpio_f(plist_id, H5FD_MPIO_COLLECTIVE_F, hdf_error) - call h5dwrite_f(& - dset_vx, H5T_NATIVE_DOUBLE,& - vx(1:nx,ic,xstart(3):xend(3)),& - data_count, hdf_error, file_space_id = filespace,& - mem_space_id = memspace, xfer_prp = plist_id) - call h5pclose_f(plist_id, hdf_error) - call h5dclose_f(dset_vx, hdf_error) - - call h5sclose_f(filespace, hdf_error) - call h5sclose_f(memspace, hdf_error) - - ! Redefine dataspaces to match size of staggered grid arrays in x - dims(1) = nxm - data_count(1) = nxm - - ! Create dataspace (file & memory) - call h5screate_simple_f(ndims, dims, filespace, hdf_error) - call h5screate_simple_f(ndims, data_count, memspace, hdf_error) - - ! vy - dsetname = trim("vy/"//frame) - call h5lexists_f(file_id, dsetname, dsetexists, hdf_error) - if (dsetexists) call h5ldelete_f(file_id, dsetname, hdf_error) - call h5dcreate_f(file_id, dsetname, H5T_NATIVE_DOUBLE, filespace, dset_vy, hdf_error) - call h5dget_space_f(dset_vy, filespace, hdf_error) - call h5sselect_hyperslab_f (filespace, H5S_SELECT_SET_F,data_offset, data_count, hdf_error) - call h5pcreate_f(H5P_DATASET_XFER_F, plist_id, hdf_error) - call h5pset_dxpl_mpio_f(plist_id, H5FD_MPIO_COLLECTIVE_F, hdf_error) - call h5dwrite_f(& - dset_vy, H5T_NATIVE_DOUBLE,& - vy(1:nxm,ic,xstart(3):xend(3)),& - data_count, hdf_error, file_space_id = filespace,& - mem_space_id = memspace, xfer_prp = plist_id) - call h5pclose_f(plist_id, hdf_error) - call h5dclose_f(dset_vy, hdf_error) - - ! vz - dsetname = trim("vz/"//frame) - call h5lexists_f(file_id, dsetname, dsetexists, hdf_error) - if (dsetexists) call h5ldelete_f(file_id, dsetname, hdf_error) - call h5dcreate_f(file_id, dsetname, H5T_NATIVE_DOUBLE, filespace, dset_vz, hdf_error) - call h5dget_space_f(dset_vz, filespace, hdf_error) - call h5sselect_hyperslab_f (filespace, H5S_SELECT_SET_F,data_offset, data_count, hdf_error) - call h5pcreate_f(H5P_DATASET_XFER_F, plist_id, hdf_error) - call h5pset_dxpl_mpio_f(plist_id, H5FD_MPIO_COLLECTIVE_F, hdf_error) - call h5dwrite_f(& - dset_vz, H5T_NATIVE_DOUBLE,& - vz(1:nxm,ic,xstart(3):xend(3)),& - data_count, hdf_error, file_space_id = filespace,& - mem_space_id = memspace, xfer_prp = plist_id) - call h5pclose_f(plist_id, hdf_error) - call h5dclose_f(dset_vz, hdf_error) - - ! temp - dsetname = trim("temp/"//frame) - call h5lexists_f(file_id, dsetname, dsetexists, hdf_error) - if (dsetexists) call h5ldelete_f(file_id, dsetname, hdf_error) - call h5dcreate_f(file_id, dsetname, H5T_NATIVE_DOUBLE, filespace, dset_temp, hdf_error) - call h5dget_space_f(dset_temp, filespace, hdf_error) - call h5sselect_hyperslab_f (filespace, H5S_SELECT_SET_F,data_offset, data_count, hdf_error) - call h5pcreate_f(H5P_DATASET_XFER_F, plist_id, hdf_error) - call h5pset_dxpl_mpio_f(plist_id, H5FD_MPIO_COLLECTIVE_F, hdf_error) - call h5dwrite_f(& - dset_temp, H5T_NATIVE_DOUBLE,& - temp(1:nxm,ic,xstart(3):xend(3)),& - data_count, hdf_error, file_space_id = filespace,& - mem_space_id = memspace, xfer_prp = plist_id) - call h5pclose_f(plist_id, hdf_error) - call h5dclose_f(dset_temp, hdf_error) - - call h5sclose_f(filespace, hdf_error) - call h5sclose_f(memspace, hdf_error) - - call h5fclose_f(file_id, hdf_error) - - end if - - if (salinity) then - !! Repeat on refined grid to save salinity - ! Select midplane index for refined grid - ic = nymr/2 + 1 - - if (ic.le.xendr(2) .and. ic.ge.xstartr(2)) then - - ! Begin HDF5 - call h5pcreate_f(H5P_FILE_ACCESS_F, plist_id, hdf_error) - call h5pset_fapl_mpio_f(plist_id, comm, info, hdf_error) - call h5fopen_f(filename, H5F_ACC_RDWR_F, file_id, hdf_error, access_prp=plist_id) - call h5pclose_f(plist_id, hdf_error) - - ! Define dimension - ndims = 2 - dims(1) = nxmr - dims(2) = nzmr - - data_count(1) = nxmr - data_count(2) = xendr(3) - xstartr(3) + 1 - - data_offset(1) = 0 - data_offset(2) = xstartr(3) - 1 - - ! Create dataspace (file & memory) - call h5screate_simple_f(ndims, dims, filespace, hdf_error) - call h5screate_simple_f(ndims, data_count, memspace, hdf_error) - - ! Create and write salinity data - dsetname = trim("sal/"//frame) - call h5lexists_f(file_id, dsetname, dsetexists, hdf_error) - if (dsetexists) call h5ldelete_f(file_id, dsetname, hdf_error) - call h5dcreate_f(file_id, dsetname, H5T_NATIVE_DOUBLE, filespace, dset_sal, hdf_error) - call h5dget_space_f(dset_sal, filespace, hdf_error) - call h5sselect_hyperslab_f (filespace, H5S_SELECT_SET_F,data_offset, data_count, hdf_error) - call h5pcreate_f(H5P_DATASET_XFER_F, plist_id, hdf_error) - call h5pset_dxpl_mpio_f(plist_id, H5FD_MPIO_COLLECTIVE_F, hdf_error) - call h5dwrite_f(& - dset_sal, H5T_NATIVE_DOUBLE,& - sal(1:nxmr,ic,xstartr(3):xendr(3)),& - data_count, hdf_error, file_space_id = filespace,& - mem_space_id = memspace, xfer_prp = plist_id) - call h5pclose_f(plist_id, hdf_error) - call h5dclose_f(dset_sal, hdf_error) - - call h5sclose_f(filespace, hdf_error) - call h5sclose_f(memspace, hdf_error) - - call h5fclose_f(file_id, hdf_error) - - end if - end if - - if (phasefield) then - !! Repeat on refined grid to save salinity - ! Select midplane index for refined grid - ic = nymr/2 + 1 - - if (ic.le.xendr(2) .and. ic.ge.xstartr(2)) then - - ! Begin HDF5 - call h5pcreate_f(H5P_FILE_ACCESS_F, plist_id, hdf_error) - call h5pset_fapl_mpio_f(plist_id, comm, info, hdf_error) - call h5fopen_f(filename, H5F_ACC_RDWR_F, file_id, hdf_error, access_prp=plist_id) - call h5pclose_f(plist_id, hdf_error) - - ! Define dimension - ndims = 2 - dims(1) = nxmr - dims(2) = nzmr - - data_count(1) = nxmr - data_count(2) = xendr(3) - xstartr(3) + 1 - - data_offset(1) = 0 - data_offset(2) = xstartr(3) - 1 - - ! Create dataspace (file & memory) - call h5screate_simple_f(ndims, dims, filespace, hdf_error) - call h5screate_simple_f(ndims, data_count, memspace, hdf_error) - - ! Create and write salinity data - dsetname = trim("phi/"//frame) - call h5lexists_f(file_id, dsetname, dsetexists, hdf_error) - if (dsetexists) call h5ldelete_f(file_id, dsetname, hdf_error) - call h5dcreate_f(file_id, dsetname, H5T_NATIVE_DOUBLE, filespace, dset_sal, hdf_error) - call h5dget_space_f(dset_sal, filespace, hdf_error) - call h5sselect_hyperslab_f (filespace, H5S_SELECT_SET_F,data_offset, data_count, hdf_error) - call h5pcreate_f(H5P_DATASET_XFER_F, plist_id, hdf_error) - call h5pset_dxpl_mpio_f(plist_id, H5FD_MPIO_COLLECTIVE_F, hdf_error) - call h5dwrite_f(& - dset_sal, H5T_NATIVE_DOUBLE,& - phi(1:nxmr,ic,xstartr(3):xendr(3)),& - data_count, hdf_error, file_space_id = filespace,& - mem_space_id = memspace, xfer_prp = plist_id) - call h5pclose_f(plist_id, hdf_error) - call h5dclose_f(dset_sal, hdf_error) - - call h5sclose_f(filespace, hdf_error) - call h5sclose_f(memspace, hdf_error) - - call h5fclose_f(file_id, hdf_error) - - end if - end if - -end subroutine Mkmov_ycut \ No newline at end of file diff --git a/src/flow_solver/MakeMovieZCut.F90 b/src/flow_solver/MakeMovieZCut.F90 deleted file mode 100644 index b062fb6f..00000000 --- a/src/flow_solver/MakeMovieZCut.F90 +++ /dev/null @@ -1,370 +0,0 @@ -!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! -! ! -! FILE: MakeMovieZCut.F90 ! -! CONTAINS: subroutine Mkmov_zcut ! -! ! -! PURPOSE: Write down zcut snapshot for flow vars ! -! ! -!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! - -subroutine Mkmov_zcut - use param - use mpih - use hdf5 - use decomp_2d, only: xstart,xend,xstartr,xendr,DECOMP_2D_COMM_CART_X - use local_arrays, only: vz,vy,vx,temp - use mgrd_arrays, only: sal, phi, phic, tempr!, vxr, vyr - implicit none - character(70) :: filename - character(30) :: dsetname - character( 5) :: frame - - integer(HID_T) :: file_id, group_id, plist_id - integer(HID_T) :: filespace, slabspace, memspace - integer(HID_T) :: dset_vx, dset_vy, dset_vz, dset_temp, dset_sal - integer(HSIZE_T), dimension(2) :: dims, data_count, data_offset - integer :: ic, hdf_error, ndims, comm, info, ierror - - logical :: fexist, dsetexists - - ! Select mid-plane - ic = nzm/2 + 1 - - ! Record frame number and filename as strings - write(frame,"(i5.5)")nint(time/tframe) - filename = trim("outputdir/flowmov/movie_zcut.h5") - - ! MPI - call MPI_CART_SUB(DECOMP_2D_COMM_CART_X, (/.true.,.false./), comm, ierror) - info = MPI_INFO_NULL - - if(ic.le.xend(3) .and. ic.ge.xstart(3)) then - - ! Begin HDF5 - call h5pcreate_f(H5P_FILE_ACCESS_F, plist_id, hdf_error) - call h5pset_fapl_mpio_f(plist_id, comm, info, hdf_error) - - ! Check if movie file already exists, and create it if not - inquire(file=filename,exist=fexist) - if (fexist) then - call h5fopen_f(filename, H5F_ACC_RDWR_F, file_id, hdf_error, access_prp=plist_id) - else - call h5fcreate_f(filename, H5F_ACC_TRUNC_F, file_id, hdf_error,access_prp=plist_id) - call h5gcreate_f(file_id, "vx", group_id, hdf_error) - call h5gclose_f(group_id, hdf_error) - call h5gcreate_f(file_id, "vy", group_id, hdf_error) - call h5gclose_f(group_id, hdf_error) - call h5gcreate_f(file_id, "vz", group_id, hdf_error) - call h5gclose_f(group_id, hdf_error) - call h5gcreate_f(file_id, "temp", group_id, hdf_error) - call h5gclose_f(group_id, hdf_error) - if (salinity) then - call h5gcreate_f(file_id, "sal", group_id, hdf_error) - call h5gclose_f(group_id, hdf_error) - ! call h5gcreate_f(file_id, "vyr", group_id, hdf_error) - ! call h5gclose_f(group_id, hdf_error) - ! call h5gcreate_f(file_id, "vxr", group_id, hdf_error) - ! call h5gclose_f(group_id, hdf_error) - end if - if (phasefield) then - call h5gcreate_f(file_id, "phi", group_id, hdf_error) - call h5gclose_f(group_id, hdf_error) - ! call h5gcreate_f(file_id, "phic", group_id, hdf_error) - ! call h5gclose_f(group_id, hdf_error) - ! call h5gcreate_f(file_id, "tempr", group_id, hdf_error) - ! call h5gclose_f(group_id, hdf_error) - end if - end if - - call h5pclose_f(plist_id, hdf_error) - - ! Define dimension - ndims = 2 - dims(1) = nx - dims(2) = nym - - data_count(1) = nx - data_count(2) = xend(2) - xstart(2) + 1 - - data_offset(1) = 0 - data_offset(2) = xstart(2) - 1 - - ! Create dataspace (file & memory) - call h5screate_simple_f(ndims, dims, filespace, hdf_error) - call h5screate_simple_f(ndims, data_count, memspace, hdf_error) - - !! Create datasets and write data to file - ! vx - dsetname = trim("vx/"//frame) - call h5lexists_f(file_id, dsetname, dsetexists, hdf_error) - if (dsetexists) call h5ldelete_f(file_id, dsetname, hdf_error) - call h5dcreate_f(file_id, dsetname, H5T_NATIVE_DOUBLE, filespace, dset_vx, hdf_error) - call h5dget_space_f(dset_vx, filespace, hdf_error) - call h5sselect_hyperslab_f (filespace, H5S_SELECT_SET_F,data_offset, data_count, hdf_error) - call h5pcreate_f(H5P_DATASET_XFER_F, plist_id, hdf_error) - call h5pset_dxpl_mpio_f(plist_id, H5FD_MPIO_COLLECTIVE_F, hdf_error) - call h5dwrite_f(& - dset_vx, H5T_NATIVE_DOUBLE,& - vx(1:nx,xstart(2):xend(2),ic),& - data_count, hdf_error, file_space_id = filespace,& - mem_space_id = memspace, xfer_prp = plist_id) - call h5pclose_f(plist_id, hdf_error) - call h5dclose_f(dset_vx, hdf_error) - - call h5sclose_f(filespace, hdf_error) - call h5sclose_f(memspace, hdf_error) - - ! Redefine dataspaces to match size of staggered grid arrays in x - dims(1) = nxm - data_count(1) = nxm - - ! Create dataspace (file & memory) - call h5screate_simple_f(ndims, dims, filespace, hdf_error) - call h5screate_simple_f(ndims, data_count, memspace, hdf_error) - - ! vy - dsetname = trim("vy/"//frame) - call h5lexists_f(file_id, dsetname, dsetexists, hdf_error) - if (dsetexists) call h5ldelete_f(file_id, dsetname, hdf_error) - call h5dcreate_f(file_id, dsetname, H5T_NATIVE_DOUBLE, filespace, dset_vy, hdf_error) - call h5dget_space_f(dset_vy, filespace, hdf_error) - call h5sselect_hyperslab_f (filespace, H5S_SELECT_SET_F,data_offset, data_count, hdf_error) - call h5pcreate_f(H5P_DATASET_XFER_F, plist_id, hdf_error) - call h5pset_dxpl_mpio_f(plist_id, H5FD_MPIO_COLLECTIVE_F, hdf_error) - call h5dwrite_f(& - dset_vy, H5T_NATIVE_DOUBLE,& - vy(1:nxm,xstart(2):xend(2),ic),& - data_count, hdf_error, file_space_id = filespace,& - mem_space_id = memspace, xfer_prp = plist_id) - call h5pclose_f(plist_id, hdf_error) - call h5dclose_f(dset_vy, hdf_error) - - ! vz - dsetname = trim("vz/"//frame) - call h5lexists_f(file_id, dsetname, dsetexists, hdf_error) - if (dsetexists) call h5ldelete_f(file_id, dsetname, hdf_error) - call h5dcreate_f(file_id, dsetname, H5T_NATIVE_DOUBLE, filespace, dset_vz, hdf_error) - call h5dget_space_f(dset_vz, filespace, hdf_error) - call h5sselect_hyperslab_f (filespace, H5S_SELECT_SET_F,data_offset, data_count, hdf_error) - call h5pcreate_f(H5P_DATASET_XFER_F, plist_id, hdf_error) - call h5pset_dxpl_mpio_f(plist_id, H5FD_MPIO_COLLECTIVE_F, hdf_error) - call h5dwrite_f(& - dset_vz, H5T_NATIVE_DOUBLE,& - vz(1:nxm,xstart(2):xend(2),ic),& - data_count, hdf_error, file_space_id = filespace,& - mem_space_id = memspace, xfer_prp = plist_id) - call h5pclose_f(plist_id, hdf_error) - call h5dclose_f(dset_vz, hdf_error) - - ! temp - dsetname = trim("temp/"//frame) - call h5lexists_f(file_id, dsetname, dsetexists, hdf_error) - if (dsetexists) call h5ldelete_f(file_id, dsetname, hdf_error) - call h5dcreate_f(file_id, dsetname, H5T_NATIVE_DOUBLE, filespace, dset_temp, hdf_error) - call h5dget_space_f(dset_temp, filespace, hdf_error) - call h5sselect_hyperslab_f (filespace, H5S_SELECT_SET_F,data_offset, data_count, hdf_error) - call h5pcreate_f(H5P_DATASET_XFER_F, plist_id, hdf_error) - call h5pset_dxpl_mpio_f(plist_id, H5FD_MPIO_COLLECTIVE_F, hdf_error) - call h5dwrite_f(& - dset_temp, H5T_NATIVE_DOUBLE,& - temp(1:nxm,xstart(2):xend(2),ic),& - data_count, hdf_error, file_space_id = filespace,& - mem_space_id = memspace, xfer_prp = plist_id) - call h5pclose_f(plist_id, hdf_error) - call h5dclose_f(dset_temp, hdf_error) - - ! if (phasefield) then - ! ! phic - ! dsetname = trim("phic/"//frame) - ! call h5lexists_f(file_id, dsetname, dsetexists, hdf_error) - ! if (dsetexists) call h5ldelete_f(file_id, dsetname, hdf_error) - ! call h5dcreate_f(file_id, dsetname, H5T_NATIVE_DOUBLE, filespace, dset_temp, hdf_error) - ! call h5dget_space_f(dset_temp, filespace, hdf_error) - ! call h5sselect_hyperslab_f (filespace, H5S_SELECT_SET_F,data_offset, data_count, hdf_error) - ! call h5pcreate_f(H5P_DATASET_XFER_F, plist_id, hdf_error) - ! call h5pset_dxpl_mpio_f(plist_id, H5FD_MPIO_COLLECTIVE_F, hdf_error) - ! call h5dwrite_f(& - ! dset_temp, H5T_NATIVE_DOUBLE,& - ! phic(1:nxm,xstart(2):xend(2),ic),& - ! data_count, hdf_error, file_space_id = filespace,& - ! mem_space_id = memspace, xfer_prp = plist_id) - ! call h5pclose_f(plist_id, hdf_error) - ! call h5dclose_f(dset_temp, hdf_error) - ! end if - - call h5sclose_f(filespace, hdf_error) - call h5sclose_f(memspace, hdf_error) - - call h5fclose_f(file_id, hdf_error) - - end if - - if (salinity) then - !! Repeat on refined grid to save salinity - ! Select midplane index for refined grid - ic = nzmr/2 + 1 - - if(ic.le.xendr(3) .and. ic.ge.xstartr(3)) then - - ! Begin HDF5 - call h5pcreate_f(H5P_FILE_ACCESS_F, plist_id, hdf_error) - call h5pset_fapl_mpio_f(plist_id, comm, info, hdf_error) - call h5fopen_f(filename, H5F_ACC_RDWR_F, file_id, hdf_error, access_prp=plist_id) - call h5pclose_f(plist_id, hdf_error) - - ! Define dimension - ndims = 2 - dims(1) = nxmr - dims(2) = nymr - - data_count(1) = nxmr - data_count(2) = xendr(2) - xstartr(2) + 1 - - data_offset(1) = 0 - data_offset(2) = xstartr(2) - 1 - - ! Create dataspace (file & memory) - call h5screate_simple_f(ndims, dims, filespace, hdf_error) - call h5screate_simple_f(ndims, data_count, memspace, hdf_error) - - ! Create and write salinity data - dsetname = trim("sal/"//frame) - call h5lexists_f(file_id, dsetname, dsetexists, hdf_error) - if (dsetexists) call h5ldelete_f(file_id, dsetname, hdf_error) - call h5dcreate_f(file_id, dsetname, H5T_NATIVE_DOUBLE, filespace, dset_sal, hdf_error) - call h5dget_space_f(dset_sal, filespace, hdf_error) - call h5sselect_hyperslab_f (filespace, H5S_SELECT_SET_F,data_offset, data_count, hdf_error) - call h5pcreate_f(H5P_DATASET_XFER_F, plist_id, hdf_error) - call h5pset_dxpl_mpio_f(plist_id, H5FD_MPIO_COLLECTIVE_F, hdf_error) - call h5dwrite_f(& - dset_sal, H5T_NATIVE_DOUBLE,& - sal(1:nxmr,xstartr(2):xendr(2),ic),& - data_count, hdf_error, file_space_id = filespace,& - mem_space_id = memspace, xfer_prp = plist_id) - call h5pclose_f(plist_id, hdf_error) - call h5dclose_f(dset_sal, hdf_error) - - ! Create and write refined velocity data - ! dsetname = trim("vyr/"//frame) - ! call h5lexists_f(file_id, dsetname, dsetexists, hdf_error) - ! if (dsetexists) call h5ldelete_f(file_id, dsetname, hdf_error) - ! call h5dcreate_f(file_id, dsetname, H5T_NATIVE_DOUBLE, filespace, dset_sal, hdf_error) - ! call h5dget_space_f(dset_sal, filespace, hdf_error) - ! call h5sselect_hyperslab_f (filespace, H5S_SELECT_SET_F,data_offset, data_count, hdf_error) - ! call h5pcreate_f(H5P_DATASET_XFER_F, plist_id, hdf_error) - ! call h5pset_dxpl_mpio_f(plist_id, H5FD_MPIO_COLLECTIVE_F, hdf_error) - ! call h5dwrite_f(& - ! dset_sal, H5T_NATIVE_DOUBLE,& - ! vyr(1:nxmr,xstartr(2):xendr(2),ic),& - ! data_count, hdf_error, file_space_id = filespace,& - ! mem_space_id = memspace, xfer_prp = plist_id) - ! call h5pclose_f(plist_id, hdf_error) - ! call h5dclose_f(dset_sal, hdf_error) - - ! call h5sclose_f(memspace, hdf_error) - ! call h5sclose_f(filespace, hdf_error) - - ! ! Redefine dataspaces to match size of vxr arrays in x - ! dims(1) = nxr - ! data_count(1) = nxr - - ! ! Create dataspace (file & memory) - ! call h5screate_simple_f(ndims, dims, filespace, hdf_error) - ! call h5screate_simple_f(ndims, data_count, memspace, hdf_error) - - ! ! Create and write refined velocity data - ! dsetname = trim("vxr/"//frame) - ! call h5lexists_f(file_id, dsetname, dsetexists, hdf_error) - ! if (dsetexists) call h5ldelete_f(file_id, dsetname, hdf_error) - ! call h5dcreate_f(file_id, dsetname, H5T_NATIVE_DOUBLE, filespace, dset_sal, hdf_error) - ! call h5dget_space_f(dset_sal, filespace, hdf_error) - ! call h5sselect_hyperslab_f (filespace, H5S_SELECT_SET_F,data_offset, data_count, hdf_error) - ! call h5pcreate_f(H5P_DATASET_XFER_F, plist_id, hdf_error) - ! call h5pset_dxpl_mpio_f(plist_id, H5FD_MPIO_COLLECTIVE_F, hdf_error) - ! call h5dwrite_f(& - ! dset_sal, H5T_NATIVE_DOUBLE,& - ! vxr(1:nxr,xstartr(2):xendr(2),ic),& - ! data_count, hdf_error, file_space_id = filespace,& - ! mem_space_id = memspace, xfer_prp = plist_id) - ! call h5pclose_f(plist_id, hdf_error) - ! call h5dclose_f(dset_sal, hdf_error) - - call h5sclose_f(memspace, hdf_error) - call h5sclose_f(filespace, hdf_error) - - call h5fclose_f(file_id, hdf_error) - - end if - end if - - if (phasefield) then - !! Repeat on refined grid to save phi - ! Select midplane index for refined grid - ic = nzmr/2 + 1 - - if(ic.le.xendr(3) .and. ic.ge.xstartr(3)) then - - ! Begin HDF5 - call h5pcreate_f(H5P_FILE_ACCESS_F, plist_id, hdf_error) - call h5pset_fapl_mpio_f(plist_id, comm, info, hdf_error) - call h5fopen_f(filename, H5F_ACC_RDWR_F, file_id, hdf_error, access_prp=plist_id) - call h5pclose_f(plist_id, hdf_error) - - ! Define dimension - ndims = 2 - dims(1) = nxmr - dims(2) = nymr - - data_count(1) = nxmr - data_count(2) = xendr(2) - xstartr(2) + 1 - - data_offset(1) = 0 - data_offset(2) = xstartr(2) - 1 - - ! Create dataspace (file & memory) - call h5screate_simple_f(ndims, dims, filespace, hdf_error) - call h5screate_simple_f(ndims, data_count, memspace, hdf_error) - - ! Create and write phase-field data - dsetname = trim("phi/"//frame) - call h5lexists_f(file_id, dsetname, dsetexists, hdf_error) - if (dsetexists) call h5ldelete_f(file_id, dsetname, hdf_error) - call h5dcreate_f(file_id, dsetname, H5T_NATIVE_DOUBLE, filespace, dset_sal, hdf_error) - call h5dget_space_f(dset_sal, filespace, hdf_error) - call h5sselect_hyperslab_f (filespace, H5S_SELECT_SET_F,data_offset, data_count, hdf_error) - call h5pcreate_f(H5P_DATASET_XFER_F, plist_id, hdf_error) - call h5pset_dxpl_mpio_f(plist_id, H5FD_MPIO_COLLECTIVE_F, hdf_error) - call h5dwrite_f(& - dset_sal, H5T_NATIVE_DOUBLE,& - phi(1:nxmr,xstartr(2):xendr(2),ic),& - data_count, hdf_error, file_space_id = filespace,& - mem_space_id = memspace, xfer_prp = plist_id) - call h5pclose_f(plist_id, hdf_error) - call h5dclose_f(dset_sal, hdf_error) - - ! ! Create and write phase-field data - ! dsetname = trim("tempr/"//frame) - ! call h5lexists_f(file_id, dsetname, dsetexists, hdf_error) - ! if (dsetexists) call h5ldelete_f(file_id, dsetname, hdf_error) - ! call h5dcreate_f(file_id, dsetname, H5T_NATIVE_DOUBLE, filespace, dset_sal, hdf_error) - ! call h5dget_space_f(dset_sal, filespace, hdf_error) - ! call h5sselect_hyperslab_f (filespace, H5S_SELECT_SET_F,data_offset, data_count, hdf_error) - ! call h5pcreate_f(H5P_DATASET_XFER_F, plist_id, hdf_error) - ! call h5pset_dxpl_mpio_f(plist_id, H5FD_MPIO_COLLECTIVE_F, hdf_error) - ! call h5dwrite_f(& - ! dset_sal, H5T_NATIVE_DOUBLE,& - ! tempr(1:nxmr,xstartr(2):xendr(2),ic),& - ! data_count, hdf_error, file_space_id = filespace,& - ! mem_space_id = memspace, xfer_prp = plist_id) - ! call h5pclose_f(plist_id, hdf_error) - ! call h5dclose_f(dset_sal, hdf_error) - - call h5sclose_f(memspace, hdf_error) - call h5sclose_f(filespace, hdf_error) - - call h5fclose_f(file_id, hdf_error) - - end if - end if - -end subroutine Mkmov_zcut \ No newline at end of file diff --git a/src/flow_solver/HdfReadContinua.F90 b/src/h5tools/HdfReadContinua.F90 similarity index 100% rename from src/flow_solver/HdfReadContinua.F90 rename to src/h5tools/HdfReadContinua.F90 diff --git a/src/flow_solver/HdfRoutines.F90 b/src/h5tools/HdfRoutines.F90 similarity index 100% rename from src/flow_solver/HdfRoutines.F90 rename to src/h5tools/HdfRoutines.F90 diff --git a/src/h5tools/MakeMovieXCut.F90 b/src/h5tools/MakeMovieXCut.F90 new file mode 100644 index 00000000..c474ad9c --- /dev/null +++ b/src/h5tools/MakeMovieXCut.F90 @@ -0,0 +1,70 @@ +!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! +! ! +! FILE: MakeMovieXCut.F90 ! +! CONTAINS: subroutine Mkmov_xcut ! +! ! +! PURPOSE: Write down xcut snapshot for flow vars ! +! ! +!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! + +subroutine Mkmov_xcut + use mpih + use hdf5 + use decomp_2d, only: xstart,xend,xstartr,xendr,DECOMP_2D_COMM_CART_X + use local_arrays, only: vz,vy,vx,temp + use mgrd_arrays, only: sal,phi + use h5_tools + implicit none + character(70) :: filename + character(4) :: varname + + integer(HID_T) :: file_id + integer :: ic, comm + + logical :: fexist + + ! Select plane - plane next to lower wall + ic = 1 + + ! Record filename as string + filename = trim("outputdir/flowmov/movie_xcut.h5") + + ! MPI + call MPI_CART_SUB(DECOMP_2D_COMM_CART_X, (/.true.,.true./), comm, ierr) + info = MPI_INFO_NULL + + call h5_open_or_create(file_id, filename, comm, fexist) + if (.not. fexist) call h5_add_slice_groups(file_id) + + varname = 'vx' + call write_H5_plane(file_id, varname, vx(max(ic,2), xstart(2):xend(2), xstart(3):xend(3)), 'x') + + varname = 'vy' + call write_H5_plane(file_id, varname, vy(ic, xstart(2):xend(2), xstart(3):xend(3)), 'x') + + varname = 'vz' + call write_H5_plane(file_id, varname, vz(ic, xstart(2):xend(2), xstart(3):xend(3)), 'x') + + varname = 'temp' + call write_H5_plane(file_id, varname, temp(ic, xstart(2):xend(2), xstart(3):xend(3)), 'x') + + call h5fclose_f(file_id, hdf_error) + + if (salinity) then + !! Repeat on refined grid to save salinity + ! Select wall plane index for refined grid + ic = 1 + + call h5_open_or_create(file_id, filename, comm, fexist) + if (.not. fexist) call h5_add_slice_groups(file_id) + + varname = 'sal' + call write_H5_plane(file_id, varname, sal(ic, xstartr(2):xendr(2), xstartr(3):xendr(3)), 'x') + + call h5fclose_f(file_id, hdf_error) + + end if + + call MpiBarrier + +end subroutine Mkmov_xcut \ No newline at end of file diff --git a/src/h5tools/MakeMovieYCut.F90 b/src/h5tools/MakeMovieYCut.F90 new file mode 100644 index 00000000..9ff28e1d --- /dev/null +++ b/src/h5tools/MakeMovieYCut.F90 @@ -0,0 +1,96 @@ +!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! +! ! +! FILE: MakeMovieYCut.F90 ! +! CONTAINS: subroutine Mkmov_ycut ! +! ! +! PURPOSE: Write down ycut snapshot for flow vars ! +! ! +!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! + +subroutine Mkmov_ycut + use param, only: nym, nymr + use mpih + use hdf5 + use decomp_2d, only: xstart,xend,xstartr,xendr,DECOMP_2D_COMM_CART_X + use local_arrays, only: vz,vy,vx,temp + use mgrd_arrays, only: sal,phi + use h5_tools + implicit none + character(70) :: filename + character(4) :: varname + + integer(HID_T) :: file_id + integer :: ic, comm + + logical :: fexist + + ! Select mid-plane + ic = nym/2 + 1 + + ! Record filename as string + filename = trim("outputdir/flowmov/movie_ycut.h5") + + ! MPI + call MPI_CART_SUB(DECOMP_2D_COMM_CART_X, (/.false.,.true./), comm, ierr) + info = MPI_INFO_NULL + + if (ic.le.xend(2) .and. ic.ge.xstart(2)) then + + call h5_open_or_create(file_id, filename, comm, fexist) + if (.not. fexist) call h5_add_slice_groups(file_id) + + varname = 'vx' + call write_H5_plane(file_id, varname, vx(1:nx, ic, xstart(3):xend(3)), 'y') + + varname = 'vy' + call write_H5_plane(file_id, varname, vy(1:nxm, ic, xstart(3):xend(3)), 'y') + + varname = 'vz' + call write_H5_plane(file_id, varname, vz(1:nxm, ic, xstart(3):xend(3)), 'y') + + varname = 'temp' + call write_H5_plane(file_id, varname, temp(1:nxm, ic, xstart(3):xend(3)), 'y') + + call h5fclose_f(file_id, hdf_error) + + end if + + if (salinity) then + !! Repeat on refined grid to save salinity + ! Select midplane index for refined grid + ic = nymr/2 + 1 + + if (ic.le.xendr(2) .and. ic.ge.xstartr(2)) then + + call h5_open_or_create(file_id, filename, comm, fexist) + if (.not. fexist) call h5_add_slice_groups(file_id) + + varname = 'sal' + call write_H5_plane(file_id, varname, sal(1:nxmr, ic, xstartr(3):xendr(3)), 'y') + + call h5fclose_f(file_id, hdf_error) + + end if + end if + + if (phasefield) then + !! Repeat on refined grid to save salinity + ! Select midplane index for refined grid + ic = nymr/2 + 1 + + if (ic.le.xendr(2) .and. ic.ge.xstartr(2)) then + + call h5_open_or_create(file_id, filename, comm, fexist) + if (.not. fexist) call h5_add_slice_groups(file_id) + + varname = 'phi' + call write_H5_plane(file_id, varname, phi(1:nxmr, ic, xstartr(3):xendr(3)), 'y') + + call h5fclose_f(file_id, hdf_error) + + end if + end if + + call MpiBarrier + +end subroutine Mkmov_ycut \ No newline at end of file diff --git a/src/h5tools/MakeMovieZCut.F90 b/src/h5tools/MakeMovieZCut.F90 new file mode 100644 index 00000000..6241befe --- /dev/null +++ b/src/h5tools/MakeMovieZCut.F90 @@ -0,0 +1,96 @@ +!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! +! ! +! FILE: MakeMovieZCut.F90 ! +! CONTAINS: subroutine Mkmov_zcut ! +! ! +! PURPOSE: Write down zcut snapshot for flow vars ! +! ! +!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! + +subroutine Mkmov_zcut + use param, only: nzm, nzmr + use mpih + use hdf5 + use decomp_2d, only: xstart,xend,xstartr,xendr,DECOMP_2D_COMM_CART_X + use local_arrays, only: vz,vy,vx,temp, pr + use mgrd_arrays, only: sal, phi, phic, tempr + use h5_tools + implicit none + character(70) :: filename + character(4) :: varname + + integer(HID_T) :: file_id + integer :: ic, comm + + logical :: fexist + + ! Select mid-plane + ic = nzm/2 + 1 + + ! Record filename as string + filename = trim("outputdir/flowmov/movie_zcut.h5") + + ! MPI + call MPI_CART_SUB(DECOMP_2D_COMM_CART_X, (/.true.,.false./), comm, ierr) + info = MPI_INFO_NULL + + if(ic.le.xend(3) .and. ic.ge.xstart(3)) then + + call h5_open_or_create(file_id, filename, comm, fexist) + if (.not. fexist) call h5_add_slice_groups(file_id) + + varname = 'vx' + call write_H5_plane(file_id, varname, vx(1:nx, xstart(2):xend(2), ic), 'z') + + varname = 'vy' + call write_H5_plane(file_id, varname, vy(1:nxm, xstart(2):xend(2), ic), 'z') + + varname = 'vz' + call write_H5_plane(file_id, varname, vz(1:nxm, xstart(2):xend(2), ic), 'z') + + varname = 'temp' + call write_H5_plane(file_id, varname, temp(1:nxm, xstart(2):xend(2), ic), 'z') + + call h5fclose_f(file_id, hdf_error) + + end if + + if (salinity) then + !! Repeat on refined grid to save salinity + ! Select midplane index for refined grid + ic = nzmr/2 + 1 + + if(ic.le.xendr(3) .and. ic.ge.xstartr(3)) then + + call h5_open_or_create(file_id, filename, comm, fexist) + if (.not. fexist) call h5_add_slice_groups(file_id) + + varname = 'sal' + call write_H5_plane(file_id, varname, sal(1:nxmr, xstartr(2):xendr(2), ic), 'z') + + call h5fclose_f(file_id, hdf_error) + + end if + end if + + if (phasefield) then + !! Repeat on refined grid to save phi + ! Select midplane index for refined grid + ic = nzmr/2 + 1 + + if(ic.le.xendr(3) .and. ic.ge.xstartr(3)) then + + call h5_open_or_create(file_id, filename, comm, fexist) + if (.not. fexist) call h5_add_slice_groups(file_id) + + varname = 'phi' + call write_H5_plane(file_id, varname, phi(1:nxmr, xstartr(2):xendr(2), ic), 'z') + + call h5fclose_f(file_id, hdf_error) + + end if + end if + + call MpiBarrier + +end subroutine Mkmov_zcut \ No newline at end of file diff --git a/src/h5tools/h5_tools.F90 b/src/h5tools/h5_tools.F90 new file mode 100644 index 00000000..ed913612 --- /dev/null +++ b/src/h5tools/h5_tools.F90 @@ -0,0 +1,129 @@ +module h5_tools + use HDF5 + use mpih, only: MPI_INFO_NULL + use param, only: salinity, phasefield, nx, nxm, nym, nzm, time, tframe, nxmr, nymr, nzmr + use decomp_2d, only: xstart, xstartr + implicit none + + integer(HID_T) :: group_id, plist_id + integer :: info, hdf_error + +contains + +subroutine h5_open_or_create(file_id, filename, comm, fexist) + integer(HID_T), intent(out) :: file_id + integer, intent(in) :: comm + character(70), intent(in) :: filename + logical, intent(out) :: fexist + + info = MPI_INFO_NULL + + call h5pcreate_f(H5P_FILE_ACCESS_F, plist_id, hdf_error) + call h5pset_fapl_mpio_f(plist_id, comm, info, hdf_error) + + inquire(file=filename, exist=fexist) + if (fexist) then + call h5fopen_f(filename, H5F_ACC_RDWR_F, file_id, hdf_error, access_prp=plist_id) + else + call h5fcreate_f(filename, H5F_ACC_TRUNC_F, file_id, hdf_error, access_prp=plist_id) + end if + + call h5pclose_f(plist_id, hdf_error) +end subroutine + +subroutine h5_add_slice_groups(file_id) + integer(HID_T), intent(inout) :: file_id + + integer :: i + character(len=4), dimension(4) :: var_list + + var_list = [character(len=4) :: "vx", "vy", "vz", "temp"] + + do i=1,size(var_list) + call h5gcreate_f(file_id, trim(var_list(i)), group_id, hdf_error) + call h5gclose_f(group_id, hdf_error) + end do + + if (salinity) then + call h5gcreate_f(file_id, "sal", group_id, hdf_error) + call h5gclose_f(group_id, hdf_error) + end if + + if (phasefield) then + call h5gcreate_f(file_id, "phi", group_id, hdf_error) + call h5gclose_f(group_id, hdf_error) + end if +end subroutine + +subroutine write_H5_plane(file_id, varname, var, axis) + integer(HID_T), intent(inout) :: file_id + character(len=4), intent(in) :: varname + real, intent(in) :: var(:,:) + character, intent(in) :: axis + + integer, parameter :: ndims = 2 + integer(HID_T) :: filespace, slabspace, memspace, dset + integer(HSIZE_T), dimension(2) :: dims, data_count, data_offset + character(len=5) :: frame + character(len=10) :: dsetname + logical :: dexist + + if (axis=='x') then + dims = [nym, nzm] + data_offset = [xstart(2) - 1, xstart(3) - 1] + else + if (axis=='y') then + dims(2) = nzm + data_offset = [0, xstart(3) - 1] + else + dims(2) = nym + data_offset = [0, xstart(2) - 1] + end if + if (trim(varname)=='vx') then + dims(1) = nx + else + dims(1) = nxm + end if + end if + if (index('phisal', trim(varname)) /= 0) then + dims = [nxmr, nzmr] + data_offset = [0, xstartr(3) - 1] + if (axis=='x') then + dims(1) = nymr + data_offset(1) = xstartr(2) - 1 + elseif (axis=='z') then + dims(2) = nymr + data_offset(2) = xstartr(2) - 1 + end if + end if + + data_count = [size(var, dim=1), size(var, dim=2)] + + ! Create dataspaces + call h5screate_simple_f(ndims, dims, filespace, hdf_error) + call h5screate_simple_f(ndims, data_count, memspace, hdf_error) + + + write(frame,"(i5.5)")nint(time/tframe) + dsetname = trim(trim(varname)//"/"//frame) + call h5lexists_f(file_id, dsetname, dexist, hdf_error) + if (dexist) call h5ldelete_f(file_id, dsetname, hdf_error) + call h5dcreate_f(file_id, dsetname, H5T_NATIVE_DOUBLE, filespace, dset, hdf_error) + call h5dget_space_f(dset, filespace, hdf_error) + call h5sselect_hyperslab_f(filespace, H5S_SELECT_SET_F,data_offset, data_count, hdf_error) + call h5pcreate_f(H5P_DATASET_XFER_F, plist_id, hdf_error) + call h5pset_dxpl_mpio_f(plist_id, H5FD_MPIO_COLLECTIVE_F, hdf_error) + call h5dwrite_f(& + dset, H5T_NATIVE_DOUBLE,& + var,& + data_count, hdf_error, file_space_id = filespace,& + mem_space_id = memspace, xfer_prp = plist_id) + call h5pclose_f(plist_id, hdf_error) + call h5dclose_f(dset, hdf_error) + + call h5sclose_f(filespace, hdf_error) + call h5sclose_f(memspace, hdf_error) + +end subroutine + +end module h5_tools \ No newline at end of file diff --git a/src/h5tools/mean_zplane.F90 b/src/h5tools/mean_zplane.F90 new file mode 100644 index 00000000..bb3aecc1 --- /dev/null +++ b/src/h5tools/mean_zplane.F90 @@ -0,0 +1,228 @@ + +subroutine mean_zplane + use param, only: nx, nzm + use mpih + use hdf5 + use decomp_2d, only: xstart,xend,xsize,xstartr,xendr,DECOMP_2D_COMM_CART_X,nrank + use local_arrays, only: vz,vy,vx,temp + use mgrd_arrays, only: sal, phi + use h5_tools + use means + implicit none + + character(70) :: filename + character(4) :: varname + real, allocatable :: uplane(:,:), vplane(:,:), wplane(:,:), Tplane(:,:) + integer :: comm, zrank, i, j, bufsize + integer(HID_T) :: file_id + logical :: fexist + + filename = trim("outputdir/flowmov/movie_zmean.h5") + + allocate(uplane(1:nx ,xstart(2):xend(2))) + allocate(vplane(1:nxm,xstart(2):xend(2))) + allocate(wplane(1:nxm,xstart(2):xend(2))) + allocate(Tplane(1:nxm,xstart(2):xend(2))) + + call MPI_CART_SUB(DECOMP_2D_COMM_CART_X, (/.false., .true./), comm, ierr) + call MPI_COMM_RANK(comm, zrank, ierr) + + call zmean(vx(1:nx ,xstart(2):xend(2),xstart(3):xend(3)), uplane, comm, zrank) + call zmean(vy(1:nxm,xstart(2):xend(2),xstart(3):xend(3)), vplane, comm, zrank) + call zmean(vz(1:nxm,xstart(2):xend(2),xstart(3):xend(3)), wplane, comm, zrank) + call zmean(temp(1:nxm,xstart(2):xend(2),xstart(3):xend(3)), Tplane, comm, zrank) + + call MPI_CART_SUB(DECOMP_2D_COMM_CART_X, (/.true., .false./), comm, ierr) + + if (zrank==0) then + + ! Open the movie file (create if it doesn't exist) + call h5_open_or_create(file_id, filename, comm, fexist) + ! If file only just created, add groups for data storage + if (.not. fexist) call h5_add_slice_groups(file_id) + + varname='vx' + call write_H5_plane(file_id, varname, uplane, 'z') + varname='vy' + call write_H5_plane(file_id, varname, vplane, 'z') + varname='vz' + call write_H5_plane(file_id, varname, wplane, 'z') + varname='temp' + call write_H5_plane(file_id, varname, Tplane, 'z') + + ! Close HDF5 file + call h5fclose_f(file_id, hdf_error) + end if + + if (allocated(uplane)) deallocate(uplane) + if (allocated(vplane)) deallocate(vplane) + if (allocated(wplane)) deallocate(wplane) + if (allocated(Tplane)) deallocate(Tplane) + + if (salinity) then + allocate(Tplane(1:nxmr,xstartr(2):xendr(2))) + + call MPI_CART_SUB(DECOMP_2D_COMM_CART_X, (/.false., .true./), comm, ierr) + call MPI_COMM_RANK(comm, zrank, ierr) + + call zmeanr(sal(1:nxmr,xstartr(2):xendr(2),xstartr(3):xendr(3)), Tplane, comm, zrank) + + call MPI_CART_SUB(DECOMP_2D_COMM_CART_X, (/.true., .false./), comm, ierr) + + if (zrank==0) then + ! Open the movie file (create if it doesn't exist) + call h5_open_or_create(file_id, filename, comm, fexist) + ! If file only just created, add groups for data storage + if (.not. fexist) call h5_add_slice_groups(file_id) + + varname='sal' + call write_H5_plane(file_id, varname, Tplane, 'z') + + ! Close HDF5 file + call h5fclose_f(file_id, hdf_error) + end if + + if (allocated(Tplane)) deallocate(Tplane) + end if + + if (phasefield) then + allocate(Tplane(1:nxmr,xstartr(2):xendr(2))) + + call MPI_CART_SUB(DECOMP_2D_COMM_CART_X, (/.false., .true./), comm, ierr) + call MPI_COMM_RANK(comm, zrank, ierr) + + call zmeanr(phi(1:nxmr,xstartr(2):xendr(2),xstartr(3):xendr(3)), Tplane, comm, zrank) + + call MPI_CART_SUB(DECOMP_2D_COMM_CART_X, (/.true., .false./), comm, ierr) + + if (zrank==0) then + ! Open the movie file (create if it doesn't exist) + call h5_open_or_create(file_id, filename, comm, fexist) + ! If file only just created, add groups for data storage + if (.not. fexist) call h5_add_slice_groups(file_id) + + varname='sal' + call write_H5_plane(file_id, varname, Tplane, 'z') + + ! Close HDF5 file + call h5fclose_f(file_id, hdf_error) + end if + + if (allocated(Tplane)) deallocate(Tplane) + end if + +end subroutine mean_zplane + +subroutine mean_yplane + use param, only: nx, nym + use mpih + use hdf5 + use decomp_2d, only: xstart,xend,xsize,xstartr,xendr,DECOMP_2D_COMM_CART_X,nrank + use local_arrays, only: vz,vy,vx,temp + use mgrd_arrays, only: sal, phi + use h5_tools + use means + implicit none + + character(70) :: filename + character(4) :: varname + real, allocatable :: uplane(:,:), vplane(:,:), wplane(:,:), Tplane(:,:) + integer :: comm, yrank, i, j, bufsize + integer(HID_T) :: file_id + logical :: fexist + + filename = trim("outputdir/flowmov/movie_ymean.h5") + + allocate(uplane(1:nx ,xstart(3):xend(3))) + allocate(vplane(1:nxm,xstart(3):xend(3))) + allocate(wplane(1:nxm,xstart(3):xend(3))) + allocate(Tplane(1:nxm,xstart(3):xend(3))) + + call MPI_CART_SUB(DECOMP_2D_COMM_CART_X, (/.true., .false./), comm, ierr) + call MPI_COMM_RANK(comm, yrank, ierr) + + call ymean(vx(1:nx ,xstart(2):xend(2),xstart(3):xend(3)), uplane, comm, yrank) + call ymean(vy(1:nxm,xstart(2):xend(2),xstart(3):xend(3)), vplane, comm, yrank) + call ymean(vz(1:nxm,xstart(2):xend(2),xstart(3):xend(3)), wplane, comm, yrank) + call ymean(temp(1:nxm,xstart(2):xend(2),xstart(3):xend(3)), Tplane, comm, yrank) + + call MPI_CART_SUB(DECOMP_2D_COMM_CART_X, (/.false., .true./), comm, ierr) + + if (yrank==0) then + + ! Open the movie file (create if it doesn't exist) + call h5_open_or_create(file_id, filename, comm, fexist) + ! If file only just created, add groups for data storage + if (.not. fexist) call h5_add_slice_groups(file_id) + + varname='vx' + call write_H5_plane(file_id, varname, uplane, 'y') + varname='vy' + call write_H5_plane(file_id, varname, vplane, 'y') + varname='vz' + call write_H5_plane(file_id, varname, wplane, 'y') + varname='temp' + call write_H5_plane(file_id, varname, Tplane, 'y') + + ! Close HDF5 file + call h5fclose_f(file_id, hdf_error) + end if + + if (allocated(uplane)) deallocate(uplane) + if (allocated(vplane)) deallocate(vplane) + if (allocated(wplane)) deallocate(wplane) + if (allocated(Tplane)) deallocate(Tplane) + + if (salinity) then + allocate(Tplane(1:nxmr,xstartr(3):xendr(3))) + + call MPI_CART_SUB(DECOMP_2D_COMM_CART_X, (/.true., .false./), comm, ierr) + call MPI_COMM_RANK(comm, yrank, ierr) + + call ymeanr(sal(1:nxmr,xstartr(2):xendr(2),xstartr(3):xendr(3)), Tplane, comm, yrank) + + call MPI_CART_SUB(DECOMP_2D_COMM_CART_X, (/.false., .true./), comm, ierr) + + if (yrank==0) then + ! Open the movie file (create if it doesn't exist) + call h5_open_or_create(file_id, filename, comm, fexist) + ! If file only just created, add groups for data storage + if (.not. fexist) call h5_add_slice_groups(file_id) + + varname='sal' + call write_H5_plane(file_id, varname, Tplane, 'y') + + ! Close HDF5 file + call h5fclose_f(file_id, hdf_error) + end if + + if (allocated(Tplane)) deallocate(Tplane) + end if + + if (phasefield) then + allocate(Tplane(1:nxmr,xstartr(3):xendr(3))) + + call MPI_CART_SUB(DECOMP_2D_COMM_CART_X, (/.true., .false./), comm, ierr) + call MPI_COMM_RANK(comm, yrank, ierr) + + call ymeanr(phi(1:nxmr,xstartr(2):xendr(2),xstartr(3):xendr(3)), Tplane, comm, yrank) + + call MPI_CART_SUB(DECOMP_2D_COMM_CART_X, (/.false., .true./), comm, ierr) + + if (yrank==0) then + ! Open the movie file (create if it doesn't exist) + call h5_open_or_create(file_id, filename, comm, fexist) + ! If file only just created, add groups for data storage + if (.not. fexist) call h5_add_slice_groups(file_id) + + varname='sal' + call write_H5_plane(file_id, varname, Tplane, 'y') + + ! Close HDF5 file + call h5fclose_f(file_id, hdf_error) + end if + + if (allocated(Tplane)) deallocate(Tplane) + end if + +end subroutine mean_yplane diff --git a/src/h5tools/means.F90 b/src/h5tools/means.F90 new file mode 100644 index 00000000..aca6a3c4 --- /dev/null +++ b/src/h5tools/means.F90 @@ -0,0 +1,157 @@ +module means + use mpih + use decomp_2d + use param, only: nym, nymr, nzm, nzmr + implicit none + +contains + +subroutine ymean(var, plane, comm, yrank) + real, intent(in) :: var(:,xstart(2):,xstart(3):) + real, intent(out) :: plane(:,xstart(3):) + integer, intent(in) :: comm, yrank + + integer :: i, j, k, bufsize, n1 + + plane(:,:) = 0.0 + + n1 = size(var, dim=1) + do i=xstart(3),xend(3) + do j=xstart(2),xend(2) + do k=1,n1 + plane(k,i) = plane(k,i) + var(k,j,i) + end do + end do + end do + + bufsize = n1*xsize(3) + if (yrank==0) then + call MPI_REDUCE( & + MPI_IN_PLACE, & + plane(1:n1,xstart(3):xend(3)), & + bufsize, MDP, MPI_SUM, 0, comm, ierr & + ) + else + call MPI_REDUCE( & + plane(1:n1,xstart(3):xend(3)), & + plane(1:n1,xstart(3):xend(3)), & + bufsize, MDP, MPI_SUM, 0, comm, ierr & + ) + end if + + plane = plane/nym + +end subroutine ymean + +subroutine zmean(var, plane, comm, zrank) + real, intent(in) :: var(:,xstart(2):,xstart(3):) + real, intent(out) :: plane(:,xstart(2):) + integer, intent(in) :: comm, zrank + + integer :: i, j, k, bufsize, n1 + + plane(:,:) = 0.0 + + n1 = size(var, dim=1) + do i=xstart(3),xend(3) + do j=xstart(2),xend(2) + do k=1,n1 + plane(k,j) = plane(k,j) + var(k,j,i) + end do + end do + end do + + bufsize = n1*xsize(2) + if (zrank==0) then + call MPI_REDUCE( & + MPI_IN_PLACE, & + plane(1:n1,xstart(2):xend(2)), & + bufsize, MDP, MPI_SUM, 0, comm, ierr & + ) + else + call MPI_REDUCE( & + plane(1:n1,xstart(2):xend(2)), & + plane(1:n1,xstart(2):xend(2)), & + bufsize, MDP, MPI_SUM, 0, comm, ierr & + ) + end if + + plane = plane/nzm + +end subroutine zmean + +subroutine ymeanr(var, plane, comm, yrank) + real, intent(in) :: var(:,xstartr(2):,xstartr(3):) + real, intent(out) :: plane(:,xstartr(3):) + integer, intent(in) :: comm, yrank + + integer :: i, j, k, bufsize, n1 + + plane(:,:) = 0.0 + + n1 = size(var, dim=1) + do i=xstartr(3),xendr(3) + do j=xstartr(2),xendr(2) + do k=1,n1 + plane(k,i) = plane(k,i) + var(k,j,i) + end do + end do + end do + + bufsize = n1*xsizer(3) + if (yrank==0) then + call MPI_REDUCE( & + MPI_IN_PLACE, & + plane(1:n1,xstartr(3):xendr(3)), & + bufsize, MDP, MPI_SUM, 0, comm, ierr & + ) + else + call MPI_REDUCE( & + plane(1:n1,xstartr(3):xendr(3)), & + plane(1:n1,xstartr(3):xendr(3)), & + bufsize, MDP, MPI_SUM, 0, comm, ierr & + ) + end if + + plane = plane/nymr + +end subroutine ymeanr + +subroutine zmeanr(var, plane, comm, zrank) + real, intent(in) :: var(:,xstartr(2):,xstartr(3):) + real, intent(out) :: plane(:,xstartr(2):) + integer, intent(in) :: comm, zrank + + integer :: i, j, k, bufsize, n1 + + plane(:,:) = 0.0 + + n1 = size(var, dim=1) + do i=xstartr(3),xendr(3) + do j=xstartr(2),xendr(2) + do k=1,n1 + plane(k,j) = plane(k,j) + var(k,j,i) + end do + end do + end do + + bufsize = n1*xsizer(2) + if (zrank==0) then + call MPI_REDUCE( & + MPI_IN_PLACE, & + plane(1:n1,xstartr(2):xendr(2)), & + bufsize, MDP, MPI_SUM, 0, comm, ierr & + ) + else + call MPI_REDUCE( & + plane(1:n1,xstartr(2):xendr(2)), & + plane(1:n1,xstartr(2):xendr(2)), & + bufsize, MDP, MPI_SUM, 0, comm, ierr & + ) + end if + + plane = plane/nzmr + +end subroutine zmeanr + +end module means \ No newline at end of file diff --git a/src/main.F90 b/src/main.F90 index 7c330b9c..e0abf0f8 100644 --- a/src/main.F90 +++ b/src/main.F90 @@ -18,7 +18,7 @@ program AFiD integer :: prow=0,pcol=0 integer :: lfactor,lfactor2 character(100) :: arg - logical :: nanexist + logical :: nanexist, write_mean_planes=.true. real,allocatable,dimension(:,:) :: dummy,dscan,dbot integer :: comm,ierror,row_id,row_coords(2),ic,jc,kc @@ -196,6 +196,10 @@ program AFiD call Mkmov_xcut call Mkmov_ycut call Mkmov_zcut + if (write_mean_planes) then + call mean_yplane + call mean_zplane + end if if (ismaster) write(*,*) "Writing 3D fields" call WriteFlowField(.false.) @@ -289,6 +293,10 @@ program AFiD call Mkmov_xcut call Mkmov_ycut call Mkmov_zcut + if (write_mean_planes) then + call mean_yplane + call mean_zplane + end if endif if(ntime.eq.1.or.mod(time,tout).lt.dt) then