diff --git a/sources/io.F90 b/sources/io.F90 index 5d528d8..066230c 100644 --- a/sources/io.F90 +++ b/sources/io.F90 @@ -44,6 +44,13 @@ module io ! subroutine interfaces ! + interface read_snapshot_parameter +#ifdef HDF5 + module procedure read_snapshot_parameter_string_h5 + module procedure read_snapshot_parameter_integer_h5 + module procedure read_snapshot_parameter_double_h5 +#endif /* HDF5 */ + end interface interface write_attribute #ifdef HDF5 module procedure write_scalar_attribute_string_h5 @@ -158,6 +165,7 @@ module io ! declare public subroutines ! public :: initialize_io, finalize_io + public :: read_snapshot_parameter public :: read_restart_snapshot, write_restart_snapshot, write_snapshot public :: restart_from_snapshot public :: next_tout @@ -747,6 +755,348 @@ module io !=============================================================================== ! #ifdef HDF5 +! +!=============================================================================== +! +! subroutine READ_SNAPSHOT_PARAMETER_STRING_H5: +! -------------------------------------------- +! +! Subroutine reads a string parameter from the restart snapshot. +! +! Arguments: +! +! pname - the parameter name; +! pvalue - the parameter value; +! iret - the success flag (the success is 0, failure otherwise); +! +!=============================================================================== +! + subroutine read_snapshot_parameter_string_h5(pname, pvalue, iret) + +! import external procedures +! + use hdf5 , only : hid_t + use hdf5 , only : H5F_ACC_RDONLY_F + use hdf5 , only : H5T_NATIVE_CHARACTER + use hdf5 , only : hid_t, hsize_t, size_t + use hdf5 , only : h5open_f, h5close_f + use hdf5 , only : h5fopen_f, h5fclose_f + use hdf5 , only : h5gopen_f, h5gclose_f + use hdf5 , only : h5aexists_by_name_f + use hdf5 , only : h5tcopy_f, h5tset_size_f + use hdf5 , only : h5aopen_by_name_f, h5aread_f, h5aclose_f + use iso_fortran_env, only : error_unit + use mpitools , only : nproc + use parameters , only : get_parameter_string, get_parameter_integer + +! local variables are not implicit by default +! + implicit none + +! subroutine arguments +! + character(len=*), intent(in) :: pname + character(len=*), intent(inout) :: pvalue + integer , intent(inout) :: iret + +! local variables +! + logical :: info + character(len=255) :: rpath = "./" + character(len=255) :: rname + integer :: nrest = 0 + integer :: np + integer(hid_t) :: fid, gid, tid, aid + integer(size_t) :: aln + integer(hsize_t) :: am(1) = 1 + +! local parameters +! + character(len=*), parameter :: loc = & + 'IO::read_snapshot_parameter_string_h5()' +! +!------------------------------------------------------------------------------- +! +! reset the success flag +! + iret = 0 + +! get the path and the number of the restart snapshot +! + call get_parameter_string ("restart_path" , rpath) + call get_parameter_integer("restart_number", nrest) + +! generate the filename of the restart snapshot +! + info = .false. + np = nproc + 1 + do while (.not. info .and. np >= 0) + np = np - 1 + write (rname, "(a,'r',i6.6,'_',i5.5,'.h5')") trim(rpath), nrest, np + inquire(file = rname, exist = info) + end do + +! procees if file exists +! + if (info) then + call h5open_f(iret) + if (iret >= 0) then + call h5fopen_f(rname, H5F_ACC_RDONLY_F, fid, iret) + if (iret >= 0) then + call h5gopen_f(fid, 'attributes', gid, iret) + if (iret >= 0) then + call h5aexists_by_name_f(gid, '.', trim(pname), info, iret) + if (info .and. iret >= 0) then + call h5aopen_by_name_f(gid, '.', trim(pname), aid, iret) + if (iret >= 0) then + aln = len(pvalue) + call h5tcopy_f(H5T_NATIVE_CHARACTER, tid, iret) + call h5tset_size_f(tid, aln, iret) + call h5aread_f(aid, tid, pvalue, am(:), iret) + call h5aclose_f(aid, iret) + end if + end if + call h5gclose_f(gid, iret) + end if + call h5fclose_f(fid, iret) + end if + call h5close_f(iret) + end if + else + write(error_unit,"('[', a, ']: ', a)") trim(loc) & + , "Snapshot " // trim(rname) // " file does not exist!" + iret = 1 + end if + +!------------------------------------------------------------------------------- +! + end subroutine read_snapshot_parameter_string_h5 +! +!=============================================================================== +! +! subroutine READ_SNAPSHOT_PARAMETER_INTEGER_H5: +! --------------------------------------------- +! +! Subroutine reads an integer parameter from the restart snapshot. +! +! Arguments: +! +! pname - the parameter name; +! pvalue - the parameter value; +! iret - the success flag (the success is 0, failure otherwise); +! +!=============================================================================== +! + subroutine read_snapshot_parameter_integer_h5(pname, pvalue, iret) + +! import external procedures +! + use hdf5 , only : hid_t + use hdf5 , only : H5F_ACC_RDONLY_F + use hdf5 , only : H5T_NATIVE_INTEGER + use hdf5 , only : hid_t, hsize_t, size_t + use hdf5 , only : h5open_f, h5close_f + use hdf5 , only : h5fopen_f, h5fclose_f + use hdf5 , only : h5gopen_f, h5gclose_f + use hdf5 , only : h5aexists_by_name_f + use hdf5 , only : h5aopen_by_name_f, h5aread_f, h5aclose_f + use iso_fortran_env, only : error_unit + use mpitools , only : nproc + use parameters , only : get_parameter_string, get_parameter_integer + +! local variables are not implicit by default +! + implicit none + +! subroutine arguments +! + character(len=*), intent(in) :: pname + integer , intent(inout) :: pvalue + integer , intent(inout) :: iret + +! local variables +! + logical :: info + character(len=255) :: rpath = "./" + character(len=255) :: rname + integer :: nrest = 0 + integer :: np + integer(hid_t) :: fid, gid, aid + integer(size_t) :: aln + integer(hsize_t) :: am(1) = 1 + +! local parameters +! + character(len=*), parameter :: loc = & + 'IO::read_snapshot_parameter_integer_h5()' +! +!------------------------------------------------------------------------------- +! +! reset the success flag +! + iret = 0 + +! get the path and the number of the restart snapshot +! + call get_parameter_string ("restart_path" , rpath) + call get_parameter_integer("restart_number", nrest) + +! generate the filename of the restart snapshot +! + info = .false. + np = nproc + 1 + do while (.not. info .and. np >= 0) + np = np - 1 + write (rname, "(a,'r',i6.6,'_',i5.5,'.h5')") trim(rpath), nrest, np + inquire(file = rname, exist = info) + end do + +! procees if file exists +! + if (info) then + call h5open_f(iret) + if (iret >= 0) then + call h5fopen_f(rname, H5F_ACC_RDONLY_F, fid, iret) + if (iret >= 0) then + call h5gopen_f(fid, 'attributes', gid, iret) + if (iret >= 0) then + call h5aexists_by_name_f(gid, '.', trim(pname), info, iret) + if (info .and. iret >= 0) then + call h5aopen_by_name_f(gid, '.', trim(pname), aid, iret) + if (iret >= 0) then + call h5aread_f(aid, H5T_NATIVE_INTEGER, pvalue, am(:), iret) + call h5aclose_f(aid, iret) + end if + end if + call h5gclose_f(gid, iret) + end if + call h5fclose_f(fid, iret) + end if + call h5close_f(iret) + end if + else + write(error_unit,"('[', a, ']: ', a)") trim(loc) & + , "Snapshot " // trim(rname) // " file does not exist!" + iret = 1 + end if + +!------------------------------------------------------------------------------- +! + end subroutine read_snapshot_parameter_integer_h5 +! +!=============================================================================== +! +! subroutine READ_SNAPSHOT_PARAMETER_DOUBLE_H5: +! -------------------------------------------- +! +! Subroutine reads a double precision real parameter from the restart +! snapshot. +! +! Arguments: +! +! pname - the parameter name; +! pvalue - the parameter value; +! iret - the success flag (the success is 0, failure otherwise); +! +!=============================================================================== +! + subroutine read_snapshot_parameter_double_h5(pname, pvalue, iret) + +! import external procedures +! + use hdf5 , only : hid_t + use hdf5 , only : H5F_ACC_RDONLY_F + use hdf5 , only : H5T_NATIVE_DOUBLE + use hdf5 , only : hid_t, hsize_t, size_t + use hdf5 , only : h5open_f, h5close_f + use hdf5 , only : h5fopen_f, h5fclose_f + use hdf5 , only : h5gopen_f, h5gclose_f + use hdf5 , only : h5aexists_by_name_f + use hdf5 , only : h5aopen_by_name_f, h5aread_f, h5aclose_f + use iso_fortran_env, only : error_unit + use mpitools , only : nproc + use parameters , only : get_parameter_string, get_parameter_integer + +! local variables are not implicit by default +! + implicit none + +! subroutine arguments +! + character(len=*), intent(in) :: pname + real(kind=8) , intent(inout) :: pvalue + integer , intent(inout) :: iret + +! local variables +! + logical :: info + character(len=255) :: rpath = "./" + character(len=255) :: rname + integer :: nrest = 0 + integer :: np + integer(hid_t) :: fid, gid, aid + integer(size_t) :: aln + integer(hsize_t) :: am(1) = 1 + +! local parameters +! + character(len=*), parameter :: loc = & + 'IO::read_snapshot_parameter_double_h5()' +! +!------------------------------------------------------------------------------- +! +! reset the success flag +! + iret = 0 + +! get the path and the number of the restart snapshot +! + call get_parameter_string ("restart_path" , rpath) + call get_parameter_integer("restart_number", nrest) + +! generate the filename of the restart snapshot +! + info = .false. + np = nproc + 1 + do while (.not. info .and. np >= 0) + np = np - 1 + write (rname, "(a,'r',i6.6,'_',i5.5,'.h5')") trim(rpath), nrest, np + inquire(file = rname, exist = info) + end do + +! procees if file exists +! + if (info) then + call h5open_f(iret) + if (iret >= 0) then + call h5fopen_f(rname, H5F_ACC_RDONLY_F, fid, iret) + if (iret >= 0) then + call h5gopen_f(fid, 'attributes', gid, iret) + if (iret >= 0) then + call h5aexists_by_name_f(gid, '.', trim(pname), info, iret) + if (info .and. iret >= 0) then + call h5aopen_by_name_f(gid, '.', trim(pname), aid, iret) + if (iret >= 0) then + call h5aread_f(aid, H5T_NATIVE_DOUBLE, pvalue, am(:), iret) + call h5aclose_f(aid, iret) + end if + end if + call h5gclose_f(gid, iret) + end if + call h5fclose_f(fid, iret) + end if + call h5close_f(iret) + end if + else + write(error_unit,"('[', a, ']: ', a)") trim(loc) & + , "Snapshot " // trim(rname) // " file does not exist!" + iret = 1 + end if + +!------------------------------------------------------------------------------- +! + end subroutine read_snapshot_parameter_double_h5 +! !=============================================================================== ! ! subroutine READ_RESTART_SNAPSHOT_H5: