DOLFIN-X
DOLFIN-X C++ interface
Namespaces | Functions
dolfinx::refinement Namespace Reference

Mesh refinement algorithms. More...

Namespaces

 PlazaRefinementND
 Implementation of the refinement method described in Plaza and Carey "Local refinement of simplicial grids based on the skeleton" (Applied Numerical Mathematics 32 (2000) 195-218)
 

Functions

mesh::Mesh refine (const mesh::Mesh &mesh, bool redistribute=true)
 Create uniformly refined mesh. More...
 
mesh::Mesh refine (const mesh::Mesh &mesh, const mesh::MeshTags< std::int8_t > &cell_markers, bool redistribute=true)
 Create locally refined mesh. More...
 
std::pair< MPI_Comm, std::map< std::int32_t, std::set< int > > > compute_edge_sharing (const mesh::Mesh &mesh)
 Compute the sharing of edges between processes. The resulting MPI_Comm is over the neighborhood of shared edges, allowing direct communication between peers. The resulting map is from local edge index to the set of neighbors (within the comm) that share that edge. More...
 
void update_logical_edgefunction (const MPI_Comm &neighbor_comm, const std::vector< std::vector< std::int32_t >> &marked_for_update, std::vector< bool > &marked_edges, const common::IndexMap &map_e)
 Transfer marked edges between processes. More...
 
std::pair< std::map< std::int32_t, std::int64_t >, Eigen::Array< double, Eigen::Dynamic, Eigen::Dynamic, Eigen::RowMajor > > create_new_vertices (const MPI_Comm &neighbor_comm, const std::map< std::int32_t, std::set< std::int32_t >> &shared_edges, const mesh::Mesh &mesh, const std::vector< bool > &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. More...
 
mesh::Mesh partition (const mesh::Mesh &old_mesh, const std::vector< std::int64_t > &cell_topology, int num_ghost_cells, const Eigen::Array< double, Eigen::Dynamic, Eigen::Dynamic, Eigen::RowMajor > &new_vertex_coordinates, bool redistribute)
 Use vertex and topology data to partition new mesh across processes. More...
 
mesh::Mesh build_local (const mesh::Mesh &old_mesh, const std::vector< std::int64_t > &cell_topology, const Eigen::Array< double, Eigen::Dynamic, Eigen::Dynamic, Eigen::RowMajor > &new_vertex_coordinates)
 Build local mesh from internal data when not running in parallel. More...
 
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 topological vertices on each process into the space of the index map. More...
 

Detailed Description

Mesh refinement algorithms.

Methods for refining meshes uniformly, or with markers, using edge bisection.

Function Documentation

◆ adjust_indices()

std::vector< std::int64_t > dolfinx::refinement::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 topological vertices on each process into the space of the index map.

Parameters
index_mapIndexMap of current mesh vertices
nNumber of new entries to be accommodated on this process
Returns
Global indices as if "n" extra values are appended on each process

◆ build_local()

mesh::Mesh dolfinx::refinement::build_local ( const mesh::Mesh old_mesh,
const std::vector< std::int64_t > &  cell_topology,
const Eigen::Array< double, Eigen::Dynamic, Eigen::Dynamic, Eigen::RowMajor > &  new_vertex_coordinates 
)

Build local mesh from internal data when not running in parallel.

Parameters
[in]old_mesh
[in]cell_topology
[in]new_vertex_coordinates
Returns
A Mesh

◆ compute_edge_sharing()

std::pair< MPI_Comm, std::map< std::int32_t, std::set< int > > > dolfinx::refinement::compute_edge_sharing ( const mesh::Mesh mesh)

Compute the sharing of edges between processes. The resulting MPI_Comm is over the neighborhood of shared edges, allowing direct communication between peers. The resulting map is from local edge index to the set of neighbors (within the comm) that share that edge.

Parameters
[in]meshMesh
Returns
pair of comm and map

◆ create_new_vertices()

std::pair< std::map< std::int32_t, std::int64_t >, Eigen::Array< double, Eigen::Dynamic, Eigen::Dynamic, Eigen::RowMajor > > dolfinx::refinement::create_new_vertices ( const MPI_Comm &  neighbor_comm,
const std::map< std::int32_t, std::set< std::int32_t >> &  shared_edges,
const mesh::Mesh mesh,
const std::vector< bool > &  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
[in]neighbor_commMPI Communicator for neighborhood
[in]shared_edges
[in]meshExisting mesh
[in]marked_edges
Returns
edge_to_new_vertex map and geometry array

◆ partition()

mesh::Mesh dolfinx::refinement::partition ( const mesh::Mesh old_mesh,
const std::vector< std::int64_t > &  cell_topology,
int  num_ghost_cells,
const Eigen::Array< double, Eigen::Dynamic, Eigen::Dynamic, Eigen::RowMajor > &  new_vertex_coordinates,
bool  redistribute 
)

Use vertex and topology data to partition new mesh across processes.

Parameters
[in]old_mesh
[in]cell_topologyTopology of cells, (vertex indices)
[in]num_ghost_cellsNumber of cells which are ghost (at end of list)
[in]new_vertex_coordinates
[in]redistributeCall graph partitioner if true
Returns
New mesh

◆ refine() [1/2]

mesh::Mesh dolfinx::refinement::refine ( const mesh::Mesh mesh,
bool  redistribute = true 
)

Create uniformly refined mesh.

Parameters
[in]meshThe mesh from which to build a refined Mesh
[in]redistributeOptional argument to redistribute the refined mesh if mesh is a distributed mesh.
Returns
A refined mesh

◆ refine() [2/2]

mesh::Mesh dolfinx::refinement::refine ( const mesh::Mesh mesh,
const mesh::MeshTags< std::int8_t > &  cell_markers,
bool  redistribute = true 
)

Create locally refined mesh.

Parameters
[in]meshThe mesh from which to build a refined Mesh
[in]cell_markersA mesh function over integers specifying which cells should be refined (value == 1) (and which should not (any other integer value)).
[in]redistributeOptional argument to redistribute the refined mesh if mesh is a distributed mesh.
Returns
A locally refined mesh

◆ update_logical_edgefunction()

void dolfinx::refinement::update_logical_edgefunction ( const MPI_Comm &  neighbor_comm,
const std::vector< std::vector< std::int32_t >> &  marked_for_update,
std::vector< bool > &  marked_edges,
const common::IndexMap map_e 
)

Transfer marked edges between processes.

Parameters
neighbor_commMPI Communicator for neighborhood
marked_for_updateLists of edges to be updates on each neighbor
marked_edgesMarked edges to be updated
map_eIndexMap for edges