Hi, thank you for the replies that gave me a better understandment of the problem. As you suggested, it is better that I explain in some more detail the situation and my requests. I am developing a MPI Fortran CFD code. The 3D CFD domain is divided in ordered cells along three computational coordinate indices. These cells are grouped in a number, say nb, of blocks. Each block has ni(b),nj(b),nk(b) cells along three coordinates for a total of ni(b)*nj(b)*nk(b) cells. The blocks are scattered between the np processes of a MPI run, in a way to get a balanced load in a run. The physical variables r (density) and p (pressure) are stored in a structure of double precision arrays for each process n in the range (0:np-1)
b=1,nb(n)
var(b)%r(1:ni(b),1:nj(b),1:nk(b))
var(b)%p(1:ni(b),1:nj(b),1:nk(b))
n=0,np
where nb(n) is the number of blocks hosted in process n.
Now, the central point is to save these arrays in hdf5 format in the most efficient way. Of course, the first choice is parallel collective writing. In hdf5 forum, they suggested me to create a compound datatype for each process and then to write in parallel an array of compound datatype scattered between processes.
In order to do this, hdf5 library needs to know the offset locations of the members of the structure. One could provide them by:
- their dimensions (but one needs that the data are contiguous in memory)
- the offset of their memory locations
and I write down a simple code that I provided in the first post of the topic, where I found that in the MPI run the memory locations of the structure members are scrambled in unpredictable manner.
Now you will find the hdf5 version of the same code below, where a compound datatype is created, populated with the 2*nb arrays in each process, and then written to disk: to my purpose, in a MPI run only the rank 0 will write data since I need to understand the memory layout behavior. The code has an input argument: if 0 or unspecified, it will use offset based on the dimension of arrays, otherwise will use memory locations offsets.
I have tested the code with intel and gnu compilers and the results are the same:
serial version (no -DMPIV preprocessing flag):
- dimension offset (0) : run went fine and file created, but values in the file are wrong
- memory location offset (1) : run went fine and file created, and values OK
parallel version (-DMPIV preprocessing flag):
- dimension offset (0) : run went fine and file created, but values in the file are wrong
- memory location offset (1) : hdf5 library complains about incorrect memory locations
major: Datatype
minor: Unable to insert object
#001: H5Tcompound.c line 458 in H5T__insert(): member extends past end of compound type
major: Datatype
minor: Unable to insert object
#000: H5Tcompound.c line 370 in H5Tinsert(): unable to insert member
So, in MPI implementation, there are no working alternatives.
module var_mod
type vt
double precision,dimension(:,:,:),allocatable :: r,p
end type vt
type(vt),target,dimension(:),allocatable :: var
end module var_mod
PROGRAM main
USE hdf5
USE ISO_C_BINDING
#ifdef MPIV
USE MPI
#endif
USE var_mod
IMPLICIT NONE
! FILES
CHARACTER(LEN=*),PARAMETER :: H5FILE_NAME = "mpivars.h5"
CHARACTER(LEN=*),PARAMETER :: DATASETNAME = "CompoundOfArrays"
CHARACTER(LEN=10) :: varname,FIRSTARG
INTEGER,PARAMETER :: NB = 4
INTEGER,PARAMETER :: RANK = 3
!----------------------------------------------------------------
INTEGER :: b,ilen,ioff,ist,i,j,k,me,np,mpierr
INTEGER(hid_t) :: file,dataset,space
INTEGER(hsize_t) :: DIM(1) = (/1/) ! Dataspace dimensions
INTEGER(SIZE_T) :: offset,type_size,type_sized ! Size of the double precision datatype
INTEGER :: hdferr
TYPE(C_PTR) :: f_ptr
INTEGER(hid_t) :: vwrite ! File datatype identifier
INTEGER(hid_t),dimension(nb) :: tidr,tidp ! /* Nested Array Datatype ID */
INTEGER(HSIZE_T),DIMENSION(3) :: tdims
integer,dimension(:),allocatable :: ni,nj,nk
!----------------------------------------------------------------
me = 0
#ifdef MPIV
call mpi_init(mpierr)
call mpi_comm_size(mpi_comm_world,np,mpierr)
call mpi_comm_rank(mpi_comm_world,me,mpierr)
#endif
if (me==0) then
ioff = 0
if (command_argument_count()/=0) then
call get_command_argument(1,FIRSTARG,ilen,ist)
read (FIRSTARG,*) ioff
end if
allocate (ni(nb))
allocate (nj(nb))
allocate (nk(nb))
allocate (var(nb))
do b = 1,NB
ni(b) = b/3+1
nj(b) = b
nk(b) = b/5+1
allocate (var(b)%r(ni(b),nj(b),nk(b)))
allocate (var(b)%p(ni(b),nj(b),nk(b)))
END DO
!
! Initialize FORTRAN interface.
!
CALL h5open_f(hdferr)
!
! Initialize the data
!
do b = 1,nb
DO k = 1,nk(b)
DO j = 1,nj(b)
DO i = 1,ni(b)
var(b)%r(i,j,k) = 10000*me+1000*b+i*100+j*10+k
var(b)%p(i,j,k) = -var(b)%r(i,j,k)
END DO
END DO
END DO
end do
!
! Create the file.
!
CALL H5Fcreate_f(H5FILE_NAME,H5F_ACC_TRUNC_F,file,hdferr)
CALL h5tget_size_f(H5T_NATIVE_DOUBLE,type_sized,hdferr)
type_size = sum(ni*nj*nk)*2*type_sized
if(ioff>0) type_size = H5OFFSETOF(C_LOC(var),C_LOC(var(nb)%p))+H5OFFSETOF(C_LOC(var(nb)%r),C_LOC(var(nb)%p))
write (*,*) 'var loc:',LOC(var)
do i = 1,nb
write (*,*) 'item:',i
write (*,*) 'var loc:',LOC(var(i))
write (*,*) 'var%r loc:',LOC(var(i)%r)
write (*,*) 'var%p loc:',LOC(var(i)%p)
end do
do i = 1,nb
write (*,*) 'item:',i
if (i>1) write (*,*) 'offset var i/i-1:',LOC(var(i))-LOC(var(i-1))
write (*,*) 'offset var%r/var:',LOC(var(i)%r)-LOC(var(i))
write (*,*) 'offset var%p/var:',LOC(var(i)%p)-LOC(var(i))
write (*,*) 'offset var%p/var%r:',LOC(var(i)%p)-LOC(var(i)%r)
end do
!
! Create the memory data type.
!
CALL H5Tcreate_f(H5T_COMPOUND_F,type_size,vwrite,hdferr)
do b = 1,nb
tdims = (/ni(b),nj(b),nk(b)/)
CALL h5tarray_create_f(H5T_NATIVE_DOUBLE,RANK,tdims,tidr(b),hdferr)
CALL h5tarray_create_f(H5T_NATIVE_DOUBLE,RANK,tdims,tidp(b),hdferr)
end do
offset = 0
do b = 1,nb
write (varname,10) b
10 format('var_r_b',I0)
if(ioff>0) offset = H5OFFSETOF(C_LOC(var),C_LOC(var(b)%r))
CALL h5tinsert_f(vwrite,varname,offset,tidr(b),hdferr)
if(ioff==0) offset = offset+type_sized*ni(b)*nj(b)*nk(b)
write (varname,11) b
11 format('var_p_b',I0)
if(ioff>0) offset = H5OFFSETOF(C_LOC(var),C_LOC(var(b)%p))
CALL h5tinsert_f(vwrite,varname,offset,tidp(b),hdferr)
if(ioff==0) offset = offset+type_sized*ni(b)*nj(b)*nk(b)
end do
!
! Create the data space.
!
!
CALL H5Screate_simple_f(1,dim,space,hdferr)
!
! Create the dataset.
!
CALL H5Dcreate_f(file,DATASETNAME,vwrite,space,dataset,hdferr)
!
! Write data to the dataset
!
f_ptr = C_LOC(var)
CALL H5Dwrite_f(dataset,vwrite,f_ptr,hdferr)
!
! Release resources
!
CALL H5Tclose_f(vwrite,hdferr)
CALL H5Sclose_f(space,hdferr)
CALL H5Dclose_f(dataset,hdferr)
CALL H5Fclose_f(file,hdferr)
do b = 1,nb
deallocate (var(b)%r,var(b)%p)
END DO
deallocate (var)
deallocate (ni,nj,nk)
end if
#ifdef MPIV
call mpi_finalize(mpierr)
#endif
END PROGRAM main
