24 template <
int dim,
int spacedim>
34 , cell_quadrature(quadrature)
35 , face_quadrature(face_quadrature)
36 , hp_capability_enabled(false)
37 , cell_update_flags(update_flags)
38 , neighbor_cell_update_flags(update_flags)
39 , face_update_flags(face_update_flags)
40 , neighbor_face_update_flags(face_update_flags)
41 , local_dof_indices(fe.n_dofs_per_cell())
42 , neighbor_dof_indices(fe.n_dofs_per_cell())
47 template <
int dim,
int spacedim>
59 , cell_quadrature(quadrature)
60 , face_quadrature(face_quadrature)
61 , hp_capability_enabled(false)
62 , cell_update_flags(update_flags)
63 , neighbor_cell_update_flags(neighbor_update_flags)
64 , face_update_flags(face_update_flags)
65 , neighbor_face_update_flags(neighbor_face_update_flags)
66 , local_dof_indices(fe.n_dofs_per_cell())
67 , neighbor_dof_indices(fe.n_dofs_per_cell())
72 template <
int dim,
int spacedim>
90 template <
int dim,
int spacedim>
104 neighbor_update_flags,
107 neighbor_face_update_flags)
112 template <
int dim,
int spacedim>
120 : mapping_collection(&mapping_collection)
121 , fe_collection(&fe_collection)
122 , cell_quadrature_collection(cell_quadrature_collection)
123 , face_quadrature_collection(face_quadrature_collection)
124 , hp_capability_enabled(true)
125 , cell_update_flags(cell_update_flags)
126 , neighbor_cell_update_flags(cell_update_flags)
127 , face_update_flags(face_update_flags)
128 , neighbor_face_update_flags(face_update_flags)
136 template <
int dim,
int spacedim>
146 : mapping_collection(&mapping_collection)
147 , fe_collection(&fe_collection)
148 , cell_quadrature_collection(cell_quadrature_collection)
149 , face_quadrature_collection(face_quadrature_collection)
150 , hp_capability_enabled(true)
151 , cell_update_flags(cell_update_flags)
152 , neighbor_cell_update_flags(neighbor_cell_update_flags)
153 , face_update_flags(face_update_flags)
154 , neighbor_face_update_flags(neighbor_face_update_flags)
162 template <
int dim,
int spacedim>
169 :
ScratchData(fe_collection.get_reference_cell_default_linear_mapping(),
171 cell_quadrature_collection,
173 face_quadrature_collection,
179 template <
int dim,
int spacedim>
188 :
ScratchData(fe_collection.get_reference_cell_default_linear_mapping(),
190 cell_quadrature_collection,
192 neighbor_cell_update_flags,
193 face_quadrature_collection,
195 neighbor_face_update_flags)
200 template <
int dim,
int spacedim>
203 : mapping(scratch.mapping)
205 , cell_quadrature(scratch.cell_quadrature)
206 , face_quadrature(scratch.face_quadrature)
207 , mapping_collection(scratch.mapping_collection)
208 , fe_collection(scratch.fe_collection)
209 , cell_quadrature_collection(scratch.cell_quadrature_collection)
210 , face_quadrature_collection(scratch.face_quadrature_collection)
211 , hp_capability_enabled(scratch.hp_capability_enabled)
212 , cell_update_flags(scratch.cell_update_flags)
213 , neighbor_cell_update_flags(scratch.neighbor_cell_update_flags)
214 , face_update_flags(scratch.face_update_flags)
215 , neighbor_face_update_flags(scratch.neighbor_face_update_flags)
216 , local_dof_indices(scratch.local_dof_indices)
217 , neighbor_dof_indices(scratch.neighbor_dof_indices)
218 , user_data_storage(scratch.user_data_storage)
219 , internal_data_storage(scratch.internal_data_storage)
224 template <
int dim,
int spacedim>
229 if (hp_capability_enabled ==
false)
232 fe_values = std::make_unique<FEValues<dim, spacedim>>(
233 *mapping, *fe, cell_quadrature, cell_update_flags);
235 fe_values->reinit(cell);
236 local_dof_indices.resize(fe_values->dofs_per_cell);
237 cell->get_dof_indices(local_dof_indices);
238 current_fe_values = fe_values.get();
244 hp_fe_values = std::make_unique<hp::FEValues<dim, spacedim>>(
247 cell_quadrature_collection,
250 hp_fe_values->
reinit(cell);
251 const auto &fe_values = hp_fe_values->get_present_fe_values();
254 (*fe_collection)[cell->active_fe_index()].n_dofs_per_cell(),
255 fe_values.dofs_per_cell);
256 local_dof_indices.resize(fe_values.dofs_per_cell);
257 cell->get_dof_indices(local_dof_indices);
259 current_fe_values = &fe_values;
266 template <
int dim,
int spacedim>
270 const unsigned int face_no)
272 if (hp_capability_enabled ==
false)
275 fe_face_values = std::make_unique<FEFaceValues<dim, spacedim>>(
276 *mapping, *fe, face_quadrature, face_update_flags);
278 fe_face_values->reinit(cell, face_no);
279 local_dof_indices.resize(fe->n_dofs_per_cell());
280 cell->get_dof_indices(local_dof_indices);
281 current_fe_values = fe_face_values.get();
282 return *fe_face_values;
286 return reinit(cell, cell, face_no);
292 template <
int dim,
int spacedim>
298 const unsigned int face_no)
300 Assert(hp_capability_enabled, ExcOnlyAvailableWithHP());
302 if (!hp_fe_face_values)
303 hp_fe_face_values = std::make_unique<hp::FEFaceValues<dim, spacedim>>(
306 face_quadrature_collection,
309 if (neighbor_cell == cell)
311 hp_fe_face_values->reinit(cell, face_no);
320 unsigned int dominated_fe_index = fe_collection->find_dominated_fe(
321 {cell->active_fe_index(), neighbor_cell->active_fe_index()});
330 hp_fe_face_values->reinit(cell,
336 const auto &fe_face_values = hp_fe_face_values->get_present_fe_values();
337 const auto &fe = (*fe_collection)[cell->active_fe_index()];
339 local_dof_indices.resize(fe.n_dofs_per_cell());
340 cell->get_dof_indices(local_dof_indices);
342 current_fe_values = &fe_face_values;
343 return fe_face_values;
348 template <
int dim,
int spacedim>
352 const unsigned int face_no,
353 const unsigned int subface_no)
357 if (hp_capability_enabled ==
false)
359 if (!fe_subface_values)
361 std::make_unique<FESubfaceValues<dim, spacedim>>(
362 *mapping, *fe, face_quadrature, face_update_flags);
363 fe_subface_values->reinit(cell, face_no, subface_no);
364 local_dof_indices.resize(fe->n_dofs_per_cell());
365 cell->get_dof_indices(local_dof_indices);
367 current_fe_values = fe_subface_values.get();
368 return *fe_subface_values;
372 return reinit(cell, cell, face_no, subface_no);
376 return reinit(cell, face_no);
381 template <
int dim,
int spacedim>
387 const unsigned int face_no,
388 const unsigned int subface_no)
390 Assert(hp_capability_enabled, ExcOnlyAvailableWithHP());
394 if (!hp_fe_subface_values)
395 hp_fe_subface_values =
396 std::make_unique<hp::FESubfaceValues<dim, spacedim>>(
399 face_quadrature_collection,
402 if (neighbor_cell == cell)
404 hp_fe_subface_values->reinit(cell, face_no, subface_no);
413 unsigned int dominated_fe_index = fe_collection->find_dominated_fe(
414 {cell->active_fe_index(), neighbor_cell->active_fe_index()});
423 hp_fe_subface_values->reinit(cell,
430 const auto &fe_subface_values =
431 hp_fe_subface_values->get_present_fe_values();
432 const auto &fe = (*fe_collection)[cell->active_fe_index()];
434 local_dof_indices.resize(fe.n_dofs_per_cell());
435 cell->get_dof_indices(local_dof_indices);
437 current_fe_values = &fe_subface_values;
438 return fe_subface_values;
442 return reinit(cell, neighbor_cell, face_no);
448 template <
int dim,
int spacedim>
452 const unsigned int face_no,
455 const unsigned int face_no_neighbor)
467 template <
int dim,
int spacedim>
471 const unsigned int face_no,
472 const unsigned int sub_face_no,
475 const unsigned int face_no_neighbor,
476 const unsigned int sub_face_no_neighbor)
478 if (hp_capability_enabled ==
false)
480 if (!interface_fe_values)
481 interface_fe_values =
482 std::make_unique<FEInterfaceValues<dim, spacedim>>(
483 *mapping, *fe, face_quadrature, face_update_flags);
485 interface_fe_values->reinit(cell,
490 sub_face_no_neighbor);
494 if (!interface_fe_values)
495 interface_fe_values =
496 std::make_unique<FEInterfaceValues<dim, spacedim>>(
499 face_quadrature_collection,
507 unsigned int dominated_fe_index = fe_collection->find_dominated_fe(
508 {cell->active_fe_index(), cell_neighbor->active_fe_index()});
517 interface_fe_values->reinit(cell,
522 sub_face_no_neighbor,
527 current_fe_values = &interface_fe_values->get_fe_face_values(0);
528 current_neighbor_fe_values = &interface_fe_values->get_fe_face_values(1);
530 local_dof_indices = interface_fe_values->get_interface_dof_indices();
531 return *interface_fe_values;
536 template <
int dim,
int spacedim>
541 if (hp_capability_enabled ==
false)
543 if (!neighbor_fe_values)
544 neighbor_fe_values = std::make_unique<FEValues<dim, spacedim>>(
545 *mapping, *fe, cell_quadrature, neighbor_cell_update_flags);
547 neighbor_fe_values->reinit(cell);
548 cell->get_dof_indices(neighbor_dof_indices);
549 current_neighbor_fe_values = neighbor_fe_values.get();
550 return *neighbor_fe_values;
554 if (!neighbor_hp_fe_values)
555 neighbor_hp_fe_values = std::make_unique<hp::FEValues<dim, spacedim>>(
558 cell_quadrature_collection,
559 neighbor_cell_update_flags);
561 neighbor_hp_fe_values->
reinit(cell);
562 const auto &neighbor_fe_values =
563 neighbor_hp_fe_values->get_present_fe_values();
566 (*fe_collection)[cell->active_fe_index()].n_dofs_per_cell(),
567 neighbor_fe_values.dofs_per_cell);
568 neighbor_dof_indices.resize(neighbor_fe_values.dofs_per_cell);
569 cell->get_dof_indices(neighbor_dof_indices);
571 current_neighbor_fe_values = &neighbor_fe_values;
572 return neighbor_fe_values;
578 template <
int dim,
int spacedim>
582 const unsigned int face_no)
584 if (hp_capability_enabled ==
false)
586 if (!neighbor_fe_face_values)
587 neighbor_fe_face_values =
588 std::make_unique<FEFaceValues<dim, spacedim>>(
589 *mapping, *fe, face_quadrature, neighbor_face_update_flags);
590 neighbor_fe_face_values->reinit(cell, face_no);
591 cell->get_dof_indices(neighbor_dof_indices);
592 current_neighbor_fe_values = neighbor_fe_face_values.get();
593 return *neighbor_fe_face_values;
597 return reinit_neighbor(cell, cell, face_no);
603 template <
int dim,
int spacedim>
609 const unsigned int face_no)
611 Assert(hp_capability_enabled, ExcOnlyAvailableWithHP());
613 if (!neighbor_hp_fe_face_values)
614 neighbor_hp_fe_face_values =
615 std::make_unique<hp::FEFaceValues<dim, spacedim>>(
618 face_quadrature_collection,
619 neighbor_face_update_flags);
621 if (neighbor_cell == cell)
623 neighbor_hp_fe_face_values->reinit(neighbor_cell, face_no);
632 unsigned int dominated_fe_index = fe_collection->find_dominated_fe(
633 {cell->active_fe_index(), neighbor_cell->active_fe_index()});
642 neighbor_hp_fe_face_values->reinit(neighbor_cell,
648 const auto &neighbor_fe_face_values =
649 neighbor_hp_fe_face_values->get_present_fe_values();
650 const auto &neighbor_fe =
651 (*fe_collection)[neighbor_cell->active_fe_index()];
653 neighbor_dof_indices.resize(neighbor_fe.n_dofs_per_cell());
654 neighbor_cell->get_dof_indices(neighbor_dof_indices);
656 current_neighbor_fe_values = &neighbor_fe_face_values;
657 return neighbor_fe_face_values;
662 template <
int dim,
int spacedim>
666 const unsigned int face_no,
667 const unsigned int subface_no)
671 if (hp_capability_enabled ==
false)
673 if (!neighbor_fe_subface_values)
674 neighbor_fe_subface_values =
675 std::make_unique<FESubfaceValues<dim, spacedim>>(
676 *mapping, *fe, face_quadrature, neighbor_face_update_flags);
677 neighbor_fe_subface_values->reinit(cell, face_no, subface_no);
678 cell->get_dof_indices(neighbor_dof_indices);
679 current_neighbor_fe_values = neighbor_fe_subface_values.get();
680 return *neighbor_fe_subface_values;
684 return reinit_neighbor(cell, cell, face_no, subface_no);
688 return reinit_neighbor(cell, face_no);
693 template <
int dim,
int spacedim>
699 const unsigned int face_no,
700 const unsigned int subface_no)
702 Assert(hp_capability_enabled, ExcOnlyAvailableWithHP());
706 if (!neighbor_hp_fe_subface_values)
707 neighbor_hp_fe_subface_values =
708 std::make_unique<hp::FESubfaceValues<dim, spacedim>>(
711 face_quadrature_collection,
712 neighbor_face_update_flags);
714 if (neighbor_cell == cell)
716 neighbor_hp_fe_subface_values->reinit(neighbor_cell,
727 unsigned int dominated_fe_index = fe_collection->find_dominated_fe(
728 {cell->active_fe_index(), neighbor_cell->active_fe_index()});
737 neighbor_hp_fe_subface_values->reinit(neighbor_cell,
744 const auto &neighbor_fe_subface_values =
745 neighbor_hp_fe_subface_values->get_present_fe_values();
746 const auto &neighbor_fe =
747 (*fe_collection)[neighbor_cell->active_fe_index()];
749 neighbor_dof_indices.resize(neighbor_fe.n_dofs_per_cell());
750 neighbor_cell->get_dof_indices(neighbor_dof_indices);
752 current_neighbor_fe_values = &neighbor_fe_subface_values;
753 return neighbor_fe_subface_values;
757 return reinit_neighbor(cell, neighbor_cell, face_no);
763 template <
int dim,
int spacedim>
767 Assert(current_fe_values !=
nullptr,
768 ExcMessage(
"You have to initialize the cache using one of the "
769 "reinit functions first!"));
770 return *current_fe_values;
775 template <
int dim,
int spacedim>
779 Assert(interface_fe_values !=
nullptr,
780 ExcMessage(
"You have to initialize the cache using one of the "
781 "reinit functions first!"));
782 return *interface_fe_values;
787 template <
int dim,
int spacedim>
791 Assert(current_neighbor_fe_values !=
nullptr,
792 ExcMessage(
"You have to initialize the cache using one of the "
793 "reinit functions first!"));
794 return *current_neighbor_fe_values;
799 template <
int dim,
int spacedim>
800 const std::vector<Point<spacedim>> &
803 return get_current_fe_values().get_quadrature_points();
808 template <
int dim,
int spacedim>
809 const std::vector<double> &
812 return get_current_fe_values().get_JxW_values();
817 template <
int dim,
int spacedim>
818 const std::vector<double> &
821 return get_current_neighbor_fe_values().get_JxW_values();
826 template <
int dim,
int spacedim>
827 const std::vector<Tensor<1, spacedim>> &
830 return get_current_fe_values().get_normal_vectors();
835 template <
int dim,
int spacedim>
836 const std::vector<Tensor<1, spacedim>> &
839 return get_current_neighbor_fe_values().get_normal_vectors();
844 template <
int dim,
int spacedim>
845 const std::vector<types::global_dof_index> &
848 return local_dof_indices;
853 template <
int dim,
int spacedim>
857 return local_dof_indices.size();
862 template <
int dim,
int spacedim>
863 const std::vector<types::global_dof_index> &
866 return neighbor_dof_indices;
871 template <
int dim,
int spacedim>
875 return neighbor_dof_indices.size();
880 template <
int dim,
int spacedim>
884 return user_data_storage;
889 template <
int dim,
int spacedim>
893 return user_data_storage;
898 template <
int dim,
int spacedim>
902 Assert(hp_capability_enabled ==
false, ExcOnlyAvailableWithoutHP());
909 template <
int dim,
int spacedim>
913 Assert(hp_capability_enabled ==
false, ExcOnlyAvailableWithoutHP());
920 template <
int dim,
int spacedim>
924 Assert(hp_capability_enabled ==
false, ExcOnlyAvailableWithoutHP());
925 return cell_quadrature;
930 template <
int dim,
int spacedim>
934 Assert(hp_capability_enabled ==
false, ExcOnlyAvailableWithoutHP());
935 return face_quadrature;
940 template <
int dim,
int spacedim>
944 Assert(hp_capability_enabled ==
true, ExcOnlyAvailableWithHP());
946 return *mapping_collection;
951 template <
int dim,
int spacedim>
955 Assert(hp_capability_enabled ==
true, ExcOnlyAvailableWithHP());
957 return *fe_collection;
962 template <
int dim,
int spacedim>
966 Assert(hp_capability_enabled ==
true, ExcOnlyAvailableWithHP());
967 return cell_quadrature_collection;
972 template <
int dim,
int spacedim>
976 Assert(hp_capability_enabled ==
true, ExcOnlyAvailableWithHP());
977 return face_quadrature_collection;
982 template <
int dim,
int spacedim>
986 return hp_capability_enabled;
991 template <
int dim,
int spacedim>
995 return cell_update_flags;
1000 template <
int dim,
int spacedim>
1004 return neighbor_cell_update_flags;
1009 template <
int dim,
int spacedim>
1013 return face_update_flags;
1018 template <
int dim,
int spacedim>
1022 return neighbor_face_update_flags;
1032 #include "scratch_data.inst"
void reinit(const Triangulation< dim, spacedim > &tria)
bool has_hp_capabilities() const
const hp::MappingCollection< dim, spacedim > & get_mapping_collection() const
const FEValuesBase< dim, spacedim > & get_current_fe_values() const
GeneralDataStorage & get_general_data_storage()
const Quadrature< dim - 1 > & get_face_quadrature() const
const std::vector< Tensor< 1, spacedim > > & get_neighbor_normal_vectors()
unsigned int n_dofs_per_cell() const
unsigned int n_neighbor_dofs_per_cell() const
const std::vector< double > & get_JxW_values() const
const Mapping< dim, spacedim > & get_mapping() const
const std::vector< double > & get_neighbor_JxW_values() const
const FiniteElement< dim, spacedim > & get_fe() const
SmartPointer< const hp::FECollection< dim, spacedim > > fe_collection
const hp::QCollection< dim > & get_cell_quadrature_collection() const
std::vector< types::global_dof_index > neighbor_dof_indices
const std::vector< types::global_dof_index > & get_neighbor_dof_indices() const
UpdateFlags get_face_update_flags() const
UpdateFlags get_neighbor_cell_update_flags() const
const hp::FECollection< dim, spacedim > & get_fe_collection() const
const std::vector< Tensor< 1, spacedim > > & get_normal_vectors() const
const FEValues< dim, spacedim > & reinit_neighbor(const typename DoFHandler< dim, spacedim >::active_cell_iterator &cell)
std::vector< types::global_dof_index > local_dof_indices
const std::vector< Point< spacedim > > & get_quadrature_points() const
const FEValuesBase< dim, spacedim > & get_current_neighbor_fe_values() const
const std::vector< types::global_dof_index > & get_local_dof_indices() const
const FEInterfaceValues< dim, spacedim > & get_current_interface_fe_values() const
UpdateFlags get_cell_update_flags() const
const Quadrature< dim > & get_cell_quadrature() const
const hp::QCollection< dim - 1 > & get_face_quadrature_collection() const
ScratchData(const Mapping< dim, spacedim > &mapping, const FiniteElement< dim, spacedim > &fe, const Quadrature< dim > &quadrature, const UpdateFlags &update_flags, const Quadrature< dim - 1 > &face_quadrature=Quadrature< dim - 1 >(), const UpdateFlags &face_update_flags=update_default)
const FEValues< dim, spacedim > & reinit(const typename DoFHandler< dim, spacedim >::active_cell_iterator &cell)
UpdateFlags get_neighbor_face_update_flags() const
#define DEAL_II_NAMESPACE_OPEN
#define DEAL_II_NAMESPACE_CLOSE
static ::ExceptionBase & ExcNotInitialized()
#define Assert(cond, exc)
#define AssertDimension(dim1, dim2)
static ::ExceptionBase & ExcMessage(std::string arg1)
const Mapping< dim, spacedim > & get_default_linear_mapping(const Triangulation< dim, spacedim > &triangulation)
void reference_cell(Triangulation< dim, spacedim > &tria, const ReferenceCell &reference_cell)
void reinit(MatrixBlock< MatrixType > &v, const BlockSparsityPattern &p)
const types::fe_index invalid_fe_index
static const unsigned int invalid_unsigned_int