deal.II version GIT relicensing-2890-gf173123fa3 2025-03-22 01:40:00+00:00
\(\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\}}\)
Loading...
Searching...
No Matches
hdf5.h
Go to the documentation of this file.
1// ------------------------------------------------------------------------
2//
3// SPDX-License-Identifier: LGPL-2.1-or-later
4// Copyright (C) 2019 - 2024 by the deal.II authors
5//
6// This file is part of the deal.II library.
7//
8// Part of the source code is dual licensed under Apache-2.0 WITH
9// LLVM-exception OR LGPL-2.1-or-later. Detailed license information
10// governing the source code and code contributions can be found in
11// LICENSE.md and CONTRIBUTING.md at the top level directory of deal.II.
12//
13// ------------------------------------------------------------------------
14
15#ifndef dealii_hdf5_h
16#define dealii_hdf5_h
17
18#include <deal.II/base/config.h>
19
20#ifdef DEAL_II_WITH_HDF5
21
23
25
26# include <hdf5.h>
27
28# include <numeric>
29
31
32// It is necessary to turn clang-format off in order to maintain the Doxygen
33// links because they are longer than 80 characters
34// clang-format off
342// clang-format on
343namespace HDF5
344{
349 {
350 protected:
355 HDF5Object(const std::string &name, const bool mpi);
356
357 public:
370 template <typename T>
371 T
372 get_attribute(const std::string &attr_name) const;
373
386 template <typename T>
387 void
388 set_attribute(const std::string &attr_name, const T value);
389
395 std::string
396 get_name() const;
397
398 protected:
404 const std::string name;
405
413 std::shared_ptr<hid_t> hdf5_reference;
414
418 const bool mpi;
419 };
420
424 class DataSet : public HDF5Object
425 {
426 friend class Group;
427
428 protected:
433 DataSet(const std::string &name, const hid_t &parent_group_id, bool mpi);
434
439 DataSet(const std::string &name,
440 const hid_t &parent_group_id,
441 const std::vector<hsize_t> &dimensions,
442 const std::shared_ptr<hid_t> &t_type,
443 const bool mpi);
444
445 public:
463 template <typename Container>
464 Container
465 read();
466
500 template <typename Container>
501 Container
502 read_selection(const std::vector<hsize_t> &coordinates);
503
504 // clang-format off
538 // clang-format on
539 template <typename Container>
540 Container
541 read_hyperslab(const std::vector<hsize_t> &offset,
542 const std::vector<hsize_t> &count);
543
576 template <typename Container>
577 Container
578 read_hyperslab(const std::vector<hsize_t> &data_dimensions,
579 const std::vector<hsize_t> &offset,
580 const std::vector<hsize_t> &stride,
581 const std::vector<hsize_t> &count,
582 const std::vector<hsize_t> &block);
583
595 template <typename number>
596 void
597 read_none();
598
617 template <typename Container>
618 void
619 write(const Container &data);
620
649 template <typename Container>
650 void
651 write_selection(const Container &data,
652 const std::vector<hsize_t> &coordinates);
653
654 // clang-format off
680 // clang-format on
681 template <typename Container>
682 void
683 write_hyperslab(const Container &data,
684 const std::vector<hsize_t> &offset,
685 const std::vector<hsize_t> &count);
686
719 template <typename Container>
720 void
721 write_hyperslab(const Container &data,
722 const std::vector<hsize_t> &data_dimensions,
723 const std::vector<hsize_t> &offset,
724 const std::vector<hsize_t> &stride,
725 const std::vector<hsize_t> &count,
726 const std::vector<hsize_t> &block);
727
745 template <typename number>
746 void
747 write_none();
748
765 bool
766 get_query_io_mode() const;
767
771 void
772 set_query_io_mode(const bool new_query_io_mode);
773
788 std::string
789 get_io_mode();
790
807 H5D_mpio_actual_io_mode_t
809
828 std::string
830
852 std::uint32_t
854
873 std::string
875
896 std::uint32_t
898
904 std::vector<hsize_t>
905 get_dimensions() const;
906
910 unsigned int
911 get_size() const;
912
916 unsigned int
917 get_rank() const;
918
919 private:
923 unsigned int rank;
924
929 std::vector<hsize_t> dimensions;
930
934 std::shared_ptr<hid_t> dataspace;
935
939 unsigned int size;
940
951
955 H5D_mpio_actual_io_mode_t io_mode;
956
963
970 };
971
975 class Group : public HDF5Object
976 {
977 protected:
982 {
986 open,
990 create
991 };
1000 Group(const std::string &name,
1001 const Group &parent_group,
1002 const bool mpi,
1003 const GroupAccessMode mode);
1004
1010 Group(const std::string &name, const bool mpi);
1011
1012 public:
1016 Group
1017 open_group(const std::string &name) const;
1018
1022 Group
1023 create_group(const std::string &name) const;
1024
1028 DataSet
1029 open_dataset(const std::string &name) const;
1030
1041 template <typename number>
1042 DataSet
1043 create_dataset(const std::string &name,
1044 const std::vector<hsize_t> &dimensions) const;
1045
1064 template <typename Container>
1065 void
1066 write_dataset(const std::string &name, const Container &data) const;
1067 };
1068
1072 class File : public Group
1073 {
1074 public:
1079 {
1083 open,
1087 create
1088 };
1089
1095 File(const std::string &name, const FileAccessMode mode);
1096
1104 File(const std::string &name,
1105 const FileAccessMode mode,
1106 const MPI_Comm mpi_communicator);
1107
1108 private:
1116 File(const std::string &name,
1117 const FileAccessMode mode,
1118 const bool mpi,
1119 const MPI_Comm mpi_communicator);
1120 };
1121
1122 namespace internal
1123 {
1135 template <typename number>
1136 std::shared_ptr<hid_t>
1138
1148 template <typename number>
1149 std::vector<hsize_t>
1150 get_container_dimensions(const std::vector<number> &data);
1151
1156 template <typename number>
1157 std::vector<hsize_t>
1159
1164 template <typename number>
1165 std::vector<hsize_t>
1167
1172 template <typename number>
1173 unsigned int
1174 get_container_size(const std::vector<number> &data);
1175
1180 template <typename number>
1181 unsigned int
1183
1188 template <typename number>
1189 unsigned int
1191
1210 template <typename Container>
1211 std::enable_if_t<
1212 std::is_same_v<Container, std::vector<typename Container::value_type>>,
1213 Container>
1214 initialize_container(const std::vector<hsize_t> &dimensions);
1215
1219 template <typename Container>
1220 std::enable_if_t<
1221 std::is_same_v<Container, Vector<typename Container::value_type>>,
1222 Container>
1223 initialize_container(const std::vector<hsize_t> &dimensions);
1224
1228 template <typename Container>
1229 std::enable_if_t<
1230 std::is_same_v<Container, FullMatrix<typename Container::value_type>>,
1231 Container>
1232 initialize_container(const std::vector<hsize_t> &dimensions);
1233
1240 inline void
1241 set_plist(hid_t &plist, const bool mpi);
1242
1251 inline void
1252 release_plist(hid_t &plist,
1253 H5D_mpio_actual_io_mode_t &io_mode,
1254 std::uint32_t &local_no_collective_cause,
1255 std::uint32_t &global_no_collective_cause,
1256 const bool mpi,
1257 const bool query_io_mode);
1258
1262 inline std::string
1263 no_collective_cause_to_string(const std::uint32_t no_collective_cause);
1264 } // namespace internal
1265
1266
1267
1268 // definitions
1269
1270 namespace internal
1271 {
1272 template <typename number>
1273 std::shared_ptr<hid_t>
1275 {
1276 static_assert(std::is_same_v<number, float> ||
1277 std::is_same_v<number, double> ||
1278 std::is_same_v<number, int> ||
1279 std::is_same_v<number, bool> ||
1280 std::is_same_v<number, unsigned int> ||
1281 std::is_same_v<number, std::complex<float>> ||
1282 std::is_same_v<number, std::complex<double>>,
1283 "The data type you are trying to get the HDF5 tag for "
1284 "is not supported by this function.");
1285
1286 if (std::is_same_v<number, float>)
1287 {
1288 return std::make_shared<hid_t>(H5T_NATIVE_FLOAT);
1289 }
1290 else if (std::is_same_v<number, double>)
1291 {
1292 return std::make_shared<hid_t>(H5T_NATIVE_DOUBLE);
1293 }
1294 else if (std::is_same_v<number, int>)
1295 {
1296 return std::make_shared<hid_t>(H5T_NATIVE_INT);
1297 }
1298 else if (std::is_same_v<number, unsigned int>)
1299 {
1300 return std::make_shared<hid_t>(H5T_NATIVE_UINT);
1301 }
1302 else if (std::is_same_v<number, bool>)
1303 {
1304 return std::make_shared<hid_t>(H5T_NATIVE_HBOOL);
1305 }
1306 else if (std::is_same_v<number, std::complex<float>>)
1307 {
1308 std::shared_ptr<hid_t> t_type =
1309 std::shared_ptr<hid_t>(new hid_t, [](hid_t *pointer) {
1310 // Release the HDF5 resource
1311 const herr_t ret = H5Tclose(*pointer);
1312 AssertNothrow(ret >= 0, ExcInternalError());
1313 (void)ret;
1314 delete pointer;
1315 });
1316
1317 *t_type = H5Tcreate(H5T_COMPOUND, sizeof(std::complex<float>));
1318 // The C++ standards committee agreed to mandate that the storage
1319 // format used for the std::complex type be binary-compatible with
1320 // the C99 type, i.e. an array T[2] with consecutive real [0] and
1321 // imaginary [1] parts.
1322 herr_t ret = H5Tinsert(*t_type, "r", 0, H5T_NATIVE_FLOAT);
1323 Assert(ret >= 0, ExcInternalError());
1324 ret = H5Tinsert(*t_type, "i", sizeof(float), H5T_NATIVE_FLOAT);
1325 Assert(ret >= 0, ExcInternalError());
1326 (void)ret;
1327 return t_type;
1328 }
1329 else if (std::is_same_v<number, std::complex<double>>)
1330 {
1331 std::shared_ptr<hid_t> t_type =
1332 std::shared_ptr<hid_t>(new hid_t, [](hid_t *pointer) {
1333 // Release the HDF5 resource
1334 const herr_t ret = H5Tclose(*pointer);
1335 AssertNothrow(ret >= 0, ExcInternalError());
1336 (void)ret;
1337 delete pointer;
1338 });
1339 *t_type = H5Tcreate(H5T_COMPOUND, sizeof(std::complex<double>));
1340 // The C++ standards committee agreed to mandate that the storage
1341 // format used for the std::complex type be binary-compatible with
1342 // the C99 type, i.e. an array T[2] with consecutive real [0] and
1343 // imaginary [1] parts.
1344 herr_t ret = H5Tinsert(*t_type, "r", 0, H5T_NATIVE_DOUBLE);
1345 Assert(ret >= 0, ExcInternalError());
1346 ret = H5Tinsert(*t_type, "i", sizeof(double), H5T_NATIVE_DOUBLE);
1347 Assert(ret >= 0, ExcInternalError());
1348 (void)ret;
1349 return t_type;
1350 }
1351
1352 // The function should not reach this point
1354 return {};
1355 }
1356
1357
1358
1359 template <typename number>
1360 std::vector<hsize_t>
1361 get_container_dimensions(const std::vector<number> &data)
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.size()};
1374 return dimensions;
1375 }
1376
1377
1378
1379 template <typename number>
1380 std::vector<hsize_t>
1382 {
1383 std::vector<hsize_t> dimensions = {data.m(), data.n()};
1384 return dimensions;
1385 }
1386
1387
1388
1389 template <typename number>
1390 unsigned int
1391 get_container_size(const std::vector<number> &data)
1392 {
1393 return static_cast<unsigned int>(data.size());
1394 }
1395
1396
1397
1398 template <typename number>
1399 unsigned int
1401 {
1402 return static_cast<unsigned int>(data.size());
1403 }
1404
1405
1406
1407 template <typename number>
1408 unsigned int
1410 {
1411 return static_cast<unsigned int>(data.m() * data.n());
1412 }
1413
1414
1415
1416 template <typename Container>
1417 std::enable_if_t<
1418 std::is_same_v<Container, std::vector<typename Container::value_type>>,
1419 Container>
1420 initialize_container(const std::vector<hsize_t> &dimensions)
1421 {
1422 return Container(std::accumulate(
1423 dimensions.begin(), dimensions.end(), 1, std::multiplies<int>()));
1424 }
1425
1426
1427
1428 template <typename Container>
1429 std::enable_if_t<
1430 std::is_same_v<Container, Vector<typename Container::value_type>>,
1431 Container>
1432 initialize_container(const std::vector<hsize_t> &dimensions)
1433 {
1434 return Container(std::accumulate(
1435 dimensions.begin(), dimensions.end(), 1, std::multiplies<int>()));
1436 }
1437
1438
1439
1440 template <typename Container>
1441 std::enable_if_t<
1442 std::is_same_v<Container, FullMatrix<typename Container::value_type>>,
1443 Container>
1444 initialize_container(const std::vector<hsize_t> &dimensions)
1445 {
1446 // If the rank is higher than 2, then remove single-dimensional entries
1447 // from the shape defined by dimensions. This is equivalent to the squeeze
1448 // function of python/numpy. For example the following code would convert
1449 // the vector {1,3,1,2} to {3,2}
1450 std::vector<hsize_t> squeezed_dimensions;
1451
1452 if (dimensions.size() > 2)
1453 {
1454 for (const auto &dimension : dimensions)
1455 {
1456 if (dimension > 1)
1457 squeezed_dimensions.push_back(dimension);
1458 }
1459 }
1460 else
1461 {
1462 squeezed_dimensions = dimensions;
1463 }
1464
1465 AssertDimension(squeezed_dimensions.size(), 2);
1466 return Container(squeezed_dimensions[0], squeezed_dimensions[1]);
1467 }
1468
1469
1470 inline void
1471 set_plist(hid_t &plist, const bool mpi)
1472 {
1473 if (mpi)
1474 {
1475# ifdef DEAL_II_WITH_MPI
1476 plist = H5Pcreate(H5P_DATASET_XFER);
1477 Assert(plist >= 0, ExcInternalError());
1478 const herr_t ret = H5Pset_dxpl_mpio(plist, H5FD_MPIO_COLLECTIVE);
1479 (void)ret;
1480 Assert(ret >= 0, ExcInternalError());
1481# else
1483# endif
1484 }
1485 else
1486 {
1487 plist = H5P_DEFAULT;
1488 }
1489
1490 (void)plist;
1491 (void)mpi;
1492 }
1493
1494
1495 inline void
1496 release_plist(hid_t &plist,
1497 H5D_mpio_actual_io_mode_t &io_mode,
1498 std::uint32_t &local_no_collective_cause,
1499 std::uint32_t &global_no_collective_cause,
1500 const bool mpi,
1501 const bool query_io_mode)
1502 {
1503 if (mpi)
1504 {
1505# ifdef DEAL_II_WITH_MPI
1506 herr_t ret;
1507 (void)ret;
1508 if (query_io_mode)
1509 {
1510 ret = H5Pget_mpio_actual_io_mode(plist, &io_mode);
1511 Assert(ret >= 0, ExcInternalError());
1512 ret =
1513 H5Pget_mpio_no_collective_cause(plist,
1514 &local_no_collective_cause,
1515 &global_no_collective_cause);
1516 Assert(ret >= 0, ExcInternalError());
1517 }
1518 ret = H5Pclose(plist);
1519 Assert(ret >= 0, ExcInternalError());
1520# else
1522# endif
1523 }
1524
1525 (void)plist;
1526 (void)io_mode;
1527 (void)local_no_collective_cause;
1528 (void)global_no_collective_cause;
1529 (void)mpi;
1530 (void)query_io_mode;
1531 }
1532
1533
1534 inline std::string
1535 no_collective_cause_to_string(const std::uint32_t no_collective_cause)
1536 {
1537 std::string message;
1538
1539 auto append_to_message = [&message](const char *p) {
1540 if (message.size() > 0)
1541 message += ", ";
1542 message += p;
1543 };
1544
1545 // The first is not a bitmask comparison, the rest are bitmask
1546 // comparisons.
1547 // https://support.hdfgroup.org/HDF5/doc/RM/RM_H5P.html#Property-GetMpioNoCollectiveCause
1548 // See H5Ppublic.h
1549 // Hex codes are used because the HDF5 Group can deprecate some of the
1550 // enum codes. For example the enum code H5D_MPIO_FILTERS is not defined
1551 // in 1.10.2 because it is possible to use compressed datasets with the
1552 // MPI/IO driver.
1553
1554 // H5D_MPIO_COLLECTIVE
1555 if (no_collective_cause == 0x00)
1556 {
1557 append_to_message("H5D_MPIO_COLLECTIVE");
1558 }
1559 // H5D_MPIO_SET_INDEPENDENT
1560 if ((no_collective_cause & 0x01) == 0x01)
1561 {
1562 append_to_message("H5D_MPIO_SET_INDEPENDENT");
1563 }
1564 // H5D_MPIO_DATATYPE_CONVERSION
1565 if ((no_collective_cause & 0x02) == 0x02)
1566 {
1567 append_to_message("H5D_MPIO_DATATYPE_CONVERSION");
1568 }
1569 // H5D_MPIO_DATA_TRANSFORMS
1570 if ((no_collective_cause & 0x04) == 0x04)
1571 {
1572 append_to_message("H5D_MPIO_DATA_TRANSFORMS");
1573 }
1574 // H5D_MPIO_NOT_SIMPLE_OR_SCALAR_DATASPACES
1575 if ((no_collective_cause & 0x10) == 0x10)
1576 {
1577 append_to_message("H5D_MPIO_NOT_SIMPLE_OR_SCALAR_DATASPACES");
1578 }
1579 // H5D_MPIO_NOT_CONTIGUOUS_OR_CHUNKED_DATASET
1580 if ((no_collective_cause & 0x20) == 0x20)
1581 {
1582 append_to_message("H5D_MPIO_NOT_CONTIGUOUS_OR_CHUNKED_DATASET");
1583 }
1584 // H5D_MPIO_FILTERS
1585 if ((no_collective_cause & 0x40) == 0x40)
1586 {
1587 append_to_message("H5D_MPIO_FILTERS");
1588 }
1589 return message;
1590 }
1591 } // namespace internal
1592
1593
1594 template <typename T>
1595 T
1596 HDF5Object::get_attribute(const std::string &attr_name) const
1597 {
1598 const std::shared_ptr<hid_t> t_type = internal::get_hdf5_datatype<T>();
1599 T value;
1600 hid_t attr;
1601 herr_t ret;
1602
1603 attr = H5Aopen(*hdf5_reference, attr_name.data(), H5P_DEFAULT);
1604 Assert(attr >= 0, ExcMessage("Error at H5Aopen"));
1605 (void)ret;
1606 ret = H5Aread(attr, *t_type, &value);
1607 Assert(ret >= 0, ExcMessage("Error at H5Aread"));
1608 (void)ret;
1609 ret = H5Aclose(attr);
1610 Assert(ret >= 0, ExcMessage("Error at H5Aclose"));
1611
1612 return value;
1613 }
1614
1615
1616
1617 template <>
1618 inline std::string
1619 HDF5Object::get_attribute(const std::string &attr_name) const
1620 {
1621 // Reads a UTF8 variable string
1622 //
1623 // code inspired from
1624 // https://support.hdfgroup.org/ftp/HDF5/examples/misc-examples/vlstratt.c
1625 //
1626 // In the case of a variable length string the user does not have to reserve
1627 // memory for string_out. H5Aread will reserve the memory and the
1628 // user has to free the memory.
1629 //
1630 // Todo:
1631 // - Use H5Dvlen_reclaim instead of free
1632
1633 char *string_out;
1634 hid_t attr;
1635 hid_t type;
1636 herr_t ret;
1637
1638 /* Create a datatype to refer to. */
1639 type = H5Tcopy(H5T_C_S1);
1640 Assert(type >= 0, ExcInternalError());
1641
1642 // Python strings are encoded in UTF8
1643 ret = H5Tset_cset(type, H5T_CSET_UTF8);
1644 Assert(type >= 0, ExcInternalError());
1645
1646 ret = H5Tset_size(type, H5T_VARIABLE);
1647 Assert(ret >= 0, ExcInternalError());
1648
1649 attr = H5Aopen(*hdf5_reference, attr_name.data(), H5P_DEFAULT);
1650 Assert(attr >= 0, ExcInternalError());
1651
1652 ret = H5Aread(attr, type, &string_out);
1653 Assert(ret >= 0, ExcInternalError());
1654
1655 std::string string_value(string_out);
1656 // The memory of the variable length string has to be freed.
1657 // H5Dvlen_reclaim could be also used
1658 free(string_out);
1659 ret = H5Tclose(type);
1660 Assert(ret >= 0, ExcInternalError());
1661
1662 ret = H5Aclose(attr);
1663 Assert(ret >= 0, ExcInternalError());
1664
1665 (void)ret;
1666 return string_value;
1667 }
1668
1669
1670
1671 template <typename T>
1672 void
1673 HDF5Object::set_attribute(const std::string &attr_name, const T value)
1674 {
1675 hid_t attr;
1676 hid_t aid;
1677 herr_t ret;
1678
1679 const std::shared_ptr<hid_t> t_type = internal::get_hdf5_datatype<T>();
1680
1681
1682 /*
1683 * Create scalar attribute.
1684 */
1685 aid = H5Screate(H5S_SCALAR);
1686 Assert(aid >= 0, ExcMessage("Error at H5Screate"));
1687 attr = H5Acreate2(*hdf5_reference,
1688 attr_name.data(),
1689 *t_type,
1690 aid,
1691 H5P_DEFAULT,
1692 H5P_DEFAULT);
1693 Assert(attr >= 0, ExcMessage("Error at H5Acreate2"));
1694
1695 /*
1696 * Write scalar attribute.
1697 */
1698 ret = H5Awrite(attr, *t_type, &value);
1699 Assert(ret >= 0, ExcMessage("Error at H5Awrite"));
1700
1701 ret = H5Sclose(aid);
1702 Assert(ret >= 0, ExcMessage("Error at H5Sclose"));
1703 ret = H5Aclose(attr);
1704 Assert(ret >= 0, ExcMessage("Error at H5Aclose"));
1705
1706 (void)ret;
1707 }
1708
1709
1710
1711 template <>
1712 inline void
1713 HDF5Object::set_attribute(const std::string &attr_name,
1714 const std::string value) // NOLINT
1715 {
1716 // Writes a UTF8 variable string
1717 //
1718 // code inspired from
1719 // https://support.hdfgroup.org/ftp/HDF5/examples/misc-examples/vlstratt.c
1720
1721 hid_t attr;
1722 hid_t aid;
1723 hid_t t_type;
1724 herr_t ret;
1725
1726 /* Create a datatype to refer to. */
1727 t_type = H5Tcopy(H5T_C_S1);
1728 Assert(t_type >= 0, ExcInternalError());
1729
1730 // Python strings are encoded in UTF8
1731 ret = H5Tset_cset(t_type, H5T_CSET_UTF8);
1732 Assert(t_type >= 0, ExcInternalError());
1733
1734 ret = H5Tset_size(t_type, H5T_VARIABLE);
1735 Assert(ret >= 0, ExcInternalError());
1736
1737 /*
1738 * Create scalar attribute.
1739 */
1740 aid = H5Screate(H5S_SCALAR);
1741 Assert(aid >= 0, ExcMessage("Error at H5Screate"));
1742 attr = H5Acreate2(
1743 *hdf5_reference, attr_name.data(), t_type, aid, H5P_DEFAULT, H5P_DEFAULT);
1744 Assert(attr >= 0, ExcMessage("Error at H5Acreate2"));
1745
1746 /*
1747 * Write scalar attribute.
1748 * In most of the cases H5Awrite and H5Dwrite take a pointer to the data.
1749 * But in the particular case of a variable length string, H5Awrite takes
1750 * the address of the pointer of the string.
1751 */
1752 const char *c_string_value = value.c_str();
1753 ret = H5Awrite(attr, t_type, &c_string_value);
1754 Assert(ret >= 0, ExcInternalError());
1755
1756 ret = H5Sclose(aid);
1757 Assert(ret >= 0, ExcMessage("Error at H5Sclose"));
1758 ret = H5Aclose(attr);
1759 Assert(ret >= 0, ExcMessage("Error at H5Aclose"));
1760
1761 (void)ret;
1762 }
1763
1764
1765
1766 template <typename Container>
1767 Container
1769 {
1770 const std::shared_ptr<hid_t> t_type =
1771 internal::get_hdf5_datatype<typename Container::value_type>();
1772 hid_t plist;
1773 herr_t ret;
1774
1775 Container data = internal::initialize_container<Container>(dimensions);
1776
1777 internal::set_plist(plist, mpi);
1778
1779 ret = H5Dread(*hdf5_reference,
1780 *t_type,
1781 H5S_ALL,
1782 H5S_ALL,
1783 plist,
1785 Assert(ret >= 0, ExcInternalError());
1786
1787
1789 io_mode,
1792 mpi,
1794
1795 (void)ret;
1796 return data;
1797 }
1798
1799
1800
1801 template <typename Container>
1802 Container
1803 DataSet::read_selection(const std::vector<hsize_t> &coordinates)
1804 {
1805 Assert(coordinates.size() % rank == 0,
1806 ExcMessage(
1807 "The dimension of coordinates has to be divisible by the rank"));
1808 const std::shared_ptr<hid_t> t_type =
1809 internal::get_hdf5_datatype<typename Container::value_type>();
1810 hid_t plist;
1811 hid_t memory_dataspace;
1812 herr_t ret;
1813
1814 std::vector<hsize_t> data_dimensions{
1815 static_cast<hsize_t>(coordinates.size() / rank)};
1816
1817 Container data = internal::initialize_container<Container>(data_dimensions);
1818
1819 memory_dataspace = H5Screate_simple(1, data_dimensions.data(), nullptr);
1820 Assert(memory_dataspace >= 0, ExcMessage("Error at H5Screate_simple"));
1821 ret = H5Sselect_elements(*dataspace,
1822 H5S_SELECT_SET,
1823 data.size(),
1824 coordinates.data());
1825 Assert(ret >= 0, ExcMessage("Error at H5Sselect_elements"));
1826
1827 internal::set_plist(plist, mpi);
1828
1829 ret = H5Dread(*hdf5_reference,
1830 *t_type,
1831 memory_dataspace,
1832 *dataspace,
1833 plist,
1835 Assert(ret >= 0, ExcMessage("Error at H5Dread"));
1836
1838 io_mode,
1841 mpi,
1843
1844 ret = H5Sclose(memory_dataspace);
1845 Assert(ret >= 0, ExcMessage("Error at H5SClose"));
1846
1847 (void)ret;
1848 return data;
1849 }
1850
1851
1852
1853 template <typename Container>
1854 Container
1855 DataSet::read_hyperslab(const std::vector<hsize_t> &offset,
1856 const std::vector<hsize_t> &count)
1857 {
1858 const std::shared_ptr<hid_t> t_type =
1859 internal::get_hdf5_datatype<typename Container::value_type>();
1860 hid_t plist;
1861 hid_t memory_dataspace;
1862 herr_t ret;
1863
1864 // In this particular overload of read_hyperslab the data_dimensions are
1865 // the same as count
1866 std::vector<hsize_t> data_dimensions = count;
1867
1868 Container data = internal::initialize_container<Container>(data_dimensions);
1869
1870 memory_dataspace =
1871 H5Screate_simple(data_dimensions.size(), data_dimensions.data(), nullptr);
1872 Assert(memory_dataspace >= 0, ExcMessage("Error at H5Screate_simple"));
1873 ret = H5Sselect_hyperslab(*dataspace,
1874 H5S_SELECT_SET,
1875 offset.data(),
1876 nullptr,
1877 count.data(),
1878 nullptr);
1879 Assert(ret >= 0, ExcMessage("Error at H5Sselect_hyperslab"));
1880
1881 internal::set_plist(plist, mpi);
1882
1883 ret = H5Dread(*hdf5_reference,
1884 *t_type,
1885 memory_dataspace,
1886 *dataspace,
1887 plist,
1889 Assert(ret >= 0, ExcMessage("Error at H5Dread"));
1890
1892 io_mode,
1895 mpi,
1897
1898 ret = H5Sclose(memory_dataspace);
1899 Assert(ret >= 0, ExcMessage("Error at H5SClose"));
1900
1901 (void)ret;
1902 return data;
1903 }
1904
1905
1906
1907 template <typename Container>
1908 Container
1909 DataSet::read_hyperslab(const std::vector<hsize_t> &data_dimensions,
1910 const std::vector<hsize_t> &offset,
1911 const std::vector<hsize_t> &stride,
1912 const std::vector<hsize_t> &count,
1913 const std::vector<hsize_t> &block)
1914 {
1915 const std::shared_ptr<hid_t> t_type =
1916 internal::get_hdf5_datatype<typename Container::value_type>();
1917 hid_t plist;
1918 hid_t memory_dataspace;
1919 herr_t ret;
1920
1921 Container data = internal::initialize_container<Container>(data_dimensions);
1922
1923 memory_dataspace =
1924 H5Screate_simple(data_dimensions.size(), data_dimensions.data(), nullptr);
1925 Assert(memory_dataspace >= 0, ExcMessage("Error at H5Screate_simple"));
1926 ret = H5Sselect_hyperslab(*dataspace,
1927 H5S_SELECT_SET,
1928 offset.data(),
1929 stride.data(),
1930 count.data(),
1931 block.data());
1932 Assert(ret >= 0, ExcMessage("Error at H5Sselect_hyperslab"));
1933
1934 internal::set_plist(plist, mpi);
1935
1936 ret = H5Dread(*hdf5_reference,
1937 *t_type,
1938 memory_dataspace,
1939 *dataspace,
1940 plist,
1942 Assert(ret >= 0, ExcMessage("Error at H5Dread"));
1943
1945 io_mode,
1948 mpi,
1950
1951 ret = H5Sclose(memory_dataspace);
1952 Assert(ret >= 0, ExcMessage("Error at H5SClose"));
1953
1954 (void)ret;
1955 return data;
1956 }
1957
1958
1959
1960 template <typename number>
1961 void
1963 {
1964 const std::shared_ptr<hid_t> t_type = internal::get_hdf5_datatype<number>();
1965 const std::vector<hsize_t> data_dimensions = {0};
1966
1967 hid_t memory_dataspace;
1968 hid_t plist;
1969 herr_t ret;
1970
1971 memory_dataspace = H5Screate_simple(1, data_dimensions.data(), nullptr);
1972 Assert(memory_dataspace >= 0, ExcMessage("Error at H5Screate_simple"));
1973 ret = H5Sselect_none(*dataspace);
1974 Assert(ret >= 0, ExcMessage("H5Sselect_none"));
1975
1976 internal::set_plist(plist, mpi);
1977
1978 // The pointer of data can safely be nullptr, see the discussion at the HDF5
1979 // forum:
1980 // https://forum.hdfgroup.org/t/parallel-i-o-does-not-support-filters-yet/884/17
1981 ret = H5Dread(
1982 *hdf5_reference, *t_type, memory_dataspace, *dataspace, plist, nullptr);
1983 Assert(ret >= 0, ExcMessage("Error at H5Dread"));
1984
1986 io_mode,
1989 mpi,
1991
1992 ret = H5Sclose(memory_dataspace);
1993 Assert(ret >= 0, ExcMessage("Error at H5SClose"));
1994
1995 (void)ret;
1996 }
1997
1998
1999
2000 template <typename Container>
2001 void
2002 DataSet::write(const Container &data)
2003 {
2005 const std::shared_ptr<hid_t> t_type =
2006 internal::get_hdf5_datatype<typename Container::value_type>();
2007 hid_t plist;
2008 herr_t ret;
2009
2010 internal::set_plist(plist, mpi);
2011
2012 ret = H5Dwrite(*hdf5_reference,
2013 *t_type,
2014 H5S_ALL,
2015 H5S_ALL,
2016 plist,
2018 Assert(ret >= 0, ExcMessage("Error at H5Dwrite"));
2019
2021 io_mode,
2024 mpi,
2026
2027 (void)ret;
2028 }
2029
2030
2031
2032 template <typename Container>
2033 void
2035 const std::vector<hsize_t> &coordinates)
2036 {
2037 AssertDimension(coordinates.size(), data.size() * rank);
2038 const std::shared_ptr<hid_t> t_type =
2039 internal::get_hdf5_datatype<typename Container::value_type>();
2040 const std::vector<hsize_t> data_dimensions =
2042
2043 hid_t memory_dataspace;
2044 hid_t plist;
2045 herr_t ret;
2046
2047
2048 memory_dataspace = H5Screate_simple(1, data_dimensions.data(), nullptr);
2049 Assert(memory_dataspace >= 0, ExcMessage("Error at H5Screate_simple"));
2050 ret = H5Sselect_elements(*dataspace,
2051 H5S_SELECT_SET,
2052 data.size(),
2053 coordinates.data());
2054 Assert(ret >= 0, ExcMessage("Error at H5Sselect_elements"));
2055
2056 internal::set_plist(plist, mpi);
2057
2058 ret = H5Dwrite(*hdf5_reference,
2059 *t_type,
2060 memory_dataspace,
2061 *dataspace,
2062 plist,
2064 Assert(ret >= 0, ExcMessage("Error at H5Dwrite"));
2065
2067 io_mode,
2070 mpi,
2072
2073 ret = H5Sclose(memory_dataspace);
2074 Assert(ret >= 0, ExcMessage("Error at H5SClose"));
2075
2076 (void)ret;
2077 }
2078
2079
2080
2081 template <typename Container>
2082 void
2084 const std::vector<hsize_t> &offset,
2085 const std::vector<hsize_t> &count)
2086 {
2087 AssertDimension(std::accumulate(count.begin(),
2088 count.end(),
2089 1,
2090 std::multiplies<unsigned int>()),
2092 const std::shared_ptr<hid_t> t_type =
2093 internal::get_hdf5_datatype<typename Container::value_type>();
2094 // In this particular overload of write_hyperslab the data_dimensions are
2095 // the same as count
2096 const std::vector<hsize_t> &data_dimensions = count;
2097
2098 hid_t memory_dataspace;
2099 hid_t plist;
2100 herr_t ret;
2101
2102 memory_dataspace =
2103 H5Screate_simple(data_dimensions.size(), data_dimensions.data(), nullptr);
2104 Assert(memory_dataspace >= 0, ExcMessage("Error at H5Screate_simple"));
2105 ret = H5Sselect_hyperslab(*dataspace,
2106 H5S_SELECT_SET,
2107 offset.data(),
2108 nullptr,
2109 count.data(),
2110 nullptr);
2111 Assert(ret >= 0, ExcMessage("Error at H5Sselect_hyperslab"));
2112
2113 internal::set_plist(plist, mpi);
2114
2115 ret = H5Dwrite(*hdf5_reference,
2116 *t_type,
2117 memory_dataspace,
2118 *dataspace,
2119 plist,
2121 Assert(ret >= 0, ExcMessage("Error at H5Dwrite"));
2122
2124 io_mode,
2127 mpi,
2129
2130 ret = H5Sclose(memory_dataspace);
2131 Assert(ret >= 0, ExcMessage("Error at H5Sclose"));
2132
2133 (void)ret;
2134 }
2135
2136
2137
2138 template <typename Container>
2139 void
2141 const std::vector<hsize_t> &data_dimensions,
2142 const std::vector<hsize_t> &offset,
2143 const std::vector<hsize_t> &stride,
2144 const std::vector<hsize_t> &count,
2145 const std::vector<hsize_t> &block)
2146 {
2147 const std::shared_ptr<hid_t> t_type =
2148 internal::get_hdf5_datatype<typename Container::value_type>();
2149
2150 hid_t memory_dataspace;
2151 hid_t plist;
2152 herr_t ret;
2153
2154 memory_dataspace =
2155 H5Screate_simple(data_dimensions.size(), data_dimensions.data(), nullptr);
2156 Assert(memory_dataspace >= 0, ExcMessage("Error at H5Screate_simple"));
2157 ret = H5Sselect_hyperslab(*dataspace,
2158 H5S_SELECT_SET,
2159 offset.data(),
2160 stride.data(),
2161 count.data(),
2162 block.data());
2163 Assert(ret >= 0, ExcMessage("Error at H5Sselect_hyperslab"));
2164
2165 internal::set_plist(plist, mpi);
2166
2167 ret = H5Dwrite(*hdf5_reference,
2168 *t_type,
2169 memory_dataspace,
2170 *dataspace,
2171 plist,
2173 Assert(ret >= 0, ExcMessage("Error at H5Dwrite"));
2174
2176 io_mode,
2179 mpi,
2181
2182 ret = H5Sclose(memory_dataspace);
2183 Assert(ret >= 0, ExcMessage("Error at H5Sclose"));
2184
2185 (void)ret;
2186 }
2187
2188
2189
2190 template <typename number>
2191 void
2193 {
2194 std::shared_ptr<hid_t> t_type = internal::get_hdf5_datatype<number>();
2195 std::vector<hsize_t> data_dimensions = {0};
2196
2197 hid_t memory_dataspace;
2198 hid_t plist;
2199 herr_t ret;
2200
2201 memory_dataspace = H5Screate_simple(1, data_dimensions.data(), nullptr);
2202 Assert(memory_dataspace >= 0, ExcMessage("Error at H5Screate_simple"));
2203 ret = H5Sselect_none(*dataspace);
2204 Assert(ret >= 0, ExcMessage("Error at H5PSselect_none"));
2205
2206 internal::set_plist(plist, mpi);
2207
2208 // The pointer of data can safely be nullptr, see the discussion at the HDF5
2209 // forum:
2210 // https://forum.hdfgroup.org/t/parallel-i-o-does-not-support-filters-yet/884/17
2211 ret = H5Dwrite(
2212 *hdf5_reference, *t_type, memory_dataspace, *dataspace, plist, nullptr);
2213 Assert(ret >= 0, ExcMessage("Error at H5Dwrite"));
2214
2216 io_mode,
2219 mpi,
2221
2222 ret = H5Sclose(memory_dataspace);
2223 Assert(ret >= 0, ExcMessage("Error at H5Sclose"));
2224
2225 (void)ret;
2226 }
2227
2228
2229
2230 template <typename number>
2231 DataSet
2232 Group::create_dataset(const std::string &name,
2233 const std::vector<hsize_t> &dimensions) const
2234 {
2235 std::shared_ptr<hid_t> t_type = internal::get_hdf5_datatype<number>();
2236 return {name, *hdf5_reference, dimensions, t_type, mpi};
2237 }
2238
2239
2240
2241 template <typename Container>
2242 void
2243 Group::write_dataset(const std::string &name, const Container &data) const
2244 {
2245 std::vector<hsize_t> dimensions = internal::get_container_dimensions(data);
2246 auto dataset =
2247 create_dataset<typename Container::value_type>(name, dimensions);
2248 dataset.write(data);
2249 }
2250} // namespace HDF5
2251
2253
2254
2255#endif // DEAL_II_WITH_HDF5
2256
2257#endif // dealii_hdf5_h
ArrayView< std::remove_reference_t< typename std::iterator_traits< Iterator >::reference >, MemorySpaceType > make_array_view(const Iterator begin, const Iterator end)
Definition array_view.h:954
unsigned int rank
Definition hdf5.h:923
unsigned int get_rank() const
Definition hdf5.cc:316
void write(const Container &data)
Definition hdf5.h:2002
bool query_io_mode
Definition hdf5.h:950
std::vector< hsize_t > get_dimensions() const
Definition hdf5.cc:200
std::vector< hsize_t > dimensions
Definition hdf5.h:929
std::uint32_t get_local_no_collective_cause_as_hdf5_type()
Definition hdf5.cc:265
void write_hyperslab(const Container &data, const std::vector< hsize_t > &offset, const std::vector< hsize_t > &count)
Definition hdf5.h:2083
unsigned int get_size() const
Definition hdf5.cc:308
void write_selection(const Container &data, const std::vector< hsize_t > &coordinates)
Definition hdf5.h:2034
void read_none()
Definition hdf5.h:1962
std::shared_ptr< hid_t > dataspace
Definition hdf5.h:934
std::string get_local_no_collective_cause()
Definition hdf5.cc:253
std::string get_io_mode()
Definition hdf5.cc:208
std::uint32_t local_no_collective_cause
Definition hdf5.h:962
std::string get_global_no_collective_cause()
Definition hdf5.cc:277
void write_none()
Definition hdf5.h:2192
Container read()
Definition hdf5.h:1768
bool get_query_io_mode() const
Definition hdf5.cc:300
std::uint32_t global_no_collective_cause
Definition hdf5.h:969
std::uint32_t get_global_no_collective_cause_as_hdf5_type()
Definition hdf5.cc:289
Container read_selection(const std::vector< hsize_t > &coordinates)
Definition hdf5.h:1803
H5D_mpio_actual_io_mode_t get_io_mode_as_hdf5_type()
Definition hdf5.cc:242
void set_query_io_mode(const bool new_query_io_mode)
Definition hdf5.cc:192
Container read_hyperslab(const std::vector< hsize_t > &offset, const std::vector< hsize_t > &count)
Definition hdf5.h:1855
unsigned int size
Definition hdf5.h:939
H5D_mpio_actual_io_mode_t io_mode
Definition hdf5.h:955
FileAccessMode
Definition hdf5.h:1079
DataSet open_dataset(const std::string &name) const
Definition hdf5.cc:383
Group open_group(const std::string &name) const
Definition hdf5.cc:367
DataSet create_dataset(const std::string &name, const std::vector< hsize_t > &dimensions) const
Definition hdf5.h:2232
Group create_group(const std::string &name) const
Definition hdf5.cc:375
GroupAccessMode
Definition hdf5.h:982
void write_dataset(const std::string &name, const Container &data) const
Definition hdf5.h:2243
std::shared_ptr< hid_t > hdf5_reference
Definition hdf5.h:413
std::string get_name() const
Definition hdf5.cc:76
T get_attribute(const std::string &attr_name) const
Definition hdf5.h:1596
const std::string name
Definition hdf5.h:404
void set_attribute(const std::string &attr_name, const T value)
Definition hdf5.h:1673
const bool mpi
Definition hdf5.h:418
#define DEAL_II_NAMESPACE_OPEN
Definition config.h:40
#define DEAL_II_NAMESPACE_CLOSE
Definition config.h:41
#define DEAL_II_ASSERT_UNREACHABLE()
static ::ExceptionBase & ExcNotImplemented()
#define Assert(cond, exc)
#define AssertDimension(dim1, dim2)
#define AssertNothrow(cond, exc)
static ::ExceptionBase & ExcInternalError()
static ::ExceptionBase & ExcMessage(std::string arg1)
#define AssertThrow(cond, exc)
std::vector< index_type > data
Definition mpi.cc:746
void release_plist(hid_t &plist, H5D_mpio_actual_io_mode_t &io_mode, std::uint32_t &local_no_collective_cause, std::uint32_t &global_no_collective_cause, const bool mpi, const bool query_io_mode)
Definition hdf5.h:1496
std::string no_collective_cause_to_string(const std::uint32_t no_collective_cause)
Definition hdf5.h:1535
std::enable_if_t< std::is_same_v< Container, std::vector< typename Container::value_type > >, Container > initialize_container(const std::vector< hsize_t > &dimensions)
Definition hdf5.h:1420
unsigned int get_container_size(const std::vector< number > &data)
Definition hdf5.h:1391
void set_plist(hid_t &plist, const bool mpi)
Definition hdf5.h:1471
std::shared_ptr< hid_t > get_hdf5_datatype()
Definition hdf5.h:1274
std::vector< hsize_t > get_container_dimensions(const std::vector< number > &data)
Definition hdf5.h:1361
Definition hdf5.h:344