24 #ifndef MOL_STATESPACE_TIME_HANDLER_H_
25 #define MOL_STATESPACE_TIME_HANDLER_H_
34 #include <deal.II/dofs/dof_handler.h>
35 #include <deal.II/dofs/dof_renumbering.h>
36 #include <deal.II/dofs/dof_tools.h>
37 #include <deal.II/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,
59 bool flux_pattern =
false,
63 dealdim>(index_setter), sparse_mkr_dynamic_(true), triangulation_(
64 triangulation), state_dof_handler_(triangulation_), state_fe_(
66 &DOpEWrapper::StaticMappingQ1<dealdim, DH>::mapping_q1), state_mesh_transfer_(
70 user_defined_dof_constr_ = NULL;
73 dealii::Triangulation<dealdim>& triangulation,
const FE<dealdim, dealdim>& state_fe,
74 dealii::Triangulation<1> & times,
75 bool flux_pattern =
false,
79 dealdim>(times, index_setter), sparse_mkr_dynamic_(true), triangulation_(
80 triangulation), state_dof_handler_(triangulation_), state_fe_(
82 &DOpEWrapper::StaticMappingQ1<dealdim, DH>::mapping_q1), state_mesh_transfer_(
86 user_defined_dof_constr_ = NULL;
90 dealii::Triangulation<dealdim>& triangulation,
92 const FE<dealdim, dealdim>& state_fe,
93 bool flux_pattern =
false,
97 dealdim>(index_setter), sparse_mkr_dynamic_(true), triangulation_(
98 triangulation), state_dof_handler_(triangulation_), state_fe_(
99 &state_fe), mapping_(&mapping), state_mesh_transfer_(NULL)
102 user_defined_dof_constr_ = NULL;
105 dealii::Triangulation<dealdim>& triangulation,
107 const FE<dealdim, dealdim>& state_fe,
108 dealii::Triangulation<1> & times,
109 bool flux_pattern =
false,
113 dealdim>(times, index_setter), sparse_mkr_dynamic_(true), triangulation_(
114 triangulation), state_dof_handler_(triangulation_), state_fe_(
115 &state_fe), mapping_(&mapping), state_mesh_transfer_(NULL)
118 user_defined_dof_constr_ = NULL;
124 state_dof_handler_.clear();
126 if (state_mesh_transfer_ != NULL)
128 delete state_mesh_transfer_;
130 if (sparsitymaker_ != NULL && sparse_mkr_dynamic_ ==
true)
132 delete sparsitymaker_;
141 const std::vector<unsigned int>& state_block_component)
146 state_dof_handler_.distribute_dofs(
GetFESystem(
"state"));
149 DoFRenumbering::component_wise(
150 static_cast<DH<dealdim, dealdim>&
>(state_dof_handler_));
152 state_dof_constraints_.clear();
153 DoFTools::make_hanging_node_constraints(
154 static_cast<DH<dealdim, dealdim>&
>(state_dof_handler_),
155 state_dof_constraints_);
157 if (GetUserDefinedDoFConstraints() != NULL
160 state_dof_handler_, state_dof_constraints_);
161 state_dof_constraints_.close();
162 state_dofs_per_block_.resize(state_n_blocks);
164 DoFTools::count_dofs_per_block(
165 static_cast<DH<dealdim, dealdim>&
>(state_dof_handler_),
166 state_dofs_per_block_, state_block_component);
168 support_points_.clear();
184 return state_dof_handler_;
202 return state_dofs_per_block_[b];
208 const std::vector<unsigned int>&
211 return state_dofs_per_block_;
217 const dealii::ConstraintMatrix&
220 return state_dof_constraints_;
228 const std::vector<VECTOR*> & local_vectors,
double t,
233 if (local_vectors.size() != 2)
235 "This function is currently not implemented for anything other than"
236 " linear interpolation of 2 DoFs.",
237 "MethodOfLine_SpaceTimeHandler::InterpolateState");
243 result = *local_vectors[0];
245 result.sadd(lambda_l, lambda_r, *local_vectors[1]);
260 const std::vector<Point<dealdim> >&
265 < std::vector<Point<dealdim> >, dealdim
267 return support_points_;
284 const FE<dealdim, dealdim>&
295 "MethodOfLines_StateSpaceTimeHandler::GetFESystem");
338 if (state_mesh_transfer_ != NULL)
340 delete state_mesh_transfer_;
341 state_mesh_transfer_ = NULL;
343 state_mesh_transfer_ =
new dealii::SolutionTransfer<dealdim, VECTOR, DH<dealdim, dealdim> >(
348 triangulation_.set_all_refine_flags();
352 GridRefinement::refine_and_coarsen_fixed_number(triangulation_,
360 GridRefinement::refine_and_coarsen_fixed_fraction(triangulation_,
368 GridRefinement::refine_and_coarsen_optimize(triangulation_,
375 "MethodOfLines_StateSpaceTimeHandler::RefineStateSpace");
377 triangulation_.prepare_coarsening_and_refinement();
378 if (state_mesh_transfer_ != NULL
380 state_mesh_transfer_->prepare_for_pure_refinement();
381 triangulation_.execute_coarsening_and_refinement();
403 VECTOR& new_values)
const
405 if (state_mesh_transfer_ != NULL
407 state_mesh_transfer_->refine_interpolate(old_values, new_values);
423 user_defined_dof_constr_ = &user_defined_dof_constr;
437 assert(sparse_mkr_dynamic_==
true);
438 if (sparsitymaker_ != NULL && sparse_mkr_dynamic_)
439 delete sparsitymaker_;
440 sparsitymaker_ = &sparsity_maker;
441 sparse_mkr_dynamic_ =
false;
446 GetSparsityMaker()
const
448 return sparsitymaker_;
450 const UserDefinedDoFConstraints<DH, dealdim>*
451 GetUserDefinedDoFConstraints()
const
453 return user_defined_dof_constr_;
455 SparsityMaker<DH, dealdim>* sparsitymaker_;
456 UserDefinedDoFConstraints<DH, dealdim>* user_defined_dof_constr_;
457 bool sparse_mkr_dynamic_;
459 dealii::Triangulation<dealdim>& triangulation_;
462 std::vector<unsigned int> state_dofs_per_block_;
464 dealii::ConstraintMatrix state_dof_constraints_;
466 const dealii::SmartPointer<const FE<dealdim, dealdim> > state_fe_;
467 const dealii::SmartPointer<
470 std::vector<Point<dealdim> > support_points_;
471 dealii::SolutionTransfer<dealdim, VECTOR, DH<dealdim, dealdim> >* state_mesh_transfer_;
const dealii::ConstraintMatrix & GetStateDoFConstraints() const
Definition: mol_statespacetimehandler.h:218
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:140
void RefineSpace(DOpEtypes::RefinementType=DOpEtypes::RefinementType::global)
Definition: mol_statespacetimehandler.h:312
MethodOfLines_StateSpaceTimeHandler(dealii::Triangulation< dealdim > &triangulation, const FE< dealdim, dealdim > &state_fe, bool flux_pattern=false, const ActiveFEIndexSetterInterface< dealdim > &index_setter=ActiveFEIndexSetterInterface< dealdim >())
Definition: mol_statespacetimehandler.h:57
void SetActiveFEIndicesState(DOpEWrapper::DoFHandler< dealdim, DH > &dof_handler)
Definition: statespacetimehandler.h:253
MethodOfLines_StateSpaceTimeHandler(dealii::Triangulation< dealdim > &triangulation, const FE< dealdim, dealdim > &state_fe, dealii::Triangulation< 1 > ×, bool flux_pattern=false, const ActiveFEIndexSetterInterface< dealdim > &index_setter=ActiveFEIndexSetterInterface< dealdim >())
Definition: mol_statespacetimehandler.h:72
Definition: dopetypes.h:52
unsigned int NewTimePointToOldTimePoint(unsigned int t) const
Definition: mol_statespacetimehandler.h:389
void RefineSpace(const RefinementContainer &ref_container)
Definition: mol_statespacetimehandler.h:331
Definition: mol_statespacetimehandler.h:53
const DOpEWrapper::Mapping< dealdim, DH > & GetMapping() const
Definition: mol_statespacetimehandler.h:191
unsigned int GetStateNDoFs(int=-1) const
Definition: mol_statespacetimehandler.h:252
Definition: timeiterator.h:62
double get_k() const
Definition: timeiterator.h:194
DOpEtypes::RefinementType GetRefType() const
Definition: refinementcontainer.cc:86
RefinementType
Definition: dopetypes.h:50
void SetSparsityMaker(SparsityMaker< DH, dealdim > &sparsity_maker)
Definition: mol_statespacetimehandler.h:435
void ComputeStateSparsityPattern(SPARSITYPATTERN &sparsity) const
Definition: mol_statespacetimehandler.h:272
Definition: refinementcontainer.h:42
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:122
virtual double GetConvergenceOrder() const
Definition: refinementcontainer.cc:76
const DOpEWrapper::DoFHandler< dealdim, DH > & GetStateDoFHandler() const
Definition: mol_statespacetimehandler.h:181
double get_right() const
Definition: timeiterator.h:184
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:420
void SpatialMeshTransferState(const VECTOR &old_values, VECTOR &new_values) const
Definition: mol_statespacetimehandler.h:402
const FE< dealdim, dealdim > & GetFESystem(std::string name) const
Definition: mol_statespacetimehandler.h:285
const std::vector< Point< dealdim > > & GetMapDoFToSupportPoints()
Definition: mol_statespacetimehandler.h:261
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:113
unsigned int GetStateDoFsPerBlock(unsigned int b, int=-1) const
Definition: mol_statespacetimehandler.h:200
MethodOfLines_StateSpaceTimeHandler(dealii::Triangulation< dealdim > &triangulation, const DOpEWrapper::Mapping< dealdim, DH > &mapping, const FE< dealdim, dealdim > &state_fe, bool flux_pattern=false, const ActiveFEIndexSetterInterface< dealdim > &index_setter=ActiveFEIndexSetterInterface< dealdim >())
Definition: mol_statespacetimehandler.h:89
bool UsesCoarsening() const
Definition: refinementcontainer.cc:94
const std::vector< unsigned int > & GetStateDoFsPerBlock(int=-1) const
Definition: mol_statespacetimehandler.h:209
double get_left() const
Definition: timeiterator.h:174
MethodOfLines_StateSpaceTimeHandler(dealii::Triangulation< dealdim > &triangulation, const DOpEWrapper::Mapping< dealdim, DH > &mapping, const FE< dealdim, dealdim > &state_fe, dealii::Triangulation< 1 > ×, bool flux_pattern=false, const ActiveFEIndexSetterInterface< dealdim > &index_setter=ActiveFEIndexSetterInterface< dealdim >())
Definition: mol_statespacetimehandler.h:104
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:227
Definition: dopetypes.h:52
virtual double GetBottomFraction() const
Definition: refinementcontainer.cc:66
std::string DOpEtypesToString(const C &)
Definition: dopetypes.h:134
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