hdf5 icon indicating copy to clipboard operation
hdf5 copied to clipboard

Allocatable components of Fortran derived types change when reading static components

Open HiSPEET opened this issue 1 year ago • 10 comments

Describe the bug A HDF5 compound type is constructed from a Fortran derived type, which contains static as well as allocatable components. Only the static components are included in the compound. When writing and reading the compound type, the allocation status of the dynamic components in the read-in data depends on the written data. In particular, the dynamic components may change their status from deallocated to allocated. This seems to happen, when they were allocated in the written data.

Expected behavior The status and contents of allocatable components in the read data should stay untouched.

Platform (please complete the following information)

  • HDF5 1.14.4.2 (MacPorts)
  • MAC OS Sonoma 14.5 (23F79)
  • GNU Fortran (MacPorts gcc12 12.3.0_4+stdlib_flag) 12.3.0
  • cmake version 3.29.3
  • serial HDF5

Additional context The attached file contains the Fortran source and the control output for a derived-type array of length 2. The latter shows that the dynamic component of the first element is deallocated. In the second element it is allocated, although the file should contain no information about the dynamic component.

hdf5_test.zip

HiSPEET avatar Jun 18 '24 13:06 HiSPEET

I'm surprised this works as reasonably as it does; I would suspect if you investigated other compiler/compiler combinations, you would encounter issues. What you would like to use instead is referred to as variable length in HDF5 terminology. I rewrote your example to use a variable length type instead. If you dump the resulting file, you get the following:

HDF5 "elements.h5" {
GROUP "/" {
   GROUP "data" {
      DATASET "elements" {
         DATATYPE  H5T_COMPOUND {
            H5T_STD_I32LE "id";
            H5T_VLEN { H5T_STD_I32LE } "nb";
         }
         DATASPACE  SIMPLE { ( 2 ) / ( 2 ) }
         DATA {
         (0): {
               100,
               ()
            },
         (1): {
               200,
               (10, 20)
            }
         }
      }
   }
}
}

and the output should be:

elem_orig%id = 100 200
elem_read%id = 100 200
elem_read( 2 )%nb = 10 20

hdf5_test.F90.gz

brtnfld avatar Jun 19 '24 16:06 brtnfld

Hi Scot,

thank you for pointing to that interesting approach. Indeed I was not aware of the variable length features of HDF5. Unfortunately, the suggested approach requires to introduce additional data structures for dynamic components and also the use of pointers in these structures, which we generally avoid. Our intended strategy is to pack dynamic components of the same type into auxiliary arrays that can be written or read in one piece.

So, the wish remains that HDF5 does not override the status of dynamic components, when reading data. The current behavior is that the dynamic components get „allocated“ if they were allocated when writing. This is naturally the case, e.g., when writing restart data. Unfortunately, this status is wrong and trying to fix it by deallocating the component leads often to a runtime error.

Best regards Joerg

Am 19.06.2024 um 18:27 schrieb Scot Breitenfeld @.***>:

I'm surprised this works as well as it does, but what you want to use instead is referred to as variable length in HDF5 terminology. I rewrote your example to use a variable length type instead. If you dump the resulting file, you get the following:

HDF5 "elements.h5" { GROUP "/" { GROUP "data" { DATASET "elements" { DATATYPE H5T_COMPOUND { H5T_STD_I32LE "id"; H5T_VLEN { H5T_STD_I32LE } "nb"; } DATASPACE SIMPLE { ( 2 ) / ( 2 ) } DATA { (0): { 100, () }, (1): { 200, (10, 20) } } } } } }

and the output should be:

elem_orig%id = 100 200 elem_read%id = 100 200 elem_read( 2 )%nb = 10 20 hdf5_test.F90.gz https://github.com/user-attachments/files/15903707/hdf5_test.F90.gz — Reply to this email directly, view it on GitHub https://github.com/HDFGroup/hdf5/issues/4577#issuecomment-2179087611, or unsubscribe https://github.com/notifications/unsubscribe-auth/AKDBWQIX6CA2OFC5GLLFAU3ZIGWPPAVCNFSM6AAAAABJQDFH5GVHI2DSMVQWIX3LMV43OSLTON2WKQ3PNVWWK3TUHMZDCNZZGA4DONRRGE. You are receiving this because you authored the thread.

HiSPEET avatar Jun 19 '24 19:06 HiSPEET

A trick to fix the corrupted allocation status in Fortran is described here:

https://stackoverflow.com/questions/64891525/reset-deallocate-nullify-a-fortran-allocatable-array-that-has-been-corrupted

This could solve our problem. But, of course I would prefer a clean solution in HDF5.

Joerg

Am 19.06.2024 um 18:27 schrieb Scot Breitenfeld @.***>:

I'm surprised this works as well as it does, but what you want to use instead is referred to as variable length in HDF5 terminology. I rewrote your example to use a variable length type instead. If you dump the resulting file, you get the following:

HDF5 "elements.h5" { GROUP "/" { GROUP "data" { DATASET "elements" { DATATYPE H5T_COMPOUND { H5T_STD_I32LE "id"; H5T_VLEN { H5T_STD_I32LE } "nb"; } DATASPACE SIMPLE { ( 2 ) / ( 2 ) } DATA { (0): { 100, () }, (1): { 200, (10, 20) } } } } } }

and the output should be:

elem_orig%id = 100 200 elem_read%id = 100 200 elem_read( 2 )%nb = 10 20 hdf5_test.F90.gz https://github.com/user-attachments/files/15903707/hdf5_test.F90.gz — Reply to this email directly, view it on GitHub https://github.com/HDFGroup/hdf5/issues/4577#issuecomment-2179087611, or unsubscribe https://github.com/notifications/unsubscribe-auth/AKDBWQIX6CA2OFC5GLLFAU3ZIGWPPAVCNFSM6AAAAABJQDFH5GVHI2DSMVQWIX3LMV43OSLTON2WKQ3PNVWWK3TUHMZDCNZZGA4DONRRGE. You are receiving this because you authored the thread.

HiSPEET avatar Jun 19 '24 20:06 HiSPEET

Unfortunately, what you are proposing is not supported in HDF5. HDF5 does not support I/O having a struct with a field that is not allocated unless that field is declared a vl type. You need to inform HDF5 what fields are written in the compound type. Currently, HDF5 expects one field to be written, but you are writing another field, which can corrupt other data in the file.

The H5Dread API does not allocate anything; it only fills the buffer being passed in, and Fortran misinterprets what C is reading and passing back. I don't think the link you provided is relevant in this case.

You can either use fixed arrays in the type, which I'm assuming you don't want to do if this is a list of nodes associated with different kinds of elements, or use an array of integers with unlimited dimension, pack all the variables in the array and store another array with element lengths.

brtnfld avatar Jun 19 '24 22:06 brtnfld

You can either use [...], or use an array of integers with unlimited dimension, pack all the variables in the array and store another array with element lengths.

This is exactly what we do. However, keeping the element type with dynamic components is essential for us.

The problem is that H5Dread fills not only to the addresses declared in the type definition but also the locations in between. This destroys the descriptors of allocatable or pointer components though not included in the type.

I see three options

  1. Clean handling in HDF5, e.g.
    • not modify memory locations that are not included in the HDF5 compound, or
    • provide the opportunity to specify fill values for those locations, or
    • provide explicit handling of allocatable components based on the standardized ISO C Fortran interface.
  2. Destroy dynamic data before calling H5Dwrite. This seems to lead H5Dread to return the state „deallocated“ .
  3. Fix the allocation status manually by calling the C function CFI_deallocate or, if this fails, another C function doing the dirty job.

The first option is preferable, but I understand that it generates work on your side.

In fortran/src/H5Df.c, line 138, the data is first read into buf_c and then transferred to the Fortran buffer buf using memcpy. This implies that the „holes“ in the compound are filled with something already contained in buf_c.

As I understand it, a solution to the problem would look like this:

  1. memcpy buf to buf_c
  2. read the compound data from file
  3. fill the declared fields in buf_c
  4. memcpy buf_c to buf

In step 3 it is important to fill only the fields that were defined with h5tinsert_f when constructing the HDF5 compound.

Joerg

Am 20.06.2024 um 00:06 schrieb Scot Breitenfeld @.***>:

Unfortunately, what you are proposing is not supported in HDF5. HDF5 does not support I/O having a struct with a field that is not allocated unless that field is declared a vl type. You need to inform HDF5 what fields are written in the compound type. Currently, HDF5 expects one field to be written, but you are writing another field, which can corrupt other data in the file.

The H5Dread API does not allocate anything; it only fills the buffer being passed in, and Fortran misinterprets what C is reading and passing back. I don't think the link you provided is relevant in this case.

You can either use fixed arrays in the type, which I'm assuming you don't want to do if this is a list of nodes associated with different kinds of elements, or use an array of integers with unlimited dimension, pack all the variables in the array and store another array with element lengths.

— Reply to this email directly, view it on GitHub https://github.com/HDFGroup/hdf5/issues/4577#issuecomment-2179509628, or unsubscribe https://github.com/notifications/unsubscribe-auth/AKDBWQNBHE6L2AGIVS5AIWLZIH6GDAVCNFSM6AAAAABJQDFH5GVHI2DSMVQWIX3LMV43OSLTON2WKQ3PNVWWK3TUHMZDCNZZGUYDSNRSHA. You are receiving this because you authored the thread.

HiSPEET avatar Jun 20 '24 06:06 HiSPEET

To illustrate the workaround with CFI_deallocate consider the appended example. The basic idea is as follows

  interface
    ! interface to CFI_deallocate from ISO Fortran binding
    integer function CFI_deallocate(a)  bind(C, name='CFI_deallocate')
      integer, allocatable :: a(:)
    end function CFI_deallocate
  end interface
  ...
  ! read static components
  ...
  ! fix dynamic components
  do i = 1, size(elem_read)
    if (allocated(elem_read(i)%nb)) then
      err = CFI_deallocate(elem_read(i)%nb)
    end if
  end do

The modified test program gives the output

check static components: differences should be zero
   elem_read%id - elem_orig%id = 0 0

check dynamic components: should be deallocated
  allocated(elem_read(1)%nb) = F
  allocated(elem_read(2)%nb) = T

dynamic components corrupted, trying to fix
  allocated(elem_read(1)%nb) = F
  allocated(elem_read(2)%nb) = F

The advantage over Fortran's deallocate is that CFI_deallocate succeeds also with corrupted array descriptors or deallocated arrays.

Nevertheless I hope that a clean solution can be provided within HDF5.

Joerg

hdf5_test.tgz

HiSPEET avatar Jun 20 '24 11:06 HiSPEET

No memory copies are happening; the code path is h5dread_ptr -> h5dread_f_c ->H5Dread. All the data (including nd) you allocate is written in the hdf5 file because you set the size in the H5Tcreate_f call. You can see this by swapping out H5T_COMPOUND_F with H5T_OPAQUE_F and forgoing the H5Tinsert. You will get the same results.

There is no place for an introspect of the data to mess with allocation states, as you've done in your example. This would have to be done on the application's side. The allocation state is an artifact of the C/Fortran interoperability, as you can see from the attached C/Fortran code, which loosely mimics what is happening in the HDF5 library. The Fortran standard may not address this scenario, so I'm unsure if it is not compiler-dependent.

ex.tar.gz

brtnfld avatar Jun 20 '24 20:06 brtnfld

Hi Scot,

All the data (including nd) you allocate is written in the hdf5 file

I believe there is still an misunderstanding: a dump of the HDF5 file produced in my example gives no indication that the dynamic component (nb) is actually written:

HDF5 "elements.h5" {
GROUP "/" {
   GROUP "data" {
      DATASET "elements" {
         DATATYPE  H5T_COMPOUND {
            H5T_STD_I32LE "id";
         }
         DATASPACE  SIMPLE { ( 2 ) / ( 2 ) }
         DATA {
         (0): {
               1
            },
         (1): {
               2
            }
         }
      }
   }
}
}

The values of elem_orig(2)%nb, namely 10 and 20, are not written. However, somehow the allocation status of elem_orig%nb is restored in elem_read%nb.

Allocatable arrays of a type that has the attribute bind(C) can be reset from the ISO-Fortran binding with CFI_deallocate. However, this excludes types with character components of length > 1 or with type-bound procedures and thus a considerable loss of functionality. A solution to the problem within HDF5 is therefore clearly preferable.

HiSPEET avatar Jun 20 '24 21:06 HiSPEET

h5dump will only report registered members in the compound type, so it does not report nb values.

brtnfld avatar Jun 20 '24 21:06 brtnfld

h5dump will only report registered members in the compound type, so it does not report nb values.

So the descriptors are written and then read, which explains the observed behavior. I suspected this and tried to prevent writing unregistered data by calling h5tpack_f, but with no effect.

HiSPEET avatar Jun 20 '24 22:06 HiSPEET