Reads marker data from file in parallel using HDF5.
| Type | Intent | Optional | Attributes | Name | ||
|---|---|---|---|---|---|---|
| class(marker_set), | intent(inout) | :: | this |
A collection of tessellation elements |
||
| integer, | intent(out) | :: | iter |
Iteration at write |
||
| real(kind=wp), | intent(out) | :: | time |
Time at write |
impure subroutine marker_set_ReadHDF5(this,iter,time) !> Reads marker data from file in parallel using HDF5. use leapUtils, only: stringtool_obj implicit none class(marker_set), intent(inout) :: this !! A collection of tessellation elements integer, intent(out) :: iter !! Iteration at write real(wp), intent(out) :: time !! Time at write ! Work variables integer, allocatable :: buffi(:) !! Integer buffer integer(leapI8), allocatable :: buffi8(:) !! Integer (kind=8) buffer real(wp), allocatable :: buffr(:) !! Real buffer integer, allocatable :: size_per_rank(:) !! Number of markers per rank integer :: max_chunk !! Max number of chunks on any MPI rank integer :: chunk_size !! Number of markers in chunk integer :: goff !! Global offset integer :: loff !! Local offset integer :: n !! Iterator type(stringtool_obj):: stringtool character(len=:), & allocatable :: filename type(hdf5_obj) :: hdf5 filename = 'HDF5/'//stringtool%RemoveExtension(this%read_file) ! Initialize HDF5 call hdf5%Initialize(this%parallel) ! Open file call hdf5%Open(filename,'R') ! Read attributes call hdf5%ReadAttributes("/","Time", time) call hdf5%ReadAttributes("/","Iter", iter) call hdf5%ReadAttributes("/","Count",this%count) ! Distribute equally the data on each mpi rank associate(mpi=>this%parallel) allocate(size_per_rank(mpi%nproc)) do n=1,mpi%nproc size_per_rank(n) = int(this%count/mpi%nproc) if (n.le.mod(this%count,mpi%nproc)) size_per_rank(n)=size_per_rank(n)+1 end do end associate ! Allocate buffers call this%parallel%Max(ceiling(size_per_rank(this%parallel%rank%mine)/real(this%R_chunk_size,wp)),max_chunk) loff = 0 chunk_size = 0 do n=1,max_chunk ! Local offset loff = loff + chunk_size ! Determine next chunk size if (loff.lt.size_per_rank(this%parallel%rank%mine)) then chunk_size = min(size_per_rank(this%parallel%rank%mine) - loff,this%R_chunk_size) ! Resize the data array to contain read Lagrangian objects call this%Resize(this%count_+chunk_size) else chunk_size = 0 end if ! Global offset goff = sum(size_per_rank(1:this%parallel%rank%mine-1)) + loff + 1 allocate(buffi8(chunk_size),buffi(chunk_size),buffr(chunk_size)) ! Read the base Lagrangian object structure if (chunk_size.eq.0) then call hdf5%Read("/","id",buffi8,offset=goff) call hdf5%Read("/","x", buffr, offset=goff) call hdf5%Read("/","y", buffr, offset=goff) call hdf5%Read("/","z", buffr, offset=goff) call hdf5%Read("/","i", buffi, offset=goff) call hdf5%Read("/","j", buffi, offset=goff) call hdf5%Read("/","k", buffi, offset=goff) call hdf5%Read("/","s", buffi, offset=goff) call hdf5%Read("/","SA",buffr, offset=goff) call hdf5%Read("/","nx",buffr, offset=goff) call hdf5%Read("/","ny",buffr, offset=goff) call hdf5%Read("/","nz",buffr, offset=goff) call hdf5%Read("/","vx",buffr, offset=goff) call hdf5%Read("/","vy",buffr, offset=goff) call hdf5%Read("/","vz",buffr, offset=goff) call hdf5%Read("/","fx",buffr, offset=goff) call hdf5%Read("/","fy",buffr, offset=goff) call hdf5%Read("/","fz",buffr, offset=goff) else select type (markers=>this%p) type is (marker_obj) call hdf5%Read("/","id",buffi8,offset=goff); markers(this%count_-chunk_size+1:this%count_)%id = buffi8 call hdf5%Read("/","x", buffr, offset=goff); markers(this%count_-chunk_size+1:this%count_)%p(1) = buffr call hdf5%Read("/","y", buffr, offset=goff); markers(this%count_-chunk_size+1:this%count_)%p(2) = buffr call hdf5%Read("/","z", buffr, offset=goff); markers(this%count_-chunk_size+1:this%count_)%p(3) = buffr call hdf5%Read("/","i", buffi, offset=goff); markers(this%count_-chunk_size+1:this%count_)%c(1) = buffi call hdf5%Read("/","j", buffi, offset=goff); markers(this%count_-chunk_size+1:this%count_)%c(2) = buffi call hdf5%Read("/","k", buffi, offset=goff); markers(this%count_-chunk_size+1:this%count_)%c(3) = buffi call hdf5%Read("/","s", buffi, offset=goff); markers(this%count_-chunk_size+1:this%count_)%s = buffi call hdf5%Read("/","SA",buffr, offset=goff); markers(this%count_-chunk_size+1:this%count_)%SA = buffr call hdf5%Read("/","nx",buffr, offset=goff); markers(this%count_-chunk_size+1:this%count_)%n(1) = buffr call hdf5%Read("/","ny",buffr, offset=goff); markers(this%count_-chunk_size+1:this%count_)%n(2) = buffr call hdf5%Read("/","nz",buffr, offset=goff); markers(this%count_-chunk_size+1:this%count_)%n(3) = buffr call hdf5%Read("/","vx",buffr, offset=goff); markers(this%count_-chunk_size+1:this%count_)%v(1) = buffr call hdf5%Read("/","vy",buffr, offset=goff); markers(this%count_-chunk_size+1:this%count_)%v(2) = buffr call hdf5%Read("/","vz",buffr, offset=goff); markers(this%count_-chunk_size+1:this%count_)%v(3) = buffr call hdf5%Read("/","fx",buffr, offset=goff); markers(this%count_-chunk_size+1:this%count_)%f(1) = buffr call hdf5%Read("/","fy",buffr, offset=goff); markers(this%count_-chunk_size+1:this%count_)%f(2) = buffr call hdf5%Read("/","fz",buffr, offset=goff); markers(this%count_-chunk_size+1:this%count_)%f(3) = buffr end select end if ! Free buffers deallocate(buffi8,buffi,buffr) ! Send markers to their mpiranks call this%Communicate() end do ! Free remaining arrays, close file and finalize HDF5 deallocate(size_per_rank) call hdf5%Close() call hdf5%Finalize() ! Set unread variables to zero select type (markers=>this%p) class is (marker_obj) do n=1,this%count_ markers(n)%pold = 0.0_wp markers(n)%vold = 0.0_wp end do end select ! Localize markers on the grid call this%Localize() return end subroutine marker_set_ReadHDF5