Reference documentation for deal.II version 9.5.0
\(\newcommand{\dealvcentcolon}{\mathrel{\mathop{:}}}\) \(\newcommand{\dealcoloneq}{\dealvcentcolon\mathrel{\mkern-1.2mu}=}\) \(\newcommand{\jump}[1]{\left[\!\left[ #1 \right]\!\right]}\) \(\newcommand{\average}[1]{\left\{\!\left\{ #1 \right\}\!\right\}}\)
Loading...
Searching...
No Matches
mg_transfer_global_coarsening.h
Go to the documentation of this file.
1// ---------------------------------------------------------------------
2//
3// Copyright (C) 2020 - 2023 by the deal.II authors
4//
5// This file is part of the deal.II library.
6//
7// The deal.II library is free software; you can use it, redistribute
8// it, and/or modify it under the terms of the GNU Lesser General
9// Public License as published by the Free Software Foundation; either
10// version 2.1 of the License, or (at your option) any later version.
11// The full text of the license can be found in the file LICENSE.md at
12// the top level directory of deal.II.
13//
14// ---------------------------------------------------------------------
15
16#ifndef dealii_mg_transfer_global_coarsening_h
17#define dealii_mg_transfer_global_coarsening_h
18
22
24
27
30
33
35
36
37
39
40// Forward declarations
41#ifndef DOXYGEN
42namespace internal
43{
44 class MGTwoLevelTransferImplementation;
45}
46
48{
49 template <int dim, int spacedim>
50 class Base;
51}
52#endif
53
54
55
60{
69 {
74 bisect,
85 };
86
91 unsigned int
93 const unsigned int degree,
94 const PolynomialCoarseningSequenceType &p_sequence);
95
101 std::vector<unsigned int>
103 const unsigned int max_degree,
104 const PolynomialCoarseningSequenceType &p_sequence);
105
116 template <int dim, int spacedim>
117 std::vector<std::shared_ptr<const Triangulation<dim, spacedim>>>
120
137 template <int dim, int spacedim>
138 std::vector<std::shared_ptr<const Triangulation<dim, spacedim>>>
142 const bool preserve_fine_triangulation,
143 const bool repartition_fine_triangulation);
144
150 template <int dim, int spacedim>
151 std::vector<std::shared_ptr<const Triangulation<dim, spacedim>>>
155 const bool repartition_fine_triangulation = false);
156
157} // namespace MGTransferGlobalCoarseningTools
158
159
163template <typename VectorType>
165{
166public:
170 virtual void
171 prolongate_and_add(VectorType &dst, const VectorType &src) const = 0;
172
176 virtual void
177 restrict_and_add(VectorType &dst, const VectorType &src) const = 0;
178
185 virtual void
186 interpolate(VectorType &dst, const VectorType &src) const = 0;
187
192 virtual void
194 const std::shared_ptr<const Utilities::MPI::Partitioner>
195 &partitioner_coarse,
196 const std::shared_ptr<const Utilities::MPI::Partitioner>
197 &partitioner_fine) = 0;
198
202 virtual std::size_t
204};
205
206
214template <typename Number>
215class MGTwoLevelTransferBase<LinearAlgebra::distributed::Vector<Number>>
216 : public Subscriptor
217{
218public:
220
224 virtual void
225 prolongate_and_add(VectorType &dst, const VectorType &src) const;
226
230 virtual void
231 restrict_and_add(VectorType &dst, const VectorType &src) const;
232
237 virtual void
238 interpolate(VectorType &dst, const VectorType &src) const = 0;
239
244 virtual void
246 const std::shared_ptr<const Utilities::MPI::Partitioner>
247 &partitioner_coarse,
248 const std::shared_ptr<const Utilities::MPI::Partitioner>
249 &partitioner_fine) = 0;
250
254 virtual std::size_t
256
257protected:
261 virtual void
264 const LinearAlgebra::distributed::Vector<Number> &src) const = 0;
265
269 virtual void
272 const LinearAlgebra::distributed::Vector<Number> &src) const = 0;
273
279 void
282
288 void
290 const VectorOperation::values op) const;
291
297 void
300
305 template <int dim, std::size_t width>
306 void
308 const std::shared_ptr<const Utilities::MPI::Partitioner>
309 &partitioner_coarse,
310 const std::shared_ptr<const Utilities::MPI::Partitioner> &partitioner_fine,
312 dim,
313 VectorizedArray<Number, width>> &constraint_info_coarse,
314 std::vector<unsigned int> & dof_indices_fine);
315
322
326 std::shared_ptr<const Utilities::MPI::Partitioner> partitioner_coarse;
327
331 std::shared_ptr<const Utilities::MPI::Partitioner> partitioner_fine;
332
341
346
351 std::shared_ptr<const Utilities::MPI::Partitioner>
353
358 std::shared_ptr<const Utilities::MPI::Partitioner> partitioner_fine_embedded;
359
365
371};
372
373
374
380template <int dim, typename VectorType>
382{
383public:
387 void
388 prolongate_and_add(VectorType &dst, const VectorType &src) const override;
389
393 void
394 restrict_and_add(VectorType &dst, const VectorType &src) const override;
395
400 void
401 interpolate(VectorType &dst, const VectorType &src) const override;
402
407 void
409 const std::shared_ptr<const Utilities::MPI::Partitioner>
410 &partitioner_coarse,
411 const std::shared_ptr<const Utilities::MPI::Partitioner> &partitioner_fine)
412 override;
413
417 std::size_t
418 memory_consumption() const override;
419};
420
421
422
429template <int dim, typename Number>
430class MGTwoLevelTransfer<dim, LinearAlgebra::distributed::Vector<Number>>
431 : public MGTwoLevelTransferBase<LinearAlgebra::distributed::Vector<Number>>
432{
434
435public:
441 void
443 const DoFHandler<dim> & dof_handler_fine,
444 const DoFHandler<dim> & dof_handler_coarse,
445 const AffineConstraints<Number> &constraint_fine =
447 const AffineConstraints<Number> &constraint_coarse =
449 const unsigned int mg_level_fine = numbers::invalid_unsigned_int,
450 const unsigned int mg_level_coarse = numbers::invalid_unsigned_int);
451
461 void
463 const DoFHandler<dim> & dof_handler_fine,
464 const DoFHandler<dim> & dof_handler_coarse,
465 const AffineConstraints<Number> &constraint_fine =
467 const AffineConstraints<Number> &constraint_coarse =
469 const unsigned int mg_level_fine = numbers::invalid_unsigned_int,
470 const unsigned int mg_level_coarse = numbers::invalid_unsigned_int);
471
485 void
486 reinit(const DoFHandler<dim> & dof_handler_fine,
487 const DoFHandler<dim> & dof_handler_coarse,
488 const AffineConstraints<Number> &constraint_fine =
490 const AffineConstraints<Number> &constraint_coarse =
492 const unsigned int mg_level_fine = numbers::invalid_unsigned_int,
493 const unsigned int mg_level_coarse = numbers::invalid_unsigned_int);
494
503 static bool
504 fast_polynomial_transfer_supported(const unsigned int fe_degree_fine,
505 const unsigned int fe_degree_coarse);
506
511 void
514 const LinearAlgebra::distributed::Vector<Number> &src) const override;
515
520 void
522 const std::shared_ptr<const Utilities::MPI::Partitioner>
523 &partitioner_coarse,
524 const std::shared_ptr<const Utilities::MPI::Partitioner> &partitioner_fine)
525 override;
526
530 std::size_t
531 memory_consumption() const override;
532
533protected:
534 void
537 const LinearAlgebra::distributed::Vector<Number> &src) const override;
538
539 void
542 const LinearAlgebra::distributed::Vector<Number> &src) const override;
543
544private:
552 struct MGTransferScheme
553 {
557 unsigned int n_coarse_cells;
558
566
574
578 unsigned int degree_coarse;
579
585 unsigned int degree_fine;
586
591
596
601
606
613 };
614
618 std::vector<MGTransferScheme> schemes;
619
626
632
636 std::vector<Number> weights; // TODO: vectorize
637
643
647 unsigned int n_components;
648
649 friend class internal::MGTwoLevelTransferImplementation;
650};
651
652
653
658template <int dim, typename VectorType>
660{
661public:
665 void
666 prolongate_and_add(VectorType &dst, const VectorType &src) const override;
667
671 void
672 restrict_and_add(VectorType &dst, const VectorType &src) const override;
673
680 void
681 interpolate(VectorType &dst, const VectorType &src) const override;
682
686 std::size_t
687 memory_consumption() const override;
688};
689
690
691
698template <int dim, typename Number>
700 LinearAlgebra::distributed::Vector<Number>>
701 : public MGTwoLevelTransferBase<LinearAlgebra::distributed::Vector<Number>>
702{
703private:
705
706public:
711 void
712 reinit(const DoFHandler<dim> & dof_handler_fine,
713 const DoFHandler<dim> & dof_handler_coarse,
714 const Mapping<dim> & mapping_fine,
715 const Mapping<dim> & mapping_coarse,
716 const AffineConstraints<Number> &constraint_fine =
718 const AffineConstraints<Number> &constraint_coarse =
720
727 void
730 const LinearAlgebra::distributed::Vector<Number> &src) const override;
731
736 void
738 const std::shared_ptr<const Utilities::MPI::Partitioner>
739 &partitioner_coarse,
740 const std::shared_ptr<const Utilities::MPI::Partitioner> &partitioner_fine)
741 override;
742
746 std::size_t
747 memory_consumption() const override;
748
749protected:
753 void
756 const LinearAlgebra::distributed::Vector<Number> &src) const override;
757
761 void
764 const LinearAlgebra::distributed::Vector<Number> &src) const override;
765
766private:
772
776 std::shared_ptr<NonMatching::MappingInfo<dim, dim, Number>> mapping_info;
777
784
788 std::unique_ptr<FiniteElement<dim>> fe_coarse;
789
794 std::vector<unsigned int> level_dof_indices_fine;
795
801 std::vector<unsigned int> level_dof_indices_fine_ptrs;
802};
803
804
805
818template <int dim, typename VectorType>
820{
821public:
825 using Number = typename VectorType::value_type;
826
836 template <typename MGTwoLevelTransferObject>
839 const std::function<void(const unsigned int, VectorType &)>
841
850 void
851 build(const std::vector<std::shared_ptr<const Utilities::MPI::Partitioner>>
853
858 void
859 build(const std::function<void(const unsigned int, VectorType &)>
861
865 void
866 prolongate(const unsigned int to_level,
867 VectorType & dst,
868 const VectorType & src) const override;
869
873 void
874 prolongate_and_add(const unsigned int to_level,
875 VectorType & dst,
876 const VectorType & src) const override;
877
881 virtual void
882 restrict_and_add(const unsigned int from_level,
883 VectorType & dst,
884 const VectorType & src) const override;
885
892 template <class InVector, int spacedim>
893 void
896 const InVector & src) const;
897
904 template <class OutVector, int spacedim>
905 void
907 OutVector & dst,
908 const MGLevelObject<VectorType> &src) const;
909
922 template <class InVector>
923 void
924 interpolate_to_mg(MGLevelObject<VectorType> &dst, const InVector &src) const;
925
932 template <class InVector, int spacedim>
933 void
936 const InVector & src) const;
937
944 std::size_t
946
950 unsigned int
951 min_level() const;
952
956 unsigned int
957 max_level() const;
958
959private:
963 void
964 initialize_dof_vector(const unsigned int level, VectorType &vector) const;
965
970
974 std::vector<std::shared_ptr<const Utilities::MPI::Partitioner>>
976};
977
978
979
985template <int dim, typename VectorType>
988 dim,
989 typename VectorType::value_type,
990 MGTransferGlobalCoarsening<dim, VectorType>>
991{
992public:
998
999protected:
1001 get_matrix_free_transfer(const unsigned int b) const override;
1002
1003private:
1008};
1009
1010
1011
1012#ifndef DOXYGEN
1013
1014/* ----------------------- Inline functions --------------------------------- */
1015
1016
1017
1018template <int dim, typename VectorType>
1019template <typename MGTwoLevelTransferObject>
1022 const std::function<void(const unsigned int, VectorType &)>
1023 &initialize_dof_vector)
1024{
1025 const unsigned int min_level = transfer.min_level();
1026 const unsigned int max_level = transfer.max_level();
1027
1028 this->transfer.resize(min_level, max_level);
1029
1030 for (unsigned int l = min_level; l <= max_level; ++l)
1031 this->transfer[l] = &const_cast<MGTwoLevelTransferBase<VectorType> &>(
1032 static_cast<const MGTwoLevelTransferBase<VectorType> &>(
1033 Utilities::get_underlying_value(transfer[l])));
1034
1035 this->build(initialize_dof_vector);
1036}
1037
1038
1039
1040template <int dim, typename VectorType>
1041void
1043 const std::vector<std::shared_ptr<const Utilities::MPI::Partitioner>>
1044 &external_partitioners)
1045{
1046 this->external_partitioners = external_partitioners;
1047
1048 if (this->external_partitioners.size() > 0)
1049 {
1050 const unsigned int min_level = transfer.min_level();
1051 const unsigned int max_level = transfer.max_level();
1052
1053 AssertDimension(this->external_partitioners.size(), transfer.n_levels());
1054
1055 for (unsigned int l = min_level + 1; l <= max_level; ++l)
1056 transfer[l]->enable_inplace_operations_if_possible(
1057 this->external_partitioners[l - 1 - min_level],
1058 this->external_partitioners[l - min_level]);
1059 }
1060}
1061
1062
1063
1064template <int dim, typename VectorType>
1065void
1067 const std::function<void(const unsigned int, VectorType &)>
1068 &initialize_dof_vector)
1069{
1070 if (initialize_dof_vector)
1071 {
1072 const unsigned int min_level = transfer.min_level();
1073 const unsigned int max_level = transfer.max_level();
1074 const unsigned int n_levels = transfer.n_levels();
1075
1076 std::vector<std::shared_ptr<const Utilities::MPI::Partitioner>>
1077 external_partitioners(n_levels);
1078
1079 for (unsigned int l = min_level; l <= max_level; ++l)
1080 {
1082 vector;
1083 initialize_dof_vector(l, vector);
1084 external_partitioners[l - min_level] = vector.get_partitioner();
1085 }
1086
1087 this->build(external_partitioners);
1088 }
1089 else
1090 {
1091 this->build();
1092 }
1093}
1094
1095
1096
1097template <int dim, typename VectorType>
1098void
1100 const unsigned int level,
1101 VectorType & vec) const
1102{
1103 AssertDimension(transfer.n_levels(), external_partitioners.size());
1104 AssertIndexRange(level, transfer.max_level() + 1);
1105
1106 const auto &external_partitioner =
1107 external_partitioners[level - transfer.min_level()];
1108
1109 if (vec.get_partitioner().get() != external_partitioner.get())
1110 vec.reinit(external_partitioner);
1111}
1112
1113
1114
1115template <int dim, typename VectorType>
1116void
1118 const unsigned int to_level,
1119 VectorType & dst,
1120 const VectorType & src) const
1121{
1122 dst = 0;
1123 prolongate_and_add(to_level, dst, src);
1124}
1125
1126
1127
1128template <int dim, typename VectorType>
1129void
1131 const unsigned int to_level,
1132 VectorType & dst,
1133 const VectorType & src) const
1134{
1135 this->transfer[to_level]->prolongate_and_add(dst, src);
1136}
1137
1138
1139
1140template <int dim, typename VectorType>
1141void
1143 const unsigned int from_level,
1144 VectorType & dst,
1145 const VectorType & src) const
1146{
1147 this->transfer[from_level]->restrict_and_add(dst, src);
1148}
1149
1150
1151
1152template <int dim, typename VectorType>
1153template <class InVector, int spacedim>
1154void
1156 const DoFHandler<dim, spacedim> &dof_handler,
1158 const InVector & src) const
1159{
1160 (void)dof_handler;
1161
1162 for (unsigned int level = dst.min_level(); level <= dst.max_level(); ++level)
1163 {
1164 initialize_dof_vector(level, dst[level]);
1165
1166 if (level == dst.max_level())
1167 dst[level].copy_locally_owned_data_from(src);
1168 else
1169 dst[level] = 0.0;
1170 }
1171}
1172
1173
1174
1175template <int dim, typename VectorType>
1176template <class OutVector, int spacedim>
1177void
1179 const DoFHandler<dim, spacedim> &dof_handler,
1180 OutVector & dst,
1181 const MGLevelObject<VectorType> &src) const
1182{
1183 (void)dof_handler;
1184
1185 dst.copy_locally_owned_data_from(src[src.max_level()]);
1186}
1187
1188
1189
1190template <int dim, typename VectorType>
1191template <class InVector>
1192void
1195 const InVector & src) const
1196{
1197 const unsigned int min_level = transfer.min_level();
1198 const unsigned int max_level = transfer.max_level();
1199
1200 AssertDimension(min_level, dst.min_level());
1201 AssertDimension(max_level, dst.max_level());
1202
1203 for (unsigned int level = min_level; level <= max_level; ++level)
1204 initialize_dof_vector(level, dst[level]);
1205
1206 dst[transfer.max_level()].copy_locally_owned_data_from(src);
1207
1208 for (unsigned int l = max_level; l > min_level; --l)
1209 this->transfer[l]->interpolate(dst[l - 1], dst[l]);
1210}
1211
1212
1213
1214template <int dim, typename VectorType>
1215template <class InVector, int spacedim>
1216void
1218 const DoFHandler<dim, spacedim> &dof_handler,
1220 const InVector & src) const
1221{
1222 (void)dof_handler;
1223
1224 this->interpolate_to_mg(dst, src);
1225}
1226
1227
1228
1229template <int dim, typename VectorType>
1230std::size_t
1232{
1233 std::size_t size = 0;
1234
1235 const unsigned int min_level = transfer.min_level();
1236 const unsigned int max_level = transfer.max_level();
1237
1238 for (unsigned int l = min_level + 1; l <= max_level; ++l)
1239 size += this->transfer[l]->memory_consumption();
1240
1241 return size;
1242}
1243
1244
1245
1246template <int dim, typename VectorType>
1247inline unsigned int
1249{
1250 return transfer.min_level();
1251}
1252
1253
1254
1255template <int dim, typename VectorType>
1256inline unsigned int
1258{
1259 return transfer.max_level();
1260}
1261
1262#endif
1263
1265
1266#endif
const std::shared_ptr< const Utilities::MPI::Partitioner > & get_partitioner() const
void resize(const unsigned int new_minlevel, const unsigned int new_maxlevel, Args &&...args)
unsigned int max_level() const
unsigned int min_level() const
unsigned int n_levels() const
const MGTransferGlobalCoarsening< dim, VectorType > & get_matrix_free_transfer(const unsigned int b) const override
MGTransferBlockGlobalCoarsening(const MGTransferGlobalCoarsening< dim, VectorType > &transfer_operator)
const MGTransferGlobalCoarsening< dim, VectorType > & transfer_operator
void copy_from_mg(const DoFHandler< dim, spacedim > &dof_handler, OutVector &dst, const MGLevelObject< VectorType > &src) const
void interpolate_to_mg(const DoFHandler< dim, spacedim > &dof_handler, MGLevelObject< VectorType > &dst, const InVector &src) const
void prolongate(const unsigned int to_level, VectorType &dst, const VectorType &src) const override
typename VectorType::value_type Number
unsigned int min_level() const
std::vector< std::shared_ptr< const Utilities::MPI::Partitioner > > external_partitioners
void copy_to_mg(const DoFHandler< dim, spacedim > &dof_handler, MGLevelObject< VectorType > &dst, const InVector &src) const
MGTransferGlobalCoarsening(const MGLevelObject< MGTwoLevelTransferObject > &transfer, const std::function< void(const unsigned int, VectorType &)> &initialize_dof_vector={})
virtual void restrict_and_add(const unsigned int from_level, VectorType &dst, const VectorType &src) const override
std::size_t memory_consumption() const
void initialize_dof_vector(const unsigned int level, VectorType &vector) const
void build(const std::vector< std::shared_ptr< const Utilities::MPI::Partitioner > > &external_partitioners={})
void build(const std::function< void(const unsigned int, VectorType &)> &initialize_dof_vector)
void interpolate_to_mg(MGLevelObject< VectorType > &dst, const InVector &src) const
void prolongate_and_add(const unsigned int to_level, VectorType &dst, const VectorType &src) const override
MGLevelObject< SmartPointer< MGTwoLevelTransferBase< VectorType > > > transfer
unsigned int max_level() const
std::shared_ptr< const Utilities::MPI::Partitioner > partitioner_fine_embedded
virtual void enable_inplace_operations_if_possible(const std::shared_ptr< const Utilities::MPI::Partitioner > &partitioner_coarse, const std::shared_ptr< const Utilities::MPI::Partitioner > &partitioner_fine)=0
void compress(LinearAlgebra::distributed::Vector< Number > &vec, const VectorOperation::values op) const
void internal_enable_inplace_operations_if_possible(const std::shared_ptr< const Utilities::MPI::Partitioner > &partitioner_coarse, const std::shared_ptr< const Utilities::MPI::Partitioner > &partitioner_fine, internal::MatrixFreeFunctions::ConstraintInfo< dim, VectorizedArray< Number, width > > &constraint_info_coarse, std::vector< unsigned int > &dof_indices_fine)
void zero_out_ghost_values(const LinearAlgebra::distributed::Vector< Number > &vec) const
virtual void restrict_and_add_internal(LinearAlgebra::distributed::Vector< Number > &dst, const LinearAlgebra::distributed::Vector< Number > &src) const =0
virtual void prolongate_and_add(VectorType &dst, const VectorType &src) const
virtual void interpolate(VectorType &dst, const VectorType &src) const =0
void update_ghost_values(const LinearAlgebra::distributed::Vector< Number > &vec) const
virtual void prolongate_and_add_internal(LinearAlgebra::distributed::Vector< Number > &dst, const LinearAlgebra::distributed::Vector< Number > &src) const =0
std::shared_ptr< const Utilities::MPI::Partitioner > partitioner_coarse_embedded
virtual void restrict_and_add(VectorType &dst, const VectorType &src) const
virtual void enable_inplace_operations_if_possible(const std::shared_ptr< const Utilities::MPI::Partitioner > &partitioner_coarse, const std::shared_ptr< const Utilities::MPI::Partitioner > &partitioner_fine)=0
virtual std::size_t memory_consumption() const =0
virtual void restrict_and_add(VectorType &dst, const VectorType &src) const =0
virtual void interpolate(VectorType &dst, const VectorType &src) const =0
virtual void prolongate_and_add(VectorType &dst, const VectorType &src) const =0
void enable_inplace_operations_if_possible(const std::shared_ptr< const Utilities::MPI::Partitioner > &partitioner_coarse, const std::shared_ptr< const Utilities::MPI::Partitioner > &partitioner_fine) override
std::shared_ptr< NonMatching::MappingInfo< dim, dim, Number > > mapping_info
internal::MatrixFreeFunctions::ConstraintInfo< dim, VectorizedArrayType > constraint_info
void prolongate_and_add_internal(LinearAlgebra::distributed::Vector< Number > &dst, const LinearAlgebra::distributed::Vector< Number > &src) const override
void interpolate(LinearAlgebra::distributed::Vector< Number > &dst, const LinearAlgebra::distributed::Vector< Number > &src) const override
void restrict_and_add_internal(LinearAlgebra::distributed::Vector< Number > &dst, const LinearAlgebra::distributed::Vector< Number > &src) const override
void reinit(const DoFHandler< dim > &dof_handler_fine, const DoFHandler< dim > &dof_handler_coarse, const Mapping< dim > &mapping_fine, const Mapping< dim > &mapping_coarse, const AffineConstraints< Number > &constraint_fine=AffineConstraints< Number >(), const AffineConstraints< Number > &constraint_coarse=AffineConstraints< Number >())
void prolongate_and_add(VectorType &dst, const VectorType &src) const override
std::size_t memory_consumption() const override
void restrict_and_add(VectorType &dst, const VectorType &src) const override
void interpolate(VectorType &dst, const VectorType &src) const override
internal::MatrixFreeFunctions::ConstraintInfo< dim, VectorizedArrayType > constraint_info_fine
void enable_inplace_operations_if_possible(const std::shared_ptr< const Utilities::MPI::Partitioner > &partitioner_coarse, const std::shared_ptr< const Utilities::MPI::Partitioner > &partitioner_fine) override
void reinit_polynomial_transfer(const DoFHandler< dim > &dof_handler_fine, const DoFHandler< dim > &dof_handler_coarse, const AffineConstraints< Number > &constraint_fine=AffineConstraints< Number >(), const AffineConstraints< Number > &constraint_coarse=AffineConstraints< Number >(), const unsigned int mg_level_fine=numbers::invalid_unsigned_int, const unsigned int mg_level_coarse=numbers::invalid_unsigned_int)
void interpolate(LinearAlgebra::distributed::Vector< Number > &dst, const LinearAlgebra::distributed::Vector< Number > &src) const override
static bool fast_polynomial_transfer_supported(const unsigned int fe_degree_fine, const unsigned int fe_degree_coarse)
void prolongate_and_add_internal(LinearAlgebra::distributed::Vector< Number > &dst, const LinearAlgebra::distributed::Vector< Number > &src) const override
internal::MatrixFreeFunctions::ConstraintInfo< dim, VectorizedArrayType > constraint_info_coarse
void reinit(const DoFHandler< dim > &dof_handler_fine, const DoFHandler< dim > &dof_handler_coarse, const AffineConstraints< Number > &constraint_fine=AffineConstraints< Number >(), const AffineConstraints< Number > &constraint_coarse=AffineConstraints< Number >(), const unsigned int mg_level_fine=numbers::invalid_unsigned_int, const unsigned int mg_level_coarse=numbers::invalid_unsigned_int)
void restrict_and_add_internal(LinearAlgebra::distributed::Vector< Number > &dst, const LinearAlgebra::distributed::Vector< Number > &src) const override
void reinit_geometric_transfer(const DoFHandler< dim > &dof_handler_fine, const DoFHandler< dim > &dof_handler_coarse, const AffineConstraints< Number > &constraint_fine=AffineConstraints< Number >(), const AffineConstraints< Number > &constraint_coarse=AffineConstraints< Number >(), const unsigned int mg_level_fine=numbers::invalid_unsigned_int, const unsigned int mg_level_coarse=numbers::invalid_unsigned_int)
void enable_inplace_operations_if_possible(const std::shared_ptr< const Utilities::MPI::Partitioner > &partitioner_coarse, const std::shared_ptr< const Utilities::MPI::Partitioner > &partitioner_fine) override
void interpolate(VectorType &dst, const VectorType &src) const override
void restrict_and_add(VectorType &dst, const VectorType &src) const override
void prolongate_and_add(VectorType &dst, const VectorType &src) const override
std::size_t memory_consumption() const override
Abstract base class for mapping classes.
Definition mapping.h:317
#define DEAL_II_NAMESPACE_OPEN
Definition config.h:472
#define DEAL_II_NAMESPACE_CLOSE
Definition config.h:473
unsigned int level
Definition grid_out.cc:4618
#define AssertDimension(dim1, dim2)
#define AssertIndexRange(index, range)
void interpolate(const DoFHandler< dim, spacedim > &dof1, const InVector &u1, const DoFHandler< dim, spacedim > &dof2, OutVector &u2)
std::vector< std::shared_ptr< const Triangulation< dim, spacedim > > > create_geometric_coarsening_sequence(const Triangulation< dim, spacedim > &tria)
unsigned int create_next_polynomial_coarsening_degree(const unsigned int degree, const PolynomialCoarseningSequenceType &p_sequence)
std::vector< unsigned int > create_polynomial_coarsening_sequence(const unsigned int max_degree, const PolynomialCoarseningSequenceType &p_sequence)
std::enable_if_t< std::is_fundamental< T >::value, std::size_t > memory_consumption(const T &t)
Tensor< 2, dim, Number > l(const Tensor< 2, dim, Number > &F, const Tensor< 2, dim, Number > &dF_dt)
T & get_underlying_value(T &p)
Definition utilities.h:1612
static const unsigned int invalid_unsigned_int
Definition types.h:213
const ::Triangulation< dim, spacedim > & tria