24 #ifndef ELEMENTDATACONTAINER_H_
25 #define ELEMENTDATACONTAINER_H_
33 #include <dofs/dof_handler.h>
34 #include <multigrid/mg_dof_handler.h>
35 #include <hp/dof_handler.h>
37 using namespace dealii;
55 template<
template<
int,
int>
class DH,
typename VECTOR,
int dim>
63 "Dummy class, this constructor should never get called.",
64 "ElementDataContainer<dealii::DoFHandler<dim> , VECTOR, dim>::ElementDataContainer"));
78 template<
typename VECTOR,
int dim>
104 template<
template<
int,
int>
class FE,
typename SPARSITYPATTERN,
int dopedim,
int dealdim>
106 UpdateFlags update_flags,
108 dopedim, dealdim>& sth,
110 typename dealii::DoFHandler<dim>::active_cell_iterator>& element,
111 const std::map<std::string,
const Vector<double>*> ¶m_values,
112 const std::map<std::string, const VECTOR*> &domain_values) :
113 edcinternal::ElementDataContainerInternal<VECTOR, dim>(param_values,
114 domain_values), _element(element), _state_fe_values(
115 sth.GetMapping(), (sth.GetFESystem(
"state")), quad,
116 update_flags), _control_fe_values(sth.GetMapping(),
117 (sth.GetFESystem(
"control")), quad, update_flags)
119 _state_index = sth.GetStateIndex();
120 if (_state_index == 1)
124 _n_q_points_per_element = quad.size();
125 _n_dofs_per_element = element[0]->get_fe().dofs_per_cell;
146 template<
template<
int,
int>
class FE,
typename SPARSITYPATTERN>
148 UpdateFlags update_flags,
152 typename dealii::DoFHandler<dim>::active_cell_iterator>& element,
153 const std::map<std::string,
const Vector<double>*> ¶m_values,
154 const std::map<std::string, const VECTOR*> &domain_values) :
155 edcinternal::ElementDataContainerInternal<VECTOR, dim>(param_values,
156 domain_values), _element(element), _state_fe_values(
157 sth.GetMapping(), (sth.GetFESystem(
"state")), quad,
158 update_flags), _control_fe_values(sth.GetMapping(),
159 (sth.GetFESystem(
"state")), quad, update_flags)
161 _state_index = sth.GetStateIndex();
162 _control_index = element.size();
163 _n_q_points_per_element = quad.size();
164 _n_dofs_per_element = element[0]->get_fe().dofs_per_cell;
183 GetNDoFsPerElement()
const;
187 GetMaterialId()
const;
189 GetNbrMaterialId(
unsigned int face)
const;
191 GetFaceBoundaryIndicator(
unsigned int face)
const;
193 GetIsAtBoundary()
const;
195 GetElementDiameter()
const;
199 GetFEValuesState()
const;
201 GetFEValuesControl()
const;
207 GetStateIndex()
const;
209 GetControlIndex()
const;
213 unsigned int _state_index;
214 unsigned int _control_index;
216 const std::vector<typename dealii::DoFHandler<dim>::active_cell_iterator> & _element;
220 unsigned int _n_q_points_per_element;
221 unsigned int _n_dofs_per_element;
237 template<
typename VECTOR,
int dim>
263 template<
template<
int,
int>
class FE,
typename SPARSITYPATTERN,
int dopedim,
int dealdim>
265 const Quadrature<dim>& quad,
266 UpdateFlags update_flags,
268 VECTOR, dopedim, dealdim>& sth
271 typename dealii::MGDoFHandler<dim>::active_cell_iterator>& element,
272 const std::map<std::string,
const Vector<double>*> ¶m_values,
273 const std::map<std::string, const VECTOR*> &domain_values)
274 : edcinternal::ElementDataContainerInternal<VECTOR, dim>(
275 param_values, domain_values), _element(element), _state_fe_values(
276 sth.GetMapping(), (sth.GetFESystem(
"state")), quad,
277 update_flags), _control_fe_values(sth.GetMapping(),
278 (sth.GetFESystem(
"control")), quad, update_flags)
280 _state_index = sth.GetStateIndex();
281 if (_state_index == 1)
285 _n_q_points_per_element = quad.size();
286 _n_dofs_per_element = element[0]->get_fe().dofs_per_cell;
311 template<
template<
int,
int>
class FE,
typename SPARSITYPATTERN>
313 UpdateFlags update_flags,
315 SPARSITYPATTERN, VECTOR, dim>& sth,
317 typename dealii::MGDoFHandler<dim>::active_cell_iterator>& element,
318 const std::map<std::string,
const Vector<double>*> ¶m_values,
319 const std::map<std::string, const VECTOR*> &domain_values)
320 : edcinternal::ElementDataContainerInternal<VECTOR, dim>(
321 param_values, domain_values), _element(element), _state_fe_values(
322 sth.GetMapping(), (sth.GetFESystem(
"state")), quad,
323 update_flags), _control_fe_values(sth.GetMapping(),
324 (sth.GetFESystem(
"state")), quad, update_flags)
326 _state_index = sth.GetStateIndex();
327 _control_index = element.size();
328 _n_q_points_per_element = quad.size();
329 _n_dofs_per_element = element[0]->get_fe().dofs_per_cell;
350 GetNDoFsPerElement()
const;
354 GetMaterialId()
const;
356 GetNbrMaterialId(
unsigned int face)
const;
358 GetFaceBoundaryIndicator(
unsigned int face)
const;
360 GetIsAtBoundary()
const;
362 GetElementDiameter()
const;
363 inline Point<dim> GetCenter()
const;
365 GetFEValuesState()
const;
367 GetFEValuesControl()
const;
373 GetStateIndex()
const;
375 GetControlIndex()
const;
379 unsigned int _state_index;
380 unsigned int _control_index;
382 const std::vector<typename dealii::MGDoFHandler<dim>::active_cell_iterator> & _element;
386 unsigned int _n_q_points_per_element;
387 unsigned int _n_dofs_per_element;
394 template<
typename VECTOR,
int dim>
420 template<
template<
int,
int>
class FE,
typename SPARSITYPATTERN,
int dopedim,
int dealdim>
422 UpdateFlags update_flags,
424 VECTOR, dopedim, dealdim>& sth,
427 const std::map<std::string,
const Vector<double>*> ¶m_values,
428 const std::map<std::string, const VECTOR*> &domain_values) :
429 edcinternal::ElementDataContainerInternal<VECTOR, dim>(param_values,
430 domain_values), _element(element), _state_hp_fe_values(
431 sth.GetMapping(), (sth.GetFESystem(
"state")), q_collection,
432 update_flags), _control_hp_fe_values(sth.GetMapping(),
433 (sth.GetFESystem(
"control")), q_collection, update_flags), _q_collection(
436 _state_index = sth.GetStateIndex();
437 if (_state_index == 1)
460 template<
template<
int,
int>
class FE,
typename SPARSITYPATTERN>
462 UpdateFlags update_flags,
467 const std::map<std::string,
const Vector<double>*> ¶m_values,
468 const std::map<std::string, const VECTOR*> &domain_values) :
469 edcinternal::ElementDataContainerInternal<VECTOR, dim>(param_values,
470 domain_values), _element(element), _state_hp_fe_values(
471 sth.GetMapping(), (sth.GetFESystem(
"state")), q_collection,
472 update_flags), _control_hp_fe_values(sth.GetMapping(),
473 (sth.GetFESystem(
"state")), q_collection, update_flags), _q_collection(
476 _state_index = sth.GetStateIndex();
477 _control_index = element.size();
497 GetNDoFsPerElement()
const;
501 GetMaterialId()
const;
503 GetNbrMaterialId(
unsigned int face)
const;
505 GetFaceBoundaryIndicator(
unsigned int face)
const;
507 GetIsAtBoundary()
const;
509 GetElementDiameter()
const;
514 GetFEValuesState()
const;
516 GetFEValuesControl()
const;
520 GetStateIndex()
const;
522 GetControlIndex()
const;
523 const std::map<std::string, const VECTOR*> &
524 GetDomainValues()
const;
532 std::vector<double>& values)
const;
539 std::vector<Vector<double> >& values)
const;
545 template<
int targetdim>
549 std::vector<Tensor<1, targetdim> >& values)
const;
554 template<
int targetdim>
558 std::vector<std::vector<Tensor<1, targetdim> > >& values)
const;
562 unsigned int _state_index;
563 unsigned int _control_index;
566 typename dealii::hp::DoFHandler<dim>::active_cell_iterator> & _element;
570 const hp::QCollection<dim>& _q_collection;
577 template<
typename VECTOR,
int dim>
581 _state_fe_values.reinit(_element[this->GetStateIndex()]);
583 if (this->GetControlIndex() < _element.size())
584 _control_fe_values.reinit(_element[this->GetControlIndex()]);
588 template<
typename VECTOR,
int dim>
592 return _n_dofs_per_element;
596 template<
typename VECTOR,
int dim>
600 return _n_q_points_per_element;
604 template<
typename VECTOR,
int dim>
608 return _element[0]->material_id();
612 template<
typename VECTOR,
int dim>
615 unsigned int face)
const
617 if (_element[0]->neighbor_index(face) != -1)
618 return _element[0]->neighbor(face)->material_id();
620 throw DOpEException(
"There is no neighbor with number " + face,
621 "ElementDataContainer::GetNbrMaterialId");
625 template<
typename VECTOR,
int dim>
628 unsigned int face)
const
630 return _element[0]->face(face)->boundary_indicator();
634 template<
typename VECTOR,
int dim>
638 return _element[0]->at_boundary();
641 template<
typename VECTOR,
int dim>
645 return _element[0]->diameter();
648 template<
typename VECTOR,
int dim>
652 return _element[0]->center();
656 template<
typename VECTOR,
int dim>
660 return _state_fe_values;
664 template<
typename VECTOR,
int dim>
668 return _control_fe_values;
673 template<
typename VECTOR,
int dim>
682 template<
typename VECTOR,
int dim>
684 ElementDataContainer<dealii::DoFHandler, VECTOR, dim>::GetControlIndex()
const
686 return _control_index;
702 template<
typename VECTOR,
int dim>
706 _state_fe_values.reinit(_element[this->GetStateIndex()]);
708 if (this->GetControlIndex() < _element.size())
709 _control_fe_values.reinit(_element[this->GetControlIndex()]);
713 template<
typename VECTOR,
int dim>
717 return _n_dofs_per_element;
721 template<
typename VECTOR,
int dim>
725 return _n_q_points_per_element;
729 template<
typename VECTOR,
int dim>
733 return _element[0]->material_id();
737 template<
typename VECTOR,
int dim>
740 unsigned int face)
const
742 if (_element[0]->neighbor_index(face) != -1)
743 return _element[0]->neighbor(face)->material_id();
745 throw DOpEException(
"There is no neighbor with number " + face,
746 "ElementDataContainer::GetNbrMaterialId");
750 template<
typename VECTOR,
int dim>
753 unsigned int face)
const
755 return _element[0]->face(face)->boundary_indicator();
759 template<
typename VECTOR,
int dim>
763 return _element[0]->at_boundary();
766 template<
typename VECTOR,
int dim>
770 return _element[0]->diameter();
773 template<
typename VECTOR,
int dim>
777 return _element[0]->center();
781 template<
typename VECTOR,
int dim>
785 return _state_fe_values;
789 template<
typename VECTOR,
int dim>
793 return _control_fe_values;
798 template<
typename VECTOR,
int dim>
807 template<
typename VECTOR,
int dim>
809 ElementDataContainer<dealii::MGDoFHandler, VECTOR, dim>::GetControlIndex()
const
811 return _control_index;
824 template<
typename VECTOR,
int dim>
828 _state_hp_fe_values.reinit(_element[this->GetStateIndex()]);
830 if (this->GetControlIndex() < _element.size())
831 _control_hp_fe_values.reinit(_element[this->GetControlIndex()]);
835 template<
typename VECTOR,
int dim>
839 return _element[0]->get_fe().dofs_per_cell;
842 template<
typename VECTOR,
int dim>
846 return (_q_collection[_element[0]->active_fe_index()]).size();
850 template<
typename VECTOR,
int dim>
854 return _element[0]->material_id();
858 template<
typename VECTOR,
int dim>
861 unsigned int face)
const
863 if (_element[0]->neighbor_index(face) != -1)
864 return _element[0]->neighbor(face)->material_id();
866 throw DOpEException(
"There is no neighbor with number" + face,
867 "HpElementDataContainer::GetNbrMaterialId");
871 template<
typename VECTOR,
int dim>
874 unsigned int face)
const
876 return _element[0]->face(face)->boundary_indicator();
881 template<
typename VECTOR,
int dim>
885 return _element[0]->at_boundary();
890 template<
typename VECTOR,
int dim>
894 return _element[0]->diameter();
898 template<
typename VECTOR,
int dim>
902 return _element[0]->center();
906 template<
typename VECTOR,
int dim>
913 template<
typename VECTOR,
int dim>
921 template<
typename VECTOR,
int dim>
928 template<
typename VECTOR,
int dim>
930 ElementDataContainer<dealii::hp::DoFHandler, VECTOR, dim>::GetControlIndex()
const
932 return _control_index;
ElementDataContainer(const Quadrature< dim > &quad, UpdateFlags update_flags, StateSpaceTimeHandler< FE, dealii::MGDoFHandler, SPARSITYPATTERN, VECTOR, dim > &sth, const std::vector< typename dealii::MGDoFHandler< dim >::active_cell_iterator > &element, const std::map< std::string, const Vector< double > * > ¶m_values, const std::map< std::string, const VECTOR * > &domain_values)
Definition: elementdatacontainer.h:312
ElementDataContainer()
Definition: elementdatacontainer.h:60
ElementDataContainer(const hp::QCollection< dim > &q_collection, UpdateFlags update_flags, SpaceTimeHandler< FE, dealii::hp::DoFHandler, SPARSITYPATTERN, VECTOR, dopedim, dealdim > &sth, const std::vector< typename DOpEWrapper::DoFHandler< dim, dealii::hp::DoFHandler >::active_cell_iterator > &element, const std::map< std::string, const Vector< double > * > ¶m_values, const std::map< std::string, const VECTOR * > &domain_values)
Definition: elementdatacontainer.h:421
Definition: elementdatacontainer.h:56
~ElementDataContainer()
Definition: elementdatacontainer.h:166
ElementDataContainer(const Quadrature< dim > &quad, UpdateFlags update_flags, StateSpaceTimeHandler< FE, dealii::DoFHandler, SPARSITYPATTERN, VECTOR, dim > &sth, const std::vector< typename dealii::DoFHandler< dim >::active_cell_iterator > &element, const std::map< std::string, const Vector< double > * > ¶m_values, const std::map< std::string, const VECTOR * > &domain_values)
Definition: elementdatacontainer.h:147
ElementDataContainer(const Quadrature< dim > &quad, UpdateFlags update_flags, SpaceTimeHandler< FE, dealii::DoFHandler, SPARSITYPATTERN, VECTOR, dopedim, dealdim > &sth, const std::vector< typename dealii::DoFHandler< dim >::active_cell_iterator > &element, const std::map< std::string, const Vector< double > * > ¶m_values, const std::map< std::string, const VECTOR * > &domain_values)
Definition: elementdatacontainer.h:105
~ElementDataContainer()
Definition: elementdatacontainer.h:333
Definition: elementdatacontainer_internal.h:46
Definition: dofhandler_wrapper.h:51
Definition: statespacetimehandler.h:59
~ElementDataContainer()
Definition: elementdatacontainer.h:479
Definition: spacetimehandler.h:71
ElementDataContainer(const Quadrature< dim > &quad, UpdateFlags update_flags, SpaceTimeHandler< FE, dealii::MGDoFHandler, SPARSITYPATTERN, VECTOR, dopedim, dealdim > &sth, const std::vector< typename dealii::MGDoFHandler< dim >::active_cell_iterator > &element, const std::map< std::string, const Vector< double > * > ¶m_values, const std::map< std::string, const VECTOR * > &domain_values)
Definition: elementdatacontainer.h:264
ElementDataContainer(const hp::QCollection< dim > &q_collection, UpdateFlags update_flags, StateSpaceTimeHandler< FE, dealii::hp::DoFHandler, SPARSITYPATTERN, VECTOR, dim > &sth, const std::vector< typename DOpEWrapper::DoFHandler< dim, dealii::hp::DoFHandler >::active_cell_iterator > &element, const std::map< std::string, const Vector< double > * > ¶m_values, const std::map< std::string, const VECTOR * > &domain_values)
Definition: elementdatacontainer.h:461
Definition: fevalues_wrapper.h:43
Definition: dopeexception.h:35
Definition: fevalues_wrapper.h:186