From a9c5de0b370502baa111357a9063dcc2f83a1e15 Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?J=C3=B8rgen=20Schartum=20Dokken?= Date: Tue, 16 Apr 2024 10:24:23 +0100 Subject: [PATCH] Improve mesh refinement documentation (#3010) * Document refinement functions * Explain parent_facet_index * Fix Garth's comments * More comments * Revert spaces * Apply suggestions from code review * Fix docs * Apply suggestions from code review Co-authored-by: Chris Richardson * Tweak formatting. * Typo. * Typos --------- Co-authored-by: Chris Richardson Co-authored-by: Jack S. Hale --- cpp/dolfinx/refinement/plaza.cpp | 16 +++++++++-- cpp/dolfinx/refinement/plaza.h | 48 +++++++++++++++++++++++++++----- cpp/dolfinx/refinement/utils.cpp | 3 -- cpp/dolfinx/refinement/utils.h | 42 +++++++++++++++------------- 4 files changed, 76 insertions(+), 33 deletions(-) diff --git a/cpp/dolfinx/refinement/plaza.cpp b/cpp/dolfinx/refinement/plaza.cpp index 89df2f56a8e..2691ed619af 100644 --- a/cpp/dolfinx/refinement/plaza.cpp +++ b/cpp/dolfinx/refinement/plaza.cpp @@ -24,13 +24,23 @@ using namespace dolfinx::refinement; namespace { //----------------------------------------------------------------------------- -// 2D version of subdivision allowing for uniform subdivision (flag) +/// 2D version of subdivision allowing for uniform subdivision (flag) +/// @param[in] indices Vector containing the +/// global indices for the original vertices and potential new vertices at each +/// edge. If an edge is not refined its corresponding entry is -1. Size +/// `num_vertices + num_edges` +/// @param[in] longest_edge Local index of the longest edge in the triangle. +/// @param[in] uniform If true, the triangle is subdivided into four similar +/// sub-triangles. +/// @returns Local indices for each sub-divived triangle std::vector get_triangles(std::span indices, const std::int32_t longest_edge, bool uniform) { - // v0 and v1 are at ends of longest_edge (e2) opposite vertex has same - // index as longest_edge + // NOTE: The assumption below is based on the UFC ordering of a triangle, i.e. + // that the N-th edge of a triangle is the edge where the N-th vertex is not + // part of the set. v0 and v1 are at ends of longest_edge (e2) opposite vertex + // has same index as longest_edge const std::int32_t v0 = (longest_edge + 1) % 3; const std::int32_t v1 = (longest_edge + 2) % 3; const std::int32_t v2 = longest_edge; diff --git a/cpp/dolfinx/refinement/plaza.h b/cpp/dolfinx/refinement/plaza.h index 5e17771177c..206cfb9bc4b 100644 --- a/cpp/dolfinx/refinement/plaza.h +++ b/cpp/dolfinx/refinement/plaza.h @@ -116,13 +116,16 @@ compute_parent_facets(std::span simplex_set) /// (cell local indexing). A flag indicates if a uniform subdivision is /// preferable in 2D. /// -/// @param[in] indices Vector indicating which edges are to be -/// split (value >=0) +/// @param[in] indices 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 /// @param[in] longest_edge Vector indicating the longest edge for each -/// triangle. For tdim=2, one entry, for tdim=3, four entries. +/// 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. /// @param[in] tdim Topological dimension (2 or 3) -/// @param[in] uniform Make a "uniform" subdivision with all triangles -/// being similar shape +/// @param[in] uniform Make a "uniform" subdivision with all triangles being +/// similar shape /// @return std::vector get_simplices(std::span indices, @@ -136,7 +139,15 @@ void enforce_rules(MPI_Comm comm, const graph::AdjacencyList& shared_edges, const mesh::Topology& topology, std::span long_edge); -/// Get the longest edge of each face (using local mesh index) +/// @brief 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 +/// +/// @param[in] mesh The mesh +/// @return A tuple (longest edge, edge ratio ok) where longest edge gives the +/// local index of the longest edge for each face. template std::pair, std::vector> face_long_edge(const mesh::Mesh& mesh) @@ -253,7 +264,30 @@ face_long_edge(const mesh::Mesh& mesh) return std::pair(std::move(long_edge), std::move(edge_ratio_ok)); } -/// Convenient interface for both uniform and marker refinement +/// @brief 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. +/// +/// @param[in] neighbor_comm Neighbourhood communciator scattering owned edges +/// to processes with ghosts +/// @param[in] marked_edges A marker for all edges on the process (local + +/// ghosts) indicating if an edge should be refined +/// @param[in] shared_edges For each local edge on a process map to ghost +/// processes +/// @param[in] mesh The mesh +/// @param[in] long_edge A map from each face to its longest edge. Index is +/// local to the process. +/// @param[in] edge_ratio_ok 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 +/// @param[in] option Option to compute additional information relating refined +/// and original mesh entities +/// @return (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 std::tuple, std::vector, std::array, std::vector, diff --git a/cpp/dolfinx/refinement/utils.cpp b/cpp/dolfinx/refinement/utils.cpp index 9e4d33b3fc7..e249329a8a6 100644 --- a/cpp/dolfinx/refinement/utils.cpp +++ b/cpp/dolfinx/refinement/utils.cpp @@ -103,9 +103,6 @@ refinement::adjust_indices(const common::IndexMap& map, std::int32_t n) { // NOTE: Is this effectively concatenating index maps? - // Add in an extra "n" indices at the end of the current local_range - // of "index_map", and adjust existing indices to match. - // Get offset for 'n' for this process const std::int64_t num_local = n; std::int64_t global_offset = 0; diff --git a/cpp/dolfinx/refinement/utils.h b/cpp/dolfinx/refinement/utils.h index e4cc6ce5fe3..7faa9664c8c 100644 --- a/cpp/dolfinx/refinement/utils.h +++ b/cpp/dolfinx/refinement/utils.h @@ -46,9 +46,12 @@ std::int64_t local_to_global(std::int32_t local_index, /// Create geometric points of new Mesh, from current Mesh and a /// edge_to_vertex map listing the new local points (midpoints of those /// edges) -/// @param mesh -/// @param local_edge_to_new_vertex -/// @return array of points +/// +/// @param mesh Current mesh +/// @param local_edge_to_new_vertex A map from a local edge to the new +/// global vertex index that will be inserted at its midpoint +/// @return (1) Array of new (flattened) mesh geometry and (2) its +/// multi-dimensional shape template std::pair, std::array> create_new_geometry( const mesh::Mesh& mesh, @@ -135,12 +138,14 @@ void update_logical_edgefunction( /// Communicate new vertices with MPI to all affected processes. /// /// @param[in] comm MPI Communicator for neighborhood -/// @param[in] shared_edges +/// @param[in] shared_edges For each local edge on a process map to ghost +/// processes /// @param[in] mesh Existing mesh -/// @param[in] marked_edges -/// @return (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. +/// @param[in] marked_edges A marker for all edges on the process (local + +/// ghosts) indicating if an edge should be refined +/// @return (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::tuple, std::vector, std::array> @@ -153,7 +158,8 @@ create_new_vertices(MPI_Comm comm, std::shared_ptr edge_index_map = mesh.topology()->index_map(1); - // Add new edge midpoints to list of vertices + // Add new edge midpoints to list of vertices. The new vertex will be owned by + // the process owning the edge. int n = 0; std::map local_edge_to_new_vertex; for (int local_i = 0; local_i < edge_index_map->size_local(); ++local_i) @@ -303,18 +309,14 @@ mesh::Mesh partition(const mesh::Mesh& old_mesh, } } -/// @todo Fix docstring. It is unclear. -/// -/// @brief Add indices to account for extra n values on this process. +/// @brief Given an index map, add "n" extra indices at the end of local range /// -/// This is a utility to help add new topological vertices on each -/// process into the space of the index map. -/// -/// @param[in] map Index map for the current mesh vertices -/// @param[in] n Number of new entries to be accommodated on this -/// process -/// @return Global indices as if "n" extra values are appended on each -/// process +/// @note The returned global indices (local and ghosts) are adjust for the +/// addition of new local indices. +/// @param[in] map Index map +/// @param[in] n Number of new local indices +/// @return New global indices for both owned and ghosted indices in input +/// index map. std::vector adjust_indices(const common::IndexMap& map, std::int32_t n);