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