Large array data leads to chunk size exhaustion

hmenke commented 2 years ago

Prerequisites

[x] Please check that a similar issue isn't already filed: https://github.com/issues?q=is%3Aissue+user%3Atriqs

Description

For very large data like lattice susceptibilities the chunk size of the HDF5 data may exceed 4 GB which is not allowed and will crash.

Chunk maximums

The maximum number of elements in a chunk is 2³²-1 which is equal to 4,294,967,295. If the number of elements in a chunk is set via H5Pset_chunk to a value greater than 2³²-1, then H5Pset_chunk will fail.

Steps to Reproduce

from triqs.gf import MeshBrZone, MeshImFreq, MeshProduct, Gf
from triqs.lattice import BravaisLattice, BrillouinZone
from h5 import HDFArchive

BL = BravaisLattice(units = [(1,0,0) , (0,1,0)])
iW_mesh = MeshImFreq(beta=4.0, S='Boson', n_max=1)
iw_mesh = MeshImFreq(beta=4.0, S='Fermion', n_max=80)
br_zone = MeshBrZone(BrillouinZone(BL), 40)
prod_mesh = MeshProduct(iW_mesh, iw_mesh, br_zone)

chi0_wnk = Gf(mesh=prod_mesh, indices=[list(range(6))]*4)

with HDFArchive('/tmp/trash.h5', 'a') as ar:
    ar['chi0_wnk'] = chi0_wnk

Expected behavior: No error.

Actual behavior:

HDF5-DIAG: Error detected in HDF5 (1.10.7) thread 1:
  #000: H5D.c line 152 in H5Dcreate2(): unable to create dataset
    major: Dataset
    minor: Unable to initialize object
  #001: H5Dint.c line 338 in H5D__create_named(): unable to create and link to dataset
    major: Dataset
    minor: Unable to initialize object
  #002: H5L.c line 1605 in H5L_link_object(): unable to create new link to object
    major: Links
    minor: Unable to initialize object
  #003: H5L.c line 1846 in H5L__create_real(): can't insert link
    major: Links
    minor: Unable to insert object
  #004: H5Gtraverse.c line 848 in H5G_traverse(): internal path traversal failed
    major: Symbol table
    minor: Object not found
  #005: H5Gtraverse.c line 624 in H5G__traverse_real(): traversal operator failed
    major: Symbol table
    minor: Callback failed
  #006: H5L.c line 1652 in H5L__link_cb(): unable to create object
    major: Links
    minor: Unable to initialize object
  #007: H5Oint.c line 2496 in H5O_obj_create(): unable to open object
    major: Object header
    minor: Can't open object
  #008: H5Doh.c line 300 in H5O__dset_create(): unable to create dataset
    major: Dataset
    minor: Unable to initialize object
  #009: H5Dint.c line 1307 in H5D__create(): unable to construct layout information
    major: Dataset
  #010: H5Dchunk.c line 892 in H5D__chunk_construct(): unable to set chunk sizes
    major: Dataset
    minor: Bad value
  #011: H5Dchunk.c line 851 in H5D__chunk_set_sizes(): chunk size must be < 4GB
    major: Dataset
    minor: Unable to initialize object
HDFArchive is in trouble with the array [...]
Traceback (most recent call last):
  File "test.py", line 16, in <module>
    ar['chi0_wnk'] = chi0_wnk
  File "/home/menke/Code/triqs/triqs_unstable/install/lib/python3.6/site-packages/h5/archive.py", line 171, in __setitem__
    for k, v in list(d.items()) : SubGroup[k] = v
  File "/home/menke/Code/triqs/triqs_unstable/install/lib/python3.6/site-packages/h5/archive.py", line 175, in __setitem__
    self._write( key, numpy.array(val,copy=1,order='C') )
  File "/home/menke/Code/triqs/triqs_unstable/install/lib/python3.6/site-packages/h5/archive_basic_layer.py", line 80, in _write
    h5.h5_write(self._group, key, val)
RuntimeError: .. Error occurred at Fri May 13 18:02:14 2022

.. Error .. calling C++ overload 
.. h5_write_bare(group g, std::string name, PyObject * ob) -> void 
.. in implementation of function _h5py.h5_write
.. C++ error was : 
Cannot create the dataset data in the group /chi0_wnk

Versions

d6081913d2d42b258dfe362d17933b3d178764d3

Formatting

Please use markdown in your issue message. A useful summary of commands can be found here.

Additional Information

Any additional information, configuration or data that might be necessary to reproduce the issue.

hmenke commented 2 years ago

Possible patch:

diff --git a/c++/h5/array_interface.cpp b/c++/h5/array_interface.cpp
index 79d4c1c..6254d6d 100644
--- a/c++/h5/array_interface.cpp
+++ b/c++/h5/array_interface.cpp
@@ -78,7 +78,10 @@ namespace h5::array_interface {
     if (compress and (v.rank() != 0)) {
       int n_dims = v.rank();
       std::vector<hsize_t> chunk_dims(n_dims);
-      for (int i = 0; i < v.rank(); ++i) chunk_dims[i] = std::max(v.slab.count[i], hsize_t{1});
+      hsize_t const dtype_size = H5Tget_size(v.ty);
+      hsize_t const chunk_max  = hsize_t{1UL << 32} - hsize_t{1}; // 2^32 - 1 = 4 GB
+      hsize_t const root_dim   = std::pow(chunk_max / dtype_size, 1. / n_dims);
+      for (int i = 0; i < v.rank(); ++i) chunk_dims[i] = std::clamp(v.slab.count[i], hsize_t{1}, root_dim);
       cparms = H5Pcreate(H5P_DATASET_CREATE);
       H5Pset_chunk(cparms, n_dims, chunk_dims.data());
       H5Pset_deflate(cparms, 1);

Wentzell commented 2 years ago

Dear @hmenke,

Thank you for pointing out the issue and for proposing a solution! This is now fixed in 33b094d.

I have opted, when possible, to keep the chunk dimensions the same as the right-most dimensions in the array.

TRIQS / h5