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 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::InterpolateControl");
334 result = *local_vectors[0];
336 result.sadd(lambda_l, lambda_r, *local_vectors[1]);
341 const std::vector<VECTOR*> & local_vectors,
double t,
346 if (local_vectors.size() != 2)
348 "This function is currently not implemented for anything other than"
349 " linear interpolation of 2 DoFs.",
350 "MethodOfLine_SpaceTimeHandler::InterpolateState");
356 result = *local_vectors[0];
358 result.sadd(lambda_l, lambda_r, *local_vectors[1]);
383 return constraints_.
n_dofs(name);
391 return constraints_.
n_dofs(
"global");
400 return constraints_.
n_dofs(
"local");
407 const std::vector<Point<dim> >&
413 return support_points_;
424 this->GetSparsityMaker()->ComputeSparsityPattern(
441 else if (name ==
"control")
448 "MethodOfLines_MultiMesh_SpaceTimeHandler::GetFESystem");
520 template<
typename NUMBER>
544 if (state_mesh_transfer_ != NULL)
546 delete state_mesh_transfer_;
547 state_mesh_transfer_ = NULL;
550 DH>(state_dof_handler_);
554 state_triangulation_.set_all_refine_flags();
559 GridRefinement::refine_and_coarsen_fixed_number(
567 GridRefinement::refine_and_coarsen_fixed_fraction(
575 GridRefinement::refine_and_coarsen_optimize(state_triangulation_,
581 this->FlagIfLeftIsNotFinest(state_triangulation_,
582 control_triangulation_);
587 "MethodOfLines_MultiMesh_SpaceTimeHandler::RefineSpace");
589 state_triangulation_.prepare_coarsening_and_refinement();
590 if (state_mesh_transfer_ != NULL)
591 state_mesh_transfer_->prepare_for_pure_refinement();
593 state_triangulation_.execute_coarsening_and_refinement();
612 if (control_mesh_transfer_ != NULL)
614 delete control_mesh_transfer_;
615 control_mesh_transfer_ = NULL;
617 #if dope_dimension == deal_II_dimension
619 DH>(control_dof_handler_);
623 control_triangulation_.set_all_refine_flags();
629 GridRefinement::refine_and_coarsen_fixed_number(
637 GridRefinement::refine_and_coarsen_fixed_fraction(
645 GridRefinement::refine_and_coarsen_optimize(control_triangulation_,
649 else if (
"finest-of-both")
651 this->FlagIfLeftIsNotFinest(control_triangulation_,
652 state_triangulation_);
657 "MethodOfLines_MultiMesh_SpaceTimeHandler::RefineSpace");
659 control_triangulation_.prepare_coarsening_and_refinement();
662 if (control_mesh_transfer_ != NULL)
663 control_mesh_transfer_->prepare_for_pure_refinement();
665 control_triangulation_.execute_coarsening_and_refinement();
687 VECTOR& new_values)
const
689 if (control_mesh_transfer_ != NULL)
690 control_mesh_transfer_->refine_interpolate(old_values, new_values);
699 VECTOR& new_values)
const
701 if (state_mesh_transfer_ != NULL)
702 state_mesh_transfer_->refine_interpolate(old_values, new_values);
717 user_defined_dof_constr_ = &constraints_maker;
731 if (sparsitymaker_ != NULL && sparse_mkr_dynamic_)
732 delete sparsitymaker_;
733 sparsitymaker_ = &sparsity_maker;
734 sparse_mkr_dynamic_ =
false;
746 state_dof_handler_.clear();
747 state_triangulation_.clear();
748 state_triangulation_.copy_triangulation(tria);
749 state_dof_handler_.initialize(state_triangulation_, *state_fe_);
751 if (state_mesh_transfer_ != NULL)
752 delete state_mesh_transfer_;
753 state_mesh_transfer_ = NULL;
766 control_dof_handler_.clear();
767 control_triangulation_.clear();
768 control_triangulation_.copy_triangulation(tria);
769 control_dof_handler_.initialize(control_triangulation_, *control_fe_);
771 if (control_mesh_transfer_ != NULL)
772 delete control_mesh_transfer_;
773 control_mesh_transfer_ = NULL;
778 GetSparsityMaker()
const
780 return sparsitymaker_;
782 const UserDefinedDoFConstraints<DH, dim>*
783 GetUserDefinedDoFConstraints()
const
785 return user_defined_dof_constr_;
789 FlagIfLeftIsNotFinest(dealii::Triangulation<dim>& left,
790 const dealii::Triangulation<dim>& right);
792 SparsityMaker<DH, dim>* sparsitymaker_;
793 UserDefinedDoFConstraints<DH, dim>* user_defined_dof_constr_;
795 dealii::Triangulation<dim> control_triangulation_;
796 dealii::Triangulation<dim>& state_triangulation_;
800 std::vector<unsigned int> control_dofs_per_block_;
801 std::vector<unsigned int> state_dofs_per_block_;
803 dealii::ConstraintMatrix control_dof_constraints_;
804 dealii::ConstraintMatrix state_dof_constraints_;
806 const dealii::SmartPointer<const FE<dim, dim>> control_fe_;
807 const dealii::SmartPointer<const FE<dim, dim>> state_fe_;
811 std::vector<Point<dim> > support_points_;
813 Constraints constraints_;
816 bool sparse_mkr_dynamic_;
820 #if dope_dimension == deal_II_dimension
827 dealii::DoFHandler, dealii::BlockSparsityPattern,
828 dealii::BlockVector<double>, dope_dimension>::ComputeControlSparsityPattern(
829 dealii::BlockSparsityPattern & sparsity)
const
831 const std::vector<unsigned int>& blocks = this->GetControlDoFsPerBlock();
832 dealii::BlockCompressedSimpleSparsityPattern csp(blocks.size(),
834 for (
unsigned int i = 0; i < blocks.size(); i++)
836 for (
unsigned int j = 0; j < blocks.size(); j++)
838 csp.block(i, j).reinit(this->GetControlDoFsPerBlock(i),
839 this->GetControlDoFsPerBlock(j));
844 dealii::DoFTools::make_sparsity_pattern(
845 static_cast<const dealii::DoFHandler<deal_II_dimension>&
>(this->GetControlDoFHandler()),
847 this->GetControlDoFConstraints().condense(csp);
848 sparsity.copy_from(csp);
856 dealii::DoFHandler, dealii::SparsityPattern, dealii::Vector<double>,
857 dope_dimension>::ComputeControlSparsityPattern(
858 dealii::SparsityPattern & sparsity)
const
860 const unsigned int total_dofs = this->GetControlNDoFs();
861 dealii::CompressedSimpleSparsityPattern csp(total_dofs, total_dofs);
863 dealii::DoFTools::make_sparsity_pattern(
864 static_cast<const dealii::DoFHandler<deal_II_dimension>&
>(this->GetControlDoFHandler()),
866 this->GetControlDoFConstraints().condense(csp);
867 sparsity.copy_from(csp);
876 dealii::hp::DoFHandler, dealii::BlockSparsityPattern,
877 dealii::BlockVector<double>, dope_dimension>::ComputeControlSparsityPattern(
878 dealii::BlockSparsityPattern & sparsity)
const
880 const std::vector<unsigned int>& blocks = this->GetControlDoFsPerBlock();
881 dealii::BlockCompressedSimpleSparsityPattern csp(blocks.size(),
883 for (
unsigned int i = 0; i < blocks.size(); i++)
885 for (
unsigned int j = 0; j < blocks.size(); j++)
887 csp.block(i, j).reinit(this->GetControlDoFsPerBlock(i),
888 this->GetControlDoFsPerBlock(j));
893 dealii::DoFTools::make_sparsity_pattern(
894 static_cast<const dealii::hp::DoFHandler<deal_II_dimension>&
>(this->GetControlDoFHandler()),
896 this->GetControlDoFConstraints().condense(csp);
897 sparsity.copy_from(csp);
905 dealii::hp::DoFHandler, dealii::SparsityPattern, dealii::Vector<double>,
906 dope_dimension>::ComputeControlSparsityPattern(
907 dealii::SparsityPattern & sparsity)
const
909 const unsigned int total_dofs = this->GetControlNDoFs();
910 dealii::CompressedSimpleSparsityPattern csp(total_dofs, total_dofs);
912 dealii::DoFTools::make_sparsity_pattern(
913 static_cast<const dealii::hp::DoFHandler<deal_II_dimension>&
>(this->GetControlDoFHandler()),
915 this->GetControlDoFConstraints().condense(csp);
916 sparsity.copy_from(csp);
918 #endif //Endof explicit instanciation
920 template<
template<
int,
int>
class FE,
template<
int,
int>
class DH,
921 typename SPARSITYPATTERN,
typename VECTOR,
int dim>
924 dim>::FlagIfLeftIsNotFinest(dealii::Triangulation<dim>& left,
925 const dealii::Triangulation<dim>& right)
927 auto element_list = GridTools::get_finest_common_cells(left, right);
928 auto element_iter = element_list.begin();
929 for (; element_iter != element_list.end(); element_iter++)
931 if (element_iter->second->has_children())
934 element_iter->first->set_refine_flag();
void ReInitTime()
Definition: spacetimehandler_base.h:86
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:52
unsigned int GetControlNDoFs(int=-1) const
Definition: mol_multimesh_spacetimehandler.h:365
void SetSparsityMaker(SparsityMaker< DH, dim > &sparsity_maker)
Definition: mol_multimesh_spacetimehandler.h:729
void RefineControlSpace(DOpEtypes::RefinementType ref_type=DOpEtypes::RefinementType::global)
Definition: mol_multimesh_spacetimehandler.h:503
void RefineSpace(const RefinementContainer &ref_container)
Definition: mol_multimesh_spacetimehandler.h:522
void RefineStateSpace(DOpEtypes::RefinementType ref_type=DOpEtypes::RefinementType::global)
Definition: mol_multimesh_spacetimehandler.h:484
unsigned int GetConstraintNDoFs(std::string name) const
Definition: mol_multimesh_spacetimehandler.h:381
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
const std::vector< unsigned int > & GetStateDoFsPerBlock(int=-1) const
Definition: mol_multimesh_spacetimehandler.h:284
ControlType
Definition: dopetypes.h:103
MethodOfLines_MultiMesh_SpaceTimeHandler(dealii::Triangulation< dim > &triangulation, const FE< dim, dim > &control_fe, const FE< dim, dim > &state_fe, dealii::Triangulation< 1 > ×, const Constraints &c, DOpEtypes::ControlType type, const ActiveFEIndexSetterInterface< dim, dim > &index_setter=ActiveFEIndexSetterInterface< dim, dim >())
Definition: mol_multimesh_spacetimehandler.h:119
void SetUserDefinedDoFConstraints(UserDefinedDoFConstraints< DH, dim > &constraints_maker)
Definition: mol_multimesh_spacetimehandler.h:714
Definition: dopetypes.h:52
void RefineSpace(DOpEtypes::RefinementType ref_type=DOpEtypes::RefinementType::global)
Definition: mol_multimesh_spacetimehandler.h:464
Definition: dopetypes.h:52
const std::vector< Point< dim > > & GetMapDoFToSupportPoints()
Definition: mol_multimesh_spacetimehandler.h:408
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:673
double get_k() const
Definition: timeiterator.h:195
DOpEtypes::RefinementType GetRefType() const
Definition: refinementcontainer.cc:95
RefinementType
Definition: dopetypes.h:50
const dealii::ConstraintMatrix & GetStateDoFConstraints() const
Definition: mol_multimesh_spacetimehandler.h:308
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:744
const std::vector< unsigned int > & GetDoFsPerBlock(std::string name) const
Definition: constraints.h:159
unsigned int GetNGlobalConstraints() const
Definition: mol_multimesh_spacetimehandler.h:389
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:340
Definition: dopetypes.h:52
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:422
Definition: constraints.h:33
void SpatialMeshTransferState(const VECTOR &old_values, VECTOR &new_values) const
Definition: mol_multimesh_spacetimehandler.h:698
virtual double GetConvergenceOrder() const
Definition: refinementcontainer.cc:85
const dealii::ConstraintMatrix & GetControlDoFConstraints() const
Definition: mol_multimesh_spacetimehandler.h:300
const std::vector< unsigned int > & GetControlDoFsPerBlock(int=-1) const
Definition: mol_multimesh_spacetimehandler.h:276
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 InterpolateControl(VECTOR &result, const std::vector< VECTOR * > &local_vectors, double t, const TimeIterator &it) const
Definition: mol_multimesh_spacetimehandler.h:318
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:435
unsigned int GetStateNDoFs(int=-1) const
Definition: mol_multimesh_spacetimehandler.h:373
void IncrementControlTicket()
Definition: spacetimehandler_base.h:446
Definition: solutiontransfer_wrapper.h:51
unsigned int GetNLocalConstraints() const
Definition: mol_multimesh_spacetimehandler.h:398
void RefineStateSpace(const RefinementContainer &ref_container)
Definition: mol_multimesh_spacetimehandler.h:537
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
MethodOfLines_MultiMesh_SpaceTimeHandler(dealii::Triangulation< dim > &triangulation, const FE< dim, dim > &control_fe, const FE< dim, dim > &state_fe, dealii::Triangulation< 1 > ×, DOpEtypes::ControlType type, const ActiveFEIndexSetterInterface< dim, dim > &index_setter=ActiveFEIndexSetterInterface< dim, dim >())
Definition: mol_multimesh_spacetimehandler.h:80
Definition: dopetypes.h:52
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:686
void RefineControlSpace(const RefinementContainer &ref_container)
Definition: mol_multimesh_spacetimehandler.h:605
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:437
void ResetControlTriangulation(const dealii::Triangulation< dim > &tria)
Definition: mol_multimesh_spacetimehandler.h:764