Skip to content

Commit

Permalink
Propagate optional usage to plaza
Browse files Browse the repository at this point in the history 
Removal of duplicate functionality in plza::refine() and plaza::compute_refinement_data()
  • Loading branch information
@schnellerhase
schnellerhase committed Jul 24, 2024
1 parent 9a09ddb commit 84a6fa4
Show file tree
Hide file tree
Showing 4 changed files with 45 additions and 155 deletions.
141 changes: 19 additions & 122 deletions cpp/dolfinx/refinement/plaza.h
View file
Original file line number Diff line number Diff line change
Expand Up @@ -454,66 +454,23 @@ compute_refinement(MPI_Comm neighbor_comm,
}
} // namespace impl

/// @brief Uniform refine, optionally redistributing and optionally
/// calculating the parent-child relationships`.
///
/// @param[in] mesh Input mesh to be refined
/// @param[in] redistribute Flag to call the mesh partitioner to
/// redistribute after refinement
/// @param[in] option Control the computation of parent facets, parent
/// cells. If an option is unselected, an empty list is returned.
/// @return 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, bool redistribute, Option option)
{
auto [cell_adj, new_coords, xshape, parent_cell, parent_facet]
= compute_refinement_data(mesh, option);

if (dolfinx::MPI::size(mesh.comm()) == 1)
{
return {mesh::create_mesh(mesh.comm(), cell_adj.array(),
mesh.geometry().cmap(), new_coords, xshape,
mesh::GhostMode::none),
std::move(parent_cell), std::move(parent_facet)};
}
else
{
std::shared_ptr<const common::IndexMap> map_c
= mesh.topology()->index_map(mesh.topology()->dim());
const int num_ghost_cells = map_c->num_ghosts();
// Check if mesh has ghost cells on any rank
// FIXME: this is not a robust test. Should be user option.
int max_ghost_cells = 0;
MPI_Allreduce(&num_ghost_cells, &max_ghost_cells, 1, MPI_INT, MPI_MAX,
mesh.comm());

// Build mesh
const mesh::GhostMode ghost_mode = max_ghost_cells == 0
? mesh::GhostMode::none
: mesh::GhostMode::shared_facet;
return {partition<T>(mesh, cell_adj, std::span(new_coords), xshape,
redistribute, ghost_mode),
std::move(parent_cell), std::move(parent_facet)};
}
}

/// @brief Refine with markers, optionally redistributing, and
/// optionally calculating the parent-child relationships.
///
/// @param[in] mesh Input mesh to be refined
/// @param[in] edges Indices of the edges that should be split by this
/// refinement
/// @param[in] edges Optional indices of the edges that should be split by this
/// refinement, if optional is not set, a uniform refinement is performend, same
/// behavior as passing a list of all indices.
/// @param[in] redistribute Flag to call the Mesh Partitioner to
/// redistribute after refinement
/// @param[in] option Control the computation of parent facets, parent
/// cells. If an option is unselected, an empty list is returned.
/// @return New 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(const mesh::Mesh<T>& mesh,
std::optional<std::span<const std::int32_t>> edges, bool redistribute,
Option option)
{
auto [cell_adj, new_vertex_coords, xshape, parent_cell, parent_facet]
= compute_refinement_data(mesh, edges, option);
Expand Down Expand Up @@ -547,71 +504,6 @@ refine(const mesh::Mesh<T>& mesh, std::span<const std::int32_t> edges,
}
}

/// @brief Refine mesh returning new mesh data.
///
/// @param[in] mesh Input mesh to be refined
/// @param[in] option Control computation of parent facets and parent
/// cells. If an option is unselected, an empty list is returned.
/// @return 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, Option option)
{
common::Timer t0("PLAZA: refine");
auto topology = mesh.topology();
assert(topology);

if (topology->cell_type() != mesh::CellType::triangle
and topology->cell_type() != mesh::CellType::tetrahedron)
{
throw std::runtime_error("Cell type not supported");
}

auto map_e = topology->index_map(1);
if (!map_e)
throw std::runtime_error("Edges must be initialised");

// Get sharing ranks for each edge
graph::AdjacencyList<int> edge_ranks = map_e->index_to_dest_ranks();

// Create unique list of ranks that share edges (owners of ghosts
// plus ranks that ghost owned indices)
std::vector<int> ranks(edge_ranks.array().begin(), edge_ranks.array().end());
std::ranges::sort(ranks);
auto [unique_end, range_end] = std::ranges::unique(ranks);
ranks.erase(unique_end, range_end);

// Convert edge_ranks from global rank to to neighbourhood ranks
std::ranges::transform(edge_ranks.array(), edge_ranks.array().begin(),
[&ranks](auto r)
{
auto it = std::ranges::lower_bound(ranks, r);
assert(it != ranks.end() and *it == r);
return std::distance(ranks.begin(), it);
});

MPI_Comm comm;
MPI_Dist_graph_create_adjacent(mesh.comm(), ranks.size(), ranks.data(),
MPI_UNWEIGHTED, ranks.size(), ranks.data(),
MPI_UNWEIGHTED, MPI_INFO_NULL, false, &comm);

const auto [long_edge, edge_ratio_ok] = impl::face_long_edge(mesh);
auto [cell_adj, new_vertex_coords, xshape, parent_cell, parent_facet]
= impl::compute_refinement(
comm,
std::vector<std::int8_t>(map_e->size_local() + map_e->num_ghosts(),
true),
edge_ranks, mesh, long_edge, edge_ratio_ok, option);
MPI_Comm_free(&comm);

return {std::move(cell_adj), std::move(new_vertex_coords), xshape,
std::move(parent_cell), std::move(parent_facet)};
}

/// Refine with markers returning new mesh data.
///
/// @param[in] mesh Input mesh to be refined
Expand All @@ -626,7 +518,8 @@ 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)
std::optional<std::span<const std::int32_t>> edges,
Option option)
{
common::Timer t0("PLAZA: refine");
auto topology = mesh.topology();
Expand Down Expand Up @@ -663,17 +556,21 @@ compute_refinement_data(const mesh::Mesh<T>& mesh,

// Get number of neighbors
std::vector<std::int8_t> marked_edges(
map_e->size_local() + map_e->num_ghosts(), false);
map_e->size_local() + map_e->num_ghosts(), !edges.has_value());
std::vector<std::vector<std::int32_t>> marked_for_update(ranks.size());
for (auto edge : edges)

if (edges.has_value())
{
if (!marked_edges[edge])
for (auto edge : edges.value())
{
marked_edges[edge] = true;
if (!marked_edges[edge])
{
marked_edges[edge] = true;

// If it is a shared edge, add all sharing neighbors to update set
for (int rank : edge_ranks.links(edge))
marked_for_update[rank].push_back(edge);
// If it is a shared edge, add all sharing neighbors to update set
for (int rank : edge_ranks.links(edge))
marked_for_update[rank].push_back(edge);
}
}
}

Expand Down
5 changes: 1 addition & 4 deletions cpp/dolfinx/refinement/refine.h
View file
Original file line number Diff line number Diff line change
Expand Up @@ -38,10 +38,7 @@ mesh::Mesh<T> refine(const mesh::Mesh<T>& mesh,
}

auto [refined_mesh, parent_cell, parent_facet]
= edges.has_value()
? plaza::refine(mesh, edges.value(), redistribute,
plaza::Option::none)
: plaza::refine(mesh, redistribute, plaza::Option::none);
= plaza::refine(mesh, edges, redistribute, plaza::Option::none);

// Report the number of refined cellse
const int D = topology->dim();
Expand Down
9 changes: 3 additions & 6 deletions python/dolfinx/mesh.py
View file
Original file line number Diff line number Diff line change
Expand Up @@ -385,12 +385,9 @@ def refine_plaza(
Returns:
Refined mesh, list of parent cell for each refine cell, and list
"""
if edges is None:
mesh1, cells, facets = _cpp.refinement.refine_plaza(mesh._cpp_object, redistribute, option)
else:
mesh1, cells, facets = _cpp.refinement.refine_plaza(
mesh._cpp_object, edges, redistribute, option
)
mesh1, cells, facets = _cpp.refinement.refine_plaza(
mesh._cpp_object, edges, redistribute, option
)
return Mesh(mesh1, mesh._ufl_domain), cells, facets


Expand Down
45 changes: 22 additions & 23 deletions python/dolfinx/wrappers/refinement.cpp
View file
Original file line number Diff line number Diff line change
Expand Up @@ -4,20 +4,25 @@
//
// SPDX-License-Identifier: LGPL-3.0-or-later

#include "array.h"
#include <concepts>
#include <optional>
#include <span>

#include <nanobind/nanobind.h>
#include <nanobind/ndarray.h>
#include <nanobind/stl/optional.h>
#include <nanobind/stl/shared_ptr.h>
#include <nanobind/stl/tuple.h>
#include <nanobind/stl/vector.h>

#include <dolfinx/mesh/Mesh.h>
#include <dolfinx/mesh/MeshTags.h>
#include <dolfinx/refinement/interval.h>
#include <dolfinx/refinement/plaza.h>
#include <dolfinx/refinement/refine.h>
#include <dolfinx/refinement/utils.h>
#include <nanobind/nanobind.h>
#include <nanobind/ndarray.h>
#include <nanobind/stl/shared_ptr.h>
#include <nanobind/stl/tuple.h>
#include <nanobind/stl/vector.h>
#include <optional>

#include "array.h"

namespace nb = nanobind;

Expand Down Expand Up @@ -74,31 +79,25 @@ void export_refinement_with_variable_mesh_type(nb::module_& m)
nb::arg("mesh"), nb::arg("edges"), nb::arg("redistribute"),
nb::arg("ghost_mode"));

m.def(
"refine_plaza",
[](const dolfinx::mesh::Mesh<T>& mesh0, bool redistribute,
dolfinx::refinement::plaza::Option option)
{
auto [mesh1, cell, facet]
= dolfinx::refinement::plaza::refine(mesh0, redistribute, option);
return std::tuple{std::move(mesh1), as_nbarray(std::move(cell)),
as_nbarray(std::move(facet))};
},
nb::arg("mesh"), nb::arg("redistribute"), nb::arg("option"));

m.def(
"refine_plaza",
[](const dolfinx::mesh::Mesh<T>& mesh0,
nb::ndarray<const std::int32_t, nb::ndim<1>, nb::c_contig> edges,
std::optional<
nb::ndarray<const std::int32_t, nb::ndim<1>, nb::c_contig>>
edges,
bool redistribute, dolfinx::refinement::plaza::Option option)
{
std::optional<std::span<const std::int32_t>> cpp_edges(std::nullopt);
if (edges.has_value())
cpp_edges.emplace(
std::span(edges.value().data(), edges.value().size()));

auto [mesh1, cell, facet] = dolfinx::refinement::plaza::refine(
mesh0, std::span<const std::int32_t>(edges.data(), edges.size()),
redistribute, option);
mesh0, cpp_edges, redistribute, option);
return std::tuple{std::move(mesh1), as_nbarray(std::move(cell)),
as_nbarray(std::move(facet))};
},
nb::arg("mesh"), nb::arg("edges"), nb::arg("redistribute"),
nb::arg("mesh"), nb::arg("edges") = nb::none(), nb::arg("redistribute"),
nb::arg("option"));
}

Expand Down

0 comments on commit 84a6fa4

Please sign in to comment.