Refinement (dolfinx::refinement
)
-
namespace refinement
Mesh refinement algorithms.
Methods for refining meshes uniformly, or with markers, using edge bisection.
Functions
-
template<std::floating_point T>
mesh::Mesh<T> refine(const mesh::Mesh<T> &mesh, bool redistribute = true) Create a uniformly refined mesh.
- Parameters:
mesh – [in] Mesh to create a new, refined mesh from
redistribute – [in] If
true
refined mesh is re-partitioned across MPI ranks.
- Returns:
Refined mesh
-
template<std::floating_point T>
mesh::Mesh<T> refine(const mesh::Mesh<T> &mesh, std::span<const std::int32_t> edges, bool redistribute = true) Create a locally refined mesh.
- Parameters:
mesh – [in] Mesh to create a new, refined mesh from.
edges – [in] Indices of the edges that should be split during refinement. mesh::compute_incident_entities can be used to compute the edges that are incident to other entities, e.g. incident to cells.
redistribute – [in] If
true
refined mesh is re-partitioned across MPI ranks.
- Returns:
Refined mesh.
-
void update_logical_edgefunction(MPI_Comm comm, const std::vector<std::vector<std::int32_t>> &marked_for_update, std::vector<std::int8_t> &marked_edges, const common::IndexMap &map)
Communicate edge markers between processes that share edges.
- Parameters:
comm – [in] MPI Communicator for neighborhood
marked_for_update – [in] Lists of edges to be updated on each neighbor.
marked_for_update[r]
is the list of edge indices that are marked by the caller and are shared with local MPI rankr
.marked_edges – [inout] Marker for each edge on the calling process
map – [in] Index map for the mesh edges
-
template<std::floating_point T>
std::tuple<std::map<std::int32_t, std::int64_t>, std::vector<T>, std::array<std::size_t, 2>> create_new_vertices(MPI_Comm comm, const graph::AdjacencyList<int> &shared_edges, const mesh::Mesh<T> &mesh, std::span<const std::int8_t> marked_edges) Add new vertex for each marked edge, and create new_vertex_coordinates and global_edge->new_vertex map.
Communicate new vertices with MPI to all affected processes.
- Parameters:
comm – [in] MPI Communicator for neighborhood
shared_edges – [in] For each local edge on a process map to ghost processes
mesh – [in] Existing mesh
marked_edges – [in] A marker for all edges on the process (local + ghosts) indicating if an edge should be refined
- Returns:
(0) map from local edge index to new vertex global index, (1) the coordinates of the new vertices (row-major storage) and (2) the shape of the new coordinates.
-
template<std::floating_point T>
mesh::Mesh<T> partition(const mesh::Mesh<T> &old_mesh, const graph::AdjacencyList<std::int64_t> &cell_topology, std::span<const T> new_coords, std::array<std::size_t, 2> xshape, bool redistribute, mesh::GhostMode ghost_mode) Use vertex and topology data to partition new mesh across processes
- Parameters:
old_mesh – [in]
cell_topology – [in] Topology of cells, (vertex indices)
new_coords – [in] New coordinates, row-major storage
xshape – [in] The shape of
new_coords
redistribute – [in] Call graph partitioner if true
ghost_mode – [in] None or shared_facet
- Returns:
New mesh
-
std::vector<std::int64_t> adjust_indices(const common::IndexMap &map, std::int32_t n)
Given an index map, add “n” extra indices at the end of local range.
Note
The returned global indices (local and ghosts) are adjust for the addition of new local indices.
- Parameters:
map – [in] Index map
n – [in] Number of new local indices
- Returns:
New global indices for both owned and ghosted indices in input index map.
-
std::array<std::vector<std::int32_t>, 2> transfer_facet_meshtag(const mesh::MeshTags<std::int32_t> &tags0, const mesh::Topology &topology1, std::span<const std::int32_t> cell, std::span<const std::int8_t> facet)
Transfer facet MeshTags from coarse mesh to refined mesh.
Note
The refined mesh must not have been redistributed during refinement
Note
GhostMode must be GhostMode.none
- Parameters:
tags0 – [in] Tags on the parent mesh
topology1 – [in] Refined mesh topology
cell – [in] Parent cell of each cell in refined mesh
facet – [in] Local facets of parent in each cell in refined mesh
- Returns:
(0) entities and (1) values on the refined topology
-
std::array<std::vector<std::int32_t>, 2> transfer_cell_meshtag(const mesh::MeshTags<std::int32_t> &tags0, const mesh::Topology &topology1, std::span<const std::int32_t> parent_cell)
Transfer cell MeshTags from coarse mesh to refined mesh.
Note
The refined mesh must not have been redistributed during refinement.
Note
GhostMode must be GhostMode.none
- Parameters:
tags0 – [in] Tags on the parent mesh
topology1 – [in] Refined mesh topology
parent_cell – [in] Parent cell of each cell in refined mesh
- Returns:
(0) entities and (1) values on the refined topology
-
namespace impl
Functions
-
std::int64_t local_to_global(std::int32_t local_index, const common::IndexMap &map)
Compute global index.
-
template<std::floating_point T>
std::pair<std::vector<T>, std::array<std::size_t, 2>> create_new_geometry(const mesh::Mesh<T> &mesh, const std::map<std::int32_t, std::int64_t> &local_edge_to_new_vertex) Create geometric points of new Mesh, from current Mesh and a edge_to_vertex map listing the new local points (midpoints of those edges)
- Parameters:
mesh – Current mesh
local_edge_to_new_vertex – A map from a local edge to the new global vertex index that will be inserted at its midpoint
- Returns:
(1) Array of new (flattened) mesh geometry and (2) its multi-dimensional shape
-
std::int64_t local_to_global(std::int32_t local_index, const common::IndexMap &map)
-
namespace plaza
Plaza mesh refinement.
Functions for the refinement method described in Plaza and Carey “Local refinement of simplicial grids based on the skeleton”, Applied Numerical Mathematics 32 (2000), 195-218.
Enums
-
enum class Option : int
Options for data to compute during mesh refinement.
Values:
-
enumerator none
No extra data
-
enumerator parent_cell
Compute list with the parent cell index for each new cell
-
enumerator parent_facet
Compute list of the cell-local facet indices in the parent cell of each facet in each new cell (or -1 if no match)
-
enumerator parent_cell_and_facet
Both cell and facet parent data
-
enumerator none
Functions
-
template<std::floating_point T>
std::tuple<mesh::Mesh<T>, std::vector<std::int32_t>, std::vector<std::int8_t>> refine(const mesh::Mesh<T> &mesh, bool redistribute, Option option) Uniform refine, optionally redistributing and optionally calculating the parent-child relationships`.
- Parameters:
mesh – [in] Input mesh to be refined
redistribute – [in] Flag to call the mesh partitioner to redistribute after refinement
option – [in] Control the computation of parent facets, parent cells. If an option is unselected, an empty list is returned.
- Returns:
Refined mesh and optional parent cell index, parent facet indices
-
template<std::floating_point T>
std::tuple<mesh::Mesh<T>, std::vector<std::int32_t>, std::vector<std::int8_t>> refine(const mesh::Mesh<T> &mesh, std::span<const std::int32_t> edges, bool redistribute, Option option) Refine with markers, optionally redistributing, and optionally calculating the parent-child relationships.
- Parameters:
mesh – [in] Input mesh to be refined
edges – [in] Indices of the edges that should be split by this refinement
redistribute – [in] Flag to call the Mesh Partitioner to redistribute after refinement
option – [in] Control the computation of parent facets, parent cells. If an option is unselected, an empty list is returned.
- Returns:
New Mesh and optional parent cell index, parent facet indices
-
template<std::floating_point T>
std::tuple<graph::AdjacencyList<std::int64_t>, std::vector<T>, std::array<std::size_t, 2>, std::vector<std::int32_t>, std::vector<std::int8_t>> compute_refinement_data(const mesh::Mesh<T> &mesh, Option option) Refine mesh returning new mesh data.
- Parameters:
mesh – [in] Input mesh to be refined
option – [in] Control computation of parent facets and parent cells. If an option is unselected, an empty list is returned.
- Returns:
New mesh data: cell topology, vertex coordinates, vertex coordinates shape, and optional parent cell index, and parent facet indices.
-
template<std::floating_point T>
std::tuple<graph::AdjacencyList<std::int64_t>, std::vector<T>, std::array<std::size_t, 2>, std::vector<std::int32_t>, std::vector<std::int8_t>> compute_refinement_data(const mesh::Mesh<T> &mesh, std::span<const std::int32_t> edges, Option option) Refine with markers returning new mesh data.
- Parameters:
mesh – [in] Input mesh to be refined
edges – [in] Indices of the edges that should be split by this refinement
option – [in] Control the computation of parent facets, parent cells. If an option is unselected, an empty list is returned.
- Returns:
New mesh data: cell topology, vertex coordinates and parent cell index, and stored parent facet indices (if requested).
-
namespace impl
Functions
-
template<int tdim>
std::vector<std::int8_t> compute_parent_facets(std::span<const std::int32_t> simplex_set) Computes the parent-child facet relationship
- Parameters:
simplex_set – - index into indices for each child cell
- Returns:
mapping from child to parent facets, using cell-local index
-
std::vector<std::int32_t> get_simplices(std::span<const std::int64_t> indices, std::span<const std::int32_t> longest_edge, int tdim, bool uniform)
Get the subdivision of an original simplex into smaller simplices, for a given set of marked edges, and the longest edge of each facet (cell local indexing). A flag indicates if a uniform subdivision is preferable in 2D.
- Parameters:
indices – [in] Vector containing the global indices for the original vertices and potential new vertices at each edge. Size (num_vertices + num_edges). If an edge is not refined its corresponding entry is -1
longest_edge – [in] Vector indicating the longest edge for each triangle in the cell. For triangular cells (2D) there is only one value, and for tetrahedra (3D) there are four values, one for each facet. The values give the local edge indices of the cell.
tdim – [in] Topological dimension (2 or 3)
uniform – [in] Make a “uniform” subdivision with all triangles being similar shape
- Returns:
-
void enforce_rules(MPI_Comm comm, const graph::AdjacencyList<int> &shared_edges, std::vector<std::int8_t> &marked_edges, const mesh::Topology &topology, std::span<const std::int32_t> long_edge)
Propagate edge markers according to rules (longest edge of each face must be marked, if any edge of face is marked)
-
template<std::floating_point T>
std::pair<std::vector<std::int32_t>, std::vector<std::int8_t>> face_long_edge(const mesh::Mesh<T> &mesh) Get the longest edge of each face (using local mesh index)
Note
Edge ratio ok returns an array in 2D, where it checks if the ratio between the shortest and longest edge of a cell is bigger than sqrt(2)/2. In 3D an empty array is returned
- Parameters:
mesh – [in] The mesh
- Returns:
A tuple (longest edge, edge ratio ok) where longest edge gives the local index of the longest edge for each face.
-
template<std::floating_point T>
std::tuple<graph::AdjacencyList<std::int64_t>, std::vector<T>, std::array<std::size_t, 2>, std::vector<std::int32_t>, std::vector<std::int8_t>> compute_refinement(MPI_Comm neighbor_comm, const std::vector<std::int8_t> &marked_edges, const graph::AdjacencyList<int> &shared_edges, const mesh::Mesh<T> &mesh, const std::vector<std::int32_t> &long_edge, const std::vector<std::int8_t> &edge_ratio_ok, plaza::Option option) Convenient interface for both uniform and marker refinement.
Note
The parent facet map gives you the map from a cell given by parent cell map to the local index (relative to the cell), e.g. the i-th entry of parent facets relates to the local facet index of the i-th entry parent cell.
- Parameters:
neighbor_comm – [in] Neighbourhood communciator scattering owned edges to processes with ghosts
marked_edges – [in] A marker for all edges on the process (local + ghosts) indicating if an edge should be refined
shared_edges – [in] For each local edge on a process map to ghost processes
mesh – [in] The mesh
long_edge – [in] A map from each face to its longest edge. Index is local to the process.
edge_ratio_ok – [in] For each face in a 2D mesh this error contains a marker indicating if the ratio between smallest and largest edge is bigger than sqrt(2)/2
option – [in] Option to compute additional information relating refined and original mesh entities
- Returns:
(0) The new mesh topology, (1) the new flattened mesh geometry, (3) Shape of the new geometry_shape, (4) Map from new cells to parent cells and (5) map from refined facets to parent facets.
-
template<int tdim>
-
enum class Option : int
-
template<std::floating_point T>