DOLFINx
DOLFINx C++ interface
xdmf_meshtags.h
1// Copyright (C) 2020 Michal Habera
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 "pugixml.hpp"
10#include "xdmf_mesh.h"
11#include "xdmf_utils.h"
12#include <dolfinx/common/MPI.h>
13#include <dolfinx/common/log.h>
14#include <dolfinx/mesh/MeshTags.h>
15#include <hdf5.h>
16#include <string>
17#include <vector>
18#include <xtl/xspan.hpp>
19
20namespace dolfinx
21{
22namespace io
23{
24namespace xdmf_meshtags
25{
26
28template <typename T>
29void add_meshtags(MPI_Comm comm, const mesh::MeshTags<T>& meshtags,
30 pugi::xml_node& xml_node, const hid_t h5_id,
31 const std::string name)
32{
33 LOG(INFO) << "XDMF: add meshtags (" << name << ")";
34 // Get mesh
35 assert(meshtags.mesh());
36 std::shared_ptr<const mesh::Mesh> mesh = meshtags.mesh();
37 const int dim = meshtags.dim();
38
39 const std::int32_t num_local_entities
40 = mesh->topology().index_map(dim)->size_local();
41
42 // Find number of tagged entities in local range
43 auto it = std::lower_bound(meshtags.indices().begin(),
44 meshtags.indices().end(), num_local_entities);
45 const int num_active_entities = std::distance(meshtags.indices().begin(), it);
46
47 const std::string path_prefix = "/MeshTags/" + name;
49 comm, xml_node, h5_id, path_prefix, mesh->topology(), mesh->geometry(),
50 dim,
51 xtl::span<const std::int32_t>(meshtags.indices().data(),
52 num_active_entities));
53
54 // Add attribute node with values
55 pugi::xml_node attribute_node = xml_node.append_child("Attribute");
56 assert(attribute_node);
57 attribute_node.append_attribute("Name") = name.c_str();
58 attribute_node.append_attribute("AttributeType") = "Scalar";
59 attribute_node.append_attribute("Center") = "Cell";
60
61 std::int64_t global_num_values = 0;
62 const std::int64_t local_num_values = num_active_entities;
63 MPI_Allreduce(&local_num_values, &global_num_values, 1, MPI_INT64_T, MPI_SUM,
64 comm);
65 const std::int64_t num_local = num_active_entities;
66 std::int64_t offset = 0;
67 MPI_Exscan(&num_local, &offset, 1, MPI_INT64_T, MPI_SUM, comm);
68 const bool use_mpi_io = (dolfinx::MPI::size(comm) > 1);
69 xdmf_utils::add_data_item(
70 attribute_node, h5_id, path_prefix + std::string("/Values"),
71 xtl::span<const T>(meshtags.values().data(), num_active_entities), offset,
72 {global_num_values, 1}, "", use_mpi_io);
73}
74
75} // namespace xdmf_meshtags
76} // namespace io
77} // namespace dolfinx
int size(MPI_Comm comm)
Return size of the group (number of processes) associated with the communicator.
Definition: MPI.cpp:83
void add_topology_data(MPI_Comm comm, pugi::xml_node &xml_node, const hid_t h5_id, const std::string path_prefix, const mesh::Topology &topology, const mesh::Geometry &geometry, int cell_dim, const xtl::span< const std::int32_t > &entities)
Add Topology xml node.
Definition: xdmf_mesh.cpp:18