From 63ccfaa86fa6fcf8f8bad4a61fcdd56ba08206d7 Mon Sep 17 00:00:00 2001
From: Pedro Costa
Date: Mon, 4 Mar 2024 12:21:13 +0100
Subject: [PATCH] Fix errors in code for HDF5 backend.
---
src/load.f90 | 10 +++++-----
1 file changed, 5 insertions(+), 5 deletions(-)
diff --git a/src/load.f90 b/src/load.f90
index 4485ab23..d2cf370f 100644
--- a/src/load.f90
+++ b/src/load.f90
@@ -14,7 +14,7 @@ module mod_load
use mod_utils, only: f_sizeof
implicit none
private
- public load_all,io_field
+ public load_all,io_field,io_field_hdf5
contains
subroutine load_all(io,filename,comm,ng,nh,lo,hi,u,v,w,p,time,istep)
!
@@ -543,6 +543,7 @@ end subroutine load_one_local
!
#if defined(_USE_HDF5)
subroutine io_field_hdf5(io,filename,varname,ng,nh,lo,hi,var,meta,x_g,y_g,z_g)
+ use hdf5
!
! collective single field data I/O using HDF5
!
@@ -575,7 +576,6 @@ subroutine io_field_hdf5(io,filename,varname,ng,nh,lo,hi,var,meta,x_g,y_g,z_g)
integer(HSIZE_T) , dimension(3) :: data_count
integer(HSSIZE_T), dimension(3) :: data_offset
integer(HSSIZE_T), dimension(3) :: halo_offset
- type(MPI_INFO) :: info = MPI_INFO_NULL
!
n(:) = hi(:)-lo(:)+1
sizes(:) = ng(:)
@@ -591,7 +591,7 @@ subroutine io_field_hdf5(io,filename,varname,ng,nh,lo,hi,var,meta,x_g,y_g,z_g)
select case(io)
case('r')
call h5pcreate_f(H5P_FILE_ACCESS_F,plist_id,ierr)
- call h5pset_fapl_mpio_f(plist_id,MPI_COMM_WORLD%MPI_VAL,info%MPI_VAL,ierr)
+ call h5pset_fapl_mpio_f(plist_id,MPI_COMM_WORLD,MPI_INFO_NULL,ierr)
call h5fopen_f(filename,H5F_ACC_RDONLY_F,file_id,ierr,access_prp=plist_id)
call h5pclose_f(plist_id,ierr)
!
@@ -618,11 +618,11 @@ subroutine io_field_hdf5(io,filename,varname,ng,nh,lo,hi,var,meta,x_g,y_g,z_g)
call h5dclose_f(dset,ierr)
call h5fclose_f(file_id,ierr)
end if
- call MPI_Bcast(meta,2,MPI_REAL_RP,0,MPI_COMM_WORLD)
+ call MPI_Bcast(meta,2,MPI_REAL_RP,0,MPI_COMM_WORLD,ierr)
case('w')
call h5screate_simple_f(ndims,dims,filespace,ierr)
call h5pcreate_f(H5P_FILE_ACCESS_F,plist_id,ierr)
- call h5pset_fapl_mpio_f(plist_id,MPI_COMM_WORLD%MPI_VAL,info%MPI_VAL,ierr)
+ call h5pset_fapl_mpio_f(plist_id,MPI_COMM_WORLD,MPI_INFO_NULL,ierr)
call h5fcreate_f(filename,H5F_ACC_TRUNC_F,file_id,ierr,access_prp=plist_id)
call h5pclose_f(plist_id,ierr)
!