DOLFINx
DOLFINx C++ interface
xdmf_utils.h
1// Copyright (C) 2012 Chris N. Richardson
2//
3// This file is part of DOLFINx (https://www.fenicsproject.org)
4//
5// SPDX-License-Identifier: LGPL-3.0-or-later
6
7#pragma once
8
9#include "HDF5Interface.h"
10#include "pugixml.hpp"
11#include <array>
12#include <dolfinx/common/utils.h>
13#include <dolfinx/mesh/cell_types.h>
14#include <filesystem>
15#include <string>
16#include <utility>
17#include <vector>
18#include <xtl/xspan.hpp>
19
20namespace pugi
21{
22class xml_node;
23} // namespace pugi
24
25namespace dolfinx
26{
27
28namespace fem
29{
30template <typename T>
31class Function;
32} // namespace fem
33
34namespace fem
35{
36class CoordinateElement;
37}
38
39namespace mesh
40{
41class Mesh;
42}
43
44namespace io::xdmf_utils
45{
46
47// Get DOLFINx cell type string from XML topology node
48// @return DOLFINx cell type and polynomial degree
49std::pair<std::string, int> get_cell_type(const pugi::xml_node& topology_node);
50
51// Return (0) HDF5 filename and (1) path in HDF5 file from a DataItem
52// node
53std::array<std::string, 2> get_hdf5_paths(const pugi::xml_node& dataitem_node);
54
55std::filesystem::path
56get_hdf5_filename(const std::filesystem::path& xdmf_filename);
57
59std::vector<std::int64_t> get_dataset_shape(const pugi::xml_node& dataset_node);
60
62std::int64_t get_num_cells(const pugi::xml_node& topology_node);
63
66std::vector<double> get_point_data_values(const fem::Function<double>& u);
67std::vector<std::complex<double>>
68get_point_data_values(const fem::Function<std::complex<double>>& u);
69
71std::vector<double> get_cell_data_values(const fem::Function<double>& u);
72std::vector<std::complex<double>>
73get_cell_data_values(const fem::Function<std::complex<double>>& u);
74
76std::string vtk_cell_type_str(mesh::CellType cell_type, int num_nodes);
77
103std::pair<std::vector<std::int32_t>, std::vector<std::int32_t>>
104distribute_entity_data(const mesh::Mesh& mesh, int entity_dim,
105 const xtl::span<const std::int64_t>& entities,
106 const xtl::span<const std::int32_t>& data);
107
109template <typename T>
110void add_data_item(pugi::xml_node& xml_node, const hid_t h5_id,
111 const std::string& h5_path, const T& x, std::int64_t offset,
112 const std::vector<std::int64_t>& shape,
113 const std::string& number_type, bool use_mpi_io)
114{
115 // Add DataItem node
116 assert(xml_node);
117 pugi::xml_node data_item_node = xml_node.append_child("DataItem");
118 assert(data_item_node);
119
120 // Add dimensions attribute
121 std::string dims;
122 for (auto d : shape)
123 dims += std::to_string(d) + std::string(" ");
124 dims.pop_back();
125 data_item_node.append_attribute("Dimensions") = dims.c_str();
126
127 // Set type for topology data (needed by XDMF to prevent default to
128 // float)
129 if (!number_type.empty())
130 data_item_node.append_attribute("NumberType") = number_type.c_str();
131
132 // Add format attribute
133 if (h5_id < 0)
134 {
135 data_item_node.append_attribute("Format") = "XML";
136 assert(shape.size() == 2);
137 std::ostringstream s;
138 s.precision(16);
139 for (std::size_t i = 0; i < (std::size_t)x.size(); ++i)
140 {
141 if ((i + 1) % shape[1] == 0 and shape[1] != 0)
142 s << x.data()[i] << std::endl;
143 else
144 s << x.data()[i] << " ";
145 }
146
147 data_item_node.append_child(pugi::node_pcdata).set_value(s.str().c_str());
148 }
149 else
150 {
151 data_item_node.append_attribute("Format") = "HDF";
152
153 // Get name of HDF5 file, including path
154 const std::filesystem::path p = HDF5Interface::get_filename(h5_id);
155 const std::filesystem::path filename = p.filename().c_str();
156
157 // Add HDF5 filename and HDF5 internal path to XML file
158 const std::string xdmf_path
159 = filename.string() + std::string(":") + h5_path;
160 data_item_node.append_child(pugi::node_pcdata).set_value(xdmf_path.c_str());
161
162 // Compute data offset and range of values
163 std::int64_t local_shape0 = x.size();
164 for (std::size_t i = 1; i < shape.size(); ++i)
165 {
166 assert(local_shape0 % shape[i] == 0);
167 local_shape0 /= shape[i];
168 }
169
170 const std::array local_range{offset, offset + local_shape0};
171 HDF5Interface::write_dataset(h5_id, h5_path, x.data(), local_range, shape,
172 use_mpi_io, false);
173
174 // Add partitioning attribute to dataset
175 // std::vector<std::size_t> partitions;
176 // std::vector<std::size_t> offset_tmp(1, offset);
177 // dolfinx::MPI::gather(comm, offset_tmp, partitions);
178 // dolfinx::MPI::broadcast(comm, partitions);
179 // HDF5Interface::add_attribute(h5_id, h5_path, "partition", partitions);
180 }
181}
182
183} // namespace io::xdmf_utils
184} // namespace dolfinx
static std::filesystem::path get_filename(hid_t handle)
Get filename.
Definition: HDF5Interface.cpp:136
static void write_dataset(const hid_t handle, const std::string &dataset_path, const T *data, const std::array< std::int64_t, 2 > &range, const std::vector< std::int64_t > &global_size, bool use_mpi_io, bool use_chunking)
Write data to existing HDF file as defined by range blocks on each process.
constexpr std::array< std::int64_t, 2 > local_range(int rank, std::int64_t N, int size)
Return local range for the calling process, partitioning the global [0, N - 1] range across all ranks...
Definition: MPI.h:82
CellType
Cell type identifier.
Definition: cell_types.h:22