24 #ifndef MOL_STATESPACE_TIME_HANDLER_H_
25 #define MOL_STATESPACE_TIME_HANDLER_H_
34 #include <dofs/dof_handler.h>
35 #include <dofs/dof_renumbering.h>
36 #include <dofs/dof_tools.h>
37 #include <lac/constraint_matrix.h>
38 #include <deal.II/numerics/solution_transfer.h>
39 #include <deal.II/grid/grid_refinement.h>
51 template<
template<
int,
int>
class FE,
template<
int,
int>
class DH,
typename SPARSITYPATTERN,
52 typename VECTOR,
int dealdim>
54 DH, SPARSITYPATTERN, VECTOR, dealdim>
58 dealii::Triangulation<dealdim>& triangulation,
const FE<dealdim, dealdim>& state_fe,
62 dealdim>(index_setter), sparse_mkr_dynamic_(true), triangulation_(
63 triangulation), state_dof_handler_(triangulation_), state_fe_(
65 &DOpEWrapper::StaticMappingQ1<dealdim, DH>::mapping_q1), state_mesh_transfer_(
69 user_defined_dof_constr_ = NULL;
72 dealii::Triangulation<dealdim>& triangulation,
const FE<dealdim, dealdim>& state_fe,
73 dealii::Triangulation<1> & times,
77 dealdim>(times, index_setter), sparse_mkr_dynamic_(true), triangulation_(
78 triangulation), state_dof_handler_(triangulation_), state_fe_(
80 &DOpEWrapper::StaticMappingQ1<dealdim, DH>::mapping_q1), state_mesh_transfer_(
84 user_defined_dof_constr_ = NULL;
88 dealii::Triangulation<dealdim>& triangulation,
90 const FE<dealdim, dealdim>& state_fe,
94 dealdim>(index_setter), sparse_mkr_dynamic_(true), triangulation_(
95 triangulation), state_dof_handler_(triangulation_), state_fe_(
96 &state_fe), mapping_(&mapping), state_mesh_transfer_(NULL)
99 user_defined_dof_constr_ = NULL;
102 dealii::Triangulation<dealdim>& triangulation,
104 const FE<dealdim, dealdim>& state_fe,
105 dealii::Triangulation<1> & times,
109 dealdim>(times, index_setter), sparse_mkr_dynamic_(true), triangulation_(
110 triangulation), state_dof_handler_(triangulation_), state_fe_(
111 &state_fe), mapping_(&mapping), state_mesh_transfer_(NULL)
114 user_defined_dof_constr_ = NULL;
120 state_dof_handler_.clear();
122 if (state_mesh_transfer_ != NULL)
124 delete state_mesh_transfer_;
126 if (sparsitymaker_ != NULL && sparse_mkr_dynamic_ ==
true)
128 delete sparsitymaker_;
137 const std::vector<unsigned int>& state_block_component)
142 state_dof_handler_.distribute_dofs(
GetFESystem(
"state"));
143 DoFRenumbering::Cuthill_McKee(
144 static_cast<DH<dealdim, dealdim>&
>(state_dof_handler_));
145 DoFRenumbering::component_wise(
146 static_cast<DH<dealdim, dealdim>&
>(state_dof_handler_));
148 state_dof_constraints_.clear();
149 DoFTools::make_hanging_node_constraints(
150 static_cast<DH<dealdim, dealdim>&
>(state_dof_handler_),
151 state_dof_constraints_);
153 if (GetUserDefinedDoFConstraints() != NULL
156 state_dof_handler_, state_dof_constraints_);
157 state_dof_constraints_.close();
158 state_dofs_per_block_.resize(state_n_blocks);
160 DoFTools::count_dofs_per_block(
161 static_cast<DH<dealdim, dealdim>&
>(state_dof_handler_),
162 state_dofs_per_block_, state_block_component);
164 support_points_.clear();
180 return state_dof_handler_;
198 return state_dofs_per_block_[b];
204 const std::vector<unsigned int>&
207 return state_dofs_per_block_;
213 const dealii::ConstraintMatrix&
216 return state_dof_constraints_;
224 const std::vector<VECTOR*> & local_vectors,
double t,
229 if (local_vectors.size() != 2)
231 "This function is currently not implemented for anything other than"
232 " linear interpolation of 2 DoFs.",
233 "MethodOfLine_SpaceTimeHandler::InterpolateState");
239 result = *local_vectors[0];
241 result.sadd(lambda_l, lambda_r, *local_vectors[1]);
256 const std::vector<Point<dealdim> >&
261 < std::vector<Point<dealdim> >, dealdim
263 return support_points_;
280 const FE<dealdim, dealdim>&
291 "MethodOfLines_StateSpaceTimeHandler::GetFESystem");
334 if (state_mesh_transfer_ != NULL)
336 delete state_mesh_transfer_;
337 state_mesh_transfer_ = NULL;
339 state_mesh_transfer_ =
new dealii::SolutionTransfer<dealdim, VECTOR, DH<dealdim, dealdim> >(
344 triangulation_.set_all_refine_flags();
348 GridRefinement::refine_and_coarsen_fixed_number(triangulation_,
356 GridRefinement::refine_and_coarsen_fixed_fraction(triangulation_,
364 GridRefinement::refine_and_coarsen_optimize(triangulation_,
371 "MethodOfLines_StateSpaceTimeHandler::RefineStateSpace");
373 triangulation_.prepare_coarsening_and_refinement();
374 if (state_mesh_transfer_ != NULL
376 state_mesh_transfer_->prepare_for_pure_refinement();
377 triangulation_.execute_coarsening_and_refinement();
399 VECTOR& new_values)
const
401 if (state_mesh_transfer_ != NULL
403 state_mesh_transfer_->refine_interpolate(old_values, new_values);
419 user_defined_dof_constr_ = &user_defined_dof_constr;
433 if (sparsitymaker_ != NULL && sparse_mkr_dynamic_)
434 delete sparsitymaker_;
435 sparsitymaker_ = &sparsity_maker;
436 sparse_mkr_dynamic_ =
false;
441 GetSparsityMaker()
const
443 return sparsitymaker_;
445 const UserDefinedDoFConstraints<DH, dealdim>*
446 GetUserDefinedDoFConstraints()
const
448 return user_defined_dof_constr_;
450 SparsityMaker<DH, dealdim>* sparsitymaker_;
451 UserDefinedDoFConstraints<DH, dealdim>* user_defined_dof_constr_;
452 bool sparse_mkr_dynamic_;
454 dealii::Triangulation<dealdim>& triangulation_;
457 std::vector<unsigned int> state_dofs_per_block_;
459 dealii::ConstraintMatrix state_dof_constraints_;
461 const dealii::SmartPointer<const FE<dealdim, dealdim> > state_fe_;
462 const dealii::SmartPointer<
465 std::vector<Point<dealdim> > support_points_;
466 dealii::SolutionTransfer<dealdim, VECTOR, DH<dealdim, dealdim> >* state_mesh_transfer_;
const dealii::ConstraintMatrix & GetStateDoFConstraints() const
Definition: mol_statespacetimehandler.h:214
void ReInitTime()
Definition: spacetimehandler_base.h:86
void RegisterMapping(const typename DOpEWrapper::Mapping< dealdim, DH > &mapping)
Definition: userdefineddofconstraints.h:76
Definition: dopetypes.h:52
void ReInit(unsigned int state_n_blocks, const std::vector< unsigned int > &state_block_component)
Definition: mol_statespacetimehandler.h:136
MethodOfLines_StateSpaceTimeHandler(dealii::Triangulation< dealdim > &triangulation, const DOpEWrapper::Mapping< dealdim, DH > &mapping, const FE< dealdim, dealdim > &state_fe, const ActiveFEIndexSetterInterface< dealdim > &index_setter=ActiveFEIndexSetterInterface< dealdim >())
Definition: mol_statespacetimehandler.h:87
void RefineSpace(DOpEtypes::RefinementType ref_type=DOpEtypes::RefinementType::global)
Definition: mol_statespacetimehandler.h:308
void SetActiveFEIndicesState(DOpEWrapper::DoFHandler< dealdim, DH > &dof_handler)
Definition: statespacetimehandler.h:253
Definition: dopetypes.h:52
unsigned int NewTimePointToOldTimePoint(unsigned int t) const
Definition: mol_statespacetimehandler.h:385
void RefineSpace(const RefinementContainer &ref_container)
Definition: mol_statespacetimehandler.h:327
Definition: mol_statespacetimehandler.h:53
const DOpEWrapper::Mapping< dealdim, DH > & GetMapping() const
Definition: mol_statespacetimehandler.h:187
unsigned int GetStateNDoFs(int=-1) const
Definition: mol_statespacetimehandler.h:248
Definition: timeiterator.h:63
double get_k() const
Definition: timeiterator.h:195
DOpEtypes::RefinementType GetRefType() const
Definition: refinementcontainer.cc:95
RefinementType
Definition: dopetypes.h:50
void SetSparsityMaker(SparsityMaker< DH, dealdim > &sparsity_maker)
Definition: mol_statespacetimehandler.h:431
void ComputeStateSparsityPattern(SPARSITYPATTERN &sparsity) const
Definition: mol_statespacetimehandler.h:268
Definition: refinementcontainer.h:42
MethodOfLines_StateSpaceTimeHandler(dealii::Triangulation< dealdim > &triangulation, const FE< dealdim, dealdim > &state_fe, const ActiveFEIndexSetterInterface< dealdim > &index_setter=ActiveFEIndexSetterInterface< dealdim >())
Definition: mol_statespacetimehandler.h:57
virtual const dealii::Vector< float > & GetLocalErrorIndicators() const
Definition: refinementcontainer.cc:46
virtual double GetTopFraction() const
Definition: refinementcontainer.cc:56
Definition: statespacetimehandler.h:59
Definition: dopetypes.h:52
virtual ~MethodOfLines_StateSpaceTimeHandler()
Definition: mol_statespacetimehandler.h:118
virtual double GetConvergenceOrder() const
Definition: refinementcontainer.cc:85
const DOpEWrapper::DoFHandler< dealdim, DH > & GetStateDoFHandler() const
Definition: mol_statespacetimehandler.h:177
double get_right() const
Definition: timeiterator.h:185
MethodOfLines_StateSpaceTimeHandler(dealii::Triangulation< dealdim > &triangulation, const DOpEWrapper::Mapping< dealdim, DH > &mapping, const FE< dealdim, dealdim > &state_fe, dealii::Triangulation< 1 > ×, const ActiveFEIndexSetterInterface< dealdim > &index_setter=ActiveFEIndexSetterInterface< dealdim >())
Definition: mol_statespacetimehandler.h:101
void MapDoFsToSupportPoints(const DOpEWrapper::Mapping< dealdim, dealii::DoFHandler > &mapping, const DOpEWrapper::DoFHandler< dealdim, dealii::DoFHandler > &dh, VECTOR &support_points)
Definition: sth_internals.h:47
void SetUserDefinedDoFConstraints(UserDefinedDoFConstraints< DH, dealdim > &user_defined_dof_constr)
Definition: mol_statespacetimehandler.h:416
void SpatialMeshTransferState(const VECTOR &old_values, VECTOR &new_values) const
Definition: mol_statespacetimehandler.h:398
const FE< dealdim, dealdim > & GetFESystem(std::string name) const
Definition: mol_statespacetimehandler.h:281
const std::vector< Point< dealdim > > & GetMapDoFToSupportPoints()
Definition: mol_statespacetimehandler.h:257
virtual void ComputeSparsityPattern(const DOpEWrapper::DoFHandler< dim, DH > &dof_handler, dealii::BlockSparsityPattern &sparsity, const dealii::ConstraintMatrix &hanging_node_constraints, const std::vector< unsigned int > &blocks) const
Definition: sparsitymaker.h:111
unsigned int GetStateDoFsPerBlock(unsigned int b, int=-1) const
Definition: mol_statespacetimehandler.h:196
bool UsesCoarsening() const
Definition: refinementcontainer.cc:103
const std::vector< unsigned int > & GetStateDoFsPerBlock(int=-1) const
Definition: mol_statespacetimehandler.h:205
double get_left() const
Definition: timeiterator.h:175
Definition: dopeexception.h:35
virtual void InterpolateState(VECTOR &result, const std::vector< VECTOR * > &local_vectors, double t, const TimeIterator &it) const
Definition: mol_statespacetimehandler.h:223
Definition: dopetypes.h:52
virtual double GetBottomFraction() const
Definition: refinementcontainer.cc:66
MethodOfLines_StateSpaceTimeHandler(dealii::Triangulation< dealdim > &triangulation, const FE< dealdim, dealdim > &state_fe, dealii::Triangulation< 1 > ×, const ActiveFEIndexSetterInterface< dealdim > &index_setter=ActiveFEIndexSetterInterface< dealdim >())
Definition: mol_statespacetimehandler.h:71
virtual void MakeStateDoFConstraints(const DOpEWrapper::DoFHandler< dealdim, DH > &dof_handler, dealii::ConstraintMatrix &dof_constraints) const
Definition: userdefineddofconstraints.h:93
void IncrementStateTicket()
Definition: spacetimehandler_base.h:437