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"
31 #include "facedatacontainer_internal.h"
32 
33 #include <sstream>
34 
35 #include <deal.II/dofs/dof_handler.h>
36 //#include <deal.II/multigrid/mg_dof_handler.h>
37 #include <deal.II/hp/dof_handler.h>
38 
39 using namespace dealii;
40 
41 namespace DOpE
42 {
49  template<template<int, int> class DH, typename VECTOR, int dim>
50  class FaceDataContainer : public fdcinternal::FaceDataContainerInternal<
51  VECTOR, dim>
52  {
53  public:
54  FaceDataContainer()
55  {
56  throw(DOpEException(
57  "Dummy class, this constructor should never get called.",
58  "ElementDataContainer<dealii::DoFHandler , VECTOR, dim>::ElementDataContainer"));
59  }
60  ;
61  };
62 
72  template<typename VECTOR, int dim>
73  class FaceDataContainer<dealii::DoFHandler, VECTOR, dim> : public fdcinternal::FaceDataContainerInternal<
74  VECTOR, dim>
75  {
76 
77  public:
99  template<template<int, int> class FE, typename SPARSITYPATTERN, int dopedim, int dealdim>
100  FaceDataContainer(const Quadrature<dim - 1>& quad,
101  UpdateFlags update_flags,
102  SpaceTimeHandler<FE, dealii::DoFHandler, SPARSITYPATTERN, VECTOR,
103  dopedim, dealdim>& sth,
104  const std::vector<
105  typename DOpEWrapper::DoFHandler<dim, dealii::DoFHandler>::active_cell_iterator>& element,
106  const std::map<std::string, const Vector<double>*> &param_values,
107  const std::map<std::string, const VECTOR*> &domain_values,
108  bool need_neighbour) :
109  fdcinternal::FaceDataContainerInternal<VECTOR, dim>(param_values,
110  domain_values, need_neighbour), element_(element), state_fe_values_(
111  sth.GetMapping(), (sth.GetFESystem("state")), quad,
112  update_flags), control_fe_values_(sth.GetMapping(),
113  (sth.GetFESystem("control")), quad, update_flags)
114  {
115  state_index_ = sth.GetStateIndex();
116  if (state_index_ == 1)
117  control_index_ = 0;
118  else
119  control_index_ = 1;
120 
121  if (need_neighbour) //so we need FEFAcevalues etc. for the neighbour too.
122  {
123  nbr_control_fe_values_ = new DOpEWrapper::FEFaceValues<dim>(
124  sth.GetMapping(), (sth.GetFESystem("control")), quad,
125  update_flags);
126  control_fe_subface_values_ =
127  new DOpEWrapper::FESubfaceValues<dim>(sth.GetMapping(),
128  (sth.GetFESystem("control")), quad, update_flags);
129  }
130  this->PrivateConstructor(quad, update_flags, sth, need_neighbour);
131  }
150  template<template<int, int> class FE, typename SPARSITYPATTERN>
151  FaceDataContainer(const Quadrature<dim - 1>& quad,
152  UpdateFlags update_flags,
153  StateSpaceTimeHandler<FE, dealii::DoFHandler, SPARSITYPATTERN,
154  VECTOR, dim>& sth,
155  const std::vector<
156  typename DOpEWrapper::DoFHandler<dim, dealii::DoFHandler>::active_cell_iterator>& element,
157  const std::map<std::string, const Vector<double>*> &param_values,
158  const std::map<std::string, const VECTOR*> &domain_values,
159  bool need_neighbour) :
160  fdcinternal::FaceDataContainerInternal<VECTOR, dim>(param_values,
161  domain_values, need_neighbour), element_(element), state_fe_values_(
162  sth.GetMapping(), (sth.GetFESystem("state")), quad,
163  update_flags), control_fe_values_(sth.GetMapping(),
164  (sth.GetFESystem("state")), quad, update_flags)
165  {
166  state_index_ = sth.GetStateIndex();
167  control_index_ = element.size();
168  n_q_points_per_element_ = quad.size();
169  n_dofs_per_element_ = element[0]->get_fe().dofs_per_cell;
170 
171  if (need_neighbour) //so we need FEFAcevalues for the neighbour too.
172  {
173  nbr_control_fe_values_ = new DOpEWrapper::FEFaceValues<dim>(
174  sth.GetMapping(), (sth.GetFESystem("state")), quad,
175  update_flags);
176  control_fe_subface_values_ =
177  new DOpEWrapper::FESubfaceValues<dim>(sth.GetMapping(),
178  (sth.GetFESystem("state")), quad, update_flags);
179  }
180  this->PrivateConstructor(quad, update_flags, sth, need_neighbour);
181  }
182 
183  ~FaceDataContainer()
184  {
185  if (nbr_state_fe_values_ != NULL)
186  {
187  delete nbr_state_fe_values_;
188  }
189  if (nbr_control_fe_values_ != NULL)
190  {
191  delete nbr_control_fe_values_;
192  }
193  if (state_fe_subface_values_ != NULL)
194  {
195  delete state_fe_subface_values_;
196  }
197  if (control_fe_subface_values_ != NULL)
198  {
199  delete control_fe_subface_values_;
200  }
201  }
202  /*********************************************/
203  /*
204  * This function reinitializes the FEFaceValues on the actual face. Should
205  * be called prior to any of the get-functions.
206  *
207  * @param face_no The 'local number' (i.e. from the perspective of the actual element) of the
208  * actual face.
209  */
210  inline void
211  ReInit(unsigned int face_no);
212 
213  /*********************************************/
214  /*
215  * This function reinits the FESubfaceValues on the actual subface. Should
216  * be called prior to any of the get-functions.
217  *
218  * @param face_no The 'local number' (i.e. from the perspective of the actual element) of the
219  * actual face.
220  * @param subface_no The 'local number' (i.e. from the perspective of the actual element) of the
221  * actual subface.
222  */
223  inline void
224  ReInit(unsigned int face_no, unsigned int subface_no);
225 
226  /*********************************************/
227  /*
228  * This function reinitializes the FE(Sub)FaceValues on the neighbor_element.
229  * This should be called prior to any of the get nbr_functions.
230  * Assumes that ReInit is called prior to this function.
231  */
232  inline void
233  ReInitNbr();
234 
235  /*********************************************/
241  inline unsigned int
242  GetNDoFsPerElement() const;
243  inline unsigned int
244  GetNbrNDoFsPerElement() const;
245  inline unsigned int
246  GetNQPoints() const;
247  inline unsigned int
248  GetNbrNQPoints() const;
249  inline unsigned int
250  GetMaterialId() const;
251  inline unsigned int
252  GetNbrMaterialId() const;
253  inline unsigned int
254  GetNbrMaterialId(unsigned int face) const;
255  inline bool
256  GetIsAtBoundary() const;
257  inline double
258  GetElementDiameter() const;
259  inline unsigned int
260  GetBoundaryIndicator() const;
261  inline const FEFaceValuesBase<dim>&
262  GetFEFaceValuesState() const;
263  inline const FEFaceValuesBase<dim>&
264  GetFEFaceValuesControl() const;
265 
266  inline const FEFaceValuesBase<dim>&
267  GetNbrFEFaceValuesState() const;
268  inline const FEFaceValuesBase<dim>&
269  GetNbrFEFaceValuesControl() const;
270 
271  private:
272  /*
273  * Helper Functions
274  */
275  unsigned int
276  GetStateIndex() const;
277  unsigned int
278  GetControlIndex() const;
282  template<class STH>
283  void
284  PrivateConstructor(const Quadrature<dim - 1>& quad,
285  UpdateFlags update_flags, STH& sth, bool need_neighbour)
286  {
287  n_q_points_per_element_ = quad.size();
288  n_dofs_per_element_ = element_[0]->get_fe().dofs_per_cell;
289 
290  if (need_neighbour) //so we need FEFAcevalues etc. for the neighbour too.
291  {
292  nbr_state_fe_values_ = new DOpEWrapper::FEFaceValues<dim>(
293  (sth.GetFESystem("state")), quad, update_flags);
294  state_fe_subface_values_ = new DOpEWrapper::FESubfaceValues<dim>(
295  (sth.GetFESystem("state")), quad, update_flags);
296  }
297  else
298  {
299  nbr_state_fe_values_ = NULL;
300  nbr_control_fe_values_ = NULL;
301  state_fe_subface_values_ = NULL;
302  control_fe_subface_values_ = NULL;
303  }
304  // These will point to the object (i.e. FaceValues or SubfaceValues) we actually use.
305  // With this, we have the same interface to the user independently of the type (i.e. face or subface)
306  state_fe_values_ptr_ = NULL;
307  control_fe_values_ptr_ = NULL;
308  nbr_state_fe_values_ptr_ = NULL;
309  nbr_control_fe_values_ptr_ = NULL;
310  }
311  /***********************************************************/
312  //"global" member data, part of every instantiation
313  unsigned int state_index_;
314  unsigned int control_index_;
315 
316  const std::vector<
317  typename DOpEWrapper::DoFHandler<dim, dealii::DoFHandler>::active_cell_iterator> & element_;
318  DOpEWrapper::FEFaceValues<dim> state_fe_values_;
319  DOpEWrapper::FEFaceValues<dim> control_fe_values_;
320 
321  DOpEWrapper::FEFaceValues<dim>* nbr_state_fe_values_;
322  DOpEWrapper::FEFaceValues<dim>* nbr_control_fe_values_;
323 
324  DOpEWrapper::FESubfaceValues<dim>* state_fe_subface_values_;
325  DOpEWrapper::FESubfaceValues<dim>* control_fe_subface_values_;
326 
327  dealii::FEFaceValuesBase<dim>* state_fe_values_ptr_;
328  dealii::FEFaceValuesBase<dim>* control_fe_values_ptr_;
329  dealii::FEFaceValuesBase<dim>* nbr_state_fe_values_ptr_;
330  dealii::FEFaceValuesBase<dim>* nbr_control_fe_values_ptr_;
331 
332  unsigned int n_q_points_per_element_;
333  unsigned int n_dofs_per_element_;
334  };
335 
336 
337 
338 
339 
340  /****************************************************/
341  /* MGDofHandler */
342 //
343 //
344 // /**
345 // * This two classes hold all the information we need in the integrator to
346 // * integrate something over a face of a element (could be a functional, a PDE, etc.).
347 // * Of particular importance: This class holds the (Sub)FaceFEValues objects.
348 // *
349 // * @template VECTOR Type of the vector we use in our computations (i.e. Vector<double> or BlockVector<double>)
350 // * @template dim 1+ the dimension of the integral we are actually interested in.//TODO 1+??
351 // */
352 //
353 // template<typename VECTOR, int dim>
354 // class FaceDataContainer<dealii::MGDoFHandler, VECTOR, dim> : public fdcinternal::FaceDataContainerInternal<
355 // VECTOR, dim>
356 // {
357 //
358 // public:
359 // /**
360 // * Constructor. Initializes the FaceFEValues objects.
361 // *
362 // * @template FE The type of Finite Element in use here.
363 // * @template SPARSITYPATTERN The corresponding Sparsitypattern to the class-template VECTOR.
364 // * @template dopedim The dimension of the control variable.
365 // * @template dealdim The dimension of the state variable.
366 // *
367 // * @param quad Reference to the quadrature-rule which we use at the moment.
368 // * @param update_flags The update flags we need to initialize the FEValues obejcts
369 // * @param sth A reference to the SpaceTimeHandler in use.
370 // * @param element A vector of element iterators through which we gain most of the needed information (like
371 // * material_ids, n_dfos, etc.)
372 // * @param param_values A std::map containing parameter data (e.g. non space dependent data). If the control
373 // * is done by parameters, it is contained in this map at the position "control".
374 // * @param domain_values A std::map containing domain data (e.g. nodal vectors for FE-Functions). If the control
375 // * is distributed, it is contained in this map at the position "control". The state may always
376 // * be found in this map at the position "state"
377 // * @param need_neighbour Describes whether we need all the GetNbr (= Get Neighbor) functions.
378 // *
379 // */
380 // template<template<int, int> class FE, typename SPARSITYPATTERN, int dopedim, int dealdim>
381 // FaceDataContainer(const Quadrature<dim - 1>& quad,
382 // UpdateFlags update_flags,
383 // SpaceTimeHandler<FE, dealii::MGDoFHandler, SPARSITYPATTERN,
384 // VECTOR, dopedim, dealdim>& sth,
385 // const std::vector<
386 // typename dealii::MGDoFHandler<dim>::active_cell_iterator>& element,
387 // const std::map<std::string, const Vector<double>*> &param_values,
388 // const std::map<std::string, const VECTOR*> &domain_values,
389 // bool need_neighbour)
390 // : fdcinternal::FaceDataContainerInternal<VECTOR, dim>(
391 // param_values, domain_values, need_neighbour), element_(element), state_fe_values_(
392 // sth.GetMapping(), (sth.GetFESystem("state")), quad,
393 // update_flags), control_fe_values_(sth.GetMapping(),
394 // (sth.GetFESystem("control")), quad, update_flags)
395 // {
396 // state_index_ = sth.GetStateIndex();
397 // if (state_index_ == 1)
398 // control_index_ = 0;
399 // else
400 // control_index_ = 1;
401 //
402 // if (need_neighbour) //so we need FEFAcevalues etc. for the neighbour too.
403 // {
404 // nbr_control_fe_values_ = new DOpEWrapper::FEFaceValues<dim>(
405 // sth.GetMapping(), (sth.GetFESystem("control")), quad,
406 // update_flags);
407 // control_fe_subface_values_ =
408 // new DOpEWrapper::FESubfaceValues<dim>(sth.GetMapping(),
409 // (sth.GetFESystem("control")), quad, update_flags);
410 // }
411 // this->PrivateConstructor(quad, update_flags, sth, need_neighbour);
412 // }
413 // /**
414 // * Constructor. Initializes the FaceFEValues objects. For PDE only
415 // *
416 // * @template SPARSITYPATTERN The corresponding Sparsitypattern to the class-template VECTOR.
417 // *
418 // * @param quad Reference to the quadrature-rule which we use at the moment.
419 // * @param update_flags The update flags we need to initialize the FEValues obejcts
420 // * @param sth A reference to the SpaceTimeHandler in use.
421 // * @param element A vector of element iterators through which we gain most of the needed information (like
422 // * material_ids, n_dfos, etc.)
423 // * @param param_values A std::map containing parameter data (e.g. non space dependent data). If the control
424 // * is done by parameters, it is contained in this map at the position "control".
425 // * @param domain_values A std::map containing domain data (e.g. nodal vectors for FE-Functions). If the control
426 // * is distributed, it is contained in this map at the position "control". The state may always
427 // * be found in this map at the position "state".
428 // * @param need_neighbour Describes whether we need all the GetNbr (= Get Neighbor) functions.
429 // *
430 // */
431 // template<template<int, int> class FE, typename SPARSITYPATTERN>
432 // FaceDataContainer(const Quadrature<dim - 1>& quad,
433 // UpdateFlags update_flags,
434 // StateSpaceTimeHandler<FE, dealii::MGDoFHandler,
435 // SPARSITYPATTERN, VECTOR, dim>& sth,
436 // const std::vector<
437 // typename dealii::MGDoFHandler<dim>::active_cell_iterator>& element,
438 // const std::map<std::string, const Vector<double>*> &param_values,
439 // const std::map<std::string, const VECTOR*> &domain_values,
440 // bool need_neighbour)
441 // : fdcinternal::FaceDataContainerInternal<VECTOR, dim>(
442 // param_values, domain_values, need_neighbour), element_(element), state_fe_values_(
443 // sth.GetMapping(), (sth.GetFESystem("state")), quad,
444 // update_flags), control_fe_values_(sth.GetMapping(),
445 // (sth.GetFESystem("state")), quad, update_flags)
446 // {
447 // state_index_ = sth.GetStateIndex();
448 // control_index_ = element.size();
449 // n_q_points_per_element_ = quad.size();
450 // n_dofs_per_element_ = element[0]->get_fe().dofs_per_cell;
451 //
452 // if (need_neighbour) //so we need FEFAcevalues for the neighbour too.
453 // {
454 // nbr_control_fe_values_ = new DOpEWrapper::FEFaceValues<dim>(
455 // sth.GetMapping(), (sth.GetFESystem("state")), quad,
456 // update_flags);
457 // control_fe_subface_values_ =
458 // new DOpEWrapper::FESubfaceValues<dim>(sth.GetMapping(),
459 // (sth.GetFESystem("state")), quad, update_flags);
460 // }
461 // this->PrivateConstructor(quad, update_flags, sth, need_neighbour);
462 // }
463 //
464 // ~FaceDataContainer()
465 // {
466 // if (nbr_state_fe_values_ != NULL)
467 // {
468 // delete nbr_state_fe_values_;
469 // }
470 // if (nbr_control_fe_values_ != NULL)
471 // {
472 // delete nbr_control_fe_values_;
473 // }
474 // if (state_fe_subface_values_ != NULL)
475 // {
476 // delete state_fe_subface_values_;
477 // }
478 // if (control_fe_subface_values_ != NULL)
479 // {
480 // delete control_fe_subface_values_;
481 // }
482 // }
483 // /*********************************************/
484 // /*
485 // * This function reinitializes the FEFaceValues on the actual face. Should
486 // * be called prior to any of the get-functions.
487 // *
488 // * @param face_no The 'local number' (i.e. from the perspective of the actual element) of the
489 // * actual face.
490 // */
491 // inline void
492 // ReInit(unsigned int face_no);
493 //
494 // /*********************************************/
495 // /*
496 // * This function reinits the FESubfaceValues on the actual subface. Should
497 // * be called prior to any of the get-functions.
498 // *
499 // * @param face_no The 'local number' (i.e. from the perspective of the actual element) of the
500 // * actual face.
501 // * @param subface_no The 'local number' (i.e. from the perspective of the actual element) of the
502 // * actual subface.
503 // */
504 // inline void
505 // ReInit(unsigned int face_no, unsigned int subface_no);
506 //
507 // /*********************************************/
508 // /*
509 // * This function reinitializes the FE(Sub)FaceValues on the neighbor_element.
510 // * This should be called prior to any of the get nbr_functions.
511 // * Assumes that ReInit is called prior to this function.
512 // */
513 // inline void
514 // ReInitNbr();
515 //
516 // /*********************************************/
517 // /**
518 // * Get functions to extract data. They all assume that ReInit
519 // * (resp. ReInitNbr for the GetNbr* functions) is executed
520 // * before calling them.
521 // */
522 // inline unsigned int
523 // GetNDoFsPerElement() const;
524 // inline unsigned int
525 // GetNbrNDoFsPerElement() const;
526 // inline unsigned int
527 // GetNQPoints() const;
528 // inline unsigned int
529 // GetNbrNQPoints() const;
530 // inline unsigned int
531 // GetMaterialId() const;
532 // inline unsigned int
533 // GetNbrMaterialId() const;
534 // inline unsigned int
535 // GetNbrMaterialId(unsigned int face) const;
536 // inline bool
537 // GetIsAtBoundary() const;
538 // inline double
539 // GetElementDiameter() const;
540 // inline unsigned int
541 // GetBoundaryIndicator() const;
542 // inline const FEFaceValuesBase<dim>&
543 // GetFEFaceValuesState() const;
544 // inline const FEFaceValuesBase<dim>&
545 // GetFEFaceValuesControl() const;
546 //
547 // inline const FEFaceValuesBase<dim>&
548 // GetNbrFEFaceValuesState() const;
549 // inline const FEFaceValuesBase<dim>&
550 // GetNbrFEFaceValuesControl() const;
551 //
552 // private:
553 // /*
554 // * Helper Functions
555 // */
556 // unsigned int
557 // GetStateIndex() const;
558 // unsigned int
559 // GetControlIndex() const;
560 // /**
561 // * This function contains common code of the constructors.
562 // */
563 // template<class STH>
564 // void
565 // PrivateConstructor(const Quadrature<dim - 1>& quad,
566 // UpdateFlags update_flags, STH& sth, bool need_neighbour)
567 // {
568 // n_q_points_per_element_ = quad.size();
569 // n_dofs_per_element_ = element_[0]->get_fe().dofs_per_cell;
570 //
571 // if (need_neighbour) //so we need FEFAcevalues etc. for the neighbour too.
572 // {
573 // nbr_state_fe_values_ = new DOpEWrapper::FEFaceValues<dim>(
574 // (sth.GetFESystem("state")), quad, update_flags);
575 // state_fe_subface_values_ = new DOpEWrapper::FESubfaceValues<dim>(
576 // (sth.GetFESystem("state")), quad, update_flags);
577 // }
578 // else
579 // {
580 // nbr_state_fe_values_ = NULL;
581 // nbr_control_fe_values_ = NULL;
582 // state_fe_subface_values_ = NULL;
583 // control_fe_subface_values_ = NULL;
584 // }
585 // // These will point to the object (i.e. FaceValues or SubfaceValues) we actually use.
586 // // With this, we have the same interface to the user independently of the type (i.e. face or subface)
587 // state_fe_values_ptr_ = NULL;
588 // control_fe_values_ptr_ = NULL;
589 // nbr_state_fe_values_ptr_ = NULL;
590 // nbr_control_fe_values_ptr_ = NULL;
591 // }
592 // /***********************************************************/
593 // //"global" member data, part of every instantiation
594 // unsigned int state_index_;
595 // unsigned int control_index_;
596 //
597 // const std::vector<
598 // typename DOpEWrapper::DoFHandler<dim, dealii::MGDoFHandler >::active_cell_iterator> & element_;
599 // DOpEWrapper::FEFaceValues<dim> state_fe_values_;
600 // DOpEWrapper::FEFaceValues<dim> control_fe_values_;
601 //
602 // DOpEWrapper::FEFaceValues<dim>* nbr_state_fe_values_;
603 // DOpEWrapper::FEFaceValues<dim>* nbr_control_fe_values_;
604 //
605 // DOpEWrapper::FESubfaceValues<dim>* state_fe_subface_values_;
606 // DOpEWrapper::FESubfaceValues<dim>* control_fe_subface_values_;
607 //
608 // dealii::FEFaceValuesBase<dim>* state_fe_values_ptr_;
609 // dealii::FEFaceValuesBase<dim>* control_fe_values_ptr_;
610 // dealii::FEFaceValuesBase<dim>* nbr_state_fe_values_ptr_;
611 // dealii::FEFaceValuesBase<dim>* nbr_control_fe_values_ptr_;
612 //
613 // unsigned int n_q_points_per_element_;
614 // unsigned int n_dofs_per_element_;
615 // };
616 //
617  /* MGDofHandler */
618  /****************************************************/
619 
620 
621  template<typename VECTOR, int dim>
622  class FaceDataContainer<dealii::hp::DoFHandler, VECTOR, dim> : public fdcinternal::FaceDataContainerInternal<
623  VECTOR, dim>
624  {
625 
626  public:
647  template<template<int, int> class FE, typename SPARSITYPATTERN, int dopedim, int dealdim>
648  FaceDataContainer<dealii::hp::DoFHandler, VECTOR, dim>(
649  const hp::QCollection<dim - 1>& q_collection,
650  UpdateFlags update_flags,
651  SpaceTimeHandler<FE, dealii::hp::DoFHandler, SPARSITYPATTERN,
652  VECTOR, dopedim, dealdim>& sth,
653  const std::vector<
654  typename DOpEWrapper::DoFHandler<dim, dealii::hp::DoFHandler>::active_cell_iterator>& element,
655  const std::map<std::string, const Vector<double>*> &param_values,
656  const std::map<std::string, const VECTOR*> &domain_values,
657  bool need_neighbour) :
658  fdcinternal::FaceDataContainerInternal<VECTOR, dim>(param_values,
659  domain_values, need_neighbour), element_(element), state_hp_fe_values_(
660  (sth.GetFESystem("state")), q_collection, update_flags), control_hp_fe_values_(
661  (sth.GetFESystem("control")), q_collection, update_flags), q_collection_(
662  q_collection)
663  {
664  state_index_ = sth.GetStateIndex();
665  if (state_index_ == 1)
666  control_index_ = 0;
667  else
668  control_index_ = 1;
669 
670  if (need_neighbour) //so we need FEFAcevalues for the neighbour too.
671  {
672  nbr_control_hp_fe_values_ = new DOpEWrapper::HpFEFaceValues<dim>(
673  sth.GetMapping(), (sth.GetFESystem("control")), q_collection,
674  update_flags);
675  control_hp_fe_subface_values_ =
676  new DOpEWrapper::HpFESubfaceValues<dim>(sth.GetMapping(),
677  (sth.GetFESystem("control")), q_collection, update_flags);
678 
679  }
680  this->PrivateConstructor(q_collection, update_flags, sth,
681  need_neighbour);
682  }
683 
704  template<template<int, int> class FE, typename SPARSITYPATTERN, int dealdim>
705  FaceDataContainer<dealii::hp::DoFHandler, VECTOR, dim>(
706  const hp::QCollection<dim - 1>& q_collection,
707  UpdateFlags update_flags,
708  StateSpaceTimeHandler<FE, dealii::hp::DoFHandler, SPARSITYPATTERN,
709  VECTOR, dealdim>& sth,
710  const std::vector<
711  typename DOpEWrapper::DoFHandler<dim, dealii::hp::DoFHandler>::active_cell_iterator>& element,
712  const std::map<std::string, const Vector<double>*> &param_values,
713  const std::map<std::string, const VECTOR*> &domain_values,
714  bool need_neighbour) :
715  fdcinternal::FaceDataContainerInternal<VECTOR, dim>(param_values,
716  domain_values, need_neighbour), element_(element), state_hp_fe_values_(
717  (sth.GetFESystem("state")), q_collection, update_flags), control_hp_fe_values_(
718  (sth.GetFESystem("state")), q_collection, update_flags), q_collection_(
719  q_collection)
720  {
721  state_index_ = sth.GetStateIndex();
722  if (state_index_ == 1)
723  control_index_ = 0;
724  else
725  control_index_ = 1;
726 
727  if (need_neighbour)
728  {
729  nbr_control_hp_fe_values_ = new DOpEWrapper::HpFEFaceValues<dim>(
730  sth.GetMapping(), (sth.GetFESystem("state")), q_collection,
731  update_flags);
732  control_hp_fe_subface_values_ =
733  new DOpEWrapper::HpFESubfaceValues<dim>(sth.GetMapping(),
734  (sth.GetFESystem("state")), q_collection, update_flags);
735  }
736  this->PrivateConstructor(q_collection, update_flags, sth,
737  need_neighbour);
738  }
739 
743  ~FaceDataContainer<dealii::hp::DoFHandler, VECTOR, dim>()
744  {
745  delete nbr_state_hp_fe_values_;
746  nbr_state_hp_fe_values_ = NULL;
747  delete nbr_control_hp_fe_values_;
748  nbr_control_hp_fe_values_ = NULL;
749  delete state_hp_fe_subface_values_;
750  state_hp_fe_subface_values_ = NULL;
751  delete control_hp_fe_subface_values_;
752  control_hp_fe_subface_values_ = NULL;
753  }
754 
755  /*********************************************/
756  /*
757  * This function reinits the FEFaceValues on the actual face. Should
758  * be called prior to any of the get-functions.
759  *
760  * @param face_no The 'local number' (i.e. from the perspective of the actual element) of the
761  * actual face.
762  */
763  inline void
764  ReInit(unsigned int face_no);
765 
766  /*********************************************/
767  /*
768  * This function reinits the FESubfaceValues on the actual subface. Should
769  * be called prior to any of the get-functions.
770  *
771  * @param face_no The 'local number' (i.e. from the perspective of the actual element) of the
772  * actual face.
773  * @param subface_no The 'local number' (i.e. from the perspective of the actual element) of the
774  * actual subface.
775  */
776  inline void
777  ReInit(unsigned int face_no, unsigned int subface_no);
778 
779  /*********************************************/
780  /*
781  * This function reinits the FEFaceValues on the neighbor_element for the
782  * case that the neighbor-element coarser or as fine as the
783  * (previously set!) actual element. This should be called prior
784  * to any of the get nbr_functions.
785  * Assumes that ReInit is called prior to this function.
786  */
787  inline void
788  ReInitNbr();
789 
790  /*********************************************/
795  inline unsigned int
796  GetNDoFsPerElement() const;
797 
798  inline unsigned int
799  GetNbrNDoFsPerElement() const;
800  inline unsigned int
801  GetNQPoints() const;
802  inline unsigned int
803  GetNbrNQPoints() const;
804  inline unsigned int
805  GetMaterialId() const;
806  inline unsigned int
807  GetNbrMaterialId() const;
808  inline unsigned int
809  GetNbrMaterialId(unsigned int face) const;
810  inline bool
811  GetIsAtBoundary() const;
812  inline unsigned int
813  GetBoundaryIndicator() const;
814  inline double
815  GetElementDiameter() const;
816 
817  inline const FEFaceValuesBase<dim>&
818  GetFEFaceValuesState() const;
819  inline const FEFaceValuesBase<dim>&
820  GetFEFaceValuesControl() const;
821 
822  inline const FEFaceValuesBase<dim>&
823  GetNbrFEFaceValuesState() const;
824  inline const FEFaceValuesBase<dim>&
825  GetNbrFEFaceValuesControl() const;
826 
827  private:
828  inline unsigned int
829  GetStateIndex() const;
830  inline unsigned int
831  GetControlIndex() const;
832 
836  template<class STH>
837  void
838  PrivateConstructor(const hp::QCollection<dim - 1>& q_collection,
839  UpdateFlags update_flags, STH& sth, bool need_neighbour)
840  {
841  if (need_neighbour) //so we need FEFAcevalues for the neighbour too.
842  {
843  nbr_state_hp_fe_values_ = new DOpEWrapper::HpFEFaceValues<dim>(
844  sth.GetMapping(), (sth.GetFESystem("state")), q_collection,
845  update_flags);
846  state_hp_fe_subface_values_ = new DOpEWrapper::HpFESubfaceValues<
847  dim>(sth.GetMapping(), (sth.GetFESystem("state")),
848  q_collection, update_flags);
849 
850  }
851  else
852  {
853  nbr_state_hp_fe_values_ = NULL;
854  nbr_control_hp_fe_values_ = NULL;
855  state_hp_fe_subface_values_ = NULL;
856  control_hp_fe_subface_values_ = NULL;
857  }
858 
859  nbr_state_hp_fe_values_ptr_ = NULL;
860  nbr_control_hp_fe_values_ptr_ = NULL;
861  state_hp_fe_values_ptr_ = NULL;
862  control_hp_fe_values_ptr_ = NULL;
863  }
864  /***********************************************************/
865  //"global" member data, part of every instantiation
866  unsigned int state_index_;
867  unsigned int control_index_;
868  const std::vector<
869  typename DOpEWrapper::DoFHandler<dim, dealii::hp::DoFHandler>::active_cell_iterator>& element_;
870 
871  DOpEWrapper::HpFEFaceValues<dim> state_hp_fe_values_;
872  DOpEWrapper::HpFEFaceValues<dim> control_hp_fe_values_;
873 
874  DOpEWrapper::HpFEFaceValues<dim>* nbr_state_hp_fe_values_;
875  DOpEWrapper::HpFEFaceValues<dim>* nbr_control_hp_fe_values_;
876 
877  DOpEWrapper::HpFESubfaceValues<dim>* state_hp_fe_subface_values_;
878  DOpEWrapper::HpFESubfaceValues<dim>* control_hp_fe_subface_values_;
879 
880  const dealii::FEFaceValuesBase<dim>* state_hp_fe_values_ptr_;
881  const dealii::FEFaceValuesBase<dim>* control_hp_fe_values_ptr_;
882  const dealii::FEFaceValuesBase<dim>* nbr_state_hp_fe_values_ptr_;
883  const dealii::FEFaceValuesBase<dim>* nbr_control_hp_fe_values_ptr_;
884 
885  const hp::QCollection<dim - 1>& q_collection_;
886  };
887 
888  /***********************************************************************/
889  /************************IMPLEMENTATION*for*DoFHandler*********************************/
890  /***********************************************************************/
891 
892  namespace {
893 //The following will only be used in debug mode; switch of warning of unused functions
894 #pragma GCC diagnostic push
895 #pragma GCC diagnostic ignored "-Wunused-function"
896 
897  template<int dim, template<int, int> class DH>
898  bool sanity_check(const
899  typename DOpEWrapper::DoFHandler<dim, DH>::active_cell_iterator & element_,
900  unsigned int face,
901  unsigned int subface)
902  {
903  const auto neighbor_child =
904  element_->neighbor_child_on_subface(
905  face, subface);
906 
907  bool ret = false;
908  if(neighbor_child->face(element_->neighbor_of_neighbor(face)) == element_->face(face)->child(subface))
909  ret = true;
910  return ret;
911  }
912 
913  template<>
914  bool sanity_check<1,dealii::hp::DoFHandler>(const
915  typename DOpEWrapper::DoFHandler<1, dealii::hp::DoFHandler>::active_cell_iterator &,
916  unsigned int,
917  unsigned int)
918  {
919  return true;
920  }
921 
922  template<>
923  bool sanity_check<1,dealii::DoFHandler>(const
924  typename DOpEWrapper::DoFHandler<1, dealii::DoFHandler>::active_cell_iterator &,
925  unsigned int,
926  unsigned int)
927  {
928  return true;
929  }
930 //Reenable warning or unused functions
931 #pragma GCC diagnostic pop
932  }
933 
934  template<typename VECTOR, int dim>
935  void
936  FaceDataContainer<dealii::DoFHandler, VECTOR, dim>::ReInit(
937  unsigned int face_no)
938  {
939  this->SetFace(face_no);
940  state_fe_values_.reinit(element_[this->GetStateIndex()], face_no);
941  state_fe_values_ptr_ = &state_fe_values_;
942  //Make sure that the Control must be initialized.
943  if (this->GetControlIndex() < element_.size())
944  {
945  control_fe_values_.reinit(element_[this->GetControlIndex()], face_no);
946  control_fe_values_ptr_ = &control_fe_values_;
947  }
948  }
949 
950  /***********************************************************************/
951 
952  template<typename VECTOR, int dim>
953  void
954  FaceDataContainer<dealii::DoFHandler, VECTOR, dim>::ReInit(
955  unsigned int face_no, unsigned int subface_no)
956  {
957  this->SetFace(face_no);
958  this->SetSubFace(subface_no);
959  state_fe_subface_values_->reinit(element_[this->GetStateIndex()], face_no,
960  subface_no);
961  state_fe_values_ptr_ = state_fe_subface_values_;
962  //Make sure that the Control must be initialized.
963  if (this->GetControlIndex() < element_.size())
964  {
965  control_fe_subface_values_->reinit(element_[this->GetControlIndex()],
966  face_no, this->GetSubFace());
967  control_fe_values_ptr_ = control_fe_subface_values_;
968  }
969  }
970  /***********************************************************************/
971 
972 
973  template<typename VECTOR, int dim>
974  void
975  FaceDataContainer<dealii::DoFHandler, VECTOR, dim>::ReInitNbr()
976  {
977  Assert(this->NeedNeighbour(), ExcInternalError());
978  Assert(
979  element_[this->GetStateIndex()]->neighbor_index(this->GetFace()) != -1,
980  TriaAccessorExceptions::ExcUnusedCellAsNeighbor())
981 
982  if (element_[this->GetStateIndex()]->neighbor(this->GetFace())->has_children())
983  {
984  //if neighbor is more refined
985  const auto neighbor_child =
986  element_[this->GetStateIndex()]->neighbor_child_on_subface(
987  this->GetFace(), this->GetSubFace());
988 
989  // some sanity checks: Check, that the face and subface match and that the neighbour child
990  // is not more refined.
991  Assert((sanity_check<dim, dealii::DoFHandler>(element_[this->GetStateIndex()],
992  this->GetFace(),
993  this->GetSubFace()) == true), ExcInternalError());
994  Assert(neighbor_child->has_children() == false, ExcInternalError());
995 
996  nbr_state_fe_values_->reinit(neighbor_child,
997  element_[this->GetStateIndex()]->neighbor_of_neighbor(
998  this->GetFace()));
999  nbr_state_fe_values_ptr_ = nbr_state_fe_values_;
1000 
1001  //Make sure that the Control must be initialized.
1002  if (this->GetControlIndex() < element_.size())
1003  {
1004  const auto control_neighbor_child =
1005  element_[this->GetControlIndex()]->neighbor_child_on_subface(
1006  this->GetFace(), this->GetSubFace());
1007 
1008  nbr_control_fe_values_->reinit(control_neighbor_child,
1009  element_[this->GetControlIndex()]->neighbor_of_neighbor(
1010  this->GetFace()));
1011  nbr_control_fe_values_ptr_ = nbr_control_fe_values_;
1012  }
1013  }
1014  else if (element_[this->GetStateIndex()]->neighbor_is_coarser(
1015  this->GetFace()))
1016  {
1017  //if the neighbour is coarser
1018  Assert(
1019  element_[this->GetStateIndex()]->neighbor(this->GetFace())->level() == element_[this->GetStateIndex()]->level()-1,
1020  ExcInternalError());
1021  const auto neighbor = element_[this->GetStateIndex()]->neighbor(
1022  this->GetFace());
1023  const std::pair<unsigned int, unsigned int> faceno_subfaceno =
1024  element_[this->GetStateIndex()]->neighbor_of_coarser_neighbor(
1025  this->GetFace());
1026  const unsigned int neighbor_face_no = faceno_subfaceno.first,
1027  neighbor_subface_no = faceno_subfaceno.second;
1028  state_fe_subface_values_->reinit(neighbor, neighbor_face_no,
1029  neighbor_subface_no);
1030  nbr_state_fe_values_ptr_ = state_fe_subface_values_;
1031  if (this->GetControlIndex() < element_.size())
1032  {
1033  const auto control_neighbor =
1034  element_[this->GetControlIndex()]->neighbor(this->GetFace());
1035  const std::pair<unsigned int, unsigned int> control_faceno_subfaceno =
1036  element_[this->GetControlIndex()]->neighbor_of_coarser_neighbor(
1037  this->GetFace());
1038  const unsigned int control_neighbor_face_no =
1039  control_faceno_subfaceno.first, control_neighbor_subface_no =
1040  control_faceno_subfaceno.second;
1041  control_fe_subface_values_->reinit(control_neighbor,
1042  control_neighbor_face_no, control_neighbor_subface_no);
1043  nbr_control_fe_values_ptr_ = control_fe_subface_values_;
1044  }
1045 
1046  }
1047  else
1048  {
1049  const auto neighbor_state = element_[this->GetStateIndex()]->neighbor(
1050  this->GetFace());
1051  // neighbor element is as much refined as the
1052  Assert(neighbor_state->level() == element_[this->GetStateIndex()]->level(),
1053  ExcInternalError());
1054  nbr_state_fe_values_->reinit(neighbor_state,
1055  element_[this->GetStateIndex()]->neighbor_of_neighbor(
1056  this->GetFace()));
1057  nbr_state_fe_values_ptr_ = nbr_state_fe_values_;
1058 
1059  //Make sure that the Control must be initialized.
1060  if (this->GetControlIndex() < element_.size())
1061  {
1062  nbr_control_fe_values_->reinit(
1063  element_[this->GetControlIndex()]->neighbor(this->GetFace()),
1064  element_[this->GetControlIndex()]->neighbor_of_neighbor(
1065  this->GetFace()));
1066  nbr_control_fe_values_ptr_ = nbr_control_fe_values_;
1067  }
1068  }
1069  }
1070  /***********************************************************************/
1071  template<typename VECTOR, int dim>
1072  unsigned int
1073  FaceDataContainer<dealii::DoFHandler, VECTOR, dim>::GetNDoFsPerElement() const
1074  {
1075  return n_dofs_per_element_;
1076  }
1077 
1078  /***********************************************************************/
1079 
1080  template<typename VECTOR, int dim>
1081  unsigned int
1082  FaceDataContainer<dealii::DoFHandler, VECTOR, dim>::GetNbrNDoFsPerElement() const
1083  {
1084  return n_dofs_per_element_;
1085  }
1086 
1087  /**********************************************/
1088  template<typename VECTOR, int dim>
1089  unsigned int
1090  FaceDataContainer<dealii::DoFHandler, VECTOR, dim>::GetNQPoints() const
1091  {
1092  return n_q_points_per_element_;
1093  }
1094 
1095  /**********************************************/
1096  template<typename VECTOR, int dim>
1097  unsigned int
1098  FaceDataContainer<dealii::DoFHandler, VECTOR, dim>::GetNbrNQPoints() const
1099  {
1100  return n_q_points_per_element_;
1101  }
1102 
1103  /**********************************************/
1104  template<typename VECTOR, int dim>
1105  unsigned int
1106  FaceDataContainer<dealii::DoFHandler, VECTOR, dim>::GetMaterialId() const
1107  {
1108  return element_[0]->material_id();
1109  }
1110 
1111  /**********************************************/
1112  template<typename VECTOR, int dim>
1113  unsigned int
1114  FaceDataContainer<dealii::DoFHandler, VECTOR, dim>::GetNbrMaterialId() const
1115  {
1116  return this->GetNbrMaterialId(this->GetFace());
1117  }
1118 
1119  /**********************************************/
1120  template<typename VECTOR, int dim>
1121  unsigned int
1122  FaceDataContainer<dealii::DoFHandler, VECTOR, dim>::GetNbrMaterialId(
1123  unsigned int face) const
1124  {
1125  if (element_[0]->neighbor_index(face) != -1)
1126  return element_[0]->neighbor(face)->material_id();
1127  else
1128  {
1129  std::stringstream out;
1130  out << "There is no neighbor with number " << face;
1131  throw DOpEException(out.str(),
1132  "FaceDataContainer::GetNbrMaterialId");
1133  }
1134  }
1135 
1136  /**********************************************/
1137  template<typename VECTOR, int dim>
1138  bool
1139  FaceDataContainer<dealii::DoFHandler, VECTOR, dim>::GetIsAtBoundary() const
1140  {
1141  return element_[0]->face(this->GetFace())->at_boundary();
1142  }
1143 
1144  /**********************************************/
1145  template<typename VECTOR, int dim>
1146  double
1147  FaceDataContainer<dealii::DoFHandler, VECTOR, dim>::GetElementDiameter() const
1148  {
1149  return element_[0]->face(this->GetFace())->diameter();
1150  }
1151 
1152  /**********************************************/
1153 
1154  template<typename VECTOR, int dim>
1155  unsigned int
1156  FaceDataContainer<dealii::DoFHandler, VECTOR, dim>::GetBoundaryIndicator() const
1157  {
1158 #if DEAL_II_VERSION_GTE(8,3,0)
1159  return element_[0]->face(this->GetFace())->boundary_id();
1160 #else
1161  return element_[0]->face(this->GetFace())->boundary_indicator();
1162 #endif
1163  }
1164 
1165  /**********************************************/
1166  template<typename VECTOR, int dim>
1167  const FEFaceValuesBase<dim>&
1168  FaceDataContainer<dealii::DoFHandler, VECTOR, dim>::GetFEFaceValuesState() const
1169  {
1170  return *state_fe_values_ptr_;
1171  }
1172 
1173  /**********************************************/
1174  template<typename VECTOR, int dim>
1175  const FEFaceValuesBase<dim>&
1176  FaceDataContainer<dealii::DoFHandler, VECTOR, dim>::GetFEFaceValuesControl() const
1177  {
1178  return *control_fe_values_ptr_;
1179  }
1180  /**********************************************/
1181  template<typename VECTOR, int dim>
1182  const FEFaceValuesBase<dim>&
1183  FaceDataContainer<dealii::DoFHandler, VECTOR, dim>::GetNbrFEFaceValuesState() const
1184  {
1185  return *nbr_state_fe_values_ptr_;
1186  }
1187 
1188  /**********************************************/
1189  template<typename VECTOR, int dim>
1190  const FEFaceValuesBase<dim>&
1191  FaceDataContainer<dealii::DoFHandler, VECTOR, dim>::GetNbrFEFaceValuesControl() const
1192  {
1193  return *nbr_control_fe_values_ptr_;
1194  }
1195 
1196  template<typename VECTOR, int dim>
1197  unsigned int
1198  FaceDataContainer<dealii::DoFHandler, VECTOR, dim>::GetStateIndex() const
1199  {
1200  return state_index_;
1201  }
1202 
1203  /***********************************************************************/
1204 
1205  template<typename VECTOR, int dim>
1206  unsigned int
1207  FaceDataContainer<dealii::DoFHandler, VECTOR, dim>::GetControlIndex() const
1208  {
1209  return control_index_;
1210  }
1211 
1212  /***********************************************************************/
1213  /************************END*OF*IMPLEMENTATION**************************/
1214  /***********************************************************************/
1215  /***********************************************************************/
1216  /*****************IMPLEMENTATION for MGDoFHandler*********************/
1217  /***********************************************************************/
1218 //
1219 //
1220 //
1221 // template<typename VECTOR, int dim>
1222 // void
1223 // FaceDataContainer<dealii::MGDoFHandler, VECTOR, dim>::ReInit(
1224 // unsigned int face_no)
1225 // {
1226 // this->SetFace(face_no);
1227 // state_fe_values_.reinit(element_[this->GetStateIndex()], face_no);
1228 // state_fe_values_ptr_ = &state_fe_values_;
1229 // //Make sure that the Control must be initialized.
1230 // if (this->GetControlIndex() < element_.size())
1231 // {
1232 // control_fe_values_.reinit(element_[this->GetControlIndex()], face_no);
1233 // control_fe_values_ptr_ = &control_fe_values_;
1234 // }
1235 // }
1236 //
1237 // /***********************************************************************/
1238 //
1239 // template<typename VECTOR, int dim>
1240 // void
1241 // FaceDataContainer<dealii::MGDoFHandler, VECTOR, dim>::ReInit(
1242 // unsigned int face_no, unsigned int subface_no)
1243 // {
1244 // this->SetFace(face_no);
1245 // this->SetSubFace(subface_no);
1246 // state_fe_subface_values_->reinit(element_[this->GetStateIndex()], face_no,
1247 // subface_no);
1248 // state_fe_values_ptr_ = state_fe_subface_values_;
1249 // //Make sure that the Control must be initialized.
1250 // if (this->GetControlIndex() < element_.size())
1251 // {
1252 // control_fe_subface_values_->reinit(element_[this->GetControlIndex()],
1253 // face_no, this->GetSubFace());
1254 // control_fe_values_ptr_ = control_fe_subface_values_;
1255 // }
1256 // }
1257 // /***********************************************************************/
1258 //
1259 // template<typename VECTOR, int dim>
1260 // void
1261 // FaceDataContainer<dealii::MGDoFHandler, VECTOR, dim>::ReInitNbr()
1262 // {
1263 // Assert(this->NeedNeighbour(), ExcInternalError());
1264 // Assert(
1265 // element_[this->GetStateIndex()]->neighbor_index(this->GetFace()) != -1,
1266 // TriaAccessorExceptions::ExcUnusedCellAsNeighbor())
1267 //
1268 // if (element_[this->GetStateIndex()]->neighbor(this->GetFace())->has_children())
1269 // {
1270 // //if neighbor is more refined
1271 // const auto neighbor_child =
1272 // element_[this->GetStateIndex()]->neighbor_child_on_subface(
1273 // this->GetFace(), this->GetSubFace());
1274 //
1275 // // some sanity checks: Check, that the face and subface match and that the neighbour child
1276 // // is not more refined.
1277 // Assert((sanity_check<dim, dealii::hp::DoFHandler>(element_[this->GetStateIndex()],
1278 // this->GetFace(),
1279 // this->GetSubFace()) == true), ExcInternalError());
1280 // Assert(neighbor_child->has_children() == false, ExcInternalError());
1281 //
1282 // nbr_state_fe_values_->reinit(neighbor_child,
1283 // element_[this->GetStateIndex()]->neighbor_of_neighbor(
1284 // this->GetFace()));
1285 // nbr_state_fe_values_ptr_ = nbr_state_fe_values_;
1286 //
1287 // //Make sure that the Control must be initialized.
1288 // if (this->GetControlIndex() < element_.size())
1289 // {
1290 // const auto control_neighbor_child =
1291 // element_[this->GetControlIndex()]->neighbor_child_on_subface(
1292 // this->GetFace(), this->GetSubFace());
1293 //
1294 // nbr_control_fe_values_->reinit(control_neighbor_child,
1295 // element_[this->GetControlIndex()]->neighbor_of_neighbor(
1296 // this->GetFace()));
1297 // nbr_control_fe_values_ptr_ = nbr_control_fe_values_;
1298 // }
1299 // }
1300 // else if (element_[this->GetStateIndex()]->neighbor_is_coarser(
1301 // this->GetFace()))
1302 // {
1303 // //if the neighbour is coarser
1304 // Assert(
1305 // element_[this->GetStateIndex()]->neighbor(this->GetFace())->level() == element_[this->GetStateIndex()]->level()-1,
1306 // ExcInternalError());
1307 // const auto neighbor = element_[this->GetStateIndex()]->neighbor(
1308 // this->GetFace());
1309 // const std::pair<unsigned int, unsigned int> faceno_subfaceno =
1310 // element_[this->GetStateIndex()]->neighbor_of_coarser_neighbor(
1311 // this->GetFace());
1312 // const unsigned int neighbor_face_no = faceno_subfaceno.first,
1313 // neighbor_subface_no = faceno_subfaceno.second;
1314 // state_fe_subface_values_->reinit(neighbor, neighbor_face_no,
1315 // neighbor_subface_no);
1316 // nbr_state_fe_values_ptr_ = state_fe_subface_values_;
1317 // if (this->GetControlIndex() < element_.size())
1318 // {
1319 // const auto control_neighbor =
1320 // element_[this->GetControlIndex()]->neighbor(this->GetFace());
1321 // const std::pair<unsigned int, unsigned int> control_faceno_subfaceno =
1322 // element_[this->GetControlIndex()]->neighbor_of_coarser_neighbor(
1323 // this->GetFace());
1324 // const unsigned int control_neighbor_face_no =
1325 // control_faceno_subfaceno.first, control_neighbor_subface_no =
1326 // control_faceno_subfaceno.second;
1327 // control_fe_subface_values_->reinit(control_neighbor,
1328 // control_neighbor_face_no, control_neighbor_subface_no);
1329 // nbr_control_fe_values_ptr_ = control_fe_subface_values_;
1330 // }
1331 //
1332 // }
1333 // else
1334 // {
1335 // const auto neighbor_state = element_[this->GetStateIndex()]->neighbor(
1336 // this->GetFace());
1337 // // neighbor element is as much refined as the
1338 // Assert(neighbor_state->level() == element_[this->GetStateIndex()]->level(),
1339 // ExcInternalError());
1340 // nbr_state_fe_values_->reinit(neighbor_state,
1341 // element_[this->GetStateIndex()]->neighbor_of_neighbor(
1342 // this->GetFace()));
1343 // nbr_state_fe_values_ptr_ = nbr_state_fe_values_;
1344 //
1345 // //Make sure that the Control must be initialized.
1346 // if (this->GetControlIndex() < element_.size())
1347 // {
1348 // nbr_control_fe_values_->reinit(
1349 // element_[this->GetControlIndex()]->neighbor(this->GetFace()),
1350 // element_[this->GetControlIndex()]->neighbor_of_neighbor(
1351 // this->GetFace()));
1352 // nbr_control_fe_values_ptr_ = nbr_control_fe_values_;
1353 // }
1354 // }
1355 // }
1356 // /***********************************************************************/
1357 // template<typename VECTOR, int dim>
1358 // unsigned int
1359 // FaceDataContainer<dealii::MGDoFHandler, VECTOR, dim>::GetNDoFsPerElement() const
1360 // {
1361 // return n_dofs_per_element_;
1362 // }
1363 //
1364 // /***********************************************************************/
1365 //
1366 // template<typename VECTOR, int dim>
1367 // unsigned int
1368 // FaceDataContainer<dealii::MGDoFHandler, VECTOR, dim>::GetNbrNDoFsPerElement() const
1369 // {
1370 // return n_dofs_per_element_;
1371 // }
1372 //
1373 // /**********************************************/
1374 // template<typename VECTOR, int dim>
1375 // unsigned int
1376 // FaceDataContainer<dealii::MGDoFHandler, VECTOR, dim>::GetNQPoints() const
1377 // {
1378 // return n_q_points_per_element_;
1379 // }
1380 //
1381 // /**********************************************/
1382 // template<typename VECTOR, int dim>
1383 // unsigned int
1384 // FaceDataContainer<dealii::MGDoFHandler, VECTOR, dim>::GetNbrNQPoints() const
1385 // {
1386 // return n_q_points_per_element_;
1387 // }
1388 //
1389 // /**********************************************/
1390 // template<typename VECTOR, int dim>
1391 // unsigned int
1392 // FaceDataContainer<dealii::MGDoFHandler, VECTOR, dim>::GetMaterialId() const
1393 // {
1394 // return element_[0]->material_id();
1395 // }
1396 //
1397 // /**********************************************/
1398 // template<typename VECTOR, int dim>
1399 // unsigned int
1400 // FaceDataContainer<dealii::MGDoFHandler, VECTOR, dim>::GetNbrMaterialId() const
1401 // {
1402 // return this->GetNbrMaterialId(this->GetFace());
1403 // }
1404 //
1405 // /**********************************************/
1406 // template<typename VECTOR, int dim>
1407 // unsigned int
1408 // FaceDataContainer<dealii::MGDoFHandler, VECTOR, dim>::GetNbrMaterialId(
1409 // unsigned int face) const
1410 // {
1411 // if (element_[0]->neighbor_index(face) != -1)
1412 // return element_[0]->neighbor(face)->material_id();
1413 // else
1414 // {
1415 // std::stringstream out;
1416 // out << "There is no neighbor with number " << face;
1417 // throw DOpEException(out.str(),
1418 // "FaceDataContainer::GetNbrMaterialId");
1419 // }
1420 // }
1421 //
1422 // /**********************************************/
1423 // template<typename VECTOR, int dim>
1424 // bool
1425 // FaceDataContainer<dealii::MGDoFHandler, VECTOR, dim>::GetIsAtBoundary() const
1426 // {
1427 // return element_[0]->face(this->GetFace())->at_boundary();
1428 // }
1429 //
1430 // /**********************************************/
1431 // template<typename VECTOR, int dim>
1432 // double
1433 // FaceDataContainer<dealii::MGDoFHandler, VECTOR, dim>::GetElementDiameter() const
1434 // {
1435 // return element_[0]->face(this->GetFace())->diameter();
1436 // }
1437 //
1438 // /**********************************************/
1439 //
1440 // template<typename VECTOR, int dim>
1441 // unsigned int
1442 // FaceDataContainer<dealii::MGDoFHandler, VECTOR, dim>::GetBoundaryIndicator() const
1443 // {
1444 // return element_[0]->face(this->GetFace())->boundary_indicator();
1445 // }
1446 //
1447 // /**********************************************/
1448 // template<typename VECTOR, int dim>
1449 // const FEFaceValuesBase<dim>&
1450 // FaceDataContainer<dealii::MGDoFHandler, VECTOR, dim>::GetFEFaceValuesState() const
1451 // {
1452 // return *state_fe_values_ptr_;
1453 // }
1454 //
1455 // /**********************************************/
1456 // template<typename VECTOR, int dim>
1457 // const FEFaceValuesBase<dim>&
1458 // FaceDataContainer<dealii::MGDoFHandler, VECTOR, dim>::GetFEFaceValuesControl() const
1459 // {
1460 // return *control_fe_values_ptr_;
1461 // }
1462 // /**********************************************/
1463 // template<typename VECTOR, int dim>
1464 // const FEFaceValuesBase<dim>&
1465 // FaceDataContainer<dealii::MGDoFHandler, VECTOR, dim>::GetNbrFEFaceValuesState() const
1466 // {
1467 // return *nbr_state_fe_values_ptr_;
1468 // }
1469 //
1470 // /**********************************************/
1471 // template<typename VECTOR, int dim>
1472 // const FEFaceValuesBase<dim>&
1473 // FaceDataContainer<dealii::MGDoFHandler, VECTOR, dim>::GetNbrFEFaceValuesControl() const
1474 // {
1475 // return *nbr_control_fe_values_ptr_;
1476 // }
1477 //
1478 // template<typename VECTOR, int dim>
1479 // unsigned int
1480 // FaceDataContainer<dealii::MGDoFHandler, VECTOR, dim>::GetStateIndex() const
1481 // {
1482 // return state_index_;
1483 // }
1484 //
1485 // /***********************************************************************/
1486 //
1487 // template<typename VECTOR, int dim>
1488 // unsigned int
1489 // FaceDataContainer<dealii::MGDoFHandler, VECTOR, dim>::GetControlIndex() const
1490 // {
1491 // return control_index_;
1492 // }
1493 //
1494  /***********************************************************************/
1495  /************************END*OF*IMPLEMENTATION**************************/
1496  /***********************************************************************/
1497  /***********************************************************************/
1498  /*****************IMPLEMENTATION for hp::DoFHandler*********************/
1499  /***********************************************************************/
1500 
1501 
1502 
1503 
1504 
1505  template<typename VECTOR, int dim>
1506  void
1507  FaceDataContainer<dealii::hp::DoFHandler, VECTOR, dim>::ReInit(
1508  unsigned int face_no)
1509  {
1510  this->SetFace(face_no);
1511  state_hp_fe_values_.reinit(element_[this->GetStateIndex()], face_no);
1512  state_hp_fe_values_ptr_ = &state_hp_fe_values_.get_present_fe_values();
1513  //Make sure that the Control must be initialized.
1514  if (this->GetControlIndex() < element_.size())
1515  {
1516  control_hp_fe_values_.reinit(element_[this->GetControlIndex()], face_no);
1517  control_hp_fe_values_ptr_ =
1518  &control_hp_fe_values_.get_present_fe_values();
1519  }
1520  }
1521 
1522  /***********************************************************************/
1523 
1524  template<typename VECTOR, int dim>
1525  void
1526  FaceDataContainer<dealii::hp::DoFHandler, VECTOR, dim>::ReInit(
1527  unsigned int face_no, unsigned int subface_no)
1528  {
1529  this->SetFace(face_no);
1530  this->SetSubFace(subface_no);
1531  state_hp_fe_subface_values_->reinit(element_[this->GetStateIndex()],
1532  this->GetFace(), this->GetSubFace());
1533  state_hp_fe_values_ptr_ =
1534  &state_hp_fe_subface_values_->get_present_fe_values();
1535  //Make sure that the Control must be initialized.
1536  if (this->GetControlIndex() < element_.size())
1537  {
1538  control_hp_fe_subface_values_->reinit(element_[this->GetControlIndex()],
1539  this->GetFace(), this->GetSubFace());
1540  control_hp_fe_values_ptr_ =
1541  &control_hp_fe_subface_values_->get_present_fe_values();
1542  }
1543  }
1544 
1545  /***********************************************************************/
1546  template<typename VECTOR, int dim>
1547  void
1548  FaceDataContainer<dealii::hp::DoFHandler, VECTOR, dim>::ReInitNbr()
1549  {
1550  Assert(this->NeedNeighbour(), ExcInternalError());
1551  Assert(
1552  element_[this->GetStateIndex()]->neighbor_index(this->GetFace()) != -1,
1553  TriaAccessorExceptions::ExcUnusedCellAsNeighbor())
1554 
1555  if (element_[this->GetStateIndex()]->neighbor(this->GetFace())->has_children())
1556  {
1557  //if neighbor is more refined
1558  const auto neighbor_child =
1559  element_[this->GetStateIndex()]->neighbor_child_on_subface(
1560  this->GetFace(), this->GetSubFace());
1561 
1562  // some sanity checks: Check, that the face and subface match and that the neighbour child
1563  // is not more refined.
1564  Assert(
1565  neighbor_child->face(element_[this->GetStateIndex()]->neighbor_of_neighbor(this->GetFace())) == element_[this->GetStateIndex()]->face(this->GetFace())->child(this->GetSubFace()),
1566  ExcInternalError());
1567  Assert(neighbor_child->has_children() == false, ExcInternalError());
1568 
1569  nbr_state_hp_fe_values_->reinit(neighbor_child,
1570  element_[this->GetStateIndex()]->neighbor_of_neighbor(
1571  this->GetFace()));
1572  nbr_state_hp_fe_values_ptr_ =
1573  &nbr_state_hp_fe_values_->get_present_fe_values();
1574 
1575  //Make sure that the Control must be initialized.
1576  if (this->GetControlIndex() < element_.size())
1577  {
1578  const auto control_neighbor_child =
1579  element_[this->GetControlIndex()]->neighbor_child_on_subface(
1580  this->GetFace(), this->GetSubFace());
1581 
1582  nbr_control_hp_fe_values_->reinit(control_neighbor_child,
1583  element_[this->GetControlIndex()]->neighbor_of_neighbor(
1584  this->GetFace()));
1585  nbr_control_hp_fe_values_ptr_ =
1586  &nbr_control_hp_fe_values_->get_present_fe_values();
1587  }
1588  }
1589  else if (element_[this->GetStateIndex()]->neighbor_is_coarser(
1590  this->GetFace()))
1591  {
1592  //if the neighbour is coarser
1593  Assert(
1594  element_[this->GetStateIndex()]->neighbor(this->GetFace())->level() == element_[this->GetStateIndex()]->level()-1,
1595  ExcInternalError());
1596  const auto neighbor = element_[this->GetStateIndex()]->neighbor(
1597  this->GetFace());
1598  const std::pair<unsigned int, unsigned int> faceno_subfaceno =
1599  element_[this->GetStateIndex()]->neighbor_of_coarser_neighbor(
1600  this->GetFace());
1601  const unsigned int neighbor_face_no = faceno_subfaceno.first,
1602  neighbor_subface_no = faceno_subfaceno.second;
1603  state_hp_fe_subface_values_->reinit(neighbor, neighbor_face_no,
1604  neighbor_subface_no);
1605  nbr_state_hp_fe_values_ptr_ =
1606  &state_hp_fe_subface_values_->get_present_fe_values();
1607  if (this->GetControlIndex() < element_.size())
1608  {
1609  const auto control_neighbor =
1610  element_[this->GetControlIndex()]->neighbor(this->GetFace());
1611  const std::pair<unsigned int, unsigned int> control_faceno_subfaceno =
1612  element_[this->GetControlIndex()]->neighbor_of_coarser_neighbor(
1613  this->GetFace());
1614  const unsigned int control_neighbor_face_no =
1615  control_faceno_subfaceno.first, control_neighbor_subface_no =
1616  control_faceno_subfaceno.second;
1617 
1618  control_hp_fe_subface_values_->reinit(control_neighbor,
1619  control_neighbor_face_no, control_neighbor_subface_no);
1620  nbr_control_hp_fe_values_ptr_ =
1621  &control_hp_fe_subface_values_->get_present_fe_values();
1622  }
1623  }
1624  else
1625  {
1626  const auto neighbor_state = element_[this->GetStateIndex()]->neighbor(
1627  this->GetFace());
1628  // neighbor element is as much refined as the
1629  Assert(neighbor_state->level() == element_[this->GetStateIndex()]->level(),
1630  ExcInternalError());
1631  nbr_state_hp_fe_values_->reinit(neighbor_state,
1632  element_[this->GetStateIndex()]->neighbor_of_neighbor(
1633  this->GetFace()));
1634  nbr_state_hp_fe_values_ptr_ =
1635  &nbr_state_hp_fe_values_->get_present_fe_values();
1636 
1637  //Make sure that the Control must be initialized.
1638  if (this->GetControlIndex() < element_.size())
1639  {
1640  nbr_control_hp_fe_values_->reinit(
1641  element_[this->GetControlIndex()]->neighbor(this->GetFace()),
1642  element_[this->GetControlIndex()]->neighbor_of_neighbor(
1643  this->GetFace()));
1644  nbr_control_hp_fe_values_ptr_ =
1645  &nbr_control_hp_fe_values_->get_present_fe_values();
1646  }
1647  }
1648  }
1649 
1650  /*********************************************/
1651 
1652  template<typename VECTOR, int dim>
1653  unsigned int
1654  FaceDataContainer<dealii::hp::DoFHandler, VECTOR, dim>::GetNDoFsPerElement() const
1655  {
1656  return element_[0]->get_fe().dofs_per_cell;
1657  }
1658 
1659  /*********************************************/
1660 
1661  template<typename VECTOR, int dim>
1662  unsigned int
1663  FaceDataContainer<dealii::hp::DoFHandler, VECTOR, dim>::GetNbrNDoFsPerElement() const
1664  {
1665  if (element_[0]->neighbor_index(this->GetFace()) != -1)
1666  return element_[0]->neighbor(this->GetFace())->get_fe().dofs_per_cell;
1667  else
1668  throw DOpEException(
1669  "There is no neighbor with number" + this->GetFace(),
1670  "HpFaceDataContainer::GetNbrNDoFsPerElement");
1671  }
1672  /*********************************************/
1673  template<typename VECTOR, int dim>
1674  unsigned int
1675  FaceDataContainer<dealii::hp::DoFHandler, VECTOR, dim>::GetNQPoints() const
1676  {
1677  return q_collection_[element_[0]->active_fe_index()].size();
1678  }
1679 
1680  /*********************************************/
1681  template<typename VECTOR, int dim>
1682  unsigned int
1683  FaceDataContainer<dealii::hp::DoFHandler, VECTOR, dim>::GetNbrNQPoints() const
1684  {
1685  if (element_[0]->neighbor_index(this->GetFace()) != -1)
1686  return q_collection_[element_[0]->neighbor(this->GetFace())->active_fe_index()].size();
1687  else
1688  throw DOpEException(
1689  "There is no neighbor with number" + this->GetFace(),
1690  "HpFaceDataContainer::GetNbrNQPoints");
1691  }
1692  /*********************************************/
1693 
1694  template<typename VECTOR, int dim>
1695  unsigned int
1696  FaceDataContainer<dealii::hp::DoFHandler, VECTOR, dim>::GetMaterialId() const
1697  {
1698  return element_[0]->material_id();
1699  }
1700  /*********************************************/
1701  template<typename VECTOR, int dim>
1702  unsigned int
1703  FaceDataContainer<dealii::hp::DoFHandler, VECTOR, dim>::GetNbrMaterialId() const
1704  {
1705  return this->GetNbrMaterialId(this->GetFace());
1706  }
1707 
1708  /*********************************************/
1709 
1710  template<typename VECTOR, int dim>
1711  unsigned int
1712  FaceDataContainer<dealii::hp::DoFHandler, VECTOR, dim>::GetNbrMaterialId(
1713  unsigned int face) const
1714  {
1715  if (element_[0]->neighbor_index(face) != -1)
1716  return element_[0]->neighbor(face)->material_id();
1717  else
1718  {
1719  std::stringstream out;
1720  out << "There is no neighbor with number " << face;
1721  throw DOpEException(out.str(),
1722  "HpFaceDataContainer::GetNbrMaterialId");
1723  }
1724  }
1725 
1726  /*********************************************/
1727 
1728  template<typename VECTOR, int dim>
1729  double
1730  FaceDataContainer<dealii::hp::DoFHandler, VECTOR, dim>::GetElementDiameter() const
1731  {
1732  return element_[0]->face(this->GetFace())->diameter();
1733  }
1734 
1735  /**********************************************/
1736 
1737  template<typename VECTOR, int dim>
1738  bool
1739  FaceDataContainer<dealii::hp::DoFHandler, VECTOR, dim>::GetIsAtBoundary() const
1740  {
1741  return element_[0]->face(this->GetFace())->at_boundary();
1742  }
1743 
1744  /**********************************************/
1745  template<typename VECTOR, int dim>
1746  unsigned int
1747  FaceDataContainer<dealii::hp::DoFHandler, VECTOR, dim>::GetBoundaryIndicator() const
1748  {
1749 #if DEAL_II_VERSION_GTE(8,3,0)
1750  return element_[0]->face(this->GetFace())->boundary_id();
1751 #else
1752  return element_[0]->face(this->GetFace())->boundary_indicator();
1753 #endif
1754  }
1755 
1756  /*********************************************/
1757  template<typename VECTOR, int dim>
1758  const FEFaceValuesBase<dim>&
1759  FaceDataContainer<dealii::hp::DoFHandler, VECTOR, dim>::GetFEFaceValuesState() const
1760  {
1761  return *state_hp_fe_values_ptr_;
1762  }
1763  /*********************************************/
1764  template<typename VECTOR, int dim>
1765  const FEFaceValuesBase<dim>&
1766  FaceDataContainer<dealii::hp::DoFHandler, VECTOR, dim>::GetFEFaceValuesControl() const
1767  {
1768  return *control_hp_fe_values_ptr_;
1769  }
1770  /*********************************************/
1771  template<typename VECTOR, int dim>
1772  const FEFaceValuesBase<dim>&
1773  FaceDataContainer<dealii::hp::DoFHandler, VECTOR, dim>::GetNbrFEFaceValuesState() const
1774  {
1775  return *nbr_state_hp_fe_values_ptr_;
1776  }
1777  /*********************************************/
1778  template<typename VECTOR, int dim>
1779  const FEFaceValuesBase<dim>&
1780  FaceDataContainer<dealii::hp::DoFHandler, VECTOR, dim>::GetNbrFEFaceValuesControl() const
1781  {
1782  return *nbr_control_hp_fe_values_ptr_;
1783  }
1784 
1785  /**********************************************/
1786  template<typename VECTOR, int dim>
1787  unsigned int
1788  FaceDataContainer<dealii::hp::DoFHandler, VECTOR, dim>::GetStateIndex() const
1789  {
1790  return state_index_;
1791  }
1792  /*********************************************/
1793  template<typename VECTOR, int dim>
1794  unsigned int
1795  FaceDataContainer<dealii::hp::DoFHandler, VECTOR, dim>::GetControlIndex() const
1796  {
1797  return control_index_;
1798  }
1799 } //end of namespace
1800 
1801 #endif /* FACEDATACONTAINER_H_ */
DOpEWrapper::FESubfaceValues
Definition: fevalues_wrapper.h:136
DOpE::FaceDataContainer< dealii::DoFHandler, VECTOR, dim >::~FaceDataContainer
~FaceDataContainer()
Definition: facedatacontainer.h:183
DOpE::fdcinternal::FaceDataContainerInternal
Definition: facedatacontainer_internal.h:46
DOpE::FaceDataContainer< dealii::DoFHandler, VECTOR, dim >::FaceDataContainer
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:151
facedatacontainer_internal.h
DOpE::FaceDataContainer::FaceDataContainer
FaceDataContainer()
Definition: facedatacontainer.h:54
DOpEWrapper::HpFEFaceValues
Definition: fevalues_wrapper.h:218
DOpE::FaceDataContainer
Definition: facedatacontainer.h:50
DOpEWrapper::DoFHandler
Definition: dofhandler_wrapper.h:51
DOpE::StateSpaceTimeHandler
Definition: statespacetimehandler.h:59
DOpE::SpaceTimeHandler
Definition: spacetimehandler.h:71
DOpEWrapper::HpFESubfaceValues
Definition: fevalues_wrapper.h:259
fevalues_wrapper.h
DOpE::FaceDataContainer< dealii::DoFHandler, VECTOR, dim >::FaceDataContainer
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:100
DOpEWrapper::FEFaceValues
Definition: fevalues_wrapper.h:88
DOpE::DOpEException
Definition: dopeexception.h:35
Generated on Mon Jan 25 2016 11:32:18 for DOpE by   doxygen 1.8.6