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 class Function;
31 } // namespace function
32 
33 namespace io::xdmf_utils
34 {
35 
36 // Get DOLFINX cell type string from XML topology node
37 // @return DOLFINX cell type and polynomial degree
38 std::pair<std::string, int> get_cell_type(const pugi::xml_node& topology_node);
39 
40 // Return (0) HDF5 filename and (1) path in HDF5 file from a DataItem
41 // node
42 std::array<std::string, 2> get_hdf5_paths(const pugi::xml_node& dataitem_node);
43 
44 std::string get_hdf5_filename(std::string xdmf_filename);
45 
47 std::vector<std::int64_t> get_dataset_shape(const pugi::xml_node& dataset_node);
48 
50 std::int64_t get_num_cells(const pugi::xml_node& topology_node);
51 
54 std::vector<PetscScalar> get_point_data_values(const function::Function& u);
55 
57 std::vector<PetscScalar> get_cell_data_values(const function::Function& u);
58 
60 std::string vtk_cell_type_str(mesh::CellType cell_type, int num_nodes);
61 
63 template <typename T>
64 void add_data_item(pugi::xml_node& xml_node, const hid_t h5_id,
65  const std::string h5_path, const T& x,
66  const std::int64_t offset,
67  const std::vector<std::int64_t> shape,
68  const std::string number_type, const bool use_mpi_io)
69 {
70  // Add DataItem node
71  assert(xml_node);
72  pugi::xml_node data_item_node = xml_node.append_child("DataItem");
73  assert(data_item_node);
74 
75  // Add dimensions attribute
76  std::string dims;
77  for (auto d : shape)
78  dims += std::to_string(d) + " ";
79  dims.pop_back();
80  data_item_node.append_attribute("Dimensions") = dims.c_str();
81 
82  // Set type for topology data (needed by XDMF to prevent default to
83  // float)
84  if (!number_type.empty())
85  data_item_node.append_attribute("NumberType") = number_type.c_str();
86 
87  // Add format attribute
88  if (h5_id < 0)
89  {
90  data_item_node.append_attribute("Format") = "XML";
91  assert(shape.size() == 2);
92  data_item_node.append_child(pugi::node_pcdata)
93  .set_value(common::container_to_string(x, " ", 16, shape[1]).c_str());
94  }
95  else
96  {
97  data_item_node.append_attribute("Format") = "HDF";
98 
99  // Get name of HDF5 file
100  const std::string hdf5_filename = HDF5Interface::get_filename(h5_id);
101  const boost::filesystem::path p(hdf5_filename);
102 
103  // Add HDF5 filename and HDF5 internal path to XML file
104  const std::string xdmf_path = p.filename().string() + ":" + h5_path;
105  data_item_node.append_child(pugi::node_pcdata).set_value(xdmf_path.c_str());
106 
107  // Compute total number of items and check for consistency with shape
108  assert(!shape.empty());
109  std::int64_t num_items_total = 1;
110  for (auto n : shape)
111  num_items_total *= n;
112 
113  // Compute data offset and range of values
114  std::int64_t local_shape0 = x.size();
115  for (std::size_t i = 1; i < shape.size(); ++i)
116  {
117  assert(local_shape0 % shape[i] == 0);
118  local_shape0 /= shape[i];
119  }
120 
121  const std::array<std::int64_t, 2> local_range
122  = {{offset, offset + local_shape0}};
123  HDF5Interface::write_dataset(h5_id, h5_path, x.data(), local_range, shape,
124  use_mpi_io, false);
125 
126  // Add partitioning attribute to dataset
127  // std::vector<std::size_t> partitions;
128  // std::vector<std::size_t> offset_tmp(1, offset);
129  // dolfinx::MPI::gather(comm, offset_tmp, partitions);
130  // dolfinx::MPI::broadcast(comm, partitions);
131  // HDF5Interface::add_attribute(h5_id, h5_path, "partition", partitions);
132  }
133 } // namespace
134 //----------------------------------------------------------------------------
135 
136 } // namespace io::xdmf_utils
137 } // 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::io::HDF5Interface::write_dataset
static void write_dataset(const hid_t file_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_mpio, bool use_chunking)
Write data to existing HDF file as defined by range blocks on each process data: data to be written,...
dolfinx::mesh::CellType
CellType
Cell type identifier.
Definition: cell_types.h:22
dolfinx::io::HDF5Interface::get_filename
static std::string get_filename(hid_t hdf5_file_handle)
Get filename.
Definition: HDF5Interface.cpp:88