DOpE
facedatacontainer.h
Go to the documentation of this file.
1 
24 #ifndef FACEDATACONTAINER_H_
25 #define FACEDATACONTAINER_H_
26 
27 #include "spacetimehandler.h"
28 #include "statespacetimehandler.h"
29 #include "fevalues_wrapper.h"
30 #include "dopeexception.h"
32 
33 #include <dofs/dof_handler.h>
34 #include <multigrid/mg_dof_handler.h>
35 #include <hp/dof_handler.h>
36 
37 using namespace dealii;
38 
39 namespace DOpE
40 {
47  template<template<int, int> class DH, typename VECTOR, int dim>
49  VECTOR, dim>
50  {
51  public:
53  {
54  throw(DOpEException(
55  "Dummy class, this constructor should never get called.",
56  "ElementDataContainer<dealii::DoFHandler , VECTOR, dim>::ElementDataContainer"));
57  }
58  ;
59  };
60 
70  template<typename VECTOR, int dim>
71  class FaceDataContainer<dealii::DoFHandler, VECTOR, dim> : public fdcinternal::FaceDataContainerInternal<
72  VECTOR, dim>
73  {
74 
75  public:
97  template<template<int, int> class FE, typename SPARSITYPATTERN, int dopedim, int dealdim>
98  FaceDataContainer(const Quadrature<dim - 1>& quad,
99  UpdateFlags update_flags,
100  SpaceTimeHandler<FE, dealii::DoFHandler, SPARSITYPATTERN, VECTOR,
101  dopedim, dealdim>& sth,
102  const std::vector<
104  const std::map<std::string, const Vector<double>*> &param_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)
112  {
113  _state_index = sth.GetStateIndex();
114  if (_state_index == 1)
115  _control_index = 0;
116  else
117  _control_index = 1;
118 
119  if (need_neighbour) //so we need FEFAcevalues etc. for the neighbour too.
120  {
121  _nbr_control_fe_values = new DOpEWrapper::FEFaceValues<dim>(
122  sth.GetMapping(), (sth.GetFESystem("control")), quad,
123  update_flags);
124  _control_fe_subface_values =
125  new DOpEWrapper::FESubfaceValues<dim>(sth.GetMapping(),
126  (sth.GetFESystem("control")), quad, update_flags);
127  }
128  this->PrivateConstructor(quad, update_flags, sth, need_neighbour);
129  }
148  template<template<int, int> class FE, typename SPARSITYPATTERN>
149  FaceDataContainer(const Quadrature<dim - 1>& quad,
150  UpdateFlags update_flags,
151  StateSpaceTimeHandler<FE, dealii::DoFHandler, SPARSITYPATTERN,
152  VECTOR, dim>& sth,
153  const std::vector<
155  const std::map<std::string, const Vector<double>*> &param_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)
163  {
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;
168 
169  if (need_neighbour) //so we need FEFAcevalues for the neighbour too.
170  {
171  _nbr_control_fe_values = new DOpEWrapper::FEFaceValues<dim>(
172  sth.GetMapping(), (sth.GetFESystem("state")), quad,
173  update_flags);
174  _control_fe_subface_values =
175  new DOpEWrapper::FESubfaceValues<dim>(sth.GetMapping(),
176  (sth.GetFESystem("state")), quad, update_flags);
177  }
178  this->PrivateConstructor(quad, update_flags, sth, need_neighbour);
179  }
180 
182  {
183  if (_nbr_state_fe_values != NULL)
184  {
185  delete _nbr_state_fe_values;
186  }
187  if (_nbr_control_fe_values != NULL)
188  {
189  delete _nbr_control_fe_values;
190  }
191  if (_state_fe_subface_values != NULL)
192  {
193  delete _state_fe_subface_values;
194  }
195  if (_control_fe_subface_values != NULL)
196  {
197  delete _control_fe_subface_values;
198  }
199  }
200  /*********************************************/
201  /*
202  * This function reinitializes the FEFaceValues on the actual face. Should
203  * be called prior to any of the get-functions.
204  *
205  * @param face_no The 'local number' (i.e. from the perspective of the actual element) of the
206  * actual face.
207  */
208  inline void
209  ReInit(unsigned int face_no);
210 
211  /*********************************************/
212  /*
213  * This function reinits the FESubfaceValues on the actual subface. Should
214  * be called prior to any of the get-functions.
215  *
216  * @param face_no The 'local number' (i.e. from the perspective of the actual element) of the
217  * actual face.
218  * @param subface_no The 'local number' (i.e. from the perspective of the actual element) of the
219  * actual subface.
220  */
221  inline void
222  ReInit(unsigned int face_no, unsigned int subface_no);
223 
224  /*********************************************/
225  /*
226  * This function reinitializes the FE(Sub)FaceValues on the neighbor_element.
227  * This should be called prior to any of the get nbr_functions.
228  * Assumes that ReInit is called prior to this function.
229  */
230  inline void
231  ReInitNbr();
232 
233  /*********************************************/
239  inline unsigned int
240  GetNDoFsPerElement() const;
241  inline unsigned int
242  GetNbrNDoFsPerElement() const;
243  inline unsigned int
244  GetNQPoints() const;
245  inline unsigned int
246  GetNbrNQPoints() const;
247  inline unsigned int
248  GetMaterialId() const;
249  inline unsigned int
250  GetNbrMaterialId() const;
251  inline unsigned int
252  GetNbrMaterialId(unsigned int face) const;
253  inline bool
254  GetIsAtBoundary() const;
255  inline double
256  GetElementDiameter() const;
257  inline unsigned int
258  GetBoundaryIndicator() const;
259  inline const FEFaceValuesBase<dim>&
260  GetFEFaceValuesState() const;
261  inline const FEFaceValuesBase<dim>&
262  GetFEFaceValuesControl() const;
263 
264  inline const FEFaceValuesBase<dim>&
265  GetNbrFEFaceValuesState() const;
266  inline const FEFaceValuesBase<dim>&
267  GetNbrFEFaceValuesControl() const;
268 
269  private:
270  /*
271  * Helper Functions
272  */
273  unsigned int
274  GetStateIndex() const;
275  unsigned int
276  GetControlIndex() const;
280  template<class STH>
281  void
282  PrivateConstructor(const Quadrature<dim - 1>& quad,
283  UpdateFlags update_flags, STH& sth, bool need_neighbour)
284  {
285  _n_q_points_per_element = quad.size();
286  _n_dofs_per_element = _element[0]->get_fe().dofs_per_cell;
287 
288  if (need_neighbour) //so we need FEFAcevalues etc. for the neighbour too.
289  {
290  _nbr_state_fe_values = new DOpEWrapper::FEFaceValues<dim>(
291  (sth.GetFESystem("state")), quad, update_flags);
292  _state_fe_subface_values = new DOpEWrapper::FESubfaceValues<dim>(
293  (sth.GetFESystem("state")), quad, update_flags);
294  }
295  else
296  {
297  _nbr_state_fe_values = NULL;
298  _nbr_control_fe_values = NULL;
299  _state_fe_subface_values = NULL;
300  _control_fe_subface_values = NULL;
301  }
302  // These will point to the object (i.e. FaceValues or SubfaceValues) we actually use.
303  // With this, we have the same interface to the user independently of the type (i.e. face or subface)
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;
308  }
309  /***********************************************************/
310  //"global" member data, part of every instantiation
311  unsigned int _state_index;
312  unsigned int _control_index;
313 
314  const std::vector<
316  DOpEWrapper::FEFaceValues<dim> _state_fe_values;
317  DOpEWrapper::FEFaceValues<dim> _control_fe_values;
318 
319  DOpEWrapper::FEFaceValues<dim>* _nbr_state_fe_values;
320  DOpEWrapper::FEFaceValues<dim>* _nbr_control_fe_values;
321 
322  DOpEWrapper::FESubfaceValues<dim>* _state_fe_subface_values;
323  DOpEWrapper::FESubfaceValues<dim>* _control_fe_subface_values;
324 
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;
329 
330  unsigned int _n_q_points_per_element;
331  unsigned int _n_dofs_per_element;
332  };
333 
334 
335 
336 
337 
338  /****************************************************/
339  /* MGDofHandler */
340 
341 
351  template<typename VECTOR, int dim>
352  class FaceDataContainer<dealii::MGDoFHandler, VECTOR, dim> : public fdcinternal::FaceDataContainerInternal<
353  VECTOR, dim>
354  {
355 
356  public:
378  template<template<int, int> class FE, typename SPARSITYPATTERN, int dopedim, int dealdim>
379  FaceDataContainer(const Quadrature<dim - 1>& quad,
380  UpdateFlags update_flags,
381  SpaceTimeHandler<FE, dealii::MGDoFHandler, SPARSITYPATTERN,
382  VECTOR, dopedim, dealdim>& sth,
383  const std::vector<
384  typename dealii::MGDoFHandler<dim>::active_cell_iterator>& element,
385  const std::map<std::string, const Vector<double>*> &param_values,
386  const std::map<std::string, const VECTOR*> &domain_values,
387  bool need_neighbour)
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)
393  {
394  _state_index = sth.GetStateIndex();
395  if (_state_index == 1)
396  _control_index = 0;
397  else
398  _control_index = 1;
399 
400  if (need_neighbour) //so we need FEFAcevalues etc. for the neighbour too.
401  {
402  _nbr_control_fe_values = new DOpEWrapper::FEFaceValues<dim>(
403  sth.GetMapping(), (sth.GetFESystem("control")), quad,
404  update_flags);
405  _control_fe_subface_values =
406  new DOpEWrapper::FESubfaceValues<dim>(sth.GetMapping(),
407  (sth.GetFESystem("control")), quad, update_flags);
408  }
409  this->PrivateConstructor(quad, update_flags, sth, need_neighbour);
410  }
429  template<template<int, int> class FE, typename SPARSITYPATTERN>
430  FaceDataContainer(const Quadrature<dim - 1>& quad,
431  UpdateFlags update_flags,
432  StateSpaceTimeHandler<FE, dealii::MGDoFHandler,
433  SPARSITYPATTERN, VECTOR, dim>& sth,
434  const std::vector<
435  typename dealii::MGDoFHandler<dim>::active_cell_iterator>& element,
436  const std::map<std::string, const Vector<double>*> &param_values,
437  const std::map<std::string, const VECTOR*> &domain_values,
438  bool need_neighbour)
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)
444  {
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;
449 
450  if (need_neighbour) //so we need FEFAcevalues for the neighbour too.
451  {
452  _nbr_control_fe_values = new DOpEWrapper::FEFaceValues<dim>(
453  sth.GetMapping(), (sth.GetFESystem("state")), quad,
454  update_flags);
455  _control_fe_subface_values =
456  new DOpEWrapper::FESubfaceValues<dim>(sth.GetMapping(),
457  (sth.GetFESystem("state")), quad, update_flags);
458  }
459  this->PrivateConstructor(quad, update_flags, sth, need_neighbour);
460  }
461 
463  {
464  if (_nbr_state_fe_values != NULL)
465  {
466  delete _nbr_state_fe_values;
467  }
468  if (_nbr_control_fe_values != NULL)
469  {
470  delete _nbr_control_fe_values;
471  }
472  if (_state_fe_subface_values != NULL)
473  {
474  delete _state_fe_subface_values;
475  }
476  if (_control_fe_subface_values != NULL)
477  {
478  delete _control_fe_subface_values;
479  }
480  }
481  /*********************************************/
482  /*
483  * This function reinitializes the FEFaceValues on the actual face. Should
484  * be called prior to any of the get-functions.
485  *
486  * @param face_no The 'local number' (i.e. from the perspective of the actual element) of the
487  * actual face.
488  */
489  inline void
490  ReInit(unsigned int face_no);
491 
492  /*********************************************/
493  /*
494  * This function reinits the FESubfaceValues on the actual subface. Should
495  * be called prior to any of the get-functions.
496  *
497  * @param face_no The 'local number' (i.e. from the perspective of the actual element) of the
498  * actual face.
499  * @param subface_no The 'local number' (i.e. from the perspective of the actual element) of the
500  * actual subface.
501  */
502  inline void
503  ReInit(unsigned int face_no, unsigned int subface_no);
504 
505  /*********************************************/
506  /*
507  * This function reinitializes the FE(Sub)FaceValues on the neighbor_element.
508  * This should be called prior to any of the get nbr_functions.
509  * Assumes that ReInit is called prior to this function.
510  */
511  inline void
512  ReInitNbr();
513 
514  /*********************************************/
520  inline unsigned int
521  GetNDoFsPerElement() const;
522  inline unsigned int
523  GetNbrNDoFsPerElement() const;
524  inline unsigned int
525  GetNQPoints() const;
526  inline unsigned int
527  GetNbrNQPoints() const;
528  inline unsigned int
529  GetMaterialId() const;
530  inline unsigned int
531  GetNbrMaterialId() const;
532  inline unsigned int
533  GetNbrMaterialId(unsigned int face) const;
534  inline bool
535  GetIsAtBoundary() const;
536  inline double
537  GetElementDiameter() const;
538  inline unsigned int
539  GetBoundaryIndicator() const;
540  inline const FEFaceValuesBase<dim>&
541  GetFEFaceValuesState() const;
542  inline const FEFaceValuesBase<dim>&
543  GetFEFaceValuesControl() const;
544 
545  inline const FEFaceValuesBase<dim>&
546  GetNbrFEFaceValuesState() const;
547  inline const FEFaceValuesBase<dim>&
548  GetNbrFEFaceValuesControl() const;
549 
550  private:
551  /*
552  * Helper Functions
553  */
554  unsigned int
555  GetStateIndex() const;
556  unsigned int
557  GetControlIndex() const;
561  template<class STH>
562  void
563  PrivateConstructor(const Quadrature<dim - 1>& quad,
564  UpdateFlags update_flags, STH& sth, bool need_neighbour)
565  {
566  _n_q_points_per_element = quad.size();
567  _n_dofs_per_element = _element[0]->get_fe().dofs_per_cell;
568 
569  if (need_neighbour) //so we need FEFAcevalues etc. for the neighbour too.
570  {
571  _nbr_state_fe_values = new DOpEWrapper::FEFaceValues<dim>(
572  (sth.GetFESystem("state")), quad, update_flags);
573  _state_fe_subface_values = new DOpEWrapper::FESubfaceValues<dim>(
574  (sth.GetFESystem("state")), quad, update_flags);
575  }
576  else
577  {
578  _nbr_state_fe_values = NULL;
579  _nbr_control_fe_values = NULL;
580  _state_fe_subface_values = NULL;
581  _control_fe_subface_values = NULL;
582  }
583  // These will point to the object (i.e. FaceValues or SubfaceValues) we actually use.
584  // With this, we have the same interface to the user independently of the type (i.e. face or subface)
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;
589  }
590  /***********************************************************/
591  //"global" member data, part of every instantiation
592  unsigned int _state_index;
593  unsigned int _control_index;
594 
595  const std::vector<
597  DOpEWrapper::FEFaceValues<dim> _state_fe_values;
598  DOpEWrapper::FEFaceValues<dim> _control_fe_values;
599 
600  DOpEWrapper::FEFaceValues<dim>* _nbr_state_fe_values;
601  DOpEWrapper::FEFaceValues<dim>* _nbr_control_fe_values;
602 
603  DOpEWrapper::FESubfaceValues<dim>* _state_fe_subface_values;
604  DOpEWrapper::FESubfaceValues<dim>* _control_fe_subface_values;
605 
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;
610 
611  unsigned int _n_q_points_per_element;
612  unsigned int _n_dofs_per_element;
613  };
614 
615  /* MGDofHandler */
616  /****************************************************/
617 
618 
619  template<typename VECTOR, int dim>
620  class FaceDataContainer<dealii::hp::DoFHandler, VECTOR, dim> : public fdcinternal::FaceDataContainerInternal<
621  VECTOR, dim>
622  {
623 
624  public:
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,
649  SpaceTimeHandler<FE, dealii::hp::DoFHandler, SPARSITYPATTERN,
650  VECTOR, dopedim, dealdim>& sth,
651  const std::vector<
653  const std::map<std::string, const Vector<double>*> &param_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(
660  q_collection)
661  {
662  _state_index = sth.GetStateIndex();
663  if (_state_index == 1)
664  _control_index = 0;
665  else
666  _control_index = 1;
667 
668  if (need_neighbour) //so we need FEFAcevalues for the neighbour too.
669  {
670  _nbr_control_hp_fe_values = new DOpEWrapper::HpFEFaceValues<dim>(
671  sth.GetMapping(), (sth.GetFESystem("control")), q_collection,
672  update_flags);
673  _control_hp_fe_subface_values =
674  new DOpEWrapper::HpFESubfaceValues<dim>(sth.GetMapping(),
675  (sth.GetFESystem("control")), q_collection, update_flags);
676 
677  }
678  this->PrivateConstructor(q_collection, update_flags, sth,
679  need_neighbour);
680  }
681 
702  template<template<int, int> class FE, typename SPARSITYPATTERN, int dealdim>
704  const hp::QCollection<dim - 1>& q_collection,
705  UpdateFlags update_flags,
706  StateSpaceTimeHandler<FE, dealii::hp::DoFHandler, SPARSITYPATTERN,
707  VECTOR, dealdim>& sth,
708  const std::vector<
710  const std::map<std::string, const Vector<double>*> &param_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(
717  q_collection)
718  {
719  _state_index = sth.GetStateIndex();
720  if (_state_index == 1)
721  _control_index = 0;
722  else
723  _control_index = 1;
724 
725  if (need_neighbour)
726  {
727  _nbr_control_hp_fe_values = new DOpEWrapper::HpFEFaceValues<dim>(
728  sth.GetMapping(), (sth.GetFESystem("state")), q_collection,
729  update_flags);
730  _control_hp_fe_subface_values =
731  new DOpEWrapper::HpFESubfaceValues<dim>(sth.GetMapping(),
732  (sth.GetFESystem("state")), q_collection, update_flags);
733  }
734  this->PrivateConstructor(q_collection, update_flags, sth,
735  need_neighbour);
736  }
737 
742  {
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;
751  }
752 
753  /*********************************************/
754  /*
755  * This function reinits the FEFaceValues on the actual face. Should
756  * be called prior to any of the get-functions.
757  *
758  * @param face_no The 'local number' (i.e. from the perspective of the actual element) of the
759  * actual face.
760  */
761  inline void
762  ReInit(unsigned int face_no);
763 
764  /*********************************************/
765  /*
766  * This function reinits the FESubfaceValues on the actual subface. Should
767  * be called prior to any of the get-functions.
768  *
769  * @param face_no The 'local number' (i.e. from the perspective of the actual element) of the
770  * actual face.
771  * @param subface_no The 'local number' (i.e. from the perspective of the actual element) of the
772  * actual subface.
773  */
774  inline void
775  ReInit(unsigned int face_no, unsigned int subface_no);
776 
777  /*********************************************/
778  /*
779  * This function reinits the FEFaceValues on the neighbor_element for the
780  * case that the neighbor-element coarser or as fine as the
781  * (previously set!) actual element. This should be called prior
782  * to any of the get nbr_functions.
783  * Assumes that ReInit is called prior to this function.
784  */
785  inline void
786  ReInitNbr();
787 
788  /*********************************************/
793  inline unsigned int
794  GetNDoFsPerElement() const;
795 
796  inline unsigned int
797  GetNbrNDoFsPerElement() const;
798  inline unsigned int
799  GetNQPoints() const;
800  inline unsigned int
801  GetNbrNQPoints() const;
802  inline unsigned int
803  GetMaterialId() const;
804  inline unsigned int
805  GetNbrMaterialId() const;
806  inline unsigned int
807  GetNbrMaterialId(unsigned int face) const;
808  inline bool
809  GetIsAtBoundary() const;
810  inline unsigned int
811  GetBoundaryIndicator() const;
812  inline double
813  GetElementDiameter() const;
814 
815  inline const FEFaceValuesBase<dim>&
816  GetFEFaceValuesState() const;
817  inline const FEFaceValuesBase<dim>&
818  GetFEFaceValuesControl() const;
819 
820  inline const FEFaceValuesBase<dim>&
821  GetNbrFEFaceValuesState() const;
822  inline const FEFaceValuesBase<dim>&
823  GetNbrFEFaceValuesControl() const;
824 
825  private:
826  inline unsigned int
827  GetStateIndex() const;
828  inline unsigned int
829  GetControlIndex() const;
830 
834  template<class STH>
835  void
836  PrivateConstructor(const hp::QCollection<dim - 1>& q_collection,
837  UpdateFlags update_flags, STH& sth, bool need_neighbour)
838  {
839  if (need_neighbour) //so we need FEFAcevalues for the neighbour too.
840  {
841  _nbr_state_hp_fe_values = new DOpEWrapper::HpFEFaceValues<dim>(
842  sth.GetMapping(), (sth.GetFESystem("state")), q_collection,
843  update_flags);
844  _state_hp_fe_subface_values = new DOpEWrapper::HpFESubfaceValues<
845  dim>(sth.GetMapping(), (sth.GetFESystem("state")),
846  q_collection, update_flags);
847 
848  }
849  else
850  {
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;
855  }
856 
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;
861  }
862  /***********************************************************/
863  //"global" member data, part of every instantiation
864  unsigned int _state_index;
865  unsigned int _control_index;
866  const std::vector<
868 
869  DOpEWrapper::HpFEFaceValues<dim> _state_hp_fe_values;
870  DOpEWrapper::HpFEFaceValues<dim> _control_hp_fe_values;
871 
872  DOpEWrapper::HpFEFaceValues<dim>* _nbr_state_hp_fe_values;
873  DOpEWrapper::HpFEFaceValues<dim>* _nbr_control_hp_fe_values;
874 
875  DOpEWrapper::HpFESubfaceValues<dim>* _state_hp_fe_subface_values;
876  DOpEWrapper::HpFESubfaceValues<dim>* _control_hp_fe_subface_values;
877 
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;
882 
883  const hp::QCollection<dim - 1>& _q_collection;
884  };
885 
886  /***********************************************************************/
887  /************************IMPLEMENTATION*for*DoFHandler*********************************/
888  /***********************************************************************/
889 
890  namespace {
891  template<int dim, template<int, int> class DH>
892  bool sanity_check(const
894  unsigned int face,
895  unsigned int subface)
896  {
897  const auto neighbor_child =
898  _element->neighbor_child_on_subface(
899  face, subface);
900 
901  bool ret = false;
902  if(neighbor_child->face(_element->neighbor_of_neighbor(face)) == _element->face(face)->child(subface))
903  ret = true;
904  return ret;
905  }
906 
907  template<>
908  bool sanity_check<1,dealii::hp::DoFHandler>(const
910  unsigned int,
911  unsigned int)
912  {
913  return true;
914  }
915 
916  template<>
917  bool sanity_check<1,dealii::DoFHandler>(const
919  unsigned int,
920  unsigned int)
921  {
922  return true;
923  }
924  }
925 
926  template<typename VECTOR, int dim>
927  void
929  unsigned int face_no)
930  {
931  this->SetFace(face_no);
932  _state_fe_values.reinit(_element[this->GetStateIndex()], face_no);
933  _state_fe_values_ptr = &_state_fe_values;
934  //Make sure that the Control must be initialized.
935  if (this->GetControlIndex() < _element.size())
936  {
937  _control_fe_values.reinit(_element[this->GetControlIndex()], face_no);
938  _control_fe_values_ptr = &_control_fe_values;
939  }
940  }
941 
942  /***********************************************************************/
943 
944  template<typename VECTOR, int dim>
945  void
947  unsigned int face_no, unsigned int subface_no)
948  {
949  this->SetFace(face_no);
950  this->SetSubFace(subface_no);
951  _state_fe_subface_values->reinit(_element[this->GetStateIndex()], face_no,
952  subface_no);
953  _state_fe_values_ptr = _state_fe_subface_values;
954  //Make sure that the Control must be initialized.
955  if (this->GetControlIndex() < _element.size())
956  {
957  _control_fe_subface_values->reinit(_element[this->GetControlIndex()],
958  face_no, this->GetSubFace());
959  _control_fe_values_ptr = _control_fe_subface_values;
960  }
961  }
962  /***********************************************************************/
963 
964 
965  template<typename VECTOR, int dim>
966  void
968  {
969  Assert(this->NeedNeighbour(), ExcInternalError());
970  Assert(
971  _element[this->GetStateIndex()]->neighbor_index(this->GetFace()) != -1,
972  TriaAccessorExceptions::ExcUnusedCellAsNeighbor())
973 
974  if (_element[this->GetStateIndex()]->neighbor(this->GetFace())->has_children())
975  {
976  //if neighbor is more refined
977  const auto neighbor_child =
978  _element[this->GetStateIndex()]->neighbor_child_on_subface(
979  this->GetFace(), this->GetSubFace());
980 
981  // some sanity checks: Check, that the face and subface match and that the neighbour child
982  // is not more refined.
983  Assert((sanity_check<dim, dealii::DoFHandler>(_element[this->GetStateIndex()],
984  this->GetFace(),
985  this->GetSubFace()) == true), ExcInternalError());
986  Assert(neighbor_child->has_children() == false, ExcInternalError());
987 
988  _nbr_state_fe_values->reinit(neighbor_child,
989  _element[this->GetStateIndex()]->neighbor_of_neighbor(
990  this->GetFace()));
991  _nbr_state_fe_values_ptr = _nbr_state_fe_values;
992 
993  //Make sure that the Control must be initialized.
994  if (this->GetControlIndex() < _element.size())
995  {
996  const auto control_neighbor_child =
997  _element[this->GetControlIndex()]->neighbor_child_on_subface(
998  this->GetFace(), this->GetSubFace());
999 
1000  _nbr_control_fe_values->reinit(control_neighbor_child,
1001  _element[this->GetControlIndex()]->neighbor_of_neighbor(
1002  this->GetFace()));
1003  _nbr_control_fe_values_ptr = _nbr_control_fe_values;
1004  }
1005  }
1006  else if (_element[this->GetStateIndex()]->neighbor_is_coarser(
1007  this->GetFace()))
1008  {
1009  //if the neighbour is coarser
1010  Assert(
1011  _element[this->GetStateIndex()]->neighbor(this->GetFace())->level() == _element[this->GetStateIndex()]->level()-1,
1012  ExcInternalError());
1013  const auto neighbor = _element[this->GetStateIndex()]->neighbor(
1014  this->GetFace());
1015  const std::pair<unsigned int, unsigned int> faceno_subfaceno =
1016  _element[this->GetStateIndex()]->neighbor_of_coarser_neighbor(
1017  this->GetFace());
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())
1024  {
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(
1029  this->GetFace());
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;
1036  }
1037 
1038  }
1039  else
1040  {
1041  const auto neighbor_state = _element[this->GetStateIndex()]->neighbor(
1042  this->GetFace());
1043  // neighbor element is as much refined as the
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(
1048  this->GetFace()));
1049  _nbr_state_fe_values_ptr = _nbr_state_fe_values;
1050 
1051  //Make sure that the Control must be initialized.
1052  if (this->GetControlIndex() < _element.size())
1053  {
1054  _nbr_control_fe_values->reinit(
1055  _element[this->GetControlIndex()]->neighbor(this->GetFace()),
1056  _element[this->GetControlIndex()]->neighbor_of_neighbor(
1057  this->GetFace()));
1058  _nbr_control_fe_values_ptr = _nbr_control_fe_values;
1059  }
1060  }
1061  }
1062  /***********************************************************************/
1063  template<typename VECTOR, int dim>
1064  unsigned int
1066  {
1067  return _n_dofs_per_element;
1068  }
1069 
1070  /***********************************************************************/
1071 
1072  template<typename VECTOR, int dim>
1073  unsigned int
1075  {
1076  return _n_dofs_per_element;
1077  }
1078 
1079  /**********************************************/
1080  template<typename VECTOR, int dim>
1081  unsigned int
1083  {
1084  return _n_q_points_per_element;
1085  }
1086 
1087  /**********************************************/
1088  template<typename VECTOR, int dim>
1089  unsigned int
1091  {
1092  return _n_q_points_per_element;
1093  }
1094 
1095  /**********************************************/
1096  template<typename VECTOR, int dim>
1097  unsigned int
1099  {
1100  return _element[0]->material_id();
1101  }
1102 
1103  /**********************************************/
1104  template<typename VECTOR, int dim>
1105  unsigned int
1107  {
1108  return this->GetNbrMaterialId(this->GetFace());
1109  }
1110 
1111  /**********************************************/
1112  template<typename VECTOR, int dim>
1113  unsigned int
1115  unsigned int face) const
1116  {
1117  if (_element[0]->neighbor_index(face) != -1)
1118  return _element[0]->neighbor(face)->material_id();
1119  else
1120  throw DOpEException("There is no neighbor with number " + face,
1121  "FaceDataContainer::GetNbrMaterialId");
1122  }
1123 
1124  /**********************************************/
1125  template<typename VECTOR, int dim>
1126  bool
1128  {
1129  return _element[0]->face(this->GetFace())->at_boundary();
1130  }
1131 
1132  /**********************************************/
1133  template<typename VECTOR, int dim>
1134  double
1136  {
1137  return _element[0]->face(this->GetFace())->diameter();
1138  }
1139 
1140  /**********************************************/
1141 
1142  template<typename VECTOR, int dim>
1143  unsigned int
1145  {
1146  return _element[0]->face(this->GetFace())->boundary_indicator();
1147  }
1148 
1149  /**********************************************/
1150  template<typename VECTOR, int dim>
1151  const FEFaceValuesBase<dim>&
1153  {
1154  return *_state_fe_values_ptr;
1155  }
1156 
1157  /**********************************************/
1158  template<typename VECTOR, int dim>
1159  const FEFaceValuesBase<dim>&
1161  {
1162  return *_control_fe_values_ptr;
1163  }
1164  /**********************************************/
1165  template<typename VECTOR, int dim>
1166  const FEFaceValuesBase<dim>&
1168  {
1169  return *_nbr_state_fe_values_ptr;
1170  }
1171 
1172  /**********************************************/
1173  template<typename VECTOR, int dim>
1174  const FEFaceValuesBase<dim>&
1176  {
1177  return *_nbr_control_fe_values_ptr;
1178  }
1179 
1180  template<typename VECTOR, int dim>
1181  unsigned int
1183  {
1184  return _state_index;
1185  }
1186 
1187  /***********************************************************************/
1188 
1189  template<typename VECTOR, int dim>
1190  unsigned int
1191  FaceDataContainer<dealii::DoFHandler, VECTOR, dim>::GetControlIndex() const
1192  {
1193  return _control_index;
1194  }
1195 
1196  /***********************************************************************/
1197  /************************END*OF*IMPLEMENTATION**************************/
1198  /***********************************************************************/
1199  /***********************************************************************/
1200  /*****************IMPLEMENTATION for MGDoFHandler*********************/
1201  /***********************************************************************/
1202 
1203 
1204 
1205  template<typename VECTOR, int dim>
1206  void
1208  unsigned int face_no)
1209  {
1210  this->SetFace(face_no);
1211  _state_fe_values.reinit(_element[this->GetStateIndex()], face_no);
1212  _state_fe_values_ptr = &_state_fe_values;
1213  //Make sure that the Control must be initialized.
1214  if (this->GetControlIndex() < _element.size())
1215  {
1216  _control_fe_values.reinit(_element[this->GetControlIndex()], face_no);
1217  _control_fe_values_ptr = &_control_fe_values;
1218  }
1219  }
1220 
1221  /***********************************************************************/
1222 
1223  template<typename VECTOR, int dim>
1224  void
1226  unsigned int face_no, unsigned int subface_no)
1227  {
1228  this->SetFace(face_no);
1229  this->SetSubFace(subface_no);
1230  _state_fe_subface_values->reinit(_element[this->GetStateIndex()], face_no,
1231  subface_no);
1232  _state_fe_values_ptr = _state_fe_subface_values;
1233  //Make sure that the Control must be initialized.
1234  if (this->GetControlIndex() < _element.size())
1235  {
1236  _control_fe_subface_values->reinit(_element[this->GetControlIndex()],
1237  face_no, this->GetSubFace());
1238  _control_fe_values_ptr = _control_fe_subface_values;
1239  }
1240  }
1241  /***********************************************************************/
1242 
1243  template<typename VECTOR, int dim>
1244  void
1246  {
1247  Assert(this->NeedNeighbour(), ExcInternalError());
1248  Assert(
1249  _element[this->GetStateIndex()]->neighbor_index(this->GetFace()) != -1,
1250  TriaAccessorExceptions::ExcUnusedCellAsNeighbor())
1251 
1252  if (_element[this->GetStateIndex()]->neighbor(this->GetFace())->has_children())
1253  {
1254  //if neighbor is more refined
1255  const auto neighbor_child =
1256  _element[this->GetStateIndex()]->neighbor_child_on_subface(
1257  this->GetFace(), this->GetSubFace());
1258 
1259  // some sanity checks: Check, that the face and subface match and that the neighbour child
1260  // is not more refined.
1261  Assert((sanity_check<dim, dealii::hp::DoFHandler>(_element[this->GetStateIndex()],
1262  this->GetFace(),
1263  this->GetSubFace()) == true), ExcInternalError());
1264  Assert(neighbor_child->has_children() == false, ExcInternalError());
1265 
1266  _nbr_state_fe_values->reinit(neighbor_child,
1267  _element[this->GetStateIndex()]->neighbor_of_neighbor(
1268  this->GetFace()));
1269  _nbr_state_fe_values_ptr = _nbr_state_fe_values;
1270 
1271  //Make sure that the Control must be initialized.
1272  if (this->GetControlIndex() < _element.size())
1273  {
1274  const auto control_neighbor_child =
1275  _element[this->GetControlIndex()]->neighbor_child_on_subface(
1276  this->GetFace(), this->GetSubFace());
1277 
1278  _nbr_control_fe_values->reinit(control_neighbor_child,
1279  _element[this->GetControlIndex()]->neighbor_of_neighbor(
1280  this->GetFace()));
1281  _nbr_control_fe_values_ptr = _nbr_control_fe_values;
1282  }
1283  }
1284  else if (_element[this->GetStateIndex()]->neighbor_is_coarser(
1285  this->GetFace()))
1286  {
1287  //if the neighbour is coarser
1288  Assert(
1289  _element[this->GetStateIndex()]->neighbor(this->GetFace())->level() == _element[this->GetStateIndex()]->level()-1,
1290  ExcInternalError());
1291  const auto neighbor = _element[this->GetStateIndex()]->neighbor(
1292  this->GetFace());
1293  const std::pair<unsigned int, unsigned int> faceno_subfaceno =
1294  _element[this->GetStateIndex()]->neighbor_of_coarser_neighbor(
1295  this->GetFace());
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())
1302  {
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(
1307  this->GetFace());
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;
1314  }
1315 
1316  }
1317  else
1318  {
1319  const auto neighbor_state = _element[this->GetStateIndex()]->neighbor(
1320  this->GetFace());
1321  // neighbor element is as much refined as the
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(
1326  this->GetFace()));
1327  _nbr_state_fe_values_ptr = _nbr_state_fe_values;
1328 
1329  //Make sure that the Control must be initialized.
1330  if (this->GetControlIndex() < _element.size())
1331  {
1332  _nbr_control_fe_values->reinit(
1333  _element[this->GetControlIndex()]->neighbor(this->GetFace()),
1334  _element[this->GetControlIndex()]->neighbor_of_neighbor(
1335  this->GetFace()));
1336  _nbr_control_fe_values_ptr = _nbr_control_fe_values;
1337  }
1338  }
1339  }
1340  /***********************************************************************/
1341  template<typename VECTOR, int dim>
1342  unsigned int
1344  {
1345  return _n_dofs_per_element;
1346  }
1347 
1348  /***********************************************************************/
1349 
1350  template<typename VECTOR, int dim>
1351  unsigned int
1353  {
1354  return _n_dofs_per_element;
1355  }
1356 
1357  /**********************************************/
1358  template<typename VECTOR, int dim>
1359  unsigned int
1361  {
1362  return _n_q_points_per_element;
1363  }
1364 
1365  /**********************************************/
1366  template<typename VECTOR, int dim>
1367  unsigned int
1369  {
1370  return _n_q_points_per_element;
1371  }
1372 
1373  /**********************************************/
1374  template<typename VECTOR, int dim>
1375  unsigned int
1377  {
1378  return _element[0]->material_id();
1379  }
1380 
1381  /**********************************************/
1382  template<typename VECTOR, int dim>
1383  unsigned int
1385  {
1386  return this->GetNbrMaterialId(this->GetFace());
1387  }
1388 
1389  /**********************************************/
1390  template<typename VECTOR, int dim>
1391  unsigned int
1393  unsigned int face) const
1394  {
1395  if (_element[0]->neighbor_index(face) != -1)
1396  return _element[0]->neighbor(face)->material_id();
1397  else
1398  throw DOpEException("There is no neighbor with number " + face,
1399  "FaceDataContainer::GetNbrMaterialId");
1400  }
1401 
1402  /**********************************************/
1403  template<typename VECTOR, int dim>
1404  bool
1406  {
1407  return _element[0]->face(this->GetFace())->at_boundary();
1408  }
1409 
1410  /**********************************************/
1411  template<typename VECTOR, int dim>
1412  double
1414  {
1415  return _element[0]->face(this->GetFace())->diameter();
1416  }
1417 
1418  /**********************************************/
1419 
1420  template<typename VECTOR, int dim>
1421  unsigned int
1423  {
1424  return _element[0]->face(this->GetFace())->boundary_indicator();
1425  }
1426 
1427  /**********************************************/
1428  template<typename VECTOR, int dim>
1429  const FEFaceValuesBase<dim>&
1431  {
1432  return *_state_fe_values_ptr;
1433  }
1434 
1435  /**********************************************/
1436  template<typename VECTOR, int dim>
1437  const FEFaceValuesBase<dim>&
1439  {
1440  return *_control_fe_values_ptr;
1441  }
1442  /**********************************************/
1443  template<typename VECTOR, int dim>
1444  const FEFaceValuesBase<dim>&
1446  {
1447  return *_nbr_state_fe_values_ptr;
1448  }
1449 
1450  /**********************************************/
1451  template<typename VECTOR, int dim>
1452  const FEFaceValuesBase<dim>&
1454  {
1455  return *_nbr_control_fe_values_ptr;
1456  }
1457 
1458  template<typename VECTOR, int dim>
1459  unsigned int
1461  {
1462  return _state_index;
1463  }
1464 
1465  /***********************************************************************/
1466 
1467  template<typename VECTOR, int dim>
1468  unsigned int
1469  FaceDataContainer<dealii::MGDoFHandler, VECTOR, dim>::GetControlIndex() const
1470  {
1471  return _control_index;
1472  }
1473 
1474  /***********************************************************************/
1475  /************************END*OF*IMPLEMENTATION**************************/
1476  /***********************************************************************/
1477  /***********************************************************************/
1478  /*****************IMPLEMENTATION for hp::DoFHandler*********************/
1479  /***********************************************************************/
1480 
1481 
1482 
1483 
1484 
1485  template<typename VECTOR, int dim>
1486  void
1488  unsigned int face_no)
1489  {
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();
1493  //Make sure that the Control must be initialized.
1494  if (this->GetControlIndex() < _element.size())
1495  {
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();
1499  }
1500  }
1501 
1502  /***********************************************************************/
1503 
1504  template<typename VECTOR, int dim>
1505  void
1507  unsigned int face_no, unsigned int subface_no)
1508  {
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();
1515  //Make sure that the Control must be initialized.
1516  if (this->GetControlIndex() < _element.size())
1517  {
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();
1522  }
1523  }
1524 
1525  /***********************************************************************/
1526  template<typename VECTOR, int dim>
1527  void
1529  {
1530  Assert(this->NeedNeighbour(), ExcInternalError());
1531  Assert(
1532  _element[this->GetStateIndex()]->neighbor_index(this->GetFace()) != -1,
1533  TriaAccessorExceptions::ExcUnusedCellAsNeighbor())
1534 
1535  if (_element[this->GetStateIndex()]->neighbor(this->GetFace())->has_children())
1536  {
1537  //if neighbor is more refined
1538  const auto neighbor_child =
1539  _element[this->GetStateIndex()]->neighbor_child_on_subface(
1540  this->GetFace(), this->GetSubFace());
1541 
1542  // some sanity checks: Check, that the face and subface match and that the neighbour child
1543  // is not more refined.
1544  Assert(
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());
1548 
1549  _nbr_state_hp_fe_values->reinit(neighbor_child,
1550  _element[this->GetStateIndex()]->neighbor_of_neighbor(
1551  this->GetFace()));
1552  _nbr_state_hp_fe_values_ptr =
1553  &_nbr_state_hp_fe_values->get_present_fe_values();
1554 
1555  //Make sure that the Control must be initialized.
1556  if (this->GetControlIndex() < _element.size())
1557  {
1558  const auto control_neighbor_child =
1559  _element[this->GetControlIndex()]->neighbor_child_on_subface(
1560  this->GetFace(), this->GetSubFace());
1561 
1562  _nbr_control_hp_fe_values->reinit(control_neighbor_child,
1563  _element[this->GetControlIndex()]->neighbor_of_neighbor(
1564  this->GetFace()));
1565  _nbr_control_hp_fe_values_ptr =
1566  &_nbr_control_hp_fe_values->get_present_fe_values();
1567  }
1568  }
1569  else if (_element[this->GetStateIndex()]->neighbor_is_coarser(
1570  this->GetFace()))
1571  {
1572  //if the neighbour is coarser
1573  Assert(
1574  _element[this->GetStateIndex()]->neighbor(this->GetFace())->level() == _element[this->GetStateIndex()]->level()-1,
1575  ExcInternalError());
1576  const auto neighbor = _element[this->GetStateIndex()]->neighbor(
1577  this->GetFace());
1578  const std::pair<unsigned int, unsigned int> faceno_subfaceno =
1579  _element[this->GetStateIndex()]->neighbor_of_coarser_neighbor(
1580  this->GetFace());
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())
1588  {
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(
1593  this->GetFace());
1594  const unsigned int control_neighbor_face_no =
1595  control_faceno_subfaceno.first, control_neighbor_subface_no =
1596  control_faceno_subfaceno.second;
1597 
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();
1602  }
1603  }
1604  else
1605  {
1606  const auto neighbor_state = _element[this->GetStateIndex()]->neighbor(
1607  this->GetFace());
1608  // neighbor element is as much refined as the
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(
1613  this->GetFace()));
1614  _nbr_state_hp_fe_values_ptr =
1615  &_nbr_state_hp_fe_values->get_present_fe_values();
1616 
1617  //Make sure that the Control must be initialized.
1618  if (this->GetControlIndex() < _element.size())
1619  {
1620  _nbr_control_hp_fe_values->reinit(
1621  _element[this->GetControlIndex()]->neighbor(this->GetFace()),
1622  _element[this->GetControlIndex()]->neighbor_of_neighbor(
1623  this->GetFace()));
1624  _nbr_control_hp_fe_values_ptr =
1625  &_nbr_control_hp_fe_values->get_present_fe_values();
1626  }
1627  }
1628  }
1629 
1630  /*********************************************/
1631 
1632  template<typename VECTOR, int dim>
1633  unsigned int
1635  {
1636  return _element[0]->get_fe().dofs_per_cell;
1637  }
1638 
1639  /*********************************************/
1640 
1641  template<typename VECTOR, int dim>
1642  unsigned int
1644  {
1645  if (_element[0]->neighbor_index(this->GetFace()) != -1)
1646  return _element[0]->neighbor(this->GetFace())->get_fe().dofs_per_cell;
1647  else
1648  throw DOpEException(
1649  "There is no neighbor with number" + this->GetFace(),
1650  "HpFaceDataContainer::GetNbrNDoFsPerElement");
1651  }
1652  /*********************************************/
1653  template<typename VECTOR, int dim>
1654  unsigned int
1656  {
1657  return _q_collection[_element[0]->active_fe_index()].size();
1658  }
1659 
1660  /*********************************************/
1661  template<typename VECTOR, int dim>
1662  unsigned int
1664  {
1665  if (_element[0]->neighbor_index(this->GetFace()) != -1)
1666  return _q_collection[_element[0]->neighbor(this->GetFace())->active_fe_index()].size();
1667  else
1668  throw DOpEException(
1669  "There is no neighbor with number" + this->GetFace(),
1670  "HpFaceDataContainer::GetNbrNQPoints");
1671  }
1672  /*********************************************/
1673 
1674  template<typename VECTOR, int dim>
1675  unsigned int
1677  {
1678  return _element[0]->material_id();
1679  }
1680  /*********************************************/
1681  template<typename VECTOR, int dim>
1682  unsigned int
1684  {
1685  return this->GetNbrMaterialId(this->GetFace());
1686  }
1687 
1688  /*********************************************/
1689 
1690  template<typename VECTOR, int dim>
1691  unsigned int
1693  unsigned int face) const
1694  {
1695  if (_element[0]->neighbor_index(face) != -1)
1696  return _element[0]->neighbor(face)->material_id();
1697  else
1698  throw DOpEException("There is no neighbor with number" + face,
1699  "HpFaceDataContainer::GetNbrMaterialId");
1700  }
1701 
1702  /*********************************************/
1703 
1704  template<typename VECTOR, int dim>
1705  double
1707  {
1708  return _element[0]->face(this->GetFace())->diameter();
1709  }
1710 
1711  /**********************************************/
1712 
1713  template<typename VECTOR, int dim>
1714  bool
1716  {
1717  return _element[0]->face(this->GetFace())->at_boundary();
1718  }
1719 
1720  /**********************************************/
1721  template<typename VECTOR, int dim>
1722  unsigned int
1724  {
1725  return _element[0]->face(this->GetFace())->boundary_indicator();
1726  }
1727 
1728  /*********************************************/
1729  template<typename VECTOR, int dim>
1730  const FEFaceValuesBase<dim>&
1732  {
1733  return *_state_hp_fe_values_ptr;
1734  }
1735  /*********************************************/
1736  template<typename VECTOR, int dim>
1737  const FEFaceValuesBase<dim>&
1739  {
1740  return *_control_hp_fe_values_ptr;
1741  }
1742  /*********************************************/
1743  template<typename VECTOR, int dim>
1744  const FEFaceValuesBase<dim>&
1746  {
1747  return *_nbr_state_hp_fe_values_ptr;
1748  }
1749  /*********************************************/
1750  template<typename VECTOR, int dim>
1751  const FEFaceValuesBase<dim>&
1753  {
1754  return *_nbr_control_hp_fe_values_ptr;
1755  }
1756 
1757  /**********************************************/
1758  template<typename VECTOR, int dim>
1759  unsigned int
1761  {
1762  return _state_index;
1763  }
1764  /*********************************************/
1765  template<typename VECTOR, int dim>
1766  unsigned int
1767  FaceDataContainer<dealii::hp::DoFHandler, VECTOR, dim>::GetControlIndex() const
1768  {
1769  return _control_index;
1770  }
1771 } //end of namespace
1772 
1773 #endif /* FACEDATACONTAINER_H_ */
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 > * > &param_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 > * > &param_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 > * > &param_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
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 > * > &param_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