From: Johannes Heinz <43043310+jh66637@users.noreply.github.com>
Date: Mon, 22 Jan 2024 08:52:34 +0000 (+0100)
Subject: Add factory functions for FERemoteEvaluationCommunicator
X-Git-Tag: relicensing~81^2~1
X-Git-Url: https://gitweb.dealii.org/cgi-bin/gitweb.cgi?a=commitdiff_plain;h=29dd2073d2b09e1ace91232f1219148f35a558df;p=dealii.git
Add factory functions for FERemoteEvaluationCommunicator
---
diff --git a/doc/news/changes/minor/20231231Heinz b/doc/news/changes/minor/20231231Heinz
new file mode 100644
index 0000000000..306b95c60f
--- /dev/null
+++ b/doc/news/changes/minor/20231231Heinz
@@ -0,0 +1,4 @@
+New: Add utility functions that setup FERemoteEvaluationCommunicator for point-to-point
+interpolation and Nitsche-type mortaring.
+
+(Johannes Heinz, 2023/12/31)
diff --git a/include/deal.II/matrix_free/fe_remote_evaluation.h b/include/deal.II/matrix_free/fe_remote_evaluation.h
index 11bf9f874a..7bf81d25f6 100644
--- a/include/deal.II/matrix_free/fe_remote_evaluation.h
+++ b/include/deal.II/matrix_free/fe_remote_evaluation.h
@@ -19,12 +19,16 @@
#include
+#include
+#include
+
#include
#include
#include
#include
+#include
#include
DEAL_II_NAMESPACE_OPEN
@@ -446,6 +450,76 @@ private:
+/**
+ * The namespace collects convenience functions to setup
+ * FERemoteEvaluationCommunicator for typical use cases.
+ */
+namespace Utilities
+{
+ /**
+ * A factory function for the FERemoteEvaluationCommunicator in the case of
+ * point-to-point interpolation.
+ *
+ * @param[in] matrix_free MatrixFree object that is used to distribute
+ * quadrature points at non-matching faces.
+ * @param[in] non_matching_faces_marked_vertices A vector of boundary face IDs
+ * that relate to non-matching interfaces. Each boundary face ID is
+ * accompanied by a lambda that marks the vertices of cells which are
+ * considered during the search of remote points (quadrature points of faces
+ * with given boundary face ID).
+ * @param[in] quad_no Quadrature number in @p matrix_free.
+ * @param[in] dof_no DoFHandler number in @p matrix_free.
+ * @param[in] tolerance Tolerance to find remote points.
+ */
+ template >
+ FERemoteEvaluationCommunicator
+ setup_remote_communicator_faces_point_to_point_interpolation(
+ const MatrixFree &matrix_free,
+ const std::vector<
+ std::pair()>>>
+ &non_matching_faces_marked_vertices,
+ const unsigned int quad_no = 0,
+ const unsigned int dof_no = 0,
+ const double tolerance = 1e-9);
+
+
+
+ /**
+ * A factory function for the FERemoteEvaluationCommunicator in the case of
+ * Nitsche-type mortaring.
+ *
+ * @param[in] matrix_free MatrixFree object that is used to distribute
+ * @param[in] non_matching_faces_marked_vertices A vector of boundary face IDs
+ * that relate to non-matching interfaces. Each boundary face ID is
+ * accompanied by a lambda that marks the vertices of cells which are
+ * considered during the process of computing intersections.
+ * @param[in] n_q_pnts_1D The number of 1D quadrature points per intersection.
+ * @param[in] dof_no DoFHandler number in @p matrix_free.
+ * distributed to each intersection (given for one coordinate direction).
+ * @param[in] nm_mapping_info In case nm_mapping_info is not a `nullptr` it is
+ * setup for cell based loops, such that it can be used in combination with
+ * FERemoteEvaluation.
+ * @param[in] tolerance Tolerance to find intersections.
+ */
+ template >
+ FERemoteEvaluationCommunicator
+ setup_remote_communicator_faces_nitsche_type_mortaring(
+ const MatrixFree &matrix_free,
+ const std::vector<
+ std::pair()>>>
+ &non_matching_faces_marked_vertices,
+ const unsigned int n_q_pnts_1D,
+ const unsigned int dof_no = 0,
+ NonMatching::MappingInfo *nm_mapping_info = nullptr,
+ const double tolerance = 1e-9);
+} // namespace Utilities
+
+
+
/**
* Class to access data in a matrix-free loop for non-matching discretizations.
* Interfaces are named with FEEvaluation in mind.
@@ -999,6 +1073,373 @@ FERemoteEvaluationCommunicator::CopyInstructions::copy_data_entries(
+namespace Utilities
+{
+ template
+ FERemoteEvaluationCommunicator
+ setup_remote_communicator_faces_point_to_point_interpolation(
+ const MatrixFree &matrix_free,
+ const std::vector<
+ std::pair()>>>
+ &non_matching_faces_marked_vertices,
+ const unsigned int quad_no,
+ const unsigned int dof_no,
+ const double tolerance)
+ {
+ const auto &dof_handler = matrix_free.get_dof_handler(dof_no);
+ const auto &tria = dof_handler.get_triangulation();
+ const auto &mapping = *matrix_free.get_mapping_info().mapping;
+
+ // Communication objects know about the communication pattern. I.e.,
+ // they know about the cells and quadrature points that have to be
+ // evaluated at remote faces. This information is given via
+ // RemotePointEvaluation. Additionally, the communication objects
+ // have to be able to match the quadrature points of the remote
+ // points (that provide exterior information) to the quadrature points
+ // defined at the interior cell. In case of point-to-point interpolation
+ // a vector of pairs with face batch Ids and the number of faces in the
+ // batch is needed. @c FERemoteCommunicationObjectEntityBatches
+ // is a container to store this information.
+ //
+ // We need multiple communication objects (one for each non-matching face
+ // ID).
+ 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,
+ // 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());
+
+ // 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.
+ const auto face_batch_range =
+ std::make_pair(matrix_free.n_inner_face_batches(),
+ matrix_free.n_inner_face_batches() +
+ matrix_free.n_boundary_face_batches());
+
+ // Iterate over all non-matching face IDs.
+ for (const auto &[nm_face, marked_vertices] :
+ non_matching_faces_marked_vertices)
+ {
+ // Construct the communication object for every face ID:
+ // 1) RemotePointEvaluation with user specified function for marked
+ // vertices.
+ auto rpe = std::make_shared>(
+ tolerance, false, 0, marked_vertices);
+
+ // 2) Face batch IDs and number of faces in batch.
+ std::vector>
+ face_batch_id_n_faces;
+
+ // Points that are searched by rpe.
+ std::vector> points;
+
+ // Temporarily setup FEFaceEvaluation to access the quadrature points at
+ // the faces on the non-matching interface.
+ FEFaceEvaluation phi(matrix_free,
+ true,
+ dof_no,
+ quad_no);
+
+ // Iterate over the boundary faces.
+ for (unsigned int bface = 0;
+ bface < face_batch_range.second - face_batch_range.first;
+ ++bface)
+ {
+ const unsigned int face = face_batch_range.first + bface;
+
+ if (matrix_free.get_boundary_id(face) == nm_face)
+ {
+ phi.reinit(face);
+
+ // If @c face is on the current side of the non-matching
+ // interface. Add the face batch ID and the number of faces in
+ // the batch to the corresponding data structure.
+ const unsigned int n_faces =
+ matrix_free.n_active_entries_per_face_batch(face);
+ face_batch_id_n_faces.emplace_back(
+ std::make_pair(face, n_faces));
+
+ // Append the quadrature points to the points we need to search
+ // for.
+ for (unsigned int v = 0; v < n_faces; ++v)
+ {
+ for (unsigned int q : phi.quadrature_point_indices())
+ {
+ const auto point = phi.quadrature_point(q);
+ Point temp;
+ for (unsigned int i = 0; i < dim; ++i)
+ temp[i] = point[i][v];
+
+ points.push_back(temp);
+ }
+ }
+
+ // 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,
+ ExcMessage(
+ "Quadrature for given face already provided."));
+
+ global_quadrature_vector[bface] =
+ Quadrature(phi.n_q_points);
+ }
+ }
+
+ // Reinit RPE and ensure all points are found.
+ rpe->reinit(points, tria, mapping);
+ Assert(rpe->all_points_found(),
+ ExcMessage("Not all remote points found."));
+
+ // Add communication object to the list of objects.
+ FERemoteCommunicationObjectEntityBatches co;
+ co.batch_id_n_entities = face_batch_id_n_faces;
+ co.rpe = rpe;
+ comm_objects.push_back(co);
+ }
+
+ // Reinit the communicator `FERemoteEvaluationCommunicator`
+ // with the communication objects.
+ FERemoteEvaluationCommunicator remote_communicator;
+
+ remote_communicator.reinit_faces(comm_objects,
+ face_batch_range,
+ global_quadrature_vector);
+
+ return remote_communicator;
+ }
+
+
+
+ template
+ FERemoteEvaluationCommunicator
+ setup_remote_communicator_faces_nitsche_type_mortaring(
+ const MatrixFree &matrix_free,
+ const std::vector<
+ std::pair()>>>
+ &non_matching_faces_marked_vertices,
+ const unsigned int n_q_pnts_1D,
+ const unsigned int dof_no,
+ NonMatching::MappingInfo *nm_mapping_info,
+ const double tolerance)
+ {
+ const auto &dof_handler = matrix_free.get_dof_handler(dof_no);
+ const auto &tria = dof_handler.get_triangulation();
+ const auto &mapping = *matrix_free.get_mapping_info().mapping;
+
+ constexpr unsigned int n_lanes = VectorizedArray::size();
+
+ std::pair face_range =
+ std::make_pair(matrix_free.n_inner_face_batches(),
+ matrix_free.n_inner_face_batches() +
+ matrix_free.n_boundary_face_batches());
+
+ std::vector> global_quadrature_vector(
+ (matrix_free.n_inner_face_batches() +
+ matrix_free.n_boundary_face_batches()) *
+ n_lanes);
+
+ // In case of Nitsche-type mortaring a vector of face indices is
+ // needed as communication object.
+ // @c FERemoteCommunicationObjectFaces is a container to store this
+ // information.
+ //
+ // We need multiple communication objects (one for each non-matching face
+ // ID).
+ std::vector> comm_objects;
+
+ // Create bounding boxes and GridTools::Cache which is needed in
+ // the following loop.
+ std::vector> local_boxes;
+ for (const auto &cell : tria.active_cell_iterators())
+ if (cell->is_locally_owned())
+ local_boxes.emplace_back(mapping.get_bounding_box(cell));
+
+ // Create r-tree of bounding boxes
+ const auto local_tree = pack_rtree(local_boxes);
+
+ // Compress r-tree to a minimal set of bounding boxes
+ std::vector>> global_bboxes(1);
+ global_bboxes[0] = extract_rtree_level(local_tree, 0);
+
+ const GridTools::Cache cache(tria, mapping);
+
+ // Iterate over all sides of the non-matching interface.
+ for (const auto &[nm_face, marked_vertices] :
+ non_matching_faces_marked_vertices)
+ {
+ // 1) compute cell face pairs
+ std::vector<
+ std::pair::cell_iterator, unsigned int>>
+ cell_face_pairs;
+
+ std::vector indices;
+
+ for (unsigned int face = face_range.first; face < face_range.second;
+ ++face)
+ {
+ if (matrix_free.get_boundary_id(face) == nm_face)
+ {
+ for (unsigned int v = 0;
+ v < matrix_free.n_active_entries_per_face_batch(face);
+ ++v)
+ {
+ const auto &[c, f] = matrix_free.get_face_iterator(face, v);
+
+ cell_face_pairs.emplace_back(std::make_pair(c, f));
+ indices.push_back(face * n_lanes + v);
+ }
+ }
+ }
+
+ // 2) Create RPE.
+ // In the Nitsche-type case we do not collect points for the setup
+ // of RemotePointEvaluation. Instead we compute intersections between
+ // the faces and setup RemotePointEvaluation with the computed
+ // intersections.
+
+ // Build intersection requests. Intersection requests
+ // correspond to vertices at faces.
+ std::vector>> intersection_requests;
+ for (const auto &[cell, f] : cell_face_pairs)
+ {
+ std::vector> vertices(cell->face(f)->n_vertices());
+ std::copy_n(mapping.get_vertices(cell, f).begin(),
+ cell->face(f)->n_vertices(),
+ vertices.begin());
+ intersection_requests.emplace_back(vertices);
+ }
+
+ // Compute intersection data with user specified function for marked
+ // vertices.
+ auto intersection_data =
+ GridTools::internal::distributed_compute_intersection_locations<
+ dim - 1>(cache,
+ intersection_requests,
+ global_bboxes,
+ marked_vertices(),
+ tolerance);
+
+ // Convert to RPE.
+ std::vector> mapped_quadratures_recv_comp;
+
+ auto rpe =
+ std::make_shared>();
+ rpe->reinit(
+ intersection_data
+ .template convert_to_distributed_compute_point_locations_internal<
+ dim>(n_q_pnts_1D, tria, mapping, &mapped_quadratures_recv_comp),
+ tria,
+ mapping);
+
+ // 3) Fill global quadrature vector.
+ for (unsigned int i = 0; i < intersection_requests.size(); ++i)
+ {
+ const auto idx = indices[i];
+
+ // We do not use a structural binding here, since with
+ // C++17 caputuring structural bindings in lambdas leads
+ // to an ill formed program.
+ const auto &cell = std::get<0>(cell_face_pairs[i]);
+ const auto &f = std::get<1>(cell_face_pairs[i]);
+
+ std::vector> q_points;
+ std::vector weights;
+ for (unsigned int ptr = intersection_data.recv_ptrs[i];
+ ptr < intersection_data.recv_ptrs[i + 1];
+ ++ptr)
+ {
+ const auto &quad = mapped_quadratures_recv_comp[ptr];
+
+ const auto &ps = quad.get_points();
+ std::transform(
+ ps.begin(),
+ ps.end(),
+ std::back_inserter(q_points),
+ [&](const Point &p) {
+ return mapping.project_real_point_to_unit_point_on_face(
+ cell, f, p);
+ });
+
+ const auto &ws = quad.get_weights();
+ weights.insert(weights.end(), ws.begin(), ws.end());
+ }
+ Quadrature quad(q_points, weights);
+
+ Assert(global_quadrature_vector[idx].size() == 0,
+ ExcMessage("Quadrature for given face already provided."));
+
+ global_quadrature_vector[idx] = quad;
+ }
+
+ // Add communication object.
+ FERemoteCommunicationObject co;
+ co.indices = indices;
+ co.rpe = rpe;
+ comm_objects.push_back(co);
+ }
+
+ // Reinit the communicator with the communication objects.
+ FERemoteEvaluationCommunicator remote_communicator;
+
+ remote_communicator.reinit_faces(
+ comm_objects,
+ std::make_pair(0, global_quadrature_vector.size()),
+ global_quadrature_vector);
+
+ if (nm_mapping_info != nullptr)
+ {
+ std::vector<
+ std::pair::cell_iterator, unsigned int>>
+ vector_face_accessors;
+ vector_face_accessors.reserve((matrix_free.n_inner_face_batches() +
+ matrix_free.n_boundary_face_batches()) *
+ n_lanes);
+
+ // fill container for inner face batches
+ unsigned int face_batch = 0;
+ for (; face_batch < matrix_free.n_inner_face_batches(); ++face_batch)
+ {
+ for (unsigned int v = 0; v < n_lanes; ++v)
+ {
+ if (v < matrix_free.n_active_entries_per_face_batch(face_batch))
+ vector_face_accessors.push_back(
+ matrix_free.get_face_iterator(face_batch, v));
+ else
+ vector_face_accessors.push_back(
+ matrix_free.get_face_iterator(face_batch, 0));
+ }
+ }
+ // and boundary face batches
+ for (; face_batch < (matrix_free.n_inner_face_batches() +
+ matrix_free.n_boundary_face_batches());
+ ++face_batch)
+ {
+ for (unsigned int v = 0; v < n_lanes; ++v)
+ {
+ if (v < matrix_free.n_active_entries_per_face_batch(face_batch))
+ vector_face_accessors.push_back(
+ matrix_free.get_face_iterator(face_batch, v));
+ else
+ vector_face_accessors.push_back(
+ matrix_free.get_face_iterator(face_batch, 0));
+ }
+ }
+
+ nm_mapping_info->reinit_faces(vector_face_accessors,
+ global_quadrature_vector);
+ }
+
+ return remote_communicator;
+ }
+} // namespace Utilities
+
+
+
template
template
FERemoteEvaluation::FERemoteEvaluation(
diff --git a/tests/matrix_free/fe_remote_evaluation_03.cc b/tests/matrix_free/fe_remote_evaluation_03.cc
new file mode 100644
index 0000000000..9d12be17e6
--- /dev/null
+++ b/tests/matrix_free/fe_remote_evaluation_03.cc
@@ -0,0 +1,308 @@
+// ---------------------------------------------------------------------
+//
+// Copyright (C) 2023 by the deal.II authors
+//
+// This file is part of the deal.II library.
+//
+// The deal.II library is free software; you can use it, redistribute
+// it, and/or modify it under the terms of the GNU Lesser General
+// Public License as published by the Free Software Foundation; either
+// version 2.1 of the License, or (at your option) any later version.
+// The full text of the license can be found in the file LICENSE.md at
+// the top level directory of deal.II.
+//
+// ---------------------------------------------------------------------
+
+// This test checks the correct setup of FERemoteEvaluationCommunicator
+// via the factory functions
+// setup_remote_communicator_faces_point_to_point_interpolation()
+// and setup_remote_communicator_faces_nitsche_type_mortaring().
+
+#include
+
+#include
+#include
+
+#include
+
+#include
+
+#include
+
+#include "../tests.h"
+
+using namespace dealii;
+
+using Number = double;
+
+using VectorType = LinearAlgebra::distributed::Vector;
+
+constexpr unsigned int dim = 2;
+
+template
+struct Coordinates : public Function
+{
+ Coordinates()
+ : Function(dim, 0.0)
+ {}
+
+ double
+ value(const Point &p, const unsigned int component) const final
+ {
+ if (component == 0)
+ return p[component];
+ return p[component];
+ }
+};
+
+
+void
+test_point_to_point(const MatrixFree &matrix_free)
+{
+ const Triangulation &tria =
+ matrix_free.get_dof_handler().get_triangulation();
+
+ std::vector()>>>
+ non_matching_faces_marked_vertices;
+
+ non_matching_faces_marked_vertices.emplace_back(std::make_pair(1, [&]() {
+ std::vector mask(tria.n_vertices(), true);
+
+ for (const auto &cell : tria.active_cell_iterators())
+ for (auto const &f : cell->face_indices())
+ if (cell->face(f)->at_boundary() && cell->face(f)->boundary_id() == 1)
+ for (const auto v : cell->vertex_indices())
+ mask[cell->vertex_index(v)] = false;
+
+ return mask;
+ }));
+
+ auto remote_communicator =
+ Utilities::setup_remote_communicator_faces_point_to_point_interpolation(
+ matrix_free, non_matching_faces_marked_vertices);
+
+ // 1) Allocate vectors and interpolate grid coordinates in DoFs.
+ VectorType src;
+ VectorType dst;
+ matrix_free.initialize_dof_vector(src);
+ matrix_free.initialize_dof_vector(dst);
+
+ VectorTools::interpolate(*matrix_free.get_mapping_info().mapping,
+ matrix_free.get_dof_handler(),
+ Coordinates(),
+ src);
+
+ // 2) Construct remote evaluator...
+ FERemoteEvaluation> phi_remote_eval(
+ remote_communicator, matrix_free.get_dof_handler());
+
+ // ...and precompute remote values.
+ phi_remote_eval.gather_evaluate(src, EvaluationFlags::values);
+
+ // 3) Access values at quadrature points
+ const auto boundary_function =
+ [&](const auto &data, auto &dst, const auto &src, const auto face_range) {
+ FEFaceEvaluation phi(matrix_free, true);
+ auto phi_remote = phi_remote_eval.get_data_accessor();
+
+ for (unsigned int face = face_range.first; face < face_range.second;
+ ++face)
+ {
+ if (matrix_free.get_boundary_id(face) == 1)
+ {
+ phi.reinit(face);
+ phi_remote.reinit(face);
+
+ for (unsigned int q : phi.quadrature_point_indices())
+ {
+ const auto error =
+ phi.quadrature_point(q) - phi_remote.get_value(q);
+
+ // only consider active entries in global error
+ for (unsigned int d = 0; d < dim; ++d)
+ for (unsigned int v = 0;
+ v < matrix_free.n_active_entries_per_face_batch(face);
+ ++v)
+ AssertThrow(std::abs(error[d][v]) < 1e-13,
+ ExcMessage("Error too large."));
+ }
+ }
+ }
+ };
+
+ matrix_free.template loop(
+ {}, {}, boundary_function, dst, src, true);
+}
+
+void
+test_nitsche_type_mortaring(const MatrixFree &matrix_free)
+{
+ const Triangulation &tria =
+ matrix_free.get_dof_handler().get_triangulation();
+
+ std::vector()>>>
+ non_matching_faces_marked_vertices;
+
+ non_matching_faces_marked_vertices.emplace_back(std::make_pair(1, [&]() {
+ std::vector mask(tria.n_vertices(), true);
+
+ for (const auto &cell : tria.active_cell_iterators())
+ for (auto const &f : cell->face_indices())
+ if (cell->face(f)->at_boundary() && cell->face(f)->boundary_id() == 1)
+ for (const auto v : cell->vertex_indices())
+ mask[cell->vertex_index(v)] = false;
+
+ return mask;
+ }));
+
+ typename NonMatching::MappingInfo::AdditionalData
+ additional_data;
+ additional_data.use_global_weights = true;
+
+ NonMatching::MappingInfo nm_mapping_info(
+ *matrix_free.get_mapping_info().mapping,
+ update_quadrature_points,
+ additional_data);
+
+ auto remote_communicator =
+ Utilities::setup_remote_communicator_faces_nitsche_type_mortaring(
+ matrix_free, non_matching_faces_marked_vertices, 4, 0, &nm_mapping_info);
+
+ // 1) Allocate vectors and interpolate grid coordinates in DoFs.
+ VectorType src;
+ VectorType dst;
+ matrix_free.initialize_dof_vector(src);
+ matrix_free.initialize_dof_vector(dst);
+
+ VectorTools::interpolate(*matrix_free.get_mapping_info().mapping,
+ matrix_free.get_dof_handler(),
+ Coordinates(),
+ src);
+
+ // 2) Construct remote evaluator...
+ FERemoteEvaluation phi_remote_eval(
+ remote_communicator, matrix_free.get_dof_handler());
+
+ // ...and precompute remote values.
+ phi_remote_eval.gather_evaluate(src, EvaluationFlags::values);
+
+ // 3) Access values at quadrature points
+ const auto boundary_function =
+ [&](const auto &data, auto &dst, const auto &src, const auto face_range) {
+ FEFacePointEvaluation phi(
+ nm_mapping_info, matrix_free.get_dof_handler().get_fe());
+ auto phi_remote = phi_remote_eval.get_data_accessor();
+
+
+ for (unsigned int face = face_range.first; face < face_range.second;
+ ++face)
+ {
+ if (matrix_free.get_boundary_id(face) == 1)
+ {
+ for (unsigned int v = 0;
+ v < matrix_free.n_active_entries_per_face_batch(face);
+ ++v)
+ {
+ constexpr unsigned int n_lanes =
+ VectorizedArray::size();
+ phi.reinit(face * n_lanes + v);
+ phi_remote.reinit(face * n_lanes + v);
+
+ for (unsigned int q : phi.quadrature_point_indices())
+ {
+ const auto point = phi.real_point(q);
+ Tensor<1, dim, Number> temp;
+ for (unsigned int i = 0; i < dim; ++i)
+ temp[i] = point[i];
+
+ const auto error = temp - phi_remote.get_value(q);
+ for (unsigned int d = 0; d < dim; ++d)
+ AssertThrow(std::abs(error[d]) < 1e-13,
+ ExcMessage("Error too large."));
+ }
+ }
+ }
+ }
+ };
+
+ matrix_free.template loop(
+ {}, {}, boundary_function, dst, src, true);
+}
+
+int
+main(int argc, char **argv)
+{
+ Utilities::MPI::MPI_InitFinalize mpi(argc, argv, 1);
+ MPILogInitAll all;
+
+ // Construct non-matching triangulation.
+ const double length = 1.0;
+ Triangulation tria_left;
+ const unsigned int subdiv_left = 11;
+ GridGenerator::subdivided_hyper_rectangle(tria_left,
+ {subdiv_left, 2 * subdiv_left},
+ {0.0, 0.0},
+ {0.525 * length, length});
+
+ for (const auto &face : tria_left.active_face_iterators())
+ if (face->at_boundary())
+ {
+ face->set_boundary_id(0);
+ if (face->center()[0] > 0.525 * length - 1e-6)
+ face->set_boundary_id(1);
+ }
+
+ Triangulation tria_right;
+ const unsigned int subdiv_right = 4;
+ GridGenerator::subdivided_hyper_rectangle(tria_right,
+ {subdiv_right, 2 * subdiv_right},
+ {0.525 * length, 0.0},
+ {length, length});
+
+ for (const auto &face : tria_right.active_face_iterators())
+ if (face->at_boundary())
+ face->set_boundary_id(0);
+
+ // Merge triangulations with tolerance 0 to ensure no vertices
+ // are merged.
+ parallel::distributed::Triangulation tria(MPI_COMM_WORLD);
+ GridGenerator::merge_triangulations(tria_left,
+ tria_right,
+ tria,
+ /*tolerance*/ 0.,
+ /*copy_manifold_ids*/ false,
+ /*copy_boundary_ids*/ true);
+
+ // Create DoFHandler
+ DoFHandler dof_handler(tria);
+ dof_handler.distribute_dofs(FESystem(FE_DGQ(3), dim));
+
+ // Constraints are not considered in this test
+ AffineConstraints constraints;
+ constraints.close();
+
+ // Setup MatrixFree
+ typename MatrixFree::AdditionalData data;
+ data.mapping_update_flags = update_values;
+ // We are not considering inner faces in this test. Therefore,
+ // data.mapping_update_flags_inner_faces is not set.
+ data.mapping_update_flags_boundary_faces =
+ update_quadrature_points | update_values;
+
+
+ MatrixFree matrix_free;
+ matrix_free.reinit(
+ MappingQ1(), dof_handler, constraints, QGauss(4), data);
+
+ // 2) test point-to-point
+ test_point_to_point(matrix_free);
+
+ deallog << "OK" << std::endl;
+
+ // 3) test Nitsche-type
+ test_nitsche_type_mortaring(matrix_free);
+
+ deallog << "OK" << std::endl;
+
+ return 0;
+}
diff --git a/tests/matrix_free/fe_remote_evaluation_03.with_p4est=true.mpirun=1.with_cgal=true.output b/tests/matrix_free/fe_remote_evaluation_03.with_p4est=true.mpirun=1.with_cgal=true.output
new file mode 100644
index 0000000000..2d64160d8c
--- /dev/null
+++ b/tests/matrix_free/fe_remote_evaluation_03.with_p4est=true.mpirun=1.with_cgal=true.output
@@ -0,0 +1,3 @@
+
+DEAL:0::OK
+DEAL:0::OK
diff --git a/tests/matrix_free/fe_remote_evaluation_03.with_p4est=true.mpirun=2.with_cgal=true.output b/tests/matrix_free/fe_remote_evaluation_03.with_p4est=true.mpirun=2.with_cgal=true.output
new file mode 100644
index 0000000000..abeeb87a92
--- /dev/null
+++ b/tests/matrix_free/fe_remote_evaluation_03.with_p4est=true.mpirun=2.with_cgal=true.output
@@ -0,0 +1,7 @@
+
+DEAL:0::OK
+DEAL:0::OK
+
+DEAL:1::OK
+DEAL:1::OK
+