Append mode error
Created by: jychoi-hpc
Describe the bug I am working on Pixie3D to write restart files with appending mode without any success. I wrote a test case to reproduce the error.
In the test code, I am trying to write a file with appending mode but getting an error.
When writing with adios2_mode_write
, it works. But, a file with adios2_mode_append
is corrupted and cannot do bpls
(no error during the runtime though).
To Reproduce Here is my test code (writer.F90)
program writer
use mpi
use adios2
implicit none
integer(kind=8), dimension(1) :: shape_dims, start_dims, count_dims
real, dimension(:), allocatable :: myArray
integer :: inx, irank, isize, ierr, i, var1_type
type(adios2_adios) :: adios
type(adios2_io) :: io
type(adios2_variable) :: var, var1
type(adios2_attribute) :: attr
type(adios2_engine) :: engine
character(len=:), allocatable :: var1_name
integer :: comm
integer :: istatus
real*8 :: totalpe;
CHARACTER(len=32) :: arg
integer :: mode
logical :: do_read = .false.
logical :: do_append = .false.
do i = 1, iargc()
call getarg(i, arg)
if (trim(arg)=='-read') do_read = .true.
if (trim(arg)=='-append') do_append = .true.
end do
! Launch MPI
comm = MPI_COMM_WORLD
call MPI_Init(ierr)
call MPI_Comm_rank(comm, irank, ierr)
call MPI_Comm_size(comm, isize, ierr)
! Application variables
inx = 10
allocate( myArray(inx) )
do i=1,inx
myArray(i) = i-1
end do
! Variable dimensions
shape_dims(1) = isize * inx
start_dims(1) = irank * inx
count_dims(1) = inx
if (do_append) then
print *, 'Mode: append'
mode = adios2_mode_append
else
print *, 'Mode: write'
mode = adios2_mode_write
endif
! Create adios handler passing the communicator, debug mode and error flag
call adios2_init(adios, comm, adios2_debug_mode_on, ierr)
! Declare an IO process configuration inside adios
call adios2_declare_io(io, adios, "group1", ierr)
! Defines a variable to be written in bp format
! call adios2_define_variable(var, io, 'totalpe', adios2_type_dp, ierr)
call adios2_define_variable(var1, io, "totalpe", adios2_type_dp, 1, &
(/ 1_8 /), (/ 0_8 /), (/ 1_8 /), &
adios2_constant_dims, ierr)
call adios2_define_variable(var1, io, "arr1", adios2_type_real, 1, &
shape_dims, start_dims, count_dims, &
adios2_constant_dims, ierr)
call adios2_define_attribute(attr, io, "attr", "this is attr", ierr)
! Open myVector_f.bp in write mode, this launches an engine
call adios2_open(engine, io, "out.bp", mode, comm, ierr)
print *, 'adios2_open:ierr', ierr
do i=1,4
call adios2_begin_step(engine, adios2_step_mode_append, ierr)
print *, 'adios2_begin_step:ierr', ierr
! Put myArray contents to bp buffer, based on var1 metadata
call adios2_put(engine, 'totalpe', 777D0, ierr)
call adios2_put(engine, var1, myArray, adios2_mode_sync, ierr)
! Closes engine and deallocates it, becomes unreachable
call adios2_end_step(engine, ierr)
print *, 'adios2_end_step:ierr', ierr
enddo
call adios2_close(engine, ierr)
print *, 'adios2_close:ierr', ierr
! Deallocates adios and calls its destructor
call adios2_finalize(adios, ierr)
if( allocated(myArray) ) deallocate(myArray)
if( allocated(var1_name) ) deallocate(var1_name)
call MPI_Finalize(ierr)
end program writer
Makefile
FC=mpif90
ADIOS2_DIR=/opt/adios2
ADIOS2_INC=$(shell $(ADIOS2_DIR)/bin/adios2-config --fortran-flags)
ADIOS2_LIB=$(shell $(ADIOS2_DIR)/bin/adios2-config --fortran-libs)
BINS=writer
all: $(BINS)
writer: writer.F90
$(FC) -g -ffree-line-length-0 $(ADIOS2_INC) -o $@ $? $(ADIOS2_LIB)
clean:
rm -rf $(BINS) *.o *.dSYM
I am running in the following step. First, I create out.bp
with adios2_mode_write
mode
$ ./writer
$ bpls -lva out.bp
File info:
of variables: 2
of attributes: 1
statistics: Min / Max
float arr1 4*{10} = 0 / 9
string attr attr = "this is attr"
double totalpe 4*{1} = 777 / 777
it works fine
Then, I am trying to do with adios2_mode_append
:
$ ./writer -append
$ bpls -lva out.bp
$ bpls -lva out.bp/
HDF5-DIAG: Error detected in HDF5 (1.10.4) thread 0:
#000: H5F.c line 509 in H5Fopen(): unable to open file
major: File accessibilty
minor: Unable to open file
#001: H5Fint.c line 1400 in H5F__open(): unable to open file
major: File accessibilty
minor: Unable to open file
#002: H5Fint.c line 1700 in H5F_open(): unable to read superblock
major: File accessibilty
minor: Read failed
It should be BP4
file, not HDF5. It looks like the file corrupted.
Any help and advice will be appreciated.
Thanks.