From d4ffb06636635d926196a4cfdcdc5ff3cc371d81 Mon Sep 17 00:00:00 2001 From: Magdalena Schreter Date: Wed, 6 Dec 2023 16:20:19 +0100 Subject: [PATCH] add AdditionalData struct to RemotePointEvaluation --- doc/news/changes/minor/20231206Schreter | 6 ++ .../base/mpi_remote_point_evaluation.h | 81 ++++++++++++++----- .../mg_transfer_global_coarsening.templates.h | 6 +- source/base/mpi_remote_point_evaluation.cc | 37 +++++++-- tests/remote_point_evaluation/mapping_01.cc | 6 +- tests/remote_point_evaluation/mapping_03.cc | 7 +- tests/remote_point_evaluation/mapping_04.cc | 6 +- .../search_adjacent_cells.cc | 14 +++- 8 files changed, 129 insertions(+), 34 deletions(-) create mode 100644 doc/news/changes/minor/20231206Schreter diff --git a/doc/news/changes/minor/20231206Schreter b/doc/news/changes/minor/20231206Schreter new file mode 100644 index 0000000000..f3c7df10c5 --- /dev/null +++ b/doc/news/changes/minor/20231206Schreter @@ -0,0 +1,6 @@ +New: Add AdditionalData struct for configuring the RemotePointEvaluation class. +The latter is passed into a new constructor. +Deprecated: The old constructor of RemotePointEvaluation class taking the +parameters value-by-value is marked as deprecated. +
+(Magdalena Schreter-Fleischhacker, Peter Munch, 2023/12/06) diff --git a/include/deal.II/base/mpi_remote_point_evaluation.h b/include/deal.II/base/mpi_remote_point_evaluation.h index 0c441c2854..8f9934e2d3 100644 --- a/include/deal.II/base/mpi_remote_point_evaluation.h +++ b/include/deal.II/base/mpi_remote_point_evaluation.h @@ -53,9 +53,62 @@ namespace Utilities class RemotePointEvaluation { public: + /** + * AdditionalData structure that can be used to tweak parameters + * of RemotePointEvaluation. + */ + struct AdditionalData + { + public: + /** + * Constructor. + */ + AdditionalData( + const double tolerance = 1e-6, + const bool enforce_unique_mapping = false, + const unsigned int rtree_level = 0, + const std::function()> &marked_vertices = {}); + + /** + * Tolerance in terms of unit cell coordinates for determining all cells + * around a point passed to RemotePointEvaluation during reinit(). + * Depending on the problem, it might be necessary to adjust the + * tolerance in order to be able to identify a cell. Floating point + * arithmetic implies that a point will, in general, not lie exactly on + * a vertex, edge, or face. + */ + double tolerance; + + /** + * Enforce unique mapping, i.e., (one-to-one) relation of points and + * cells. + */ + bool enforce_unique_mapping; + + /** + * RTree level to be used during the construction of the bounding boxes. + */ + unsigned int rtree_level; + + /** + * Function that marks relevant vertices to make search of active cells + * around point more efficient. + */ + std::function()> marked_vertices; + }; + /** * Constructor. * + * @param additional_data Configure options for RemotePointEvaluation. + */ + RemotePointEvaluation( + const AdditionalData &additional_data = AdditionalData()); + + /** + * Constructor. This constructor is deprecated. Use the other constructor + * taking AdditionalData instead. + * * @param tolerance Tolerance in terms of unit cell coordinates for * determining all cells around a point passed to the class during * reinit(). Depending on the problem, it might be necessary to adjust @@ -67,9 +120,13 @@ namespace Utilities * @param rtree_level RTree level to be used during the construction of the bounding boxes. * @param marked_vertices Function that marks relevant vertices to make search * of active cells around point more efficient. + * + * @deprecated */ + DEAL_II_DEPRECATED_EARLY_WITH_COMMENT( + "Use the constructor with AdditionalData struct.") RemotePointEvaluation( - const double tolerance = 1e-6, + const double tolerance, const bool enforce_unique_mapping = false, const unsigned int rtree_level = 0, const std::function()> &marked_vertices = {}); @@ -312,27 +369,9 @@ namespace Utilities private: /** - * Tolerance to be used while determining the surrounding cells of a - * point. - */ - const double tolerance; - - /** - * Enforce unique mapping, i.e., (one-to-one) relation of points and - * cells. - */ - const bool enforce_unique_mapping; - - /** - * RTree level to be used during the construction of the bounding boxes. - */ - const unsigned int rtree_level; - - /** - * Function that marks relevant vertices to make search of active cells - * around point more efficient. + * Additional data with basic settings. */ - const std::function()> marked_vertices; + const AdditionalData additional_data; /** * Storage for the status of the triangulation signal. diff --git a/include/deal.II/multigrid/mg_transfer_global_coarsening.templates.h b/include/deal.II/multigrid/mg_transfer_global_coarsening.templates.h index 65680d32e8..323e9c27cb 100644 --- a/include/deal.II/multigrid/mg_transfer_global_coarsening.templates.h +++ b/include/deal.II/multigrid/mg_transfer_global_coarsening.templates.h @@ -4991,7 +4991,11 @@ template MGTwoLevelTransferNonNested>:: MGTwoLevelTransferNonNested(const AdditionalData &data) : additional_data(data) - , rpe(data.tolerance, false, data.rtree_level, {}) + , rpe(typename Utilities::MPI::RemotePointEvaluation::AdditionalData( + data.tolerance, + false, + data.rtree_level, + {})) {} template diff --git a/source/base/mpi_remote_point_evaluation.cc b/source/base/mpi_remote_point_evaluation.cc index e0c52cceff..27232e1bdf 100644 --- a/source/base/mpi_remote_point_evaluation.cc +++ b/source/base/mpi_remote_point_evaluation.cc @@ -36,7 +36,7 @@ namespace Utilities namespace MPI { template - RemotePointEvaluation::RemotePointEvaluation( + RemotePointEvaluation::AdditionalData::AdditionalData( const double tolerance, const bool enforce_unique_mapping, const unsigned int rtree_level, @@ -45,6 +45,29 @@ namespace Utilities , enforce_unique_mapping(enforce_unique_mapping) , rtree_level(rtree_level) , marked_vertices(marked_vertices) + {} + + + + template + RemotePointEvaluation::RemotePointEvaluation( + const AdditionalData &additional_data) + : additional_data(additional_data) + , ready_flag(false) + {} + + + + template + RemotePointEvaluation::RemotePointEvaluation( + const double tolerance, + const bool enforce_unique_mapping, + const unsigned int rtree_level, + const std::function()> &marked_vertices) + : additional_data(tolerance, + enforce_unique_mapping, + rtree_level, + marked_vertices) , ready_flag(false) {} @@ -88,7 +111,8 @@ namespace Utilities // compress r-tree to a minimal set of bounding boxes std::vector>> global_bboxes(1); - global_bboxes[0] = extract_rtree_level(local_tree, rtree_level); + global_bboxes[0] = + extract_rtree_level(local_tree, additional_data.rtree_level); const GridTools::Cache cache(tria, mapping); @@ -97,10 +121,11 @@ namespace Utilities cache, points, global_bboxes, - marked_vertices ? marked_vertices() : std::vector(), - tolerance, + additional_data.marked_vertices ? additional_data.marked_vertices() : + std::vector(), + additional_data.tolerance, true, - enforce_unique_mapping); + additional_data.enforce_unique_mapping); this->reinit(data, tria, mapping); #endif @@ -185,7 +210,7 @@ namespace Utilities all_points_found_flag = std::get<0>(n_owning_processes_global) > 0; } - Assert(enforce_unique_mapping == false || unique_mapping, + Assert(additional_data.enforce_unique_mapping == false || unique_mapping, ExcInternalError()); cell_data = std::make_unique(tria); diff --git a/tests/remote_point_evaluation/mapping_01.cc b/tests/remote_point_evaluation/mapping_01.cc index 9a65055679..454c7fa4ff 100644 --- a/tests/remote_point_evaluation/mapping_01.cc +++ b/tests/remote_point_evaluation/mapping_01.cc @@ -58,7 +58,11 @@ test(const bool enforce_unique_map) for (unsigned int i = 0; i <= 4; ++i) evaluation_points.emplace_back(i * 0.25, j * 0.25); - Utilities::MPI::RemotePointEvaluation eval(1e-6, enforce_unique_map); + typename Utilities::MPI::RemotePointEvaluation::AdditionalData + additional_data; + additional_data.enforce_unique_mapping = enforce_unique_map; + additional_data.tolerance = 1e-6; + Utilities::MPI::RemotePointEvaluation eval(additional_data); const auto result_avg = VectorTools::point_values<1>(mapping, diff --git a/tests/remote_point_evaluation/mapping_03.cc b/tests/remote_point_evaluation/mapping_03.cc index 589530a332..50a553d614 100644 --- a/tests/remote_point_evaluation/mapping_03.cc +++ b/tests/remote_point_evaluation/mapping_03.cc @@ -55,7 +55,12 @@ test(const bool enforce_unique_map) // should be assigned to rank 0 for unique mapping evaluation_points.emplace_back(0.5, 0.5 + std::pow(10.0, -i)); - Utilities::MPI::RemotePointEvaluation eval(1e-6, enforce_unique_map); + typename Utilities::MPI::RemotePointEvaluation::AdditionalData + additional_data; + additional_data.enforce_unique_mapping = enforce_unique_map; + additional_data.tolerance = 1e-6; + + Utilities::MPI::RemotePointEvaluation eval(additional_data); const auto result_avg = VectorTools::point_values<1>(mapping, diff --git a/tests/remote_point_evaluation/mapping_04.cc b/tests/remote_point_evaluation/mapping_04.cc index c791aa78dc..df321285a4 100644 --- a/tests/remote_point_evaluation/mapping_04.cc +++ b/tests/remote_point_evaluation/mapping_04.cc @@ -42,7 +42,11 @@ do_test(const bool enforce_unique_map, MappingQ1 mapping; - Utilities::MPI::RemotePointEvaluation eval(1.e-6, enforce_unique_map); + typename Utilities::MPI::RemotePointEvaluation::AdditionalData + additional_data; + additional_data.enforce_unique_mapping = enforce_unique_map; + + Utilities::MPI::RemotePointEvaluation eval(additional_data); eval.reinit(evaluation_points, tria, mapping); diff --git a/tests/remote_point_evaluation/search_adjacent_cells.cc b/tests/remote_point_evaluation/search_adjacent_cells.cc index 8e95d38684..458fe1b7ba 100644 --- a/tests/remote_point_evaluation/search_adjacent_cells.cc +++ b/tests/remote_point_evaluation/search_adjacent_cells.cc @@ -158,7 +158,10 @@ test(const unsigned int mapping_degree, } // initialize RPE without any marked points - dealii::Utilities::MPI::RemotePointEvaluation rpe(tolerance); + typename Utilities::MPI::RemotePointEvaluation::AdditionalData + additional_data; + additional_data.tolerance = tolerance; + dealii::Utilities::MPI::RemotePointEvaluation rpe(additional_data); rpe.reinit(points, tria, mapping); unsigned int n_points_not_found_rpe = 0; @@ -182,8 +185,13 @@ test(const unsigned int mapping_degree, // initialize RPE with all points marked std::vector marked_vertices(tria.n_vertices(), true); - dealii::Utilities::MPI::RemotePointEvaluation rpe2( - tolerance, false, 0, [marked_vertices]() { return marked_vertices; }); + + typename Utilities::MPI::RemotePointEvaluation::AdditionalData + additional_data2(tolerance, false, 0, [marked_vertices]() { + return marked_vertices; + }); + + dealii::Utilities::MPI::RemotePointEvaluation rpe2(additional_data2); rpe2.reinit(points, tria, mapping); -- 2.39.5