24 #ifndef _MOL_MULTIMESH_SPACE_TIME_HANDLER_H_
25 #define _MOL_MULTIMESH_SPACE_TIME_HANDLER_H_
35 #include <dofs/dof_handler.h>
36 #include <dofs/dof_renumbering.h>
37 #include <dofs/dof_tools.h>
38 #include <hp/mapping_collection.h>
39 #include <lac/constraint_matrix.h>
40 #include <deal.II/grid/grid_refinement.h>
55 template<
template<
int,
int>
class FE,
template<
int,
int>
class DH,
56 typename SPARSITYPATTERN,
typename VECTOR,
int dim>
58 DH, SPARSITYPATTERN, VECTOR, dim, dim>
62 dealii::Triangulation<dim>& triangulation,
63 const FE<dim, dim>& control_fe,
const FE<dim, dim>& state_fe,
68 index_setter), _state_triangulation(triangulation), _control_dof_handler(
69 _control_triangulation), _state_dof_handler(
70 _state_triangulation), _control_fe(&control_fe), _state_fe(
72 DOpEWrapper::StaticMappingQ1<dim, DH>::mapping_q1), _constraints(), _control_mesh_transfer(
73 NULL), _state_mesh_transfer(NULL), _sparse_mkr_dynamic(true)
75 _control_triangulation.copy_triangulation(_state_triangulation);
77 _user_defined_dof_constr = NULL;
81 dealii::Triangulation<dim>& triangulation,
82 const FE<dim, dim>& control_fe,
const FE<dim, dim>& state_fe,
87 type, index_setter), _state_triangulation(triangulation), _control_dof_handler(
88 _control_triangulation), _state_dof_handler(
89 _state_triangulation), _control_fe(&control_fe), _state_fe(
91 DOpEWrapper::StaticMappingQ1<dim, DH>::mapping_q1), _constraints(), _control_mesh_transfer(
92 NULL), _state_mesh_transfer(NULL), _sparse_mkr_dynamic(true)
94 _control_triangulation.copy_triangulation(_state_triangulation);
96 _user_defined_dof_constr = NULL;
100 dealii::Triangulation<dim>& triangulation,
101 const FE<dim, dim>& control_fe,
const FE<dim, dim>& state_fe,
106 index_setter), _state_triangulation(triangulation), _control_dof_handler(
107 _control_triangulation), _state_dof_handler(
108 _state_triangulation), _control_fe(&control_fe), _state_fe(
109 &state_fe), _mapping(
110 DOpEWrapper::StaticMappingQ1<dim, DH>::mapping_q1), _constraints(
111 c), _control_mesh_transfer(NULL), _state_mesh_transfer(NULL), _sparse_mkr_dynamic(
114 _control_triangulation.copy_triangulation(_state_triangulation);
116 _user_defined_dof_constr = NULL;
120 dealii::Triangulation<dim>& triangulation,
121 const FE<dim, dim>& control_fe,
const FE<dim, dim>& state_fe,
122 const dealii::Triangulation<1> & times,
const Constraints& c,
127 type, index_setter), _state_triangulation(triangulation), _control_dof_handler(
128 _control_triangulation), _state_dof_handler(
129 _state_triangulation), _control_fe(&control_fe), _state_fe(
130 &state_fe), _mapping(
131 DOpEWrapper::StaticMappingQ1<dim, DH>::mapping_q1), _constraints(
132 c), _control_mesh_transfer(NULL), _state_mesh_transfer(NULL), _sparse_mkr_dynamic(
135 _control_triangulation.copy_triangulation(_state_triangulation);
137 _user_defined_dof_constr = NULL;
142 _control_dof_handler.clear();
144 _state_dof_handler.clear();
146 if (_control_mesh_transfer != NULL)
148 delete _control_mesh_transfer;
150 if (_state_mesh_transfer != NULL)
152 delete _state_mesh_transfer;
154 if (_sparsitymaker != NULL && _sparse_mkr_dynamic ==
true)
156 delete _sparsitymaker;
165 const std::vector<unsigned int>& control_block_component,
166 unsigned int state_n_blocks,
167 const std::vector<unsigned int>& state_block_component)
170 _control_dof_handler);
171 _control_dof_handler.distribute_dofs(*_control_fe);
172 DoFRenumbering::component_wise(
static_cast<DH<dim, dim>&
>(_control_dof_handler));
174 _control_dof_constraints.clear();
175 DoFTools::make_hanging_node_constraints(
176 static_cast<DH<dim, dim>&
>(_control_dof_handler), _control_dof_constraints);
177 if (GetUserDefinedDoFConstraints() != NULL)
179 _control_dof_handler, _control_dof_constraints);
180 _control_dof_constraints.close();
182 _control_dofs_per_block.resize(control_n_blocks);
184 DoFTools::count_dofs_per_block(
185 static_cast<DH<dim, dim>&
>(_control_dof_handler), _control_dofs_per_block,
186 control_block_component);
191 _state_dof_handler.distribute_dofs(
GetFESystem(
"state"));
192 DoFRenumbering::component_wise(
static_cast<DH<dim, dim>&
>(_state_dof_handler));
194 _state_dof_constraints.clear();
195 DoFTools::make_hanging_node_constraints(
196 static_cast<DH<dim, dim>&
>(_state_dof_handler), _state_dof_constraints);
198 if (GetUserDefinedDoFConstraints() != NULL)
200 _state_dof_handler, _state_dof_constraints);
201 _state_dof_constraints.close();
203 _state_dofs_per_block.resize(state_n_blocks);
204 DoFTools::count_dofs_per_block(
static_cast<DH<dim, dim>&
>(_state_dof_handler),
205 _state_dofs_per_block, state_block_component);
207 _support_points.clear();
209 _constraints.
ReInit(_control_dofs_per_block);
227 return _control_dof_handler;
236 return _state_dof_handler;
254 return _control_dofs_per_block[b];
262 return _state_dofs_per_block[b];
275 const std::vector<unsigned int>&
278 return _control_dofs_per_block;
283 const std::vector<unsigned int>&
286 return _state_dofs_per_block;
291 const std::vector<unsigned int>&
299 const dealii::ConstraintMatrix&
302 return _control_dof_constraints;
307 const dealii::ConstraintMatrix&
310 return _state_dof_constraints;
319 const std::vector<VECTOR*> & local_vectors,
double t,
324 if (local_vectors.size() != 2)
326 "This function is currently not implemented for anything other than"
327 " linear interpolation of 2 DoFs.",
328 "MethodOfLine_SpaceTimeHandler::InterpolateState");
334 result = *local_vectors[0];
336 result.sadd(lambda_l, lambda_r, *local_vectors[1]);
361 return _constraints.
n_dofs(name);
369 return _constraints.
n_dofs(
"global");
378 return _constraints.
n_dofs(
"local");
385 const std::vector<Point<dim> >&
391 return _support_points;
402 this->GetSparsityMaker()->ComputeSparsityPattern(
419 else if (name ==
"control")
426 "MethodOfLines_MultiMesh_SpaceTimeHandler::GetFESystem");
498 template<
typename NUMBER>
522 if (_state_mesh_transfer != NULL)
524 delete _state_mesh_transfer;
525 _state_mesh_transfer = NULL;
528 DH>(_state_dof_handler);
532 _state_triangulation.set_all_refine_flags();
537 GridRefinement::refine_and_coarsen_fixed_number(
545 GridRefinement::refine_and_coarsen_fixed_fraction(
553 GridRefinement::refine_and_coarsen_optimize(_state_triangulation,
559 this->FlagIfLeftIsNotFinest(_state_triangulation,
560 _control_triangulation);
565 "MethodOfLines_MultiMesh_SpaceTimeHandler::RefineSpace");
567 _state_triangulation.prepare_coarsening_and_refinement();
568 if (_state_mesh_transfer != NULL)
569 _state_mesh_transfer->prepare_for_pure_refinement();
571 _state_triangulation.execute_coarsening_and_refinement();
590 if (_control_mesh_transfer != NULL)
592 delete _control_mesh_transfer;
593 _control_mesh_transfer = NULL;
595 #if dope_dimension == deal_II_dimension
597 DH>(_control_dof_handler);
601 _control_triangulation.set_all_refine_flags();
607 GridRefinement::refine_and_coarsen_fixed_number(
615 GridRefinement::refine_and_coarsen_fixed_fraction(
623 GridRefinement::refine_and_coarsen_optimize(_control_triangulation,
627 else if (
"finest-of-both")
629 this->FlagIfLeftIsNotFinest(_control_triangulation,
630 _state_triangulation);
635 "MethodOfLines_MultiMesh_SpaceTimeHandler::RefineSpace");
637 _control_triangulation.prepare_coarsening_and_refinement();
640 if (_control_mesh_transfer != NULL)
641 _control_mesh_transfer->prepare_for_pure_refinement();
643 _control_triangulation.execute_coarsening_and_refinement();
665 VECTOR& new_values)
const
667 if (_control_mesh_transfer != NULL)
668 _control_mesh_transfer->refine_interpolate(old_values, new_values);
677 VECTOR& new_values)
const
679 if (_state_mesh_transfer != NULL)
680 _state_mesh_transfer->refine_interpolate(old_values, new_values);
695 _user_defined_dof_constr = &constraints_maker;
709 if (_sparsitymaker != NULL && _sparse_mkr_dynamic)
710 delete _sparsitymaker;
711 _sparsitymaker = &sparsity_maker;
712 _sparse_mkr_dynamic =
false;
724 _state_dof_handler.clear();
725 _state_triangulation.clear();
726 _state_triangulation.copy_triangulation(tria);
727 _state_dof_handler.initialize(_state_triangulation, *_state_fe);
729 if (_state_mesh_transfer != NULL)
730 delete _state_mesh_transfer;
731 _state_mesh_transfer = NULL;
744 _control_dof_handler.clear();
745 _control_triangulation.clear();
746 _control_triangulation.copy_triangulation(tria);
747 _control_dof_handler.initialize(_control_triangulation, *_control_fe);
749 if (_control_mesh_transfer != NULL)
750 delete _control_mesh_transfer;
751 _control_mesh_transfer = NULL;
756 GetSparsityMaker()
const
758 return _sparsitymaker;
760 const UserDefinedDoFConstraints<DH, dim>*
761 GetUserDefinedDoFConstraints()
const
763 return _user_defined_dof_constr;
767 FlagIfLeftIsNotFinest(dealii::Triangulation<dim>& left,
768 const dealii::Triangulation<dim>& right);
770 SparsityMaker<DH, dim>* _sparsitymaker;
771 UserDefinedDoFConstraints<DH, dim>* _user_defined_dof_constr;
773 dealii::Triangulation<dim> _control_triangulation;
774 dealii::Triangulation<dim>& _state_triangulation;
778 std::vector<unsigned int> _control_dofs_per_block;
779 std::vector<unsigned int> _state_dofs_per_block;
781 dealii::ConstraintMatrix _control_dof_constraints;
782 dealii::ConstraintMatrix _state_dof_constraints;
784 const dealii::SmartPointer<const FE<dim, dim>> _control_fe;
785 const dealii::SmartPointer<const FE<dim, dim>> _state_fe;
789 std::vector<Point<dim> > _support_points;
791 Constraints _constraints;
794 bool _sparse_mkr_dynamic;
798 #if dope_dimension == deal_II_dimension
805 dealii::DoFHandler, dealii::BlockSparsityPattern,
806 dealii::BlockVector<double>, dope_dimension>::ComputeControlSparsityPattern(
807 dealii::BlockSparsityPattern & sparsity)
const
809 const std::vector<unsigned int>& blocks = this->GetControlDoFsPerBlock();
810 dealii::BlockCompressedSimpleSparsityPattern csp(blocks.size(),
812 for (
unsigned int i = 0; i < blocks.size(); i++)
814 for (
unsigned int j = 0; j < blocks.size(); j++)
816 csp.block(i, j).reinit(this->GetControlDoFsPerBlock(i),
817 this->GetControlDoFsPerBlock(j));
822 dealii::DoFTools::make_sparsity_pattern(
823 static_cast<const dealii::DoFHandler<deal_II_dimension>&
>(this->GetControlDoFHandler()),
825 this->GetControlDoFConstraints().condense(csp);
826 sparsity.copy_from(csp);
834 dealii::DoFHandler, dealii::SparsityPattern, dealii::Vector<double>,
835 dope_dimension>::ComputeControlSparsityPattern(
836 dealii::SparsityPattern & sparsity)
const
838 const unsigned int total_dofs = this->GetControlNDoFs();
839 dealii::CompressedSimpleSparsityPattern csp(total_dofs, total_dofs);
841 dealii::DoFTools::make_sparsity_pattern(
842 static_cast<const dealii::DoFHandler<deal_II_dimension>&
>(this->GetControlDoFHandler()),
844 this->GetControlDoFConstraints().condense(csp);
845 sparsity.copy_from(csp);
854 dealii::hp::DoFHandler, dealii::BlockSparsityPattern,
855 dealii::BlockVector<double>, dope_dimension>::ComputeControlSparsityPattern(
856 dealii::BlockSparsityPattern & sparsity)
const
858 const std::vector<unsigned int>& blocks = this->GetControlDoFsPerBlock();
859 dealii::BlockCompressedSimpleSparsityPattern csp(blocks.size(),
861 for (
unsigned int i = 0; i < blocks.size(); i++)
863 for (
unsigned int j = 0; j < blocks.size(); j++)
865 csp.block(i, j).reinit(this->GetControlDoFsPerBlock(i),
866 this->GetControlDoFsPerBlock(j));
871 dealii::DoFTools::make_sparsity_pattern(
872 static_cast<const dealii::hp::DoFHandler<deal_II_dimension>&
>(this->GetControlDoFHandler()),
874 this->GetControlDoFConstraints().condense(csp);
875 sparsity.copy_from(csp);
883 dealii::hp::DoFHandler, dealii::SparsityPattern, dealii::Vector<double>,
884 dope_dimension>::ComputeControlSparsityPattern(
885 dealii::SparsityPattern & sparsity)
const
887 const unsigned int total_dofs = this->GetControlNDoFs();
888 dealii::CompressedSimpleSparsityPattern csp(total_dofs, total_dofs);
890 dealii::DoFTools::make_sparsity_pattern(
891 static_cast<const dealii::hp::DoFHandler<deal_II_dimension>&
>(this->GetControlDoFHandler()),
893 this->GetControlDoFConstraints().condense(csp);
894 sparsity.copy_from(csp);
896 #endif //Endof explicit instanciation
898 template<
template<
int,
int>
class FE,
template<
int,
int>
class DH,
899 typename SPARSITYPATTERN,
typename VECTOR,
int dim>
902 dim>::FlagIfLeftIsNotFinest(dealii::Triangulation<dim>& left,
903 const dealii::Triangulation<dim>& right)
905 auto element_list = GridTools::get_finest_common_cells(left, right);
906 auto element_iter = element_list.begin();
907 for (; element_iter != element_list.end(); element_iter++)
909 if (element_iter->second->has_children())
912 element_iter->first->set_refine_flag();
void ReInitTime()
Definition: spacetimehandler_base.h:82
void SetActiveFEIndicesState(DOpEWrapper::DoFHandler< dealdim, DH > &dof_handler)
Definition: spacetimehandler.h:327
Definition: mapping_wrapper.h:48
void ReInit(unsigned int control_n_blocks, const std::vector< unsigned int > &control_block_component, unsigned int state_n_blocks, const std::vector< unsigned int > &state_block_component)
Definition: mol_multimesh_spacetimehandler.h:164
void RegisterMapping(const typename DOpEWrapper::Mapping< dealdim, DH > &mapping)
Definition: userdefineddofconstraints.h:76
Definition: dopetypes.h:50
void SetSparsityMaker(SparsityMaker< DH, dim > &sparsity_maker)
Definition: mol_multimesh_spacetimehandler.h:707
void RefineControlSpace(DOpEtypes::RefinementType ref_type=DOpEtypes::RefinementType::global)
Definition: mol_multimesh_spacetimehandler.h:481
void RefineSpace(const RefinementContainer &ref_container)
Definition: mol_multimesh_spacetimehandler.h:500
void RefineStateSpace(DOpEtypes::RefinementType ref_type=DOpEtypes::RefinementType::global)
Definition: mol_multimesh_spacetimehandler.h:462
unsigned int GetConstraintNDoFs(std::string name) const
Definition: mol_multimesh_spacetimehandler.h:359
const DOpEWrapper::DoFHandler< dim, DH > & GetStateDoFHandler() const
Definition: mol_multimesh_spacetimehandler.h:233
void SetActiveFEIndicesControl(DOpEWrapper::DoFHandler< dopedim, DH > &dof_handler)
Definition: spacetimehandler.h:348
MethodOfLines_MultiMesh_SpaceTimeHandler(dealii::Triangulation< dim > &triangulation, const FE< dim, dim > &control_fe, const FE< dim, dim > &state_fe, const dealii::Triangulation< 1 > ×, const Constraints &c, DOpEtypes::ControlType type, const ActiveFEIndexSetterInterface< dim, dim > &index_setter=ActiveFEIndexSetterInterface< dim, dim >())
Definition: mol_multimesh_spacetimehandler.h:119
const std::vector< unsigned int > & GetStateDoFsPerBlock(int=-1) const
Definition: mol_multimesh_spacetimehandler.h:284
MethodOfLines_MultiMesh_SpaceTimeHandler(dealii::Triangulation< dim > &triangulation, const FE< dim, dim > &control_fe, const FE< dim, dim > &state_fe, const dealii::Triangulation< 1 > ×, DOpEtypes::ControlType type, const ActiveFEIndexSetterInterface< dim, dim > &index_setter=ActiveFEIndexSetterInterface< dim, dim >())
Definition: mol_multimesh_spacetimehandler.h:80
ControlType
Definition: dopetypes.h:102
void SetUserDefinedDoFConstraints(UserDefinedDoFConstraints< DH, dim > &constraints_maker)
Definition: mol_multimesh_spacetimehandler.h:692
Definition: dopetypes.h:50
void RefineSpace(DOpEtypes::RefinementType ref_type=DOpEtypes::RefinementType::global)
Definition: mol_multimesh_spacetimehandler.h:442
Definition: dopetypes.h:50
const std::vector< Point< dim > > & GetMapDoFToSupportPoints()
Definition: mol_multimesh_spacetimehandler.h:386
void ComputeControlSparsityPattern(SPARSITYPATTERN &sparsity) const
MethodOfLines_MultiMesh_SpaceTimeHandler(dealii::Triangulation< dim > &triangulation, const FE< dim, dim > &control_fe, const FE< dim, dim > &state_fe, DOpEtypes::ControlType type, const ActiveFEIndexSetterInterface< dim, dim > &index_setter=ActiveFEIndexSetterInterface< dim, dim >())
Definition: mol_multimesh_spacetimehandler.h:61
Definition: timeiterator.h:63
unsigned int NewTimePointToOldTimePoint(unsigned int t) const
Definition: mol_multimesh_spacetimehandler.h:651
double get_k() const
Definition: timeiterator.h:195
DOpEtypes::RefinementType GetRefType() const
Definition: refinementcontainer.cc:95
RefinementType
Definition: dopetypes.h:48
const dealii::ConstraintMatrix & GetStateDoFConstraints() const
Definition: mol_multimesh_spacetimehandler.h:308
unsigned int GetControlNDoFs() const
Definition: mol_multimesh_spacetimehandler.h:343
const DOpEWrapper::DoFHandler< dim, DH > & GetControlDoFHandler() const
Definition: mol_multimesh_spacetimehandler.h:224
Definition: refinementcontainer.h:42
void ResetStateTriangulation(const dealii::Triangulation< dim > &tria)
Definition: mol_multimesh_spacetimehandler.h:722
const std::vector< unsigned int > & GetDoFsPerBlock(std::string name) const
Definition: constraints.h:159
unsigned int GetNGlobalConstraints() const
Definition: mol_multimesh_spacetimehandler.h:367
Definition: active_fe_index_setter_interface.h:39
unsigned int GetStateDoFsPerBlock(unsigned int b, int=-1) const
Definition: mol_multimesh_spacetimehandler.h:260
virtual const dealii::Vector< float > & GetLocalErrorIndicators() const
Definition: refinementcontainer.cc:46
virtual double GetTopFraction() const
Definition: refinementcontainer.cc:56
virtual void InterpolateState(VECTOR &result, const std::vector< VECTOR * > &local_vectors, double t, const TimeIterator &it) const
Definition: mol_multimesh_spacetimehandler.h:318
Definition: dopetypes.h:50
unsigned int GetConstraintDoFsPerBlock(std::string name, unsigned int b) const
Definition: mol_multimesh_spacetimehandler.h:268
void ComputeStateSparsityPattern(SPARSITYPATTERN &sparsity) const
Definition: mol_multimesh_spacetimehandler.h:400
Definition: constraints.h:33
void SpatialMeshTransferState(const VECTOR &old_values, VECTOR &new_values) const
Definition: mol_multimesh_spacetimehandler.h:676
const std::vector< unsigned int > & GetControlDoFsPerBlock() const
Definition: mol_multimesh_spacetimehandler.h:276
virtual double GetConvergenceOrder() const
Definition: refinementcontainer.cc:85
const dealii::ConstraintMatrix & GetControlDoFConstraints() const
Definition: mol_multimesh_spacetimehandler.h:300
double get_right() const
Definition: timeiterator.h:185
Definition: spacetimehandler.h:71
void MapDoFsToSupportPoints(const DOpEWrapper::Mapping< dealdim, dealii::DoFHandler > &mapping, const DOpEWrapper::DoFHandler< dealdim, dealii::DoFHandler > &dh, VECTOR &support_points)
Definition: sth_internals.h:47
virtual void MakeControlDoFConstraints(const DOpEWrapper::DoFHandler< dopedim, DH > &dof_handler, dealii::ConstraintMatrix &dof_constraints) const
Definition: userdefineddofconstraints.h:103
const std::vector< unsigned int > & GetConstraintDoFsPerBlock(std::string name) const
Definition: mol_multimesh_spacetimehandler.h:292
void ReInit(std::vector< unsigned int > &control_dofs_per_block)
Definition: constraints.h:108
~MethodOfLines_MultiMesh_SpaceTimeHandler()
Definition: mol_multimesh_spacetimehandler.h:140
unsigned int GetControlDoFsPerBlock(unsigned int b, int=-1) const
Definition: mol_multimesh_spacetimehandler.h:252
Definition: mol_multimesh_spacetimehandler.h:57
Definition: sparsitymaker.h:50
bool UsesCoarsening() const
Definition: refinementcontainer.cc:103
double get_left() const
Definition: timeiterator.h:175
const FE< dim, dim > & GetFESystem(std::string name) const
Definition: mol_multimesh_spacetimehandler.h:413
unsigned int GetStateNDoFs(int=-1) const
Definition: mol_multimesh_spacetimehandler.h:351
void IncrementControlTicket()
Definition: spacetimehandler_base.h:382
Definition: solutiontransfer_wrapper.h:51
unsigned int GetNLocalConstraints() const
Definition: mol_multimesh_spacetimehandler.h:376
void RefineStateSpace(const RefinementContainer &ref_container)
Definition: mol_multimesh_spacetimehandler.h:515
Definition: dopeexception.h:35
MethodOfLines_MultiMesh_SpaceTimeHandler(dealii::Triangulation< dim > &triangulation, const FE< dim, dim > &control_fe, const FE< dim, dim > &state_fe, const Constraints &c, DOpEtypes::ControlType type, const ActiveFEIndexSetterInterface< dim, dim > &index_setter=ActiveFEIndexSetterInterface< dim, dim >())
Definition: mol_multimesh_spacetimehandler.h:99
Definition: dopetypes.h:50
const DOpEWrapper::Mapping< dim, DH > & GetMapping() const
Definition: mol_multimesh_spacetimehandler.h:243
void SpatialMeshTransferControl(const VECTOR &old_values, VECTOR &new_values) const
Definition: mol_multimesh_spacetimehandler.h:664
void RefineControlSpace(const RefinementContainer &ref_container)
Definition: mol_multimesh_spacetimehandler.h:583
virtual double GetBottomFraction() const
Definition: refinementcontainer.cc:66
virtual void MakeStateDoFConstraints(const DOpEWrapper::DoFHandler< dealdim, DH > &dof_handler, dealii::ConstraintMatrix &dof_constraints) const
Definition: userdefineddofconstraints.h:93
unsigned int n_dofs(std::string name) const
Definition: constraints.h:141
void IncrementStateTicket()
Definition: spacetimehandler_base.h:373
void ResetControlTriangulation(const dealii::Triangulation< dim > &tria)
Definition: mol_multimesh_spacetimehandler.h:742