From: Marco Feder Date: Sun, 20 Aug 2023 17:12:47 +0000 (+0000) Subject: Expose rpe parameters through AdditionalData X-Git-Tag: relicensing~556^2~2 X-Git-Url: https://gitweb.dealii.org/cgi-bin/gitweb.cgi?a=commitdiff_plain;h=4cda1db952cea09203f10a122d7c360e5708d4b8;p=dealii.git Expose rpe parameters through AdditionalData --- diff --git a/include/deal.II/multigrid/mg_transfer_global_coarsening.h b/include/deal.II/multigrid/mg_transfer_global_coarsening.h index 11d7181354..301849a824 100644 --- a/include/deal.II/multigrid/mg_transfer_global_coarsening.h +++ b/include/deal.II/multigrid/mg_transfer_global_coarsening.h @@ -706,6 +706,30 @@ private: using VectorizedArrayType = VectorizedArray; public: + /** + * AdditionalData structure for construction arguments needed by + * RemotePointEvaluation. Default values are the same as the ones in + * RemotePointEvaluation. + */ + struct AdditionalData + { + AdditionalData(const double tol = 1e-6, + const bool enf_unique_mapping = false, + const unsigned int rtree_l = 0, + const std::function()> &marked_verts = {}) + : tolerance(tol) + , enforce_unique_mapping(enf_unique_mapping) + , rtree_level(rtree_l) + , marked_vertices(marked_verts) + {} + double tolerance; + bool enforce_unique_mapping; + unsigned int rtree_level; + std::function()> marked_vertices; + }; + + MGTwoLevelTransferNonNested(const AdditionalData &data = AdditionalData()); + /** * Set up transfer operator between the given DoFHandler objects ( * @p dof_handler_fine and @p dof_handler_coarse). @@ -788,7 +812,7 @@ private: * Object to evaluate shape functions on one mesh on visited support points of * the other mesh. */ - Utilities::MPI::RemotePointEvaluation rpe; + std::shared_ptr> rpe; /** * MappingInfo object needed as Mapping argument by FEPointEvaluation. 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 dfb99460fc..7dd71601ee 100644 --- a/include/deal.II/multigrid/mg_transfer_global_coarsening.templates.h +++ b/include/deal.II/multigrid/mg_transfer_global_coarsening.templates.h @@ -4100,6 +4100,18 @@ namespace internal } // namespace internal + +template +MGTwoLevelTransferNonNested>:: + MGTwoLevelTransferNonNested(const AdditionalData &data) +{ + rpe = std::make_shared>( + data.tolerance, + data.enforce_unique_mapping, + data.rtree_level, + data.marked_vertices); +} + template void MGTwoLevelTransferNonNested>:: @@ -4170,13 +4182,13 @@ MGTwoLevelTransferNonNested>:: this->partitioner_fine->global_to_local(global_dof_indices[i]); // hand points over to RPE - rpe.reinit(points, dof_handler_coarse.get_triangulation(), mapping_coarse); + rpe->reinit(points, dof_handler_coarse.get_triangulation(), mapping_coarse); // set up MappingInfo for easier data access - mapping_info = internal::fill_mapping_info(rpe); + mapping_info = internal::fill_mapping_info(*rpe); // set up constraints - const auto &cell_data = rpe.get_cell_data(); + const auto &cell_data = rpe->get_cell_data(); constraint_info.reinit(dof_handler_coarse, cell_data.cells.size(), @@ -4185,7 +4197,7 @@ MGTwoLevelTransferNonNested>:: for (unsigned int i = 0; i < cell_data.cells.size(); ++i) { typename DoFHandler::active_cell_iterator cell( - &rpe.get_triangulation(), + &rpe->get_triangulation(), cell_data.cells[i].first, cell_data.cells[i].second, &dof_handler_coarse); @@ -4295,18 +4307,18 @@ MGTwoLevelTransferNonNested>:: } }; - rpe.template evaluate_and_process(evaluation_point_results, - buffer, - evaluation_function); + rpe->template evaluate_and_process(evaluation_point_results, + buffer, + evaluation_function); // Weight operator in case some points are owned by multiple cells. - if (rpe.is_map_unique() == false) + if (rpe->is_map_unique() == false) { const auto evaluation_point_results_temp = evaluation_point_results; evaluation_point_results.clear(); - evaluation_point_results.reserve(rpe.get_point_ptrs().size() - 1); + evaluation_point_results.reserve(rpe->get_point_ptrs().size() - 1); - const auto &ptr = rpe.get_point_ptrs(); + const auto &ptr = rpe->get_point_ptrs(); for (unsigned int i = 0; i < ptr.size() - 1; ++i) { @@ -4388,7 +4400,7 @@ MGTwoLevelTransferNonNested>:: std::vector evaluation_point_results; std::vector buffer; - evaluation_point_results.resize(rpe.get_point_ptrs().size() - 1); + evaluation_point_results.resize(rpe->get_point_ptrs().size() - 1); for (unsigned int j = 0; j < evaluation_point_results.size(); ++j) { @@ -4423,9 +4435,9 @@ MGTwoLevelTransferNonNested>:: } // Weight operator in case some points are owned by multiple cells. - if (rpe.is_map_unique() == false) + if (rpe->is_map_unique() == false) { - const auto &ptr = rpe.get_point_ptrs(); + const auto &ptr = rpe->get_point_ptrs(); for (unsigned int i = 0; i < ptr.size() - 1; ++i) { @@ -4470,9 +4482,9 @@ MGTwoLevelTransferNonNested>:: } }; - rpe.template process_and_evaluate(evaluation_point_results, - buffer, - evaluation_function); + rpe->template process_and_evaluate(evaluation_point_results, + buffer, + evaluation_function); }