16#ifndef dealii_matrix_free_fe_remote_evaluation_h
17#define dealii_matrix_free_fe_remote_evaluation_h
41 template <
int dim,
int n_components,
typename value_type_>
44 using value_type =
typename internal::FEPointEvaluation::
45 EvaluatorTypeTraits<dim, dim, n_components, value_type_>::value_type;
48 EvaluatorTypeTraits<dim, dim, n_components, value_type_>
::
83 const unsigned int face_number)
const;
104 std::vector<unsigned int>
ptrs;
112 template <
int dim,
int n_components,
typename value_type_>
205 std::shared_ptr<Utilities::MPI::RemotePointEvaluation<dim>>
rpe;
218 std::vector<std::pair<unsigned int, unsigned int>>
234 std::shared_ptr<Utilities::MPI::RemotePointEvaluation<dim>>
rpe;
246 std::vector<unsigned int>
262 std::shared_ptr<Utilities::MPI::RemotePointEvaluation<dim>>
rpe;
269 std::pair<typename Triangulation<dim>::cell_iterator,
unsigned int>>
278 std::pair<typename Triangulation<dim>::cell_iterator,
unsigned int>>
308 const std::pair<unsigned int, unsigned int> &
face_range,
317 template <
typename Iterator>
329 template <
int n_components,
337 const VectorType &src,
339 const unsigned int first_selected_component,
358 std::vector<std::variant<FERemoteCommunicationObjectEntityBatches<dim>,
374 template <
typename T1,
typename T2>
378 const std::vector<T2> &src,
379 const std::vector<unsigned int> &indices);
385 template <
typename T1,
typename T2>
390 const std::vector<T2> &src,
392 unsigned int>> &cell_face_nos);
398 template <
typename T1,
typename T2>
402 const std::vector<T2> &src,
403 const std::vector<std::pair<unsigned int, unsigned int>>
404 &batch_id_n_entities);
410 template <
typename T1, std::
size_t n_lanes>
413 const unsigned int v,
419 template <
typename T1,
int rank_, std::
size_t n_lanes,
int dim_>
422 const unsigned int v,
428 template <
typename T1,
438 const unsigned int v,
445 template <
typename T1,
typename T2>
485 const unsigned int quad_no = 0,
486 const unsigned int dof_no = 0,
487 const double tolerance = 1e-9);
519 const unsigned int dof_no = 0,
521 const double tolerance = 1e-9);
537template <
int dim,
int n_components,
typename value_type>
554 template <
typename MeshType>
569 template <
typename VectorType>
654 const unsigned int face_number)
const
668 return ptrs[face_index];
678 template <
int dim,
int n_components,
typename value_type_>
685 , data_offset(
numbers::invalid_unsigned_int)
688 template <
int dim,
int n_components,
typename value_type_>
693 const unsigned int q)
const
698 return data.values[data_offset +
q];
701 template <
int dim,
int n_components,
typename value_type_>
710 return data.gradients[data_offset +
q];
713 template <
int dim,
int n_components,
typename value_type_>
716 const unsigned int index)
718 data_offset = view.get_shift(
index);
721 template <
int dim,
int n_components,
typename value_type_>
735std::vector<std::pair<unsigned int, unsigned int>>
739 return batch_id_n_entities;
743std::vector<unsigned int>
750std::vector<std::pair<typename Triangulation<dim>::cell_iterator,
unsigned int>>
753 return cell_face_nos;
767 communication_objects.
clear();
769 communication_objects.push_back(
co);
778 view.ptrs.resize(n_cells + 1);
781 for (
unsigned int face = 0; face < n_cells; ++face)
791 const std::pair<unsigned int, unsigned int> &
face_range,
795 communication_objects.
clear();
797 communication_objects.push_back(
co);
805 view.ptrs.resize(n_faces + 1);
808 for (
unsigned int face = 0; face < n_faces; ++face)
813template <
typename Iterator>
821 communication_objects.
clear();
823 communication_objects.push_back(
co);
836 unsigned int n_faces = 0;
839 cell_ptrs[cell->active_cell_index()] = n_faces;
840 n_faces += cell->n_faces();
847 for (
const auto &f : cell->face_indices())
849 const unsigned int face_index =
850 cell_ptrs[cell->active_cell_index()] + f;
860template <
int n_components,
868 const VectorType &src,
870 const unsigned int first_selected_component,
873 const bool has_ghost_elements = src.has_ghost_elements();
875 if (has_ghost_elements ==
false)
876 src.update_ghost_values();
884 [&](
const auto &
obj) {
885 CopyInstructions::copy_data(
888 VectorTools::point_values<n_components>(
890 obj.get_communication_object_pntrs());
898 [&](
const auto &
obj) {
899 CopyInstructions::copy_data(
902 VectorTools::point_gradients<n_components>(
904 obj.get_communication_object_pntrs());
912 if (has_ghost_elements ==
false)
913 src.zero_out_ghost_values();
924template <
typename T1,
typename T2>
929 const std::vector<T2> &src,
930 const std::vector<unsigned int> &indices)
935 for (
const auto idx : indices)
948template <
typename T1,
typename T2>
953 const std::vector<T2> &src,
955 unsigned int>> &cell_face_nos)
960 for (
const auto &[cell, f] : cell_face_nos)
962 for (
unsigned int j =
view.
get_shift(cell->active_cell_index(), f);
975template <
typename T1,
typename T2>
980 const std::vector<T2> &src,
981 const std::vector<std::pair<unsigned int, unsigned int>> &batch_id_n_entities)
988 for (
unsigned int v = 0; v <
n_entries; ++v)
996 copy_data_entries(dst[
j], v, src[c]);
1002template <
typename T1, std::
size_t n_lanes>
1006 const unsigned int v,
1015template <
typename T1,
int rank_, std::
size_t n_lanes,
int dim_>
1019 const unsigned int v,
1024 if constexpr (
rank_ == 1)
1026 for (
unsigned int i = 0; i < dim_; ++i)
1031 for (
unsigned int i = 0; i <
rank_; ++i)
1032 for (
unsigned int j = 0;
j < dim_; ++
j)
1033 dst[i][
j][v] = src[i][
j];
1038template <
typename T1,
1040 std::size_t n_lanes,
1048 const unsigned int v,
1051 if constexpr (
rank_ == 1)
1054 copy_data(dst[i], v, src[i]);
1058 for (
unsigned int i = 0; i <
rank_; ++i)
1060 dst[i][
j][v] = src[i][
j];
1065template <
typename T1,
typename T2>
1074 "copy_data_entries() not implemented for given arguments."));
1081 template <
int dim,
typename Number,
typename VectorizedArrayType>
1088 const unsigned int quad_no,
1089 const unsigned int dof_no,
1090 const double tolerance)
1092 const auto &dof_handler = matrix_free.get_dof_handler(
dof_no);
1093 const auto &tria = dof_handler.get_triangulation();
1094 const auto &mapping = *matrix_free.get_mapping_info().mapping;
1109 std::vector<FERemoteCommunicationObjectEntityBatches<dim>>
comm_objects;
1122 std::make_pair(matrix_free.n_inner_face_batches(),
1123 matrix_free.n_inner_face_batches() +
1124 matrix_free.n_boundary_face_batches());
1127 for (
const auto &[
nm_face, marked_vertices] :
1133 auto rpe = std::make_shared<Utilities::MPI::RemotePointEvaluation<dim>>(
1134 tolerance,
false, 0, marked_vertices);
1137 std::vector<std::pair<unsigned int, unsigned int>>
1141 std::vector<Point<dim>> points;
1151 for (
unsigned int bface = 0;
1157 if (matrix_free.get_boundary_id(face) ==
nm_face)
1164 const unsigned int n_faces =
1165 matrix_free.n_active_entries_per_face_batch(face);
1170 for (
unsigned int v = 0; v < n_faces; ++v)
1172 for (
unsigned int q :
phi.quadrature_point_indices())
1174 const auto point =
phi.quadrature_point(
q);
1176 for (
unsigned int i = 0; i < dim; ++i)
1177 temp[i] = point[i][v];
1179 points.push_back(temp);
1188 "Quadrature for given face already provided."));
1195 rpe->reinit(points, tria, mapping);
1196 Assert(rpe->all_points_found(),
1225 template <
int dim,
typename Number,
typename VectorizedArrayType>
1233 const unsigned int dof_no,
1235 const double tolerance)
1237 const auto &dof_handler = matrix_free.get_dof_handler(
dof_no);
1238 const auto &tria = dof_handler.get_triangulation();
1239 const auto &mapping = *matrix_free.get_mapping_info().mapping;
1243 std::pair<unsigned int, unsigned int>
face_range =
1244 std::make_pair(matrix_free.n_inner_face_batches(),
1245 matrix_free.n_inner_face_batches() +
1246 matrix_free.n_boundary_face_batches());
1249 (matrix_free.n_inner_face_batches() +
1250 matrix_free.n_boundary_face_batches()) *
1260 std::vector<FERemoteCommunicationObject<dim>>
comm_objects;
1265 for (
const auto &cell : tria.active_cell_iterators())
1266 if (cell->is_locally_owned())
1267 local_boxes.emplace_back(mapping.get_bounding_box(cell));
1273 std::vector<std::vector<BoundingBox<dim>>>
global_bboxes(1);
1279 for (
const auto &[
nm_face, marked_vertices] :
1284 std::pair<typename Triangulation<dim>::cell_iterator,
unsigned int>>
1287 std::vector<unsigned int> indices;
1292 if (matrix_free.get_boundary_id(face) ==
nm_face)
1294 for (
unsigned int v = 0;
1295 v < matrix_free.n_active_entries_per_face_batch(face);
1298 const auto &[c, f] = matrix_free.get_face_iterator(face, v);
1301 indices.push_back(face * n_lanes + v);
1317 std::vector<Point<dim>> vertices(cell->face(f)->n_vertices());
1318 std::copy_n(mapping.get_vertices(cell, f).begin(),
1319 cell->face(f)->n_vertices(),
1338 std::make_shared<Utilities::MPI::RemotePointEvaluation<dim>>();
1341 .
template convert_to_distributed_compute_point_locations_internal<
1349 const auto idx = indices[i];
1357 std::vector<
Point<dim - 1>> q_points;
1358 std::vector<double> weights;
1365 const auto &
ps = quad.get_points();
1369 std::back_inserter(q_points),
1371 return mapping.project_real_point_to_unit_point_on_face(
1375 const auto &
ws = quad.get_weights();
1376 weights.insert(weights.end(),
ws.begin(),
ws.end());
1381 ExcMessage(
"Quadrature for given face already provided."));
1388 co.indices = indices;
1401 [](
const auto &
q) { return q.size(); });
1411 std::pair<typename DoFHandler<dim>::cell_iterator,
unsigned int>>
1414 matrix_free.n_boundary_face_batches()) *
1421 for (
unsigned int v = 0; v < n_lanes; ++v)
1423 if (v < matrix_free.n_active_entries_per_face_batch(
face_batch))
1425 matrix_free.get_face_iterator(
face_batch, v));
1428 matrix_free.get_face_iterator(
face_batch, 0));
1432 for (;
face_batch < (matrix_free.n_inner_face_batches() +
1433 matrix_free.n_boundary_face_batches());
1436 for (
unsigned int v = 0; v < n_lanes; ++v)
1438 if (v < matrix_free.n_active_entries_per_face_batch(
face_batch))
1440 matrix_free.get_face_iterator(
face_batch, v));
1443 matrix_free.get_face_iterator(
face_batch, 0));
1457template <
int dim,
int n_components,
typename value_type>
1458template <
typename MeshType>
1462 const unsigned int first_selected_component,
1465 , first_selected_component(first_selected_component)
1466 , evaluation_flags(evaluation_flags)
1471template <
int dim,
int n_components,
typename value_type>
1472template <
typename VectorType>
1475 const VectorType &src,
1485 first_selected_component,
1488 else if (dof_handler)
1494 first_selected_component,
1501template <
int dim,
int n_components,
typename value_type>
1510template <
int dim,
int n_components,
typename value_type>
1518template <
int dim,
int n_components,
typename value_type>
1523 this->dof_handler = &dof_handler;
static void copy_data(const internal::PrecomputedEvaluationDataView &view, AlignedVector< T1 > &dst, const std::vector< T2 > &src, const std::vector< unsigned int > &indices)
static void copy_data_entries(VectorizedArray< T1, n_lanes > &dst, const unsigned int v, const T1 &src)
const internal::PrecomputedEvaluationDataView & get_view() const
void reinit_faces(const std::vector< FERemoteCommunicationObjectEntityBatches< dim > > &comm_objects, const std::pair< unsigned int, unsigned int > &face_batch_range, const std::vector< unsigned int > &quadrature_sizes)
void update_ghost_values(PrecomputedEvaluationDataType &dst, const MeshType &mesh, const VectorType &src, const EvaluationFlags::EvaluationFlags eval_flags, const unsigned int first_selected_component, const VectorTools::EvaluationFlags::EvaluationFlags vec_flags) const
internal::PrecomputedEvaluationDataView view
std::vector< std::variant< FERemoteCommunicationObjectEntityBatches< dim >, FERemoteCommunicationObject< dim >, FERemoteCommunicationObjectTwoLevel< dim > > > communication_objects
FERemoteEvaluation(const FERemoteEvaluationCommunicator< dim > &comm, const MeshType &mesh, const unsigned int first_selected_component=0, const VectorTools::EvaluationFlags::EvaluationFlags evaluation_flags=VectorTools::EvaluationFlags::avg)
ObserverPointer< const DoFHandler< dim > > dof_handler
ObserverPointer< const FERemoteEvaluationCommunicator< dim > > comm
void set_mesh(const Triangulation< dim > &tria)
internal::PrecomputedEvaluationData< dim, n_components, value_type > data
ObserverPointer< const Triangulation< dim > > tria
const VectorTools::EvaluationFlags::EvaluationFlags evaluation_flags
const unsigned int first_selected_component
internal::PrecomputedEvaluationDataAccessor< dim, n_components, value_type > get_data_accessor() const
void gather_evaluate(const VectorType &src, const EvaluationFlags::EvaluationFlags flags)
std::conditional_t< rank_==1, std::array< Number, dim >, std::array< Tensor< rank_ - 1, dim, Number >, dim > > values
static constexpr std::size_t size()
const gradient_type get_gradient(const unsigned int q) const
PrecomputedEvaluationDataAccessor(const PrecomputedEvaluationData< dim, n_components, value_type_ > &data, const PrecomputedEvaluationDataView &view)
const PrecomputedEvaluationDataView & view
typename PrecomputedEvaluationData< dim, n_components, value_type_ >::value_type value_type
const PrecomputedEvaluationData< dim, n_components, value_type_ > & data
const value_type get_value(const unsigned int q) const
void reinit(const unsigned int index)
typename PrecomputedEvaluationData< dim, n_components, value_type_ >::gradient_type gradient_type
#define DEAL_II_NAMESPACE_OPEN
#define DEAL_II_NAMESPACE_CLOSE
#define DEAL_II_NOT_IMPLEMENTED()
static ::ExceptionBase & ExcNotImplemented()
#define Assert(cond, exc)
#define AssertDimension(dim1, dim2)
#define AssertIndexRange(index, range)
static ::ExceptionBase & ExcInternalError()
static ::ExceptionBase & ExcMessage(std::string arg1)
std::vector< index_type > data
EvaluationFlags
The EvaluationFlags enum.
FERemoteEvaluationCommunicator< dim > compute_remote_communicator_faces_nitsche_type_mortaring(const MatrixFree< dim, Number, VectorizedArrayType > &matrix_free, const std::vector< std::pair< types::boundary_id, std::function< std::vector< bool >()> > > &non_matching_faces_marked_vertices, const unsigned int n_q_pnts_1D, const unsigned int dof_no=0, NonMatching::MappingInfo< dim, dim, Number > *nm_mapping_info=nullptr, const double tolerance=1e-9)
FERemoteEvaluationCommunicator< dim > compute_remote_communicator_faces_point_to_point_interpolation(const MatrixFree< dim, Number, VectorizedArrayType > &matrix_free, const std::vector< std::pair< types::boundary_id, std::function< std::vector< bool >()> > > &non_matching_faces_marked_vertices, const unsigned int quad_no=0, const unsigned int dof_no=0, const double tolerance=1e-9)
constexpr unsigned int invalid_unsigned_int
std::vector< BoundingBox< boost::geometry::dimension< typename Rtree::indexable_type >::value > > extract_rtree_level(const Rtree &tree, const unsigned int level)
RTree< typename LeafTypeIterator::value_type, IndexType, IndexableGetter > pack_rtree(const LeafTypeIterator &begin, const LeafTypeIterator &end)
std::vector< std::pair< unsigned int, unsigned int > > batch_id_n_entities
std::vector< std::pair< unsigned int, unsigned int > > get_communication_object_pntrs() const
std::shared_ptr< Utilities::MPI::RemotePointEvaluation< dim > > rpe
std::shared_ptr< Utilities::MPI::RemotePointEvaluation< dim > > rpe
std::vector< std::pair< typename Triangulation< dim >::cell_iterator, unsigned int > > cell_face_nos
std::vector< std::pair< typename Triangulation< dim >::cell_iterator, unsigned int > > get_communication_object_pntrs() const
std::shared_ptr< Utilities::MPI::RemotePointEvaluation< dim > > rpe
std::vector< unsigned int > indices
std::vector< unsigned int > get_communication_object_pntrs() const
unsigned int get_shift(const unsigned int index) const
unsigned int size() const
std::vector< unsigned int > ptrs_ptrs
std::vector< unsigned int > ptrs
AlignedVector< value_type > values
typename internal::FEPointEvaluation::EvaluatorTypeTraits< dim, dim, n_components, value_type_ >::value_type value_type
typename internal::FEPointEvaluation::EvaluatorTypeTraits< dim, dim, n_components, value_type_ >::real_gradient_type gradient_type
AlignedVector< gradient_type > gradients