From c1f6d130462c9b44b930a36480d6cce163122910 Mon Sep 17 00:00:00 2001 From: Johannes Heinz <43043310+jh66637@users.noreply.github.com> Date: Thu, 16 May 2024 20:00:31 +0200 Subject: [PATCH] remove Quadrature from FERemoteEvaluationCommunicator::reinit_faces() --- .../matrix_free/fe_remote_evaluation.h | 58 ++++++++++--------- tests/matrix_free/fe_remote_evaluation_01.cc | 23 ++++---- tests/matrix_free/fe_remote_evaluation_02.cc | 19 +++--- 3 files changed, 50 insertions(+), 50 deletions(-) diff --git a/include/deal.II/matrix_free/fe_remote_evaluation.h b/include/deal.II/matrix_free/fe_remote_evaluation.h index 930dd9ef17..14de198ad4 100644 --- a/include/deal.II/matrix_free/fe_remote_evaluation.h +++ b/include/deal.II/matrix_free/fe_remote_evaluation.h @@ -294,7 +294,7 @@ public: reinit_faces(const std::vector> &comm_objects, const std::pair &face_batch_range, - const std::vector> &quadrature_vector); + const std::vector &quadrature_sizes); /** * This function stores given communication objects @@ -305,7 +305,7 @@ public: reinit_faces( const std::vector> &comm_objects, const std::pair &face_range, - const std::vector> &quadrature_vector); + const std::vector &quadrature_sizes); /** * This function stores given communication objects @@ -317,8 +317,8 @@ public: void reinit_faces( const std::vector> &comm_objects, - const IteratorRange &cell_iterator_range, - const std::vector>> &quadrature_vector); + const IteratorRange &cell_iterator_range, + const std::vector> &quadrature_sizes); @@ -760,7 +760,7 @@ FERemoteEvaluationCommunicator::reinit_faces( const std::vector> &comm_objects, const std::pair &face_batch_range, - const std::vector> &quadrature_vector) + const std::vector &quadrature_sizes) { // erase type by converting to the base object communication_objects.clear(); @@ -768,7 +768,7 @@ FERemoteEvaluationCommunicator::reinit_faces( communication_objects.push_back(co); // fetch points and update communication patterns - const unsigned int n_cells = quadrature_vector.size(); + const unsigned int n_cells = quadrature_sizes.size(); AssertDimension(n_cells, face_batch_range.second - face_batch_range.first); // construct view: @@ -779,7 +779,7 @@ FERemoteEvaluationCommunicator::reinit_faces( view.ptrs[0] = 0; for (unsigned int face = 0; face < n_cells; ++face) { - view.ptrs[face + 1] = view.ptrs[face] + quadrature_vector[face].size(); + view.ptrs[face + 1] = view.ptrs[face] + quadrature_sizes[face]; } } @@ -788,14 +788,14 @@ void FERemoteEvaluationCommunicator::reinit_faces( const std::vector> &comm_objects, const std::pair &face_range, - const std::vector> &quadrature_vector) + const std::vector &quadrature_sizes) { // erase type communication_objects.clear(); for (const auto &co : comm_objects) communication_objects.push_back(co); - const unsigned int n_faces = quadrature_vector.size(); + const unsigned int n_faces = quadrature_sizes.size(); AssertDimension(n_faces, face_range.second - face_range.first); // construct view: @@ -805,7 +805,7 @@ FERemoteEvaluationCommunicator::reinit_faces( view.ptrs[0] = 0; for (unsigned int face = 0; face < n_faces; ++face) - view.ptrs[face + 1] = view.ptrs[face] + quadrature_vector[face].size(); + view.ptrs[face + 1] = view.ptrs[face] + quadrature_sizes[face]; } template @@ -813,15 +813,15 @@ template void FERemoteEvaluationCommunicator::reinit_faces( const std::vector> &comm_objects, - const IteratorRange &cell_iterator_range, - const std::vector>> &quadrature_vector) + const IteratorRange &cell_iterator_range, + const std::vector> &quadrature_sizes) { // erase type communication_objects.clear(); for (const auto &co : comm_objects) communication_objects.push_back(co); - const unsigned int n_cells = quadrature_vector.size(); + const unsigned int n_cells = quadrature_sizes.size(); AssertDimension(n_cells, std::distance(cell_iterator_range.begin(), cell_iterator_range.end())); @@ -850,7 +850,7 @@ FERemoteEvaluationCommunicator::reinit_faces( face_ptrs[face_index + 1] = face_ptrs[face_index] + - quadrature_vector[cell->active_cell_index()][f].size(); + quadrature_sizes[cell->active_cell_index()][f]; } } } @@ -1108,12 +1108,12 @@ namespace Utilities std::vector> comm_objects; // Additionally to the communication objects we need a vector - // that stores quadrature rules for every face batch. - // The quadrature can be empty in case of non non-matching faces, + // that stores quadrature rule sizes for every face batch. + // The quadrature can have size zero in case of non non-matching faces, // i.e. boundary faces. Internally this information is needed to correctly // access values over multiple communication objects. - std::vector> global_quadrature_vector( - matrix_free.n_boundary_face_batches()); + std::vector global_quadrature_sizes( + matrix_free.n_boundary_face_batches(), numbers::invalid_unsigned_int); // Get the range of face batches we have to look at during construction of // the communication objects. We only have to look at boundary faces. @@ -1180,15 +1180,14 @@ namespace Utilities } } - // Insert a quadrature rule of correct size into the global - // quadrature vector. First check that each face is only - // considered once. - Assert(global_quadrature_vector[bface].size() == 0, + // Insert the quadrature size into the global vector. + // First check that each face is only considered once. + Assert(global_quadrature_sizes[bface] == + numbers::invalid_unsigned_int, ExcMessage( "Quadrature for given face already provided.")); - global_quadrature_vector[bface] = - Quadrature(phi.n_q_points); + global_quadrature_sizes[bface] = phi.n_q_points; } } @@ -1210,7 +1209,7 @@ namespace Utilities remote_communicator.reinit_faces(comm_objects, face_batch_range, - global_quadrature_vector); + global_quadrature_sizes); return remote_communicator; } @@ -1388,10 +1387,17 @@ namespace Utilities // Reinit the communicator with the communication objects. FERemoteEvaluationCommunicator remote_communicator; + std::vector global_quadrature_sizes( + global_quadrature_vector.size()); + std::transform(global_quadrature_vector.cbegin(), + global_quadrature_vector.cend(), + global_quadrature_sizes.begin(), + [](const auto &q) { return q.size(); }); + remote_communicator.reinit_faces( comm_objects, std::make_pair(0, global_quadrature_vector.size()), - global_quadrature_vector); + global_quadrature_sizes); if (nm_mapping_info != nullptr) { diff --git a/tests/matrix_free/fe_remote_evaluation_01.cc b/tests/matrix_free/fe_remote_evaluation_01.cc index 039d689241..d22e3b8406 100644 --- a/tests/matrix_free/fe_remote_evaluation_01.cc +++ b/tests/matrix_free/fe_remote_evaluation_01.cc @@ -73,7 +73,7 @@ construct_comm_for_face_batches(const MatrixFree &matrix_free) // Fill a quadrature vector to keep track of the quadrature sizes on each // face - std::vector> quadrature_vector( + std::vector quadrature_sizes( matrix_free.n_boundary_face_batches()); // Points that are searched by rpe. @@ -108,8 +108,8 @@ construct_comm_for_face_batches(const MatrixFree &matrix_free) } } - // append quadrature of correct size - quadrature_vector[bface] = Quadrature(phi.n_q_points); + // append quadrature size + quadrature_sizes[bface] = phi.n_q_points; } // use rpe to search for stored points @@ -127,7 +127,7 @@ construct_comm_for_face_batches(const MatrixFree &matrix_free) // Renit the communicator `FERemoteEvaluationCommunicator` // with the communication objects. FERemoteEvaluationCommunicator remote_communicator; - remote_communicator.reinit_faces({co}, face_batch_range, quadrature_vector); + remote_communicator.reinit_faces({co}, face_batch_range, quadrature_sizes); return remote_communicator; } @@ -212,14 +212,12 @@ construct_comm_for_cell_face_nos(const MatrixFree &matrix_free) matrix_free.n_inner_face_batches() + matrix_free.n_boundary_face_batches()); - // Fill a quadrature vector to keep track of the quadrature sizes on each - // cell, face pair - std::vector>> quadrature_vector; + // Fill a vector of quadrature sizes for each cell, face pair + std::vector> quadrature_sizes; for (const auto &cell : matrix_free.get_dof_handler() .get_triangulation() .active_cell_iterators()) - quadrature_vector.emplace_back( - std::vector>(cell->n_faces())); + quadrature_sizes.emplace_back(std::vector(cell->n_faces())); // Points that are searched by rpe. std::vector> points; @@ -248,9 +246,8 @@ construct_comm_for_cell_face_nos(const MatrixFree &matrix_free) points.push_back(temp); } - // append quadrature of correct size - quadrature_vector[cell->active_cell_index()][f] = - Quadrature(phi.n_q_points); + // append correct quadrature size + quadrature_sizes[cell->active_cell_index()][f] = phi.n_q_points; } } @@ -272,7 +269,7 @@ construct_comm_for_cell_face_nos(const MatrixFree &matrix_free) remote_communicator.reinit_faces( {co}, matrix_free.get_dof_handler().get_triangulation().active_cell_iterators(), - quadrature_vector); + quadrature_sizes); return remote_communicator; } diff --git a/tests/matrix_free/fe_remote_evaluation_02.cc b/tests/matrix_free/fe_remote_evaluation_02.cc index cfb50937be..e5a4b86b51 100644 --- a/tests/matrix_free/fe_remote_evaluation_02.cc +++ b/tests/matrix_free/fe_remote_evaluation_02.cc @@ -71,9 +71,8 @@ construct_comm_for_face_batches(const MatrixFree &matrix_free) matrix_free.n_inner_face_batches() + matrix_free.n_boundary_face_batches()); - // Fill a quadrature vector to keep track of the quadrature sizes on each - // face - std::vector> quadrature_vector( + // Fill a vector to keep track of the quadrature sizes on each face + std::vector quadrature_sizes( matrix_free.n_boundary_face_batches()); // Points that are searched by rpe. @@ -109,7 +108,7 @@ construct_comm_for_face_batches(const MatrixFree &matrix_free) } // append quadrature of correct size - quadrature_vector[bface] = Quadrature(phi.n_q_points); + quadrature_sizes[bface] = phi.n_q_points; } // use rpe to search for stored points @@ -127,7 +126,7 @@ construct_comm_for_face_batches(const MatrixFree &matrix_free) // Renit the communicator `FERemoteEvaluationCommunicator` // with the communication objects. FERemoteEvaluationCommunicator remote_communicator; - remote_communicator.reinit_faces({co}, face_batch_range, quadrature_vector); + remote_communicator.reinit_faces({co}, face_batch_range, quadrature_sizes); return remote_communicator; } @@ -216,12 +215,11 @@ construct_comm_for_cell_face_nos(const MatrixFree &matrix_free) // Fill a quadrature vector to keep track of the quadrature sizes on each // cell, face pair - std::vector>> quadrature_vector; + std::vector> quadrature_sizes; for (const auto &cell : matrix_free.get_dof_handler() .get_triangulation() .active_cell_iterators()) - quadrature_vector.emplace_back( - std::vector>(cell->n_faces())); + quadrature_sizes.emplace_back(std::vector(cell->n_faces())); // Points that are searched by rpe. std::vector> points; @@ -251,8 +249,7 @@ construct_comm_for_cell_face_nos(const MatrixFree &matrix_free) } // append quadrature of correct size - quadrature_vector[cell->active_cell_index()][f] = - Quadrature(phi.n_q_points); + quadrature_sizes[cell->active_cell_index()][f] = phi.n_q_points; } } @@ -274,7 +271,7 @@ construct_comm_for_cell_face_nos(const MatrixFree &matrix_free) remote_communicator.reinit_faces( {co}, matrix_free.get_dof_handler().get_triangulation().active_cell_iterators(), - quadrature_vector); + quadrature_sizes); return remote_communicator; } -- 2.39.5