DOLFIN-X
DOLFIN-X C++ interface
ParallelRefinement.h
1 // Copyright (C) 2012-2018 Chris 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 <Eigen/Dense>
10 #include <cstdint>
11 #include <dolfinx/common/MPI.h>
12 #include <map>
13 #include <set>
14 #include <vector>
15 
16 namespace dolfinx
17 {
18 
19 namespace mesh
20 {
21 class Mesh;
22 template <typename T>
23 class MeshTags;
24 } // namespace mesh
25 
26 namespace common
27 {
28 class IndexMap;
29 }
30 
31 namespace refinement
32 {
34 
39 
41 {
42 public:
44  ParallelRefinement(const mesh::Mesh& mesh);
45 
47  ParallelRefinement(const ParallelRefinement& p) = delete;
48 
51 
54 
57  const std::vector<bool>& marked_edges() const;
58 
62  bool mark(std::int32_t edge_index);
63 
65  void mark_all();
66 
70  void mark(const mesh::MeshTags<std::int8_t>& refinement_marker);
71 
74 
79  std::map<std::int32_t, std::int64_t> create_new_vertices();
80 
85  mesh::Mesh partition(const std::vector<std::int64_t>& cell_topology,
86  bool redistribute) const;
87 
91  mesh::Mesh build_local(const std::vector<std::int64_t>& cell_topology) const;
92 
100  static std::vector<std::int64_t>
101  adjust_indices(const std::shared_ptr<const common::IndexMap>& index_map,
102  std::int32_t n);
103 
104 private:
105  // Mesh
106  const mesh::Mesh& _mesh;
107 
108  // New storage for all coordinates when creating new vertices
109  Eigen::Array<double, Eigen::Dynamic, 3, Eigen::RowMajor>
110  _new_vertex_coordinates;
111 
112  // Management of marked edges
113  std::vector<bool> _marked_edges;
114 
115  // Temporary storage for edges that have been recently marked (global
116  // index)
117  std::vector<std::vector<std::int64_t>> _marked_for_update;
118 
119  // Shared edges between processes
120  std::map<std::int32_t, std::set<std::int32_t>> _shared_edges;
121 
122  // Neighbourhood communicator
123  MPI_Comm _neighbour_comm;
124 };
125 } // namespace refinement
126 } // namespace dolfinx
dolfinx::refinement::ParallelRefinement::build_local
mesh::Mesh build_local(const std::vector< std::int64_t > &cell_topology) const
Build local mesh from internal data when not running in parallel.
Definition: ParallelRefinement.cpp:330
dolfinx::refinement::ParallelRefinement::mark_all
void mark_all()
Mark all edges in mesh.
Definition: ParallelRefinement.cpp:163
dolfinx::refinement::ParallelRefinement::partition
mesh::Mesh partition(const std::vector< std::int64_t > &cell_topology, bool redistribute) const
Use vertex and topology data to partition new mesh across processes.
Definition: ParallelRefinement.cpp:350
dolfinx::refinement::ParallelRefinement::mark
bool mark(std::int32_t edge_index)
Mark edge by index.
Definition: ParallelRefinement.cpp:140
dolfinx::refinement::ParallelRefinement::marked_edges
const std::vector< bool > & marked_edges() const
Return markers for all edges.
Definition: ParallelRefinement.cpp:135
dolfinx::refinement::ParallelRefinement::~ParallelRefinement
~ParallelRefinement()
Destructor.
Definition: ParallelRefinement.cpp:133
dolfinx::refinement::ParallelRefinement::update_logical_edgefunction
void update_logical_edgefunction()
Transfer marked edges between processes.
Definition: ParallelRefinement.cpp:193
dolfinx::refinement::ParallelRefinement
Data structure and methods for refining meshes in parallel.
Definition: ParallelRefinement.h:40
dolfinx::mesh::MeshTags
A MeshTags are used to associate mesh entities with values. The entity index (local to process) ident...
Definition: Form.h:35
dolfinx::mesh::Mesh
A Mesh consists of a set of connected and numbered mesh topological entities, and geometry data.
Definition: Mesh.h:46
dolfinx::refinement::ParallelRefinement::create_new_vertices
std::map< std::int32_t, std::int64_t > create_new_vertices()
Add new vertex for each marked edge, and create new_vertex_coordinates and global_edge->new_vertex ma...
Definition: ParallelRefinement.cpp:219
dolfinx::refinement::ParallelRefinement::adjust_indices
static std::vector< std::int64_t > adjust_indices(const std::shared_ptr< const common::IndexMap > &index_map, std::int32_t n)
Adjust indices to account for extra n values on each process This is a utility to help add new topolo...
Definition: ParallelRefinement.cpp:302
dolfinx::refinement::ParallelRefinement::ParallelRefinement
ParallelRefinement(const mesh::Mesh &mesh)
Constructor.
Definition: ParallelRefinement.cpp:89
dolfinx::refinement::ParallelRefinement::operator=
ParallelRefinement & operator=(const ParallelRefinement &p)=delete
Disable copy assignment.