Reference documentation for deal.II version Git f0d1c24e5f 2021-10-18 08:09:39 -0400
\(\newcommand{\dealvcentcolon}{\mathrel{\mathop{:}}}\) \(\newcommand{\dealcoloneq}{\dealvcentcolon\mathrel{\mkern-1.2mu}=}\) \(\newcommand{\jump}[1]{\left[\!\left[ #1 \right]\!\right]}\) \(\newcommand{\average}[1]{\left\{\!\left\{ #1 \right\}\!\right\}}\)
hdf5.h
Go to the documentation of this file.
1 // ---------------------------------------------------------------------
2 //
3 // Copyright (C) 2018 - 2021 by the deal.II authors
4 //
5 // This file is part of the deal.II library.
6 //
7 // The deal.II library is free software; you can use it, redistribute
8 // it, and/or modify it under the terms of the GNU Lesser General
9 // Public License as published by the Free Software Foundation; either
10 // version 2.1 of the License, or (at your option) any later version.
11 // The full text of the license can be found in the file LICENSE.md at
12 // the top level directory of deal.II.
13 //
14 // ---------------------------------------------------------------------
15 
16 #ifndef dealii_hdf5_h
17 #define dealii_hdf5_h
18 
19 #include <deal.II/base/config.h>
20 
21 #ifdef DEAL_II_WITH_HDF5
22 
23 # include <deal.II/base/array_view.h>
24 
25 # include <deal.II/lac/full_matrix.h>
26 
27 # include <hdf5.h>
28 
29 # include <numeric>
30 
32 
33 // It is necessary to turn clang-format off in order to maintain the Doxygen
34 // links because they are longer than 80 characters
35 // clang-format off
344 // clang-format on
345 namespace HDF5
346 {
351  {
352  protected:
357  HDF5Object(const std::string &name, const bool mpi);
358 
359  public:
372  template <typename T>
373  T
374  get_attribute(const std::string &attr_name) const;
375 
388  template <typename T>
389  void
390  set_attribute(const std::string &attr_name, const T value);
391 
397  std::string
398  get_name() const;
399 
400  protected:
406  const std::string name;
407 
415  std::shared_ptr<hid_t> hdf5_reference;
416 
420  const bool mpi;
421  };
422 
426  class DataSet : public HDF5Object
427  {
428  friend class Group;
429 
430  protected:
435  DataSet(const std::string &name, const hid_t &parent_group_id, bool mpi);
436 
441  DataSet(const std::string & name,
442  const hid_t & parent_group_id,
443  const std::vector<hsize_t> & dimensions,
444  const std::shared_ptr<hid_t> &t_type,
445  const bool mpi);
446 
447  public:
465  template <typename Container>
466  Container
467  read();
468 
502  template <typename Container>
503  Container
504  read_selection(const std::vector<hsize_t> &coordinates);
505 
506  // clang-format off
540  // clang-format on
541  template <typename Container>
542  Container
543  read_hyperslab(const std::vector<hsize_t> &offset,
544  const std::vector<hsize_t> &count);
545 
578  template <typename Container>
579  Container
580  read_hyperslab(const std::vector<hsize_t> &data_dimensions,
581  const std::vector<hsize_t> &offset,
582  const std::vector<hsize_t> &stride,
583  const std::vector<hsize_t> &count,
584  const std::vector<hsize_t> &block);
585 
597  template <typename number>
598  void
599  read_none();
600 
619  template <typename Container>
620  void
621  write(const Container &data);
622 
651  template <typename Container>
652  void
653  write_selection(const Container & data,
654  const std::vector<hsize_t> &coordinates);
655 
656  // clang-format off
682  // clang-format on
683  template <typename Container>
684  void
685  write_hyperslab(const Container & data,
686  const std::vector<hsize_t> &offset,
687  const std::vector<hsize_t> &count);
688 
721  template <typename Container>
722  void
723  write_hyperslab(const Container & data,
724  const std::vector<hsize_t> &data_dimensions,
725  const std::vector<hsize_t> &offset,
726  const std::vector<hsize_t> &stride,
727  const std::vector<hsize_t> &count,
728  const std::vector<hsize_t> &block);
729 
747  template <typename number>
748  void
749  write_none();
750 
767  bool
768  get_query_io_mode() const;
769 
773  void
774  set_query_io_mode(const bool new_query_io_mode);
775 
790  std::string
791  get_io_mode();
792 
809  H5D_mpio_actual_io_mode_t
810  get_io_mode_as_hdf5_type();
811 
830  std::string
831  get_local_no_collective_cause();
832 
853  uint32_t
854  get_local_no_collective_cause_as_hdf5_type();
855 
874  std::string
875  get_global_no_collective_cause();
876 
897  uint32_t
898  get_global_no_collective_cause_as_hdf5_type();
899 
905  std::vector<hsize_t>
906  get_dimensions() const;
907 
911  unsigned int
912  get_size() const;
913 
917  unsigned int
918  get_rank() const;
919 
920  private:
924  unsigned int rank;
925 
930  std::vector<hsize_t> dimensions;
931 
935  std::shared_ptr<hid_t> dataspace;
936 
940  unsigned int size;
941 
952 
956  H5D_mpio_actual_io_mode_t io_mode;
957 
964 
971  };
972 
976  class Group : public HDF5Object
977  {
978  protected:
982  enum class GroupAccessMode
983  {
987  open,
991  create
992  };
1001  Group(const std::string & name,
1002  const Group & parent_group,
1003  const bool mpi,
1004  const GroupAccessMode mode);
1005 
1011  Group(const std::string &name, const bool mpi);
1012 
1013  public:
1017  Group
1018  open_group(const std::string &name) const;
1019 
1023  Group
1024  create_group(const std::string &name) const;
1025 
1029  DataSet
1030  open_dataset(const std::string &name) const;
1031 
1042  template <typename number>
1043  DataSet
1044  create_dataset(const std::string & name,
1045  const std::vector<hsize_t> &dimensions) const;
1046 
1065  template <typename Container>
1066  void
1067  write_dataset(const std::string &name, const Container &data) const;
1068  };
1069 
1073  class File : public Group
1074  {
1075  public:
1079  enum class FileAccessMode
1080  {
1084  open,
1088  create
1089  };
1090 
1096  File(const std::string &name, const FileAccessMode mode);
1097 
1105  File(const std::string & name,
1106  const FileAccessMode mode,
1107  const MPI_Comm & mpi_communicator);
1108 
1109  private:
1117  File(const std::string & name,
1118  const FileAccessMode mode,
1119  const bool mpi,
1120  const MPI_Comm & mpi_communicator);
1121  };
1122 
1123  namespace internal
1124  {
1135  template <typename number>
1136  std::shared_ptr<hid_t>
1138 
1147  template <typename number>
1148  std::vector<hsize_t>
1149  get_container_dimensions(const std::vector<number> &data);
1150 
1154  template <typename number>
1155  std::vector<hsize_t>
1157 
1161  template <typename number>
1162  std::vector<hsize_t>
1164 
1168  template <typename number>
1169  unsigned int
1170  get_container_size(const std::vector<number> &data);
1171 
1175  template <typename number>
1176  unsigned int
1177  get_container_size(const Vector<number> &data);
1178 
1182  template <typename number>
1183  unsigned int
1185 
1203  template <typename Container>
1204  typename std::enable_if<
1205  std::is_same<Container,
1206  std::vector<typename Container::value_type>>::value,
1207  Container>::type
1208  initialize_container(const std::vector<hsize_t> &dimensions);
1209 
1212  template <typename Container>
1213  typename std::enable_if<
1214  std::is_same<Container, Vector<typename Container::value_type>>::value,
1215  Container>::type
1216  initialize_container(const std::vector<hsize_t> &dimensions);
1217 
1220  template <typename Container>
1221  typename std::enable_if<
1222  std::is_same<Container,
1224  Container>::type
1225  initialize_container(const std::vector<hsize_t> &dimensions);
1226 
1232  inline void
1233  set_plist(hid_t &plist, const bool mpi);
1234 
1242  inline void
1243  release_plist(hid_t & plist,
1244  H5D_mpio_actual_io_mode_t &io_mode,
1245  uint32_t & local_no_collective_cause,
1246  uint32_t & global_no_collective_cause,
1247  const bool mpi,
1248  const bool query_io_mode);
1249 
1252  inline std::string
1253  no_collective_cause_to_string(const uint32_t no_collective_cause);
1254  } // namespace internal
1255 
1256 
1257 
1258  // definitions
1259 
1260  namespace internal
1261  {
1262  template <typename number>
1263  std::shared_ptr<hid_t>
1265  {
1266  static_assert(std::is_same<number, float>::value ||
1267  std::is_same<number, double>::value ||
1268  std::is_same<number, int>::value ||
1269  std::is_same<number, bool>::value ||
1270  std::is_same<number, unsigned int>::value ||
1271  std::is_same<number, std::complex<float>>::value ||
1272  std::is_same<number, std::complex<double>>::value,
1273  "The data type you are trying to get the HDF5 tag for "
1274  "is not supported by this function.");
1275 
1276  if (std::is_same<number, float>::value)
1277  {
1278  return std::make_shared<hid_t>(H5T_NATIVE_FLOAT);
1279  }
1280  else if (std::is_same<number, double>::value)
1281  {
1282  return std::make_shared<hid_t>(H5T_NATIVE_DOUBLE);
1283  }
1284  else if (std::is_same<number, int>::value)
1285  {
1286  return std::make_shared<hid_t>(H5T_NATIVE_INT);
1287  }
1288  else if (std::is_same<number, unsigned int>::value)
1289  {
1290  return std::make_shared<hid_t>(H5T_NATIVE_UINT);
1291  }
1292  else if (std::is_same<number, bool>::value)
1293  {
1294  return std::make_shared<hid_t>(H5T_NATIVE_HBOOL);
1295  }
1296  else if (std::is_same<number, std::complex<float>>::value)
1297  {
1298  std::shared_ptr<hid_t> t_type =
1299  std::shared_ptr<hid_t>(new hid_t, [](hid_t *pointer) {
1300  // Release the HDF5 resource
1301  const herr_t ret = H5Tclose(*pointer);
1302  AssertNothrow(ret >= 0, ExcInternalError());
1303  (void)ret;
1304  delete pointer;
1305  });
1306 
1307  *t_type = H5Tcreate(H5T_COMPOUND, sizeof(std::complex<float>));
1308  // The C++ standards committee agreed to mandate that the storage
1309  // format used for the std::complex type be binary-compatible with
1310  // the C99 type, i.e. an array T[2] with consecutive real [0] and
1311  // imaginary [1] parts.
1312  herr_t ret = H5Tinsert(*t_type, "r", 0, H5T_NATIVE_FLOAT);
1313  Assert(ret >= 0, ExcInternalError());
1314  ret = H5Tinsert(*t_type, "i", sizeof(float), H5T_NATIVE_FLOAT);
1315  Assert(ret >= 0, ExcInternalError());
1316  (void)ret;
1317  return t_type;
1318  }
1319  else if (std::is_same<number, std::complex<double>>::value)
1320  {
1321  std::shared_ptr<hid_t> t_type =
1322  std::shared_ptr<hid_t>(new hid_t, [](hid_t *pointer) {
1323  // Release the HDF5 resource
1324  const herr_t ret = H5Tclose(*pointer);
1325  AssertNothrow(ret >= 0, ExcInternalError());
1326  (void)ret;
1327  delete pointer;
1328  });
1329  *t_type = H5Tcreate(H5T_COMPOUND, sizeof(std::complex<double>));
1330  // The C++ standards committee agreed to mandate that the storage
1331  // format used for the std::complex type be binary-compatible with
1332  // the C99 type, i.e. an array T[2] with consecutive real [0] and
1333  // imaginary [1] parts.
1334  herr_t ret = H5Tinsert(*t_type, "r", 0, H5T_NATIVE_DOUBLE);
1335  Assert(ret >= 0, ExcInternalError());
1336  ret = H5Tinsert(*t_type, "i", sizeof(double), H5T_NATIVE_DOUBLE);
1337  Assert(ret >= 0, ExcInternalError());
1338  (void)ret;
1339  return t_type;
1340  }
1341 
1342  // The function should not reach this point
1343  Assert(false, ExcInternalError());
1344  return {};
1345  }
1346 
1347 
1348 
1349  template <typename number>
1350  std::vector<hsize_t>
1351  get_container_dimensions(const std::vector<number> &data)
1352  {
1353  std::vector<hsize_t> dimensions = {data.size()};
1354  return dimensions;
1355  }
1356 
1357 
1358 
1359  template <typename number>
1360  std::vector<hsize_t>
1362  {
1363  std::vector<hsize_t> dimensions = {data.size()};
1364  return dimensions;
1365  }
1366 
1367 
1368 
1369  template <typename number>
1370  std::vector<hsize_t>
1372  {
1373  std::vector<hsize_t> dimensions = {data.m(), data.n()};
1374  return dimensions;
1375  }
1376 
1377 
1378 
1379  template <typename number>
1380  unsigned int
1381  get_container_size(const std::vector<number> &data)
1382  {
1383  return static_cast<unsigned int>(data.size());
1384  }
1385 
1386 
1387 
1388  template <typename number>
1389  unsigned int
1391  {
1392  return static_cast<unsigned int>(data.size());
1393  }
1394 
1395 
1396 
1397  template <typename number>
1398  unsigned int
1400  {
1401  return static_cast<unsigned int>(data.m() * data.n());
1402  }
1403 
1404 
1405 
1406  template <typename Container>
1407  typename std::enable_if<
1408  std::is_same<Container,
1409  std::vector<typename Container::value_type>>::value,
1410  Container>::type
1411  initialize_container(const std::vector<hsize_t> &dimensions)
1412  {
1413  return Container(std::accumulate(
1414  dimensions.begin(), dimensions.end(), 1, std::multiplies<int>()));
1415  }
1416 
1417 
1418 
1419  template <typename Container>
1420  typename std::enable_if<
1421  std::is_same<Container, Vector<typename Container::value_type>>::value,
1422  Container>::type
1423  initialize_container(const std::vector<hsize_t> &dimensions)
1424  {
1425  return Container(std::accumulate(
1426  dimensions.begin(), dimensions.end(), 1, std::multiplies<int>()));
1427  }
1428 
1429 
1430 
1431  template <typename Container>
1432  typename std::enable_if<
1433  std::is_same<Container,
1435  Container>::type
1436  initialize_container(const std::vector<hsize_t> &dimensions)
1437  {
1438  // If the rank is higher than 2, then remove single-dimensional entries
1439  // from the shape defined by dimensions. This is equivalent to the squeeze
1440  // function of python/numpy. For example the following code would convert
1441  // the vector {1,3,1,2} to {3,2}
1442  std::vector<hsize_t> squeezed_dimensions;
1443 
1444  if (dimensions.size() > 2)
1445  {
1446  for (const auto &dimension : dimensions)
1447  {
1448  if (dimension > 1)
1449  squeezed_dimensions.push_back(dimension);
1450  }
1451  }
1452  else
1453  {
1454  squeezed_dimensions = dimensions;
1455  }
1456 
1457  AssertDimension(squeezed_dimensions.size(), 2);
1458  return Container(squeezed_dimensions[0], squeezed_dimensions[1]);
1459  }
1460 
1461 
1462  inline void
1463  set_plist(hid_t &plist, const bool mpi)
1464  {
1465  if (mpi)
1466  {
1467 # ifdef DEAL_II_WITH_MPI
1468  plist = H5Pcreate(H5P_DATASET_XFER);
1469  Assert(plist >= 0, ExcInternalError());
1470  const herr_t ret = H5Pset_dxpl_mpio(plist, H5FD_MPIO_COLLECTIVE);
1471  (void)ret;
1472  Assert(ret >= 0, ExcInternalError());
1473 # else
1474  AssertThrow(false, ExcNotImplemented());
1475 # endif
1476  }
1477  else
1478  {
1479  plist = H5P_DEFAULT;
1480  }
1481 
1482  (void)plist;
1483  (void)mpi;
1484  }
1485 
1486 
1487  inline void
1488  release_plist(hid_t & plist,
1489  H5D_mpio_actual_io_mode_t &io_mode,
1490  uint32_t & local_no_collective_cause,
1491  uint32_t & global_no_collective_cause,
1492  const bool mpi,
1493  const bool query_io_mode)
1494  {
1495  if (mpi)
1496  {
1497 # ifdef DEAL_II_WITH_MPI
1498  herr_t ret;
1499  (void)ret;
1500  if (query_io_mode)
1501  {
1502  ret = H5Pget_mpio_actual_io_mode(plist, &io_mode);
1503  Assert(ret >= 0, ExcInternalError());
1504  ret =
1505  H5Pget_mpio_no_collective_cause(plist,
1506  &local_no_collective_cause,
1507  &global_no_collective_cause);
1508  Assert(ret >= 0, ExcInternalError());
1509  }
1510  ret = H5Pclose(plist);
1511  Assert(ret >= 0, ExcInternalError());
1512 # else
1513  AssertThrow(false, ExcNotImplemented());
1514 # endif
1515  }
1516 
1517  (void)plist;
1518  (void)io_mode;
1519  (void)local_no_collective_cause;
1520  (void)global_no_collective_cause;
1521  (void)mpi;
1522  (void)query_io_mode;
1523  }
1524 
1525 
1526  inline std::string
1527  no_collective_cause_to_string(const uint32_t no_collective_cause)
1528  {
1529  std::string message;
1530 
1531  auto append_to_message = [&message](const char *p) {
1532  if (message.length() > 0)
1533  message += ", ";
1534  message += p;
1535  };
1536 
1537  // The first is not a bitmask comparison, the rest are bitmask
1538  // comparisons.
1539  // https://support.hdfgroup.org/HDF5/doc/RM/RM_H5P.html#Property-GetMpioNoCollectiveCause
1540  // See H5Ppublic.h
1541  // Hex codes are used because the HDF5 Group can deprecate some of the
1542  // enum codes. For example the enum code H5D_MPIO_FILTERS is not defined
1543  // in 1.10.2 because it is possible to use compressed datasets with the
1544  // MPI/IO driver.
1545 
1546  // H5D_MPIO_COLLECTIVE
1547  if (no_collective_cause == 0x00)
1548  {
1549  append_to_message("H5D_MPIO_COLLECTIVE");
1550  }
1551  // H5D_MPIO_SET_INDEPENDENT
1552  if ((no_collective_cause & 0x01) == 0x01)
1553  {
1554  append_to_message("H5D_MPIO_SET_INDEPENDENT");
1555  }
1556  // H5D_MPIO_DATATYPE_CONVERSION
1557  if ((no_collective_cause & 0x02) == 0x02)
1558  {
1559  append_to_message("H5D_MPIO_DATATYPE_CONVERSION");
1560  }
1561  // H5D_MPIO_DATA_TRANSFORMS
1562  if ((no_collective_cause & 0x04) == 0x04)
1563  {
1564  append_to_message("H5D_MPIO_DATA_TRANSFORMS");
1565  }
1566  // H5D_MPIO_NOT_SIMPLE_OR_SCALAR_DATASPACES
1567  if ((no_collective_cause & 0x10) == 0x10)
1568  {
1569  append_to_message("H5D_MPIO_NOT_SIMPLE_OR_SCALAR_DATASPACES");
1570  }
1571  // H5D_MPIO_NOT_CONTIGUOUS_OR_CHUNKED_DATASET
1572  if ((no_collective_cause & 0x20) == 0x20)
1573  {
1574  append_to_message("H5D_MPIO_NOT_CONTIGUOUS_OR_CHUNKED_DATASET");
1575  }
1576  // H5D_MPIO_FILTERS
1577  if ((no_collective_cause & 0x40) == 0x40)
1578  {
1579  append_to_message("H5D_MPIO_FILTERS");
1580  }
1581  return message;
1582  }
1583  } // namespace internal
1584 
1585 
1586  template <typename T>
1587  T
1588  HDF5Object::get_attribute(const std::string &attr_name) const
1589  {
1590  const std::shared_ptr<hid_t> t_type = internal::get_hdf5_datatype<T>();
1591  T value;
1592  hid_t attr;
1593  herr_t ret;
1594 
1595  attr = H5Aopen(*hdf5_reference, attr_name.data(), H5P_DEFAULT);
1596  Assert(attr >= 0, ExcMessage("Error at H5Aopen"));
1597  (void)ret;
1598  ret = H5Aread(attr, *t_type, &value);
1599  Assert(ret >= 0, ExcMessage("Error at H5Aread"));
1600  (void)ret;
1601  ret = H5Aclose(attr);
1602  Assert(ret >= 0, ExcMessage("Error at H5Aclose"));
1603 
1604  return value;
1605  }
1606 
1607 
1608 
1609  template <>
1610  inline std::string
1611  HDF5Object::get_attribute(const std::string &attr_name) const
1612  {
1613  // Reads a UTF8 variable string
1614  //
1615  // code inspired from
1616  // https://support.hdfgroup.org/ftp/HDF5/examples/misc-examples/vlstratt.c
1617  //
1618  // In the case of a variable length string the user does not have to reserve
1619  // memory for string_out. H5Aread will reserve the memory and the
1620  // user has to free the memory.
1621  //
1622  // Todo:
1623  // - Use H5Dvlen_reclaim instead of free
1624 
1625  char * string_out;
1626  hid_t attr;
1627  hid_t type;
1628  herr_t ret;
1629 
1630  /* Create a datatype to refer to. */
1631  type = H5Tcopy(H5T_C_S1);
1632  Assert(type >= 0, ExcInternalError());
1633 
1634  // Python strings are encoded in UTF8
1635  ret = H5Tset_cset(type, H5T_CSET_UTF8);
1636  Assert(type >= 0, ExcInternalError());
1637 
1638  ret = H5Tset_size(type, H5T_VARIABLE);
1639  Assert(ret >= 0, ExcInternalError());
1640 
1641  attr = H5Aopen(*hdf5_reference, attr_name.data(), H5P_DEFAULT);
1642  Assert(attr >= 0, ExcInternalError());
1643 
1644  ret = H5Aread(attr, type, &string_out);
1645  Assert(ret >= 0, ExcInternalError());
1646 
1647  std::string string_value(string_out);
1648  // The memory of the variable length string has to be freed.
1649  // H5Dvlen_reclaim could be also used
1650  free(string_out);
1651  ret = H5Tclose(type);
1652  Assert(ret >= 0, ExcInternalError());
1653 
1654  ret = H5Aclose(attr);
1655  Assert(ret >= 0, ExcInternalError());
1656 
1657  (void)ret;
1658  return string_value;
1659  }
1660 
1661 
1662 
1663  template <typename T>
1664  void
1665  HDF5Object::set_attribute(const std::string &attr_name, const T value)
1666  {
1667  hid_t attr;
1668  hid_t aid;
1669  herr_t ret;
1670 
1671  const std::shared_ptr<hid_t> t_type = internal::get_hdf5_datatype<T>();
1672 
1673 
1674  /*
1675  * Create scalar attribute.
1676  */
1677  aid = H5Screate(H5S_SCALAR);
1678  Assert(aid >= 0, ExcMessage("Error at H5Screate"));
1679  attr = H5Acreate2(*hdf5_reference,
1680  attr_name.data(),
1681  *t_type,
1682  aid,
1683  H5P_DEFAULT,
1684  H5P_DEFAULT);
1685  Assert(attr >= 0, ExcMessage("Error at H5Acreate2"));
1686 
1687  /*
1688  * Write scalar attribute.
1689  */
1690  ret = H5Awrite(attr, *t_type, &value);
1691  Assert(ret >= 0, ExcMessage("Error at H5Awrite"));
1692 
1693  ret = H5Sclose(aid);
1694  Assert(ret >= 0, ExcMessage("Error at H5Sclose"));
1695  ret = H5Aclose(attr);
1696  Assert(ret >= 0, ExcMessage("Error at H5Aclose"));
1697 
1698  (void)ret;
1699  }
1700 
1701 
1702 
1703  template <>
1704  inline void
1705  HDF5Object::set_attribute(const std::string &attr_name,
1706  const std::string value) // NOLINT
1707  {
1708  // Writes a UTF8 variable string
1709  //
1710  // code inspired from
1711  // https://support.hdfgroup.org/ftp/HDF5/examples/misc-examples/vlstratt.c
1712 
1713  hid_t attr;
1714  hid_t aid;
1715  hid_t t_type;
1716  herr_t ret;
1717 
1718  /* Create a datatype to refer to. */
1719  t_type = H5Tcopy(H5T_C_S1);
1720  Assert(t_type >= 0, ExcInternalError());
1721 
1722  // Python strings are encoded in UTF8
1723  ret = H5Tset_cset(t_type, H5T_CSET_UTF8);
1724  Assert(t_type >= 0, ExcInternalError());
1725 
1726  ret = H5Tset_size(t_type, H5T_VARIABLE);
1727  Assert(ret >= 0, ExcInternalError());
1728 
1729  /*
1730  * Create scalar attribute.
1731  */
1732  aid = H5Screate(H5S_SCALAR);
1733  Assert(aid >= 0, ExcMessage("Error at H5Screate"));
1734  attr = H5Acreate2(
1735  *hdf5_reference, attr_name.data(), t_type, aid, H5P_DEFAULT, H5P_DEFAULT);
1736  Assert(attr >= 0, ExcMessage("Error at H5Acreate2"));
1737 
1738  /*
1739  * Write scalar attribute.
1740  * In most of the cases H5Awrite and H5Dwrite take a pointer to the data.
1741  * But in the particular case of a variable length string, H5Awrite takes
1742  * the address of the pointer of the string.
1743  */
1744  const char *c_string_value = value.c_str();
1745  ret = H5Awrite(attr, t_type, &c_string_value);
1746  Assert(ret >= 0, ExcInternalError());
1747 
1748  ret = H5Sclose(aid);
1749  Assert(ret >= 0, ExcMessage("Error at H5Sclose"));
1750  ret = H5Aclose(attr);
1751  Assert(ret >= 0, ExcMessage("Error at H5Aclose"));
1752 
1753  (void)ret;
1754  }
1755 
1756 
1757 
1758  template <typename Container>
1759  Container
1761  {
1762  const std::shared_ptr<hid_t> t_type =
1763  internal::get_hdf5_datatype<typename Container::value_type>();
1764  hid_t plist;
1765  herr_t ret;
1766 
1767  Container data = internal::initialize_container<Container>(dimensions);
1768 
1769  internal::set_plist(plist, mpi);
1770 
1771  ret = H5Dread(*hdf5_reference,
1772  *t_type,
1773  H5S_ALL,
1774  H5S_ALL,
1775  plist,
1776  make_array_view(data).data());
1777  Assert(ret >= 0, ExcInternalError());
1778 
1779 
1781  io_mode,
1782  local_no_collective_cause,
1783  global_no_collective_cause,
1784  mpi,
1785  query_io_mode);
1786 
1787  (void)ret;
1788  return data;
1789  }
1790 
1791 
1792 
1793  template <typename Container>
1794  Container
1795  DataSet::read_selection(const std::vector<hsize_t> &coordinates)
1796  {
1797  Assert(coordinates.size() % rank == 0,
1798  ExcMessage(
1799  "The dimension of coordinates has to be divisible by the rank"));
1800  const std::shared_ptr<hid_t> t_type =
1801  internal::get_hdf5_datatype<typename Container::value_type>();
1802  hid_t plist;
1803  hid_t memory_dataspace;
1804  herr_t ret;
1805 
1806  std::vector<hsize_t> data_dimensions{
1807  static_cast<hsize_t>(coordinates.size() / rank)};
1808 
1809  Container data = internal::initialize_container<Container>(data_dimensions);
1810 
1811  memory_dataspace = H5Screate_simple(1, data_dimensions.data(), nullptr);
1812  Assert(memory_dataspace >= 0, ExcMessage("Error at H5Screate_simple"));
1813  ret = H5Sselect_elements(*dataspace,
1814  H5S_SELECT_SET,
1815  data.size(),
1816  coordinates.data());
1817  Assert(ret >= 0, ExcMessage("Error at H5Sselect_elements"));
1818 
1819  internal::set_plist(plist, mpi);
1820 
1821  ret = H5Dread(*hdf5_reference,
1822  *t_type,
1823  memory_dataspace,
1824  *dataspace,
1825  plist,
1826  make_array_view(data).data());
1827  Assert(ret >= 0, ExcMessage("Error at H5Dread"));
1828 
1830  io_mode,
1831  local_no_collective_cause,
1832  global_no_collective_cause,
1833  mpi,
1834  query_io_mode);
1835 
1836  ret = H5Sclose(memory_dataspace);
1837  Assert(ret >= 0, ExcMessage("Error at H5SClose"));
1838 
1839  (void)ret;
1840  return data;
1841  }
1842 
1843 
1844 
1845  template <typename Container>
1846  Container
1847  DataSet::read_hyperslab(const std::vector<hsize_t> &offset,
1848  const std::vector<hsize_t> &count)
1849  {
1850  const std::shared_ptr<hid_t> t_type =
1851  internal::get_hdf5_datatype<typename Container::value_type>();
1852  hid_t plist;
1853  hid_t memory_dataspace;
1854  herr_t ret;
1855 
1856  // In this particular overload of read_hyperslab the data_dimensions are
1857  // the same as count
1858  std::vector<hsize_t> data_dimensions = count;
1859 
1860  Container data = internal::initialize_container<Container>(data_dimensions);
1861 
1862  memory_dataspace =
1863  H5Screate_simple(data_dimensions.size(), data_dimensions.data(), nullptr);
1864  Assert(memory_dataspace >= 0, ExcMessage("Error at H5Screate_simple"));
1865  ret = H5Sselect_hyperslab(*dataspace,
1866  H5S_SELECT_SET,
1867  offset.data(),
1868  nullptr,
1869  count.data(),
1870  nullptr);
1871  Assert(ret >= 0, ExcMessage("Error at H5Sselect_hyperslab"));
1872 
1873  internal::set_plist(plist, mpi);
1874 
1875  ret = H5Dread(*hdf5_reference,
1876  *t_type,
1877  memory_dataspace,
1878  *dataspace,
1879  plist,
1880  make_array_view(data).data());
1881  Assert(ret >= 0, ExcMessage("Error at H5Dread"));
1882 
1884  io_mode,
1885  local_no_collective_cause,
1886  global_no_collective_cause,
1887  mpi,
1888  query_io_mode);
1889 
1890  ret = H5Sclose(memory_dataspace);
1891  Assert(ret >= 0, ExcMessage("Error at H5SClose"));
1892 
1893  (void)ret;
1894  return data;
1895  }
1896 
1897 
1898 
1899  template <typename Container>
1900  Container
1901  DataSet::read_hyperslab(const std::vector<hsize_t> &data_dimensions,
1902  const std::vector<hsize_t> &offset,
1903  const std::vector<hsize_t> &stride,
1904  const std::vector<hsize_t> &count,
1905  const std::vector<hsize_t> &block)
1906  {
1907  const std::shared_ptr<hid_t> t_type =
1908  internal::get_hdf5_datatype<typename Container::value_type>();
1909  hid_t plist;
1910  hid_t memory_dataspace;
1911  herr_t ret;
1912 
1913  Container data = internal::initialize_container<Container>(data_dimensions);
1914 
1915  memory_dataspace =
1916  H5Screate_simple(data_dimensions.size(), data_dimensions.data(), nullptr);
1917  Assert(memory_dataspace >= 0, ExcMessage("Error at H5Screate_simple"));
1918  ret = H5Sselect_hyperslab(*dataspace,
1919  H5S_SELECT_SET,
1920  offset.data(),
1921  stride.data(),
1922  count.data(),
1923  block.data());
1924  Assert(ret >= 0, ExcMessage("Error at H5Sselect_hyperslab"));
1925 
1926  internal::set_plist(plist, mpi);
1927 
1928  ret = H5Dread(*hdf5_reference,
1929  *t_type,
1930  memory_dataspace,
1931  *dataspace,
1932  plist,
1933  make_array_view(data).data());
1934  Assert(ret >= 0, ExcMessage("Error at H5Dread"));
1935 
1937  io_mode,
1938  local_no_collective_cause,
1939  global_no_collective_cause,
1940  mpi,
1941  query_io_mode);
1942 
1943  ret = H5Sclose(memory_dataspace);
1944  Assert(ret >= 0, ExcMessage("Error at H5SClose"));
1945 
1946  (void)ret;
1947  return data;
1948  }
1949 
1950 
1951 
1952  template <typename number>
1953  void
1955  {
1956  const std::shared_ptr<hid_t> t_type = internal::get_hdf5_datatype<number>();
1957  const std::vector<hsize_t> data_dimensions = {0};
1958 
1959  hid_t memory_dataspace;
1960  hid_t plist;
1961  herr_t ret;
1962 
1963  memory_dataspace = H5Screate_simple(1, data_dimensions.data(), nullptr);
1964  Assert(memory_dataspace >= 0, ExcMessage("Error at H5Screate_simple"));
1965  ret = H5Sselect_none(*dataspace);
1966  Assert(ret >= 0, ExcMessage("H5Sselect_none"));
1967 
1968  internal::set_plist(plist, mpi);
1969 
1970  // The pointer of data can safely be nullptr, see the discussion at the HDF5
1971  // forum:
1972  // https://forum.hdfgroup.org/t/parallel-i-o-does-not-support-filters-yet/884/17
1973  ret = H5Dread(
1974  *hdf5_reference, *t_type, memory_dataspace, *dataspace, plist, nullptr);
1975  Assert(ret >= 0, ExcMessage("Error at H5Dread"));
1976 
1978  io_mode,
1979  local_no_collective_cause,
1980  global_no_collective_cause,
1981  mpi,
1982  query_io_mode);
1983 
1984  ret = H5Sclose(memory_dataspace);
1985  Assert(ret >= 0, ExcMessage("Error at H5SClose"));
1986 
1987  (void)ret;
1988  }
1989 
1990 
1991 
1992  template <typename Container>
1993  void
1994  DataSet::write(const Container &data)
1995  {
1997  const std::shared_ptr<hid_t> t_type =
1998  internal::get_hdf5_datatype<typename Container::value_type>();
1999  hid_t plist;
2000  herr_t ret;
2001 
2002  internal::set_plist(plist, mpi);
2003 
2004  ret = H5Dwrite(*hdf5_reference,
2005  *t_type,
2006  H5S_ALL,
2007  H5S_ALL,
2008  plist,
2009  make_array_view(data).data());
2010  Assert(ret >= 0, ExcMessage("Error at H5Dwrite"));
2011 
2013  io_mode,
2014  local_no_collective_cause,
2015  global_no_collective_cause,
2016  mpi,
2017  query_io_mode);
2018 
2019  (void)ret;
2020  }
2021 
2022 
2023 
2024  template <typename Container>
2025  void
2026  DataSet::write_selection(const Container & data,
2027  const std::vector<hsize_t> &coordinates)
2028  {
2029  AssertDimension(coordinates.size(), data.size() * rank);
2030  const std::shared_ptr<hid_t> t_type =
2031  internal::get_hdf5_datatype<typename Container::value_type>();
2032  const std::vector<hsize_t> data_dimensions =
2034 
2035  hid_t memory_dataspace;
2036  hid_t plist;
2037  herr_t ret;
2038 
2039 
2040  memory_dataspace = H5Screate_simple(1, data_dimensions.data(), nullptr);
2041  Assert(memory_dataspace >= 0, ExcMessage("Error at H5Screate_simple"));
2042  ret = H5Sselect_elements(*dataspace,
2043  H5S_SELECT_SET,
2044  data.size(),
2045  coordinates.data());
2046  Assert(ret >= 0, ExcMessage("Error at H5Sselect_elements"));
2047 
2048  internal::set_plist(plist, mpi);
2049 
2050  ret = H5Dwrite(*hdf5_reference,
2051  *t_type,
2052  memory_dataspace,
2053  *dataspace,
2054  plist,
2055  make_array_view(data).data());
2056  Assert(ret >= 0, ExcMessage("Error at H5Dwrite"));
2057 
2059  io_mode,
2060  local_no_collective_cause,
2061  global_no_collective_cause,
2062  mpi,
2063  query_io_mode);
2064 
2065  ret = H5Sclose(memory_dataspace);
2066  Assert(ret >= 0, ExcMessage("Error at H5SClose"));
2067 
2068  (void)ret;
2069  }
2070 
2071 
2072 
2073  template <typename Container>
2074  void
2075  DataSet::write_hyperslab(const Container & data,
2076  const std::vector<hsize_t> &offset,
2077  const std::vector<hsize_t> &count)
2078  {
2079  AssertDimension(std::accumulate(count.begin(),
2080  count.end(),
2081  1,
2082  std::multiplies<unsigned int>()),
2084  const std::shared_ptr<hid_t> t_type =
2085  internal::get_hdf5_datatype<typename Container::value_type>();
2086  // In this particular overload of write_hyperslab the data_dimensions are
2087  // the same as count
2088  const std::vector<hsize_t> &data_dimensions = count;
2089 
2090  hid_t memory_dataspace;
2091  hid_t plist;
2092  herr_t ret;
2093 
2094  memory_dataspace =
2095  H5Screate_simple(data_dimensions.size(), data_dimensions.data(), nullptr);
2096  Assert(memory_dataspace >= 0, ExcMessage("Error at H5Screate_simple"));
2097  ret = H5Sselect_hyperslab(*dataspace,
2098  H5S_SELECT_SET,
2099  offset.data(),
2100  nullptr,
2101  count.data(),
2102  nullptr);
2103  Assert(ret >= 0, ExcMessage("Error at H5Sselect_hyperslab"));
2104 
2105  internal::set_plist(plist, mpi);
2106 
2107  ret = H5Dwrite(*hdf5_reference,
2108  *t_type,
2109  memory_dataspace,
2110  *dataspace,
2111  plist,
2112  make_array_view(data).data());
2113  Assert(ret >= 0, ExcMessage("Error at H5Dwrite"));
2114 
2116  io_mode,
2117  local_no_collective_cause,
2118  global_no_collective_cause,
2119  mpi,
2120  query_io_mode);
2121 
2122  ret = H5Sclose(memory_dataspace);
2123  Assert(ret >= 0, ExcMessage("Error at H5Sclose"));
2124 
2125  (void)ret;
2126  }
2127 
2128 
2129 
2130  template <typename Container>
2131  void
2132  DataSet::write_hyperslab(const Container & data,
2133  const std::vector<hsize_t> &data_dimensions,
2134  const std::vector<hsize_t> &offset,
2135  const std::vector<hsize_t> &stride,
2136  const std::vector<hsize_t> &count,
2137  const std::vector<hsize_t> &block)
2138  {
2139  const std::shared_ptr<hid_t> t_type =
2140  internal::get_hdf5_datatype<typename Container::value_type>();
2141 
2142  hid_t memory_dataspace;
2143  hid_t plist;
2144  herr_t ret;
2145 
2146  memory_dataspace =
2147  H5Screate_simple(data_dimensions.size(), data_dimensions.data(), nullptr);
2148  Assert(memory_dataspace >= 0, ExcMessage("Error at H5Screate_simple"));
2149  ret = H5Sselect_hyperslab(*dataspace,
2150  H5S_SELECT_SET,
2151  offset.data(),
2152  stride.data(),
2153  count.data(),
2154  block.data());
2155  Assert(ret >= 0, ExcMessage("Error at H5Sselect_hyperslab"));
2156 
2157  internal::set_plist(plist, mpi);
2158 
2159  ret = H5Dwrite(*hdf5_reference,
2160  *t_type,
2161  memory_dataspace,
2162  *dataspace,
2163  plist,
2164  make_array_view(data).data());
2165  Assert(ret >= 0, ExcMessage("Error at H5Dwrite"));
2166 
2168  io_mode,
2169  local_no_collective_cause,
2170  global_no_collective_cause,
2171  mpi,
2172  query_io_mode);
2173 
2174  ret = H5Sclose(memory_dataspace);
2175  Assert(ret >= 0, ExcMessage("Error at H5Sclose"));
2176 
2177  (void)ret;
2178  }
2179 
2180 
2181 
2182  template <typename number>
2183  void
2185  {
2186  std::shared_ptr<hid_t> t_type = internal::get_hdf5_datatype<number>();
2187  std::vector<hsize_t> data_dimensions = {0};
2188 
2189  hid_t memory_dataspace;
2190  hid_t plist;
2191  herr_t ret;
2192 
2193  memory_dataspace = H5Screate_simple(1, data_dimensions.data(), nullptr);
2194  Assert(memory_dataspace >= 0, ExcMessage("Error at H5Screate_simple"));
2195  ret = H5Sselect_none(*dataspace);
2196  Assert(ret >= 0, ExcMessage("Error at H5PSselect_none"));
2197 
2198  internal::set_plist(plist, mpi);
2199 
2200  // The pointer of data can safely be nullptr, see the discussion at the HDF5
2201  // forum:
2202  // https://forum.hdfgroup.org/t/parallel-i-o-does-not-support-filters-yet/884/17
2203  ret = H5Dwrite(
2204  *hdf5_reference, *t_type, memory_dataspace, *dataspace, plist, nullptr);
2205  Assert(ret >= 0, ExcMessage("Error at H5Dwrite"));
2206 
2208  io_mode,
2209  local_no_collective_cause,
2210  global_no_collective_cause,
2211  mpi,
2212  query_io_mode);
2213 
2214  ret = H5Sclose(memory_dataspace);
2215  Assert(ret >= 0, ExcMessage("Error at H5Sclose"));
2216 
2217  (void)ret;
2218  }
2219 
2220 
2221 
2222  template <typename number>
2223  DataSet
2224  Group::create_dataset(const std::string & name,
2225  const std::vector<hsize_t> &dimensions) const
2226  {
2227  std::shared_ptr<hid_t> t_type = internal::get_hdf5_datatype<number>();
2228  return {name, *hdf5_reference, dimensions, t_type, mpi};
2229  }
2230 
2231 
2232 
2233  template <typename Container>
2234  void
2235  Group::write_dataset(const std::string &name, const Container &data) const
2236  {
2237  std::vector<hsize_t> dimensions = internal::get_container_dimensions(data);
2238  auto dataset =
2239  create_dataset<typename Container::value_type>(name, dimensions);
2240  dataset.write(data);
2241  }
2242 } // namespace HDF5
2243 
2245 
2246 
2247 #endif // DEAL_II_WITH_HDF5
2248 
2249 #endif // dealii_hdf5_h
size_type m() const
#define AssertNothrow(cond, exc)
Definition: exceptions.h:1524
void write_none()
Definition: hdf5.h:2184
#define AssertDimension(dim1, dim2)
Definition: exceptions.h:1655
Container read_selection(const std::vector< hsize_t > &coordinates)
Definition: hdf5.h:1795
GroupAccessMode
Definition: hdf5.h:982
std::string get_name() const
Definition: hdf5.cc:82
void read_none()
Definition: hdf5.h:1954
T get_attribute(const std::string &attr_name) const
Definition: hdf5.h:1588
std::vector< hsize_t > get_container_dimensions(const std::vector< number > &data)
Definition: hdf5.h:1351
Container read()
Definition: hdf5.h:1760
unsigned int get_container_size(const std::vector< number > &data)
Definition: hdf5.h:1381
HDF5Object(const std::string &name, const bool mpi)
Definition: hdf5.cc:74
std::shared_ptr< hid_t > hdf5_reference
Definition: hdf5.h:415
#define AssertThrow(cond, exc)
Definition: exceptions.h:1571
FileAccessMode
Definition: hdf5.h:1079
unsigned int size
Definition: hdf5.h:940
void write_selection(const Container &data, const std::vector< hsize_t > &coordinates)
Definition: hdf5.h:2026
void write(const Container &data)
Definition: hdf5.h:1994
bool query_io_mode
Definition: hdf5.h:951
static ::ExceptionBase & ExcMessage(std::string arg1)
size_type n() const
static const char T
#define Assert(cond, exc)
Definition: exceptions.h:1461
void set_attribute(const std::string &attr_name, const T value)
Definition: hdf5.h:1665
void release_plist(hid_t &plist, H5D_mpio_actual_io_mode_t &io_mode, uint32_t &local_no_collective_cause, uint32_t &global_no_collective_cause, const bool mpi, const bool query_io_mode)
Definition: hdf5.h:1488
ArrayView< typename std::remove_reference< typename std::iterator_traits< Iterator >::reference >::type, MemorySpaceType > make_array_view(const Iterator begin, const Iterator end)
Definition: array_view.h:698
DataSet create_dataset(const std::string &name, const std::vector< hsize_t > &dimensions) const
Definition: hdf5.h:2224
const std::string name
Definition: hdf5.h:406
#define DEAL_II_NAMESPACE_CLOSE
Definition: config.h:401
std::string no_collective_cause_to_string(const uint32_t no_collective_cause)
Definition: hdf5.h:1527
int create_group(const MPI_Comm &comm, const MPI_Group &group, const int tag, MPI_Comm *new_comm)
Definition: mpi.cc:181
void write_hyperslab(const Container &data, const std::vector< hsize_t > &offset, const std::vector< hsize_t > &count)
Definition: hdf5.h:2075
std::vector< hsize_t > dimensions
Definition: hdf5.h:930
unsigned int get_container_size(const FullMatrix< number > &data)
Definition: hdf5.h:1399
unsigned int rank
Definition: hdf5.h:924
std::shared_ptr< hid_t > dataspace
Definition: hdf5.h:935
void set_plist(hid_t &plist, const bool mpi)
Definition: hdf5.h:1463
void write_dataset(const std::string &name, const Container &data) const
Definition: hdf5.h:2235
Definition: hdf5.h:345
#define DEAL_II_NAMESPACE_OPEN
Definition: config.h:400
std::vector< hsize_t > get_container_dimensions(const FullMatrix< number > &data)
Definition: hdf5.h:1371
uint32_t local_no_collective_cause
Definition: hdf5.h:963
uint32_t global_no_collective_cause
Definition: hdf5.h:970
size_type size() const
static ::ExceptionBase & ExcNotImplemented()
std::shared_ptr< hid_t > get_hdf5_datatype()
Definition: hdf5.h:1264
void free(T *&pointer)
Definition: cuda.h:97
const bool mpi
Definition: hdf5.h:420
H5D_mpio_actual_io_mode_t io_mode
Definition: hdf5.h:956
std::enable_if< std::is_same< Container, FullMatrix< typename Container::value_type > >::value, Container >::type initialize_container(const std::vector< hsize_t > &dimensions)
Definition: hdf5.h:1436
Container read_hyperslab(const std::vector< hsize_t > &offset, const std::vector< hsize_t > &count)
Definition: hdf5.h:1847
static ::ExceptionBase & ExcInternalError()