10 #include "xdmf_utils.h"
11 #include <Eigen/Dense>
12 #include <boost/algorithm/string.hpp>
13 #include <boost/lexical_cast.hpp>
14 #include <dolfinx/common/IndexMap.h>
15 #include <dolfinx/mesh/Geometry.h>
16 #include <dolfinx/mesh/Mesh.h>
17 #include <dolfinx/mesh/Topology.h>
25 std::vector<T>
get_dataset(MPI_Comm comm,
const pugi::xml_node& dataset_node,
27 std::array<std::int64_t, 2> range = {{0, 0}})
38 pugi::xml_attribute format_attr = dataset_node.attribute(
"Format");
43 const std::vector shape_xml = xdmf_utils::get_dataset_shape(dataset_node);
45 const std::string format = format_attr.as_string();
46 std::vector<T> data_vector;
54 pugi::xml_node data_node = dataset_node.first_child();
56 std::string data_str = data_node.value();
59 std::vector<boost::iterator_range<std::string::iterator>> data_vector_str;
60 boost::split(data_vector_str, data_str, boost::is_any_of(
" \n"));
63 data_vector.reserve(data_vector_str.size());
64 for (
auto& v : data_vector_str)
66 if (v.begin() != v.end())
67 data_vector.push_back(
68 boost::lexical_cast<T>(boost::copy_range<std::string>(v)));
72 else if (format ==
"HDF")
75 auto paths = xdmf_utils::get_hdf5_paths(dataset_node);
78 const std::vector shape_hdf5
83 assert(!shape_hdf5.empty());
84 assert(shape_hdf5[0] != 0);
93 if (range[0] == 0 and range[1] == 0)
95 if (shape_xml == shape_hdf5)
100 else if (!shape_xml.empty() and shape_hdf5.size() == 1)
103 std::int64_t d = std::accumulate(shape_xml.begin(), shape_xml.end(), 1,
104 std::multiplies<std::int64_t>());
107 if (d * shape_xml[0] != shape_hdf5[0])
109 throw std::runtime_error(
"Data size in XDMF/XML and size of HDF5 "
110 "dataset are inconsistent");
121 throw std::runtime_error(
122 "This combination of array shapes in XDMF and HDF5 "
128 data_vector = HDF5Interface::read_dataset<T>(h5_id, paths[1], range);
131 throw std::runtime_error(
"Storage format \"" + format +
"\" is unknown");
134 if (shape_xml.empty())
136 std::int64_t size = 1;
137 for (
auto dim : shape_xml)
140 std::int64_t size_global = 0;
141 const std::int64_t size_local = data_vector.size();
142 MPI_Allreduce(&size_local, &size_global, 1, MPI_INT64_T, MPI_SUM, comm);
143 if (size != size_global)
145 throw std::runtime_error(
146 "Data sizes in attribute and size of data read are inconsistent");
static std::array< std::int64_t, 2 > local_range(int process, std::int64_t N, int size)
Return local range for given process, splitting [0, N - 1] into size() portions of almost equal size.
Definition: MPI.cpp:105
static int rank(MPI_Comm comm)
Return process rank for the communicator.
Definition: MPI.cpp:79
static int size(MPI_Comm comm)
Return size of the group (number of processes) associated with the communicator.
Definition: MPI.cpp:87
static std::vector< std::int64_t > get_dataset_shape(const hid_t handle, const std::string &dataset_path)
Get dataset shape (size of each dimension)
Definition: HDF5Interface.cpp:196
Low-level methods for reading XDMF files.
Definition: xdmf_read.h:21
std::vector< T > get_dataset(MPI_Comm comm, const pugi::xml_node &dataset_node, const hid_t h5_id, std::array< std::int64_t, 2 > range={{0, 0}})
Return data associated with a data set node.
Definition: xdmf_read.h:25