DOLFIN-X
DOLFIN-X 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 "utils.h"
12 #include <array>
13 #include <dolfinx/common/array2d.h>
14 #include <dolfinx/common/span.hpp>
15 #include <dolfinx/common/utils.h>
16 #include <dolfinx/mesh/cell_types.h>
17 #include <string>
18 #include <utility>
19 #include <vector>
20 
21 namespace pugi
22 {
23 class xml_node;
24 } // namespace pugi
25 
26 namespace dolfinx
27 {
28 
29 namespace fem
30 {
31 template <typename T>
32 class Function;
33 } // namespace fem
34 
35 namespace fem
36 {
37 class CoordinateElement;
38 }
39 
40 namespace mesh
41 {
42 class Mesh;
43 }
44 
45 namespace io::xdmf_utils
46 {
47 
48 // Get DOLFINX cell type string from XML topology node
49 // @return DOLFINX cell type and polynomial degree
50 std::pair<std::string, int> get_cell_type(const pugi::xml_node& topology_node);
51 
52 // Return (0) HDF5 filename and (1) path in HDF5 file from a DataItem
53 // node
54 std::array<std::string, 2> get_hdf5_paths(const pugi::xml_node& dataitem_node);
55 
56 std::string get_hdf5_filename(std::string xdmf_filename);
57 
59 std::vector<std::int64_t> get_dataset_shape(const pugi::xml_node& dataset_node);
60 
62 std::int64_t get_num_cells(const pugi::xml_node& topology_node);
63 
66 std::vector<double> get_point_data_values(const fem::Function<double>& u);
67 std::vector<std::complex<double>>
68 get_point_data_values(const fem::Function<std::complex<double>>& u);
69 
71 std::vector<double> get_cell_data_values(const fem::Function<double>& u);
72 std::vector<std::complex<double>>
73 get_cell_data_values(const fem::Function<std::complex<double>>& u);
74 
76 std::string vtk_cell_type_str(mesh::CellType cell_type, int num_nodes);
77 
97 std::pair<array2d<std::int32_t>, std::vector<std::int32_t>>
98 extract_local_entities(const mesh::Mesh& mesh, int entity_dim,
99  const array2d<std::int64_t>& entities,
100  const tcb::span<const std::int32_t>& values);
101 
103 template <typename T>
104 void add_data_item(pugi::xml_node& xml_node, const hid_t h5_id,
105  const std::string h5_path, const T& x,
106  const std::int64_t offset,
107  const std::vector<std::int64_t> shape,
108  const std::string number_type, const bool use_mpi_io)
109 {
110  // Add DataItem node
111  assert(xml_node);
112  pugi::xml_node data_item_node = xml_node.append_child("DataItem");
113  assert(data_item_node);
114 
115  // Add dimensions attribute
116  std::string dims;
117  for (auto d : shape)
118  dims += std::to_string(d) + " ";
119  dims.pop_back();
120  data_item_node.append_attribute("Dimensions") = dims.c_str();
121 
122  // Set type for topology data (needed by XDMF to prevent default to
123  // float)
124  if (!number_type.empty())
125  data_item_node.append_attribute("NumberType") = number_type.c_str();
126 
127  // Add format attribute
128  if (h5_id < 0)
129  {
130  data_item_node.append_attribute("Format") = "XML";
131  assert(shape.size() == 2);
132  data_item_node.append_child(pugi::node_pcdata)
133  .set_value(common::container_to_string(x, 16, shape[1]).c_str());
134  }
135  else
136  {
137  data_item_node.append_attribute("Format") = "HDF";
138 
139  // Get name of HDF5 file, including path
140  const std::string hdf5_filename = HDF5Interface::get_filename(h5_id);
141  const std::string filename = dolfinx::io::get_filename(hdf5_filename);
142 
143  // Add HDF5 filename and HDF5 internal path to XML file
144  const std::string xdmf_path = filename + ":" + h5_path;
145  data_item_node.append_child(pugi::node_pcdata).set_value(xdmf_path.c_str());
146 
147  // Compute total number of items and check for consistency with shape
148  assert(!shape.empty());
149  std::int64_t num_items_total = 1;
150  for (auto n : shape)
151  num_items_total *= n;
152 
153  // Compute data offset and range of values
154  std::int64_t local_shape0 = x.size();
155  for (std::size_t i = 1; i < shape.size(); ++i)
156  {
157  assert(local_shape0 % shape[i] == 0);
158  local_shape0 /= shape[i];
159  }
160 
161  const std::array local_range{offset, offset + local_shape0};
162  HDF5Interface::write_dataset(h5_id, h5_path, x.data(), local_range, shape,
163  use_mpi_io, false);
164 
165  // Add partitioning attribute to dataset
166  // std::vector<std::size_t> partitions;
167  // std::vector<std::size_t> offset_tmp(1, offset);
168  // dolfinx::MPI::gather(comm, offset_tmp, partitions);
169  // dolfinx::MPI::broadcast(comm, partitions);
170  // HDF5Interface::add_attribute(h5_id, h5_path, "partition", partitions);
171  }
172 }
173 
174 } // namespace io::xdmf_utils
175 } // namespace dolfinx
static std::string get_filename(hid_t handle)
Get filename.
Definition: HDF5Interface.cpp:128
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.
std::string container_to_string(const T &x, const int precision, const int linebreak)
Convert a container to string.
Definition: utils.h:62
std::string get_filename(const std::string &fullname)
Get filename from a fully qualified path and filename.
Definition: utils.cpp:13
CellType
Cell type identifier.
Definition: cell_types.h:22