24 #ifndef FACEDATACONTAINER_H_
25 #define FACEDATACONTAINER_H_
33 #include <dofs/dof_handler.h>
34 #include <multigrid/mg_dof_handler.h>
35 #include <hp/dof_handler.h>
37 using namespace dealii;
47 template<
template<
int,
int>
class DH,
typename VECTOR,
int dim>
55 "Dummy class, this constructor should never get called.",
56 "ElementDataContainer<dealii::DoFHandler , VECTOR, dim>::ElementDataContainer"));
70 template<
typename VECTOR,
int dim>
97 template<
template<
int,
int>
class FE,
typename SPARSITYPATTERN,
int dopedim,
int dealdim>
99 UpdateFlags update_flags,
101 dopedim, dealdim>& sth,
104 const std::map<std::string,
const Vector<double>*> ¶m_values,
105 const std::map<std::string, const VECTOR*> &domain_values,
106 bool need_neighbour) :
107 fdcinternal::FaceDataContainerInternal<VECTOR, dim>(param_values,
108 domain_values, need_neighbour), _element(element), _state_fe_values(
109 sth.GetMapping(), (sth.GetFESystem(
"state")), quad,
110 update_flags), _control_fe_values(sth.GetMapping(),
111 (sth.GetFESystem(
"control")), quad, update_flags)
113 _state_index = sth.GetStateIndex();
114 if (_state_index == 1)
122 sth.GetMapping(), (sth.GetFESystem(
"control")), quad,
124 _control_fe_subface_values =
126 (sth.GetFESystem(
"control")), quad, update_flags);
128 this->PrivateConstructor(quad, update_flags, sth, need_neighbour);
148 template<
template<
int,
int>
class FE,
typename SPARSITYPATTERN>
150 UpdateFlags update_flags,
155 const std::map<std::string,
const Vector<double>*> ¶m_values,
156 const std::map<std::string, const VECTOR*> &domain_values,
157 bool need_neighbour) :
158 fdcinternal::FaceDataContainerInternal<VECTOR, dim>(param_values,
159 domain_values, need_neighbour), _element(element), _state_fe_values(
160 sth.GetMapping(), (sth.GetFESystem(
"state")), quad,
161 update_flags), _control_fe_values(sth.GetMapping(),
162 (sth.GetFESystem(
"state")), quad, update_flags)
164 _state_index = sth.GetStateIndex();
165 _control_index = element.size();
166 _n_q_points_per_element = quad.size();
167 _n_dofs_per_element = element[0]->get_fe().dofs_per_cell;
172 sth.GetMapping(), (sth.GetFESystem(
"state")), quad,
174 _control_fe_subface_values =
176 (sth.GetFESystem(
"state")), quad, update_flags);
178 this->PrivateConstructor(quad, update_flags, sth, need_neighbour);
183 if (_nbr_state_fe_values != NULL)
185 delete _nbr_state_fe_values;
187 if (_nbr_control_fe_values != NULL)
189 delete _nbr_control_fe_values;
191 if (_state_fe_subface_values != NULL)
193 delete _state_fe_subface_values;
195 if (_control_fe_subface_values != NULL)
197 delete _control_fe_subface_values;
209 ReInit(
unsigned int face_no);
222 ReInit(
unsigned int face_no,
unsigned int subface_no);
240 GetNDoFsPerElement()
const;
242 GetNbrNDoFsPerElement()
const;
246 GetNbrNQPoints()
const;
248 GetMaterialId()
const;
250 GetNbrMaterialId()
const;
252 GetNbrMaterialId(
unsigned int face)
const;
254 GetIsAtBoundary()
const;
256 GetElementDiameter()
const;
258 GetBoundaryIndicator()
const;
259 inline const FEFaceValuesBase<dim>&
260 GetFEFaceValuesState()
const;
261 inline const FEFaceValuesBase<dim>&
262 GetFEFaceValuesControl()
const;
264 inline const FEFaceValuesBase<dim>&
265 GetNbrFEFaceValuesState()
const;
266 inline const FEFaceValuesBase<dim>&
267 GetNbrFEFaceValuesControl()
const;
274 GetStateIndex()
const;
276 GetControlIndex()
const;
282 PrivateConstructor(
const Quadrature<dim - 1>& quad,
283 UpdateFlags update_flags, STH& sth,
bool need_neighbour)
285 _n_q_points_per_element = quad.size();
286 _n_dofs_per_element = _element[0]->get_fe().dofs_per_cell;
291 (sth.GetFESystem(
"state")), quad, update_flags);
293 (sth.GetFESystem(
"state")), quad, update_flags);
297 _nbr_state_fe_values = NULL;
298 _nbr_control_fe_values = NULL;
299 _state_fe_subface_values = NULL;
300 _control_fe_subface_values = NULL;
304 _state_fe_values_ptr = NULL;
305 _control_fe_values_ptr = NULL;
306 _nbr_state_fe_values_ptr = NULL;
307 _nbr_control_fe_values_ptr = NULL;
311 unsigned int _state_index;
312 unsigned int _control_index;
325 dealii::FEFaceValuesBase<dim>* _state_fe_values_ptr;
326 dealii::FEFaceValuesBase<dim>* _control_fe_values_ptr;
327 dealii::FEFaceValuesBase<dim>* _nbr_state_fe_values_ptr;
328 dealii::FEFaceValuesBase<dim>* _nbr_control_fe_values_ptr;
330 unsigned int _n_q_points_per_element;
331 unsigned int _n_dofs_per_element;
351 template<
typename VECTOR,
int dim>
378 template<
template<
int,
int>
class FE,
typename SPARSITYPATTERN,
int dopedim,
int dealdim>
380 UpdateFlags update_flags,
382 VECTOR, dopedim, dealdim>& sth,
384 typename dealii::MGDoFHandler<dim>::active_cell_iterator>& element,
385 const std::map<std::string,
const Vector<double>*> ¶m_values,
386 const std::map<std::string, const VECTOR*> &domain_values,
388 : fdcinternal::FaceDataContainerInternal<VECTOR, dim>(
389 param_values, domain_values, need_neighbour), _element(element), _state_fe_values(
390 sth.GetMapping(), (sth.GetFESystem(
"state")), quad,
391 update_flags), _control_fe_values(sth.GetMapping(),
392 (sth.GetFESystem(
"control")), quad, update_flags)
394 _state_index = sth.GetStateIndex();
395 if (_state_index == 1)
403 sth.GetMapping(), (sth.GetFESystem(
"control")), quad,
405 _control_fe_subface_values =
407 (sth.GetFESystem(
"control")), quad, update_flags);
409 this->PrivateConstructor(quad, update_flags, sth, need_neighbour);
429 template<
template<
int,
int>
class FE,
typename SPARSITYPATTERN>
431 UpdateFlags update_flags,
433 SPARSITYPATTERN, VECTOR, dim>& sth,
435 typename dealii::MGDoFHandler<dim>::active_cell_iterator>& element,
436 const std::map<std::string,
const Vector<double>*> ¶m_values,
437 const std::map<std::string, const VECTOR*> &domain_values,
439 : fdcinternal::FaceDataContainerInternal<VECTOR, dim>(
440 param_values, domain_values, need_neighbour), _element(element), _state_fe_values(
441 sth.GetMapping(), (sth.GetFESystem(
"state")), quad,
442 update_flags), _control_fe_values(sth.GetMapping(),
443 (sth.GetFESystem(
"state")), quad, update_flags)
445 _state_index = sth.GetStateIndex();
446 _control_index = element.size();
447 _n_q_points_per_element = quad.size();
448 _n_dofs_per_element = element[0]->get_fe().dofs_per_cell;
453 sth.GetMapping(), (sth.GetFESystem(
"state")), quad,
455 _control_fe_subface_values =
457 (sth.GetFESystem(
"state")), quad, update_flags);
459 this->PrivateConstructor(quad, update_flags, sth, need_neighbour);
464 if (_nbr_state_fe_values != NULL)
466 delete _nbr_state_fe_values;
468 if (_nbr_control_fe_values != NULL)
470 delete _nbr_control_fe_values;
472 if (_state_fe_subface_values != NULL)
474 delete _state_fe_subface_values;
476 if (_control_fe_subface_values != NULL)
478 delete _control_fe_subface_values;
490 ReInit(
unsigned int face_no);
503 ReInit(
unsigned int face_no,
unsigned int subface_no);
521 GetNDoFsPerElement()
const;
523 GetNbrNDoFsPerElement()
const;
527 GetNbrNQPoints()
const;
529 GetMaterialId()
const;
531 GetNbrMaterialId()
const;
533 GetNbrMaterialId(
unsigned int face)
const;
535 GetIsAtBoundary()
const;
537 GetElementDiameter()
const;
539 GetBoundaryIndicator()
const;
540 inline const FEFaceValuesBase<dim>&
541 GetFEFaceValuesState()
const;
542 inline const FEFaceValuesBase<dim>&
543 GetFEFaceValuesControl()
const;
545 inline const FEFaceValuesBase<dim>&
546 GetNbrFEFaceValuesState()
const;
547 inline const FEFaceValuesBase<dim>&
548 GetNbrFEFaceValuesControl()
const;
555 GetStateIndex()
const;
557 GetControlIndex()
const;
563 PrivateConstructor(
const Quadrature<dim - 1>& quad,
564 UpdateFlags update_flags, STH& sth,
bool need_neighbour)
566 _n_q_points_per_element = quad.size();
567 _n_dofs_per_element = _element[0]->get_fe().dofs_per_cell;
572 (sth.GetFESystem(
"state")), quad, update_flags);
574 (sth.GetFESystem(
"state")), quad, update_flags);
578 _nbr_state_fe_values = NULL;
579 _nbr_control_fe_values = NULL;
580 _state_fe_subface_values = NULL;
581 _control_fe_subface_values = NULL;
585 _state_fe_values_ptr = NULL;
586 _control_fe_values_ptr = NULL;
587 _nbr_state_fe_values_ptr = NULL;
588 _nbr_control_fe_values_ptr = NULL;
592 unsigned int _state_index;
593 unsigned int _control_index;
606 dealii::FEFaceValuesBase<dim>* _state_fe_values_ptr;
607 dealii::FEFaceValuesBase<dim>* _control_fe_values_ptr;
608 dealii::FEFaceValuesBase<dim>* _nbr_state_fe_values_ptr;
609 dealii::FEFaceValuesBase<dim>* _nbr_control_fe_values_ptr;
611 unsigned int _n_q_points_per_element;
612 unsigned int _n_dofs_per_element;
619 template<
typename VECTOR,
int dim>
645 template<
template<
int,
int>
class FE,
typename SPARSITYPATTERN,
int dopedim,
int dealdim>
647 const hp::QCollection<dim - 1>& q_collection,
648 UpdateFlags update_flags,
650 VECTOR, dopedim, dealdim>& sth,
653 const std::map<std::string, const Vector<double>*> ¶m_values,
654 const std::map<std::string, const VECTOR*> &domain_values,
655 bool need_neighbour) :
657 domain_values, need_neighbour), _element(element), _state_hp_fe_values(
658 (sth.GetFESystem(
"state")), q_collection, update_flags), _control_hp_fe_values(
659 (sth.GetFESystem(
"control")), q_collection, update_flags), _q_collection(
662 _state_index = sth.GetStateIndex();
663 if (_state_index == 1)
671 sth.GetMapping(), (sth.GetFESystem(
"control")), q_collection,
673 _control_hp_fe_subface_values =
675 (sth.GetFESystem(
"control")), q_collection, update_flags);
678 this->PrivateConstructor(q_collection, update_flags, sth,
702 template<
template<
int,
int>
class FE,
typename SPARSITYPATTERN,
int dealdim>
704 const hp::QCollection<dim - 1>& q_collection,
705 UpdateFlags update_flags,
707 VECTOR, dealdim>& sth,
710 const std::map<std::string, const Vector<double>*> ¶m_values,
711 const std::map<std::string, const VECTOR*> &domain_values,
712 bool need_neighbour) :
714 domain_values, need_neighbour), _element(element), _state_hp_fe_values(
715 (sth.GetFESystem(
"state")), q_collection, update_flags), _control_hp_fe_values(
716 (sth.GetFESystem(
"state")), q_collection, update_flags), _q_collection(
719 _state_index = sth.GetStateIndex();
720 if (_state_index == 1)
728 sth.GetMapping(), (sth.GetFESystem(
"state")), q_collection,
730 _control_hp_fe_subface_values =
732 (sth.GetFESystem(
"state")), q_collection, update_flags);
734 this->PrivateConstructor(q_collection, update_flags, sth,
743 delete _nbr_state_hp_fe_values;
744 _nbr_state_hp_fe_values = NULL;
745 delete _nbr_control_hp_fe_values;
746 _nbr_control_hp_fe_values = NULL;
747 delete _state_hp_fe_subface_values;
748 _state_hp_fe_subface_values = NULL;
749 delete _control_hp_fe_subface_values;
750 _control_hp_fe_subface_values = NULL;
762 ReInit(
unsigned int face_no);
775 ReInit(
unsigned int face_no,
unsigned int subface_no);
794 GetNDoFsPerElement()
const;
797 GetNbrNDoFsPerElement()
const;
801 GetNbrNQPoints()
const;
803 GetMaterialId()
const;
805 GetNbrMaterialId()
const;
807 GetNbrMaterialId(
unsigned int face)
const;
809 GetIsAtBoundary()
const;
811 GetBoundaryIndicator()
const;
813 GetElementDiameter()
const;
815 inline const FEFaceValuesBase<dim>&
816 GetFEFaceValuesState()
const;
817 inline const FEFaceValuesBase<dim>&
818 GetFEFaceValuesControl()
const;
820 inline const FEFaceValuesBase<dim>&
821 GetNbrFEFaceValuesState()
const;
822 inline const FEFaceValuesBase<dim>&
823 GetNbrFEFaceValuesControl()
const;
827 GetStateIndex()
const;
829 GetControlIndex()
const;
836 PrivateConstructor(
const hp::QCollection<dim - 1>& q_collection,
837 UpdateFlags update_flags, STH& sth,
bool need_neighbour)
842 sth.GetMapping(), (sth.GetFESystem(
"state")), q_collection,
845 dim>(sth.GetMapping(), (sth.GetFESystem(
"state")),
846 q_collection, update_flags);
851 _nbr_state_hp_fe_values = NULL;
852 _nbr_control_hp_fe_values = NULL;
853 _state_hp_fe_subface_values = NULL;
854 _control_hp_fe_subface_values = NULL;
857 _nbr_state_hp_fe_values_ptr = NULL;
858 _nbr_control_hp_fe_values_ptr = NULL;
859 _state_hp_fe_values_ptr = NULL;
860 _control_hp_fe_values_ptr = NULL;
864 unsigned int _state_index;
865 unsigned int _control_index;
878 const dealii::FEFaceValuesBase<dim>* _state_hp_fe_values_ptr;
879 const dealii::FEFaceValuesBase<dim>* _control_hp_fe_values_ptr;
880 const dealii::FEFaceValuesBase<dim>* _nbr_state_hp_fe_values_ptr;
881 const dealii::FEFaceValuesBase<dim>* _nbr_control_hp_fe_values_ptr;
883 const hp::QCollection<dim - 1>& _q_collection;
891 template<
int dim,
template<
int,
int>
class DH>
892 bool sanity_check(
const
895 unsigned int subface)
897 const auto neighbor_child =
898 _element->neighbor_child_on_subface(
902 if(neighbor_child->face(_element->neighbor_of_neighbor(face)) == _element->face(face)->child(subface))
908 bool sanity_check<1,dealii::hp::DoFHandler>(
const
917 bool sanity_check<1,dealii::DoFHandler>(
const
926 template<
typename VECTOR,
int dim>
929 unsigned int face_no)
931 this->SetFace(face_no);
932 _state_fe_values.reinit(_element[this->GetStateIndex()], face_no);
933 _state_fe_values_ptr = &_state_fe_values;
935 if (this->GetControlIndex() < _element.size())
937 _control_fe_values.reinit(_element[this->GetControlIndex()], face_no);
938 _control_fe_values_ptr = &_control_fe_values;
944 template<
typename VECTOR,
int dim>
947 unsigned int face_no,
unsigned int subface_no)
949 this->SetFace(face_no);
950 this->SetSubFace(subface_no);
951 _state_fe_subface_values->reinit(_element[this->GetStateIndex()], face_no,
953 _state_fe_values_ptr = _state_fe_subface_values;
955 if (this->GetControlIndex() < _element.size())
957 _control_fe_subface_values->reinit(_element[this->GetControlIndex()],
958 face_no, this->GetSubFace());
959 _control_fe_values_ptr = _control_fe_subface_values;
965 template<
typename VECTOR,
int dim>
969 Assert(this->NeedNeighbour(), ExcInternalError());
971 _element[this->GetStateIndex()]->neighbor_index(this->GetFace()) != -1,
972 TriaAccessorExceptions::ExcUnusedCellAsNeighbor())
974 if (_element[this->GetStateIndex()]->neighbor(this->GetFace())->has_children())
977 const auto neighbor_child =
978 _element[this->GetStateIndex()]->neighbor_child_on_subface(
979 this->GetFace(), this->GetSubFace());
983 Assert((sanity_check<dim, dealii::DoFHandler>(_element[this->GetStateIndex()],
985 this->GetSubFace()) ==
true), ExcInternalError());
986 Assert(neighbor_child->has_children() ==
false, ExcInternalError());
988 _nbr_state_fe_values->reinit(neighbor_child,
989 _element[this->GetStateIndex()]->neighbor_of_neighbor(
991 _nbr_state_fe_values_ptr = _nbr_state_fe_values;
994 if (this->GetControlIndex() < _element.size())
996 const auto control_neighbor_child =
997 _element[this->GetControlIndex()]->neighbor_child_on_subface(
998 this->GetFace(), this->GetSubFace());
1000 _nbr_control_fe_values->reinit(control_neighbor_child,
1001 _element[this->GetControlIndex()]->neighbor_of_neighbor(
1003 _nbr_control_fe_values_ptr = _nbr_control_fe_values;
1006 else if (_element[this->GetStateIndex()]->neighbor_is_coarser(
1011 _element[this->GetStateIndex()]->neighbor(this->GetFace())->level() == _element[this->GetStateIndex()]->level()-1,
1012 ExcInternalError());
1013 const auto neighbor = _element[this->GetStateIndex()]->neighbor(
1015 const std::pair<unsigned int, unsigned int> faceno_subfaceno =
1016 _element[this->GetStateIndex()]->neighbor_of_coarser_neighbor(
1018 const unsigned int neighbor_face_no = faceno_subfaceno.first,
1019 neighbor_subface_no = faceno_subfaceno.second;
1020 _state_fe_subface_values->reinit(neighbor, neighbor_face_no,
1021 neighbor_subface_no);
1022 _nbr_state_fe_values_ptr = _state_fe_subface_values;
1023 if (this->GetControlIndex() < _element.size())
1025 const auto control_neighbor =
1026 _element[this->GetControlIndex()]->neighbor(this->GetFace());
1027 const std::pair<unsigned int, unsigned int> control_faceno_subfaceno =
1028 _element[this->GetControlIndex()]->neighbor_of_coarser_neighbor(
1030 const unsigned int control_neighbor_face_no =
1031 control_faceno_subfaceno.first, control_neighbor_subface_no =
1032 control_faceno_subfaceno.second;
1033 _control_fe_subface_values->reinit(control_neighbor,
1034 control_neighbor_face_no, control_neighbor_subface_no);
1035 _nbr_control_fe_values_ptr = _control_fe_subface_values;
1041 const auto neighbor_state = _element[this->GetStateIndex()]->neighbor(
1044 Assert(neighbor_state->level() == _element[this->GetStateIndex()]->level(),
1045 ExcInternalError());
1046 _nbr_state_fe_values->reinit(neighbor_state,
1047 _element[this->GetStateIndex()]->neighbor_of_neighbor(
1049 _nbr_state_fe_values_ptr = _nbr_state_fe_values;
1052 if (this->GetControlIndex() < _element.size())
1054 _nbr_control_fe_values->reinit(
1055 _element[this->GetControlIndex()]->neighbor(this->GetFace()),
1056 _element[this->GetControlIndex()]->neighbor_of_neighbor(
1058 _nbr_control_fe_values_ptr = _nbr_control_fe_values;
1063 template<
typename VECTOR,
int dim>
1067 return _n_dofs_per_element;
1072 template<
typename VECTOR,
int dim>
1076 return _n_dofs_per_element;
1080 template<
typename VECTOR,
int dim>
1084 return _n_q_points_per_element;
1088 template<
typename VECTOR,
int dim>
1092 return _n_q_points_per_element;
1096 template<
typename VECTOR,
int dim>
1100 return _element[0]->material_id();
1104 template<
typename VECTOR,
int dim>
1108 return this->GetNbrMaterialId(this->GetFace());
1112 template<
typename VECTOR,
int dim>
1115 unsigned int face)
const
1117 if (_element[0]->neighbor_index(face) != -1)
1118 return _element[0]->neighbor(face)->material_id();
1120 throw DOpEException(
"There is no neighbor with number " + face,
1121 "FaceDataContainer::GetNbrMaterialId");
1125 template<
typename VECTOR,
int dim>
1129 return _element[0]->face(this->GetFace())->at_boundary();
1133 template<
typename VECTOR,
int dim>
1137 return _element[0]->face(this->GetFace())->diameter();
1142 template<
typename VECTOR,
int dim>
1146 return _element[0]->face(this->GetFace())->boundary_indicator();
1150 template<
typename VECTOR,
int dim>
1151 const FEFaceValuesBase<dim>&
1154 return *_state_fe_values_ptr;
1158 template<
typename VECTOR,
int dim>
1159 const FEFaceValuesBase<dim>&
1162 return *_control_fe_values_ptr;
1165 template<
typename VECTOR,
int dim>
1166 const FEFaceValuesBase<dim>&
1169 return *_nbr_state_fe_values_ptr;
1173 template<
typename VECTOR,
int dim>
1174 const FEFaceValuesBase<dim>&
1177 return *_nbr_control_fe_values_ptr;
1180 template<
typename VECTOR,
int dim>
1184 return _state_index;
1189 template<
typename VECTOR,
int dim>
1191 FaceDataContainer<dealii::DoFHandler, VECTOR, dim>::GetControlIndex()
const
1193 return _control_index;
1205 template<
typename VECTOR,
int dim>
1208 unsigned int face_no)
1210 this->SetFace(face_no);
1211 _state_fe_values.reinit(_element[this->GetStateIndex()], face_no);
1212 _state_fe_values_ptr = &_state_fe_values;
1214 if (this->GetControlIndex() < _element.size())
1216 _control_fe_values.reinit(_element[this->GetControlIndex()], face_no);
1217 _control_fe_values_ptr = &_control_fe_values;
1223 template<
typename VECTOR,
int dim>
1226 unsigned int face_no,
unsigned int subface_no)
1228 this->SetFace(face_no);
1229 this->SetSubFace(subface_no);
1230 _state_fe_subface_values->reinit(_element[this->GetStateIndex()], face_no,
1232 _state_fe_values_ptr = _state_fe_subface_values;
1234 if (this->GetControlIndex() < _element.size())
1236 _control_fe_subface_values->reinit(_element[this->GetControlIndex()],
1237 face_no, this->GetSubFace());
1238 _control_fe_values_ptr = _control_fe_subface_values;
1243 template<
typename VECTOR,
int dim>
1247 Assert(this->NeedNeighbour(), ExcInternalError());
1249 _element[this->GetStateIndex()]->neighbor_index(this->GetFace()) != -1,
1250 TriaAccessorExceptions::ExcUnusedCellAsNeighbor())
1252 if (_element[this->GetStateIndex()]->neighbor(this->GetFace())->has_children())
1255 const auto neighbor_child =
1256 _element[this->GetStateIndex()]->neighbor_child_on_subface(
1257 this->GetFace(), this->GetSubFace());
1261 Assert((sanity_check<dim, dealii::hp::DoFHandler>(_element[this->GetStateIndex()],
1263 this->GetSubFace()) ==
true), ExcInternalError());
1264 Assert(neighbor_child->has_children() ==
false, ExcInternalError());
1266 _nbr_state_fe_values->reinit(neighbor_child,
1267 _element[this->GetStateIndex()]->neighbor_of_neighbor(
1269 _nbr_state_fe_values_ptr = _nbr_state_fe_values;
1272 if (this->GetControlIndex() < _element.size())
1274 const auto control_neighbor_child =
1275 _element[this->GetControlIndex()]->neighbor_child_on_subface(
1276 this->GetFace(), this->GetSubFace());
1278 _nbr_control_fe_values->reinit(control_neighbor_child,
1279 _element[this->GetControlIndex()]->neighbor_of_neighbor(
1281 _nbr_control_fe_values_ptr = _nbr_control_fe_values;
1284 else if (_element[this->GetStateIndex()]->neighbor_is_coarser(
1289 _element[this->GetStateIndex()]->neighbor(this->GetFace())->level() == _element[this->GetStateIndex()]->level()-1,
1290 ExcInternalError());
1291 const auto neighbor = _element[this->GetStateIndex()]->neighbor(
1293 const std::pair<unsigned int, unsigned int> faceno_subfaceno =
1294 _element[this->GetStateIndex()]->neighbor_of_coarser_neighbor(
1296 const unsigned int neighbor_face_no = faceno_subfaceno.first,
1297 neighbor_subface_no = faceno_subfaceno.second;
1298 _state_fe_subface_values->reinit(neighbor, neighbor_face_no,
1299 neighbor_subface_no);
1300 _nbr_state_fe_values_ptr = _state_fe_subface_values;
1301 if (this->GetControlIndex() < _element.size())
1303 const auto control_neighbor =
1304 _element[this->GetControlIndex()]->neighbor(this->GetFace());
1305 const std::pair<unsigned int, unsigned int> control_faceno_subfaceno =
1306 _element[this->GetControlIndex()]->neighbor_of_coarser_neighbor(
1308 const unsigned int control_neighbor_face_no =
1309 control_faceno_subfaceno.first, control_neighbor_subface_no =
1310 control_faceno_subfaceno.second;
1311 _control_fe_subface_values->reinit(control_neighbor,
1312 control_neighbor_face_no, control_neighbor_subface_no);
1313 _nbr_control_fe_values_ptr = _control_fe_subface_values;
1319 const auto neighbor_state = _element[this->GetStateIndex()]->neighbor(
1322 Assert(neighbor_state->level() == _element[this->GetStateIndex()]->level(),
1323 ExcInternalError());
1324 _nbr_state_fe_values->reinit(neighbor_state,
1325 _element[this->GetStateIndex()]->neighbor_of_neighbor(
1327 _nbr_state_fe_values_ptr = _nbr_state_fe_values;
1330 if (this->GetControlIndex() < _element.size())
1332 _nbr_control_fe_values->reinit(
1333 _element[this->GetControlIndex()]->neighbor(this->GetFace()),
1334 _element[this->GetControlIndex()]->neighbor_of_neighbor(
1336 _nbr_control_fe_values_ptr = _nbr_control_fe_values;
1341 template<
typename VECTOR,
int dim>
1345 return _n_dofs_per_element;
1350 template<
typename VECTOR,
int dim>
1354 return _n_dofs_per_element;
1358 template<
typename VECTOR,
int dim>
1362 return _n_q_points_per_element;
1366 template<
typename VECTOR,
int dim>
1370 return _n_q_points_per_element;
1374 template<
typename VECTOR,
int dim>
1378 return _element[0]->material_id();
1382 template<
typename VECTOR,
int dim>
1386 return this->GetNbrMaterialId(this->GetFace());
1390 template<
typename VECTOR,
int dim>
1393 unsigned int face)
const
1395 if (_element[0]->neighbor_index(face) != -1)
1396 return _element[0]->neighbor(face)->material_id();
1398 throw DOpEException(
"There is no neighbor with number " + face,
1399 "FaceDataContainer::GetNbrMaterialId");
1403 template<
typename VECTOR,
int dim>
1407 return _element[0]->face(this->GetFace())->at_boundary();
1411 template<
typename VECTOR,
int dim>
1415 return _element[0]->face(this->GetFace())->diameter();
1420 template<
typename VECTOR,
int dim>
1424 return _element[0]->face(this->GetFace())->boundary_indicator();
1428 template<
typename VECTOR,
int dim>
1429 const FEFaceValuesBase<dim>&
1432 return *_state_fe_values_ptr;
1436 template<
typename VECTOR,
int dim>
1437 const FEFaceValuesBase<dim>&
1440 return *_control_fe_values_ptr;
1443 template<
typename VECTOR,
int dim>
1444 const FEFaceValuesBase<dim>&
1447 return *_nbr_state_fe_values_ptr;
1451 template<
typename VECTOR,
int dim>
1452 const FEFaceValuesBase<dim>&
1455 return *_nbr_control_fe_values_ptr;
1458 template<
typename VECTOR,
int dim>
1462 return _state_index;
1467 template<
typename VECTOR,
int dim>
1469 FaceDataContainer<dealii::MGDoFHandler, VECTOR, dim>::GetControlIndex()
const
1471 return _control_index;
1485 template<
typename VECTOR,
int dim>
1488 unsigned int face_no)
1490 this->SetFace(face_no);
1491 _state_hp_fe_values.reinit(_element[this->GetStateIndex()], face_no);
1492 _state_hp_fe_values_ptr = &_state_hp_fe_values.get_present_fe_values();
1494 if (this->GetControlIndex() < _element.size())
1496 _control_hp_fe_values.reinit(_element[this->GetControlIndex()], face_no);
1497 _control_hp_fe_values_ptr =
1498 &_control_hp_fe_values.get_present_fe_values();
1504 template<
typename VECTOR,
int dim>
1507 unsigned int face_no,
unsigned int subface_no)
1509 this->SetFace(face_no);
1510 this->SetSubFace(subface_no);
1511 _state_hp_fe_subface_values->reinit(_element[this->GetStateIndex()],
1512 this->GetFace(), this->GetSubFace());
1513 _state_hp_fe_values_ptr =
1514 &_state_hp_fe_subface_values->get_present_fe_values();
1516 if (this->GetControlIndex() < _element.size())
1518 _control_hp_fe_subface_values->reinit(_element[this->GetControlIndex()],
1519 this->GetFace(), this->GetSubFace());
1520 _control_hp_fe_values_ptr =
1521 &_control_hp_fe_subface_values->get_present_fe_values();
1526 template<
typename VECTOR,
int dim>
1530 Assert(this->NeedNeighbour(), ExcInternalError());
1532 _element[this->GetStateIndex()]->neighbor_index(this->GetFace()) != -1,
1533 TriaAccessorExceptions::ExcUnusedCellAsNeighbor())
1535 if (_element[this->GetStateIndex()]->neighbor(this->GetFace())->has_children())
1538 const auto neighbor_child =
1539 _element[this->GetStateIndex()]->neighbor_child_on_subface(
1540 this->GetFace(), this->GetSubFace());
1545 neighbor_child->face(_element[this->GetStateIndex()]->neighbor_of_neighbor(this->GetFace())) == _element[this->GetStateIndex()]->face(this->GetFace())->child(this->GetSubFace()),
1546 ExcInternalError());
1547 Assert(neighbor_child->has_children() ==
false, ExcInternalError());
1549 _nbr_state_hp_fe_values->reinit(neighbor_child,
1550 _element[this->GetStateIndex()]->neighbor_of_neighbor(
1552 _nbr_state_hp_fe_values_ptr =
1553 &_nbr_state_hp_fe_values->get_present_fe_values();
1556 if (this->GetControlIndex() < _element.size())
1558 const auto control_neighbor_child =
1559 _element[this->GetControlIndex()]->neighbor_child_on_subface(
1560 this->GetFace(), this->GetSubFace());
1562 _nbr_control_hp_fe_values->reinit(control_neighbor_child,
1563 _element[this->GetControlIndex()]->neighbor_of_neighbor(
1565 _nbr_control_hp_fe_values_ptr =
1566 &_nbr_control_hp_fe_values->get_present_fe_values();
1569 else if (_element[this->GetStateIndex()]->neighbor_is_coarser(
1574 _element[this->GetStateIndex()]->neighbor(this->GetFace())->level() == _element[this->GetStateIndex()]->level()-1,
1575 ExcInternalError());
1576 const auto neighbor = _element[this->GetStateIndex()]->neighbor(
1578 const std::pair<unsigned int, unsigned int> faceno_subfaceno =
1579 _element[this->GetStateIndex()]->neighbor_of_coarser_neighbor(
1581 const unsigned int neighbor_face_no = faceno_subfaceno.first,
1582 neighbor_subface_no = faceno_subfaceno.second;
1583 _state_hp_fe_subface_values->reinit(neighbor, neighbor_face_no,
1584 neighbor_subface_no);
1585 _nbr_state_hp_fe_values_ptr =
1586 &_state_hp_fe_subface_values->get_present_fe_values();
1587 if (this->GetControlIndex() < _element.size())
1589 const auto control_neighbor =
1590 _element[this->GetControlIndex()]->neighbor(this->GetFace());
1591 const std::pair<unsigned int, unsigned int> control_faceno_subfaceno =
1592 _element[this->GetControlIndex()]->neighbor_of_coarser_neighbor(
1594 const unsigned int control_neighbor_face_no =
1595 control_faceno_subfaceno.first, control_neighbor_subface_no =
1596 control_faceno_subfaceno.second;
1598 _control_hp_fe_subface_values->reinit(control_neighbor,
1599 control_neighbor_face_no, control_neighbor_subface_no);
1600 _nbr_control_hp_fe_values_ptr =
1601 &_control_hp_fe_subface_values->get_present_fe_values();
1606 const auto neighbor_state = _element[this->GetStateIndex()]->neighbor(
1609 Assert(neighbor_state->level() == _element[this->GetStateIndex()]->level(),
1610 ExcInternalError());
1611 _nbr_state_hp_fe_values->reinit(neighbor_state,
1612 _element[this->GetStateIndex()]->neighbor_of_neighbor(
1614 _nbr_state_hp_fe_values_ptr =
1615 &_nbr_state_hp_fe_values->get_present_fe_values();
1618 if (this->GetControlIndex() < _element.size())
1620 _nbr_control_hp_fe_values->reinit(
1621 _element[this->GetControlIndex()]->neighbor(this->GetFace()),
1622 _element[this->GetControlIndex()]->neighbor_of_neighbor(
1624 _nbr_control_hp_fe_values_ptr =
1625 &_nbr_control_hp_fe_values->get_present_fe_values();
1632 template<
typename VECTOR,
int dim>
1636 return _element[0]->get_fe().dofs_per_cell;
1641 template<
typename VECTOR,
int dim>
1645 if (_element[0]->neighbor_index(this->GetFace()) != -1)
1646 return _element[0]->neighbor(this->GetFace())->get_fe().dofs_per_cell;
1649 "There is no neighbor with number" + this->GetFace(),
1650 "HpFaceDataContainer::GetNbrNDoFsPerElement");
1653 template<
typename VECTOR,
int dim>
1657 return _q_collection[_element[0]->active_fe_index()].size();
1661 template<
typename VECTOR,
int dim>
1665 if (_element[0]->neighbor_index(this->GetFace()) != -1)
1666 return _q_collection[_element[0]->neighbor(this->GetFace())->active_fe_index()].size();
1669 "There is no neighbor with number" + this->GetFace(),
1670 "HpFaceDataContainer::GetNbrNQPoints");
1674 template<
typename VECTOR,
int dim>
1678 return _element[0]->material_id();
1681 template<
typename VECTOR,
int dim>
1685 return this->GetNbrMaterialId(this->GetFace());
1690 template<
typename VECTOR,
int dim>
1693 unsigned int face)
const
1695 if (_element[0]->neighbor_index(face) != -1)
1696 return _element[0]->neighbor(face)->material_id();
1698 throw DOpEException(
"There is no neighbor with number" + face,
1699 "HpFaceDataContainer::GetNbrMaterialId");
1704 template<
typename VECTOR,
int dim>
1708 return _element[0]->face(this->GetFace())->diameter();
1713 template<
typename VECTOR,
int dim>
1717 return _element[0]->face(this->GetFace())->at_boundary();
1721 template<
typename VECTOR,
int dim>
1725 return _element[0]->face(this->GetFace())->boundary_indicator();
1729 template<
typename VECTOR,
int dim>
1730 const FEFaceValuesBase<dim>&
1733 return *_state_hp_fe_values_ptr;
1736 template<
typename VECTOR,
int dim>
1737 const FEFaceValuesBase<dim>&
1740 return *_control_hp_fe_values_ptr;
1743 template<
typename VECTOR,
int dim>
1744 const FEFaceValuesBase<dim>&
1747 return *_nbr_state_hp_fe_values_ptr;
1750 template<
typename VECTOR,
int dim>
1751 const FEFaceValuesBase<dim>&
1754 return *_nbr_control_hp_fe_values_ptr;
1758 template<
typename VECTOR,
int dim>
1762 return _state_index;
1765 template<
typename VECTOR,
int dim>
1767 FaceDataContainer<dealii::hp::DoFHandler, VECTOR, dim>::GetControlIndex()
const
1769 return _control_index;
Definition: fevalues_wrapper.h:136
FaceDataContainer(const Quadrature< dim-1 > &quad, UpdateFlags update_flags, StateSpaceTimeHandler< FE, dealii::MGDoFHandler, SPARSITYPATTERN, VECTOR, dim > &sth, const std::vector< typename dealii::MGDoFHandler< dim >::active_cell_iterator > &element, const std::map< std::string, const Vector< double > * > ¶m_values, const std::map< std::string, const VECTOR * > &domain_values, bool need_neighbour)
Definition: facedatacontainer.h:430
~FaceDataContainer()
Definition: facedatacontainer.h:181
Definition: facedatacontainer_internal.h:46
FaceDataContainer(const Quadrature< dim-1 > &quad, UpdateFlags update_flags, StateSpaceTimeHandler< FE, dealii::DoFHandler, SPARSITYPATTERN, VECTOR, dim > &sth, const std::vector< typename DOpEWrapper::DoFHandler< dim, dealii::DoFHandler >::active_cell_iterator > &element, const std::map< std::string, const Vector< double > * > ¶m_values, const std::map< std::string, const VECTOR * > &domain_values, bool need_neighbour)
Definition: facedatacontainer.h:149
FaceDataContainer(const Quadrature< dim-1 > &quad, UpdateFlags update_flags, SpaceTimeHandler< FE, dealii::MGDoFHandler, SPARSITYPATTERN, VECTOR, dopedim, dealdim > &sth, const std::vector< typename dealii::MGDoFHandler< dim >::active_cell_iterator > &element, const std::map< std::string, const Vector< double > * > ¶m_values, const std::map< std::string, const VECTOR * > &domain_values, bool need_neighbour)
Definition: facedatacontainer.h:379
FaceDataContainer()
Definition: facedatacontainer.h:52
Definition: fevalues_wrapper.h:218
Definition: facedatacontainer.h:48
Definition: dofhandler_wrapper.h:51
Definition: statespacetimehandler.h:59
Definition: spacetimehandler.h:71
Definition: fevalues_wrapper.h:259
Definition: facedatacontainer.h:620
FaceDataContainer(const Quadrature< dim-1 > &quad, UpdateFlags update_flags, SpaceTimeHandler< FE, dealii::DoFHandler, SPARSITYPATTERN, VECTOR, dopedim, dealdim > &sth, const std::vector< typename DOpEWrapper::DoFHandler< dim, dealii::DoFHandler >::active_cell_iterator > &element, const std::map< std::string, const Vector< double > * > ¶m_values, const std::map< std::string, const VECTOR * > &domain_values, bool need_neighbour)
Definition: facedatacontainer.h:98
Definition: fevalues_wrapper.h:88
Definition: dopeexception.h:35
~FaceDataContainer()
Definition: facedatacontainer.h:462