dune-pdelab  2.5-dev
onestep/localassembler.hh
Go to the documentation of this file.
1 #ifndef DUNE_PDELAB_GRIDOPERATOR_ONESTEP_LOCALASSEMBLER_HH
2 #define DUNE_PDELAB_GRIDOPERATOR_ONESTEP_LOCALASSEMBLER_HH
3 
4 #include <dune/typetree/typetree.hh>
5 
11 
13 
15 
17 
18 namespace Dune{
19  namespace PDELab{
20 
27  template<typename GO, typename LA0, typename LA1>
30  typename GO::Traits::MatrixBackend,
31  typename GO::Traits::TrialGridFunctionSpaceConstraints,
32  typename GO::Traits::TestGridFunctionSpaceConstraints>
33  {
34  public:
35 
37  typedef LA0 LocalAssemblerDT0;
38  typedef LA1 LocalAssemblerDT1;
39 
41 
44  <typename GO::Traits::MatrixBackend,
45  typename GO::Traits::TrialGridFunctionSpaceConstraints,
46  typename GO::Traits::TestGridFunctionSpaceConstraints> Base;
47 
54 
55  typedef typename LA1::LocalPatternAssemblerEngine LocalExplicitPatternAssemblerEngine;
58 
65 
66  void static_checks(){
67  static_assert((std::is_same<typename LA0::Traits::Jacobian::Pattern,
68  typename LA1::Traits::Jacobian::Pattern>::value),
69  "Received two local assemblers which are non-compatible "
70  "due to different matrix pattern types");
71  static_assert((std::is_same<typename LA0::Traits::Jacobian,
72  typename LA1::Traits::Jacobian>::value),
73  "Received two local assemblers which are non-compatible "
74  "due to different jacobian types");
75  static_assert((std::is_same<typename LA0::Traits::Solution,
76  typename LA1::Traits::Solution>::value),
77  "Received two local assemblers which are non-compatible "
78  "due to different solution vector types");
79  static_assert((std::is_same<typename LA0::Traits::Residual,
80  typename LA1::Traits::Residual>::value),
81  "Received two local assemblers which are non-compatible "
82  "due to different residual vector types");
83  }
84 
86  typedef typename Traits::RangeField Real;
87 
90 
92  OneStepLocalAssembler (LA0 & la0_, LA1 & la1_, typename Traits::Residual & const_residual_)
93  : Base(la0_.trialConstraints(),la0_.testConstraints()),
94  la0(la0_), la1(la1_),
95  const_residual(const_residual_),
96  time(0.0), dt_mode(MultiplyOperator0ByDT), stage(0),
97  pattern_engine(*this), prestage_engine(*this), residual_engine(*this), jacobian_engine(*this),
98  explicit_jacobian_residual_engine(*this)
99  { static_checks(); }
100 
104  void preStep(Real time_, Real dt_, int stages_){
105  time = time_;
106  dt = dt_;
107 
108  // This switch decides which term will be multiplied with dt
109  if(dt_mode == DivideOperator1ByDT){
110  dt_factor0 = 1.0;
111  dt_factor1 = 1.0 / dt;
112  }
113  else if(dt_mode == MultiplyOperator0ByDT){
114  dt_factor0 = dt;
115  dt_factor1 = 1.0;
116  }
117  else if(dt_mode == DoNotAssembleDT){
118  dt_factor0 = 1.0;
119  dt_factor1 = 1.0;
120  }
121  else{
122  DUNE_THROW(Dune::Exception,"Unknown mode for assembling of time step size!");
123  }
124 
125  la0.preStep(time_,dt_, stages_);
126  la1.preStep(time_,dt_, stages_);
127  }
128 
130  void setMethod(const OneStepParameters & method_){
131  osp_method = & method_;
132  }
133 
135  void setStage(int stage_){
136  stage = stage_;
137  }
138 
140 
145  dt_mode = dt_mode_;
146  }
147 
149  Real timeAtStage(int stage_){
150  return time+osp_method->d(stage_)*dt;
151  }
152 
154  Real timeAtStage(){
155  return time+osp_method->d(stage)*dt;
156  }
157 
158  void setWeight(const Real weight){
159  la0.setWeight(weight);
160  la1.setWeight(weight);
161  }
162 
165 
168  LocalPatternAssemblerEngine & localPatternAssemblerEngine
170  {
171  pattern_engine.setPattern(p);
172  return pattern_engine;
173  }
174 
177  LocalPreStageAssemblerEngine & localPreStageAssemblerEngine
178  (const std::vector<typename Traits::Solution*> & x)
179  {
180  prestage_engine.setSolutions(x);
181  prestage_engine.setConstResidual(const_residual);
182  return prestage_engine;
183  }
184 
187  LocalResidualAssemblerEngine & localResidualAssemblerEngine
188  (typename Traits::Residual & r, const typename Traits::Solution & x)
189  {
190  residual_engine.setSolution(x);
191  residual_engine.setConstResidual(const_residual);
192  residual_engine.setResidual(r);
193  return residual_engine;
194  }
195 
198  LocalJacobianAssemblerEngine & localJacobianAssemblerEngine
199  (typename Traits::Jacobian & a, const typename Traits::Solution & x)
200  {
201  jacobian_engine.setSolution(x);
202  jacobian_engine.setJacobian(a);
203  return jacobian_engine;
204  }
205 
208  LocalExplicitPatternAssemblerEngine & localExplicitPatternAssemblerEngine
209  (typename Traits::MatrixPattern & p)
210  {
211  return la1.localPatternAssemblerEngine(p);
212  }
213 
217  (typename Traits::Jacobian & a,
218  typename Traits::Residual & r0, typename Traits::Residual & r1,
219  const std::vector<typename Traits::Solution*> & x)
220  {
221  // Init pre stage engine
222  prestage_engine.setSolutions( x );
223  prestage_engine.setConstResiduals(r0,r1);
224  explicit_jacobian_residual_engine.setLocalPreStageEngine(prestage_engine);
225 
226  // Init jacobian engine
227  explicit_jacobian_residual_engine.setLocalJacobianEngine
228  (la1.localJacobianAssemblerEngine(a,*(x[stage])));
229 
230  return explicit_jacobian_residual_engine;
231  }
232 
234 
235  private:
236 
240  LA0 & la0;
241  LA1 & la1;
243 
246  const OneStepParameters * osp_method;
247 
249  typename Traits::Residual & const_residual;
250 
252  Real time;
253 
255  Real dt;
256 
271  Real dt_factor0, dt_factor1;
272 
275  DTAssemblingMode dt_mode;
276 
278  int stage;
279 
282  LocalPatternAssemblerEngine pattern_engine;
283  LocalPreStageAssemblerEngine prestage_engine;
284  LocalResidualAssemblerEngine residual_engine;
285  LocalJacobianAssemblerEngine jacobian_engine;
286  LocalExplicitJacobianResidualAssemblerEngine explicit_jacobian_residual_engine;
288  };
289 
290  }
291 }
292 #endif // DUNE_PDELAB_GRIDOPERATOR_ONESTEP_LOCALASSEMBLER_HH
OneStepLocalResidualAssemblerEngine< OneStepLocalAssembler > LocalResidualAssemblerEngine
Definition: onestep/localassembler.hh:52
Definition: onestep/localassembler.hh:139
LocalResidualAssemblerEngine & localResidualAssemblerEngine(typename Traits::Residual &r, const typename Traits::Solution &x)
Definition: onestep/localassembler.hh:188
OneStepExplicitLocalJacobianResidualAssemblerEngine< OneStepLocalAssembler > LocalExplicitJacobianResidualAssemblerEngine
Definition: onestep/localassembler.hh:57
void setJacobian(Jacobian &jacobian_)
Definition: onestep/jacobianengine.hh:77
Traits::RangeField Real
The local operators type for real numbers e.g. time.
Definition: onestep/localassembler.hh:86
const GO::Traits::TestGridFunctionSpaceConstraints & testConstraints() const
get the constraints on the test grid function space
Definition: assemblerutilities.hh:204
void setWeight(const Real weight)
Definition: onestep/localassembler.hh:158
void setConstResidual(Residual &const_residual_)
Definition: prestageengine.hh:106
LocalPreStageAssemblerEngine & localPreStageAssemblerEngine(const std::vector< typename Traits::Solution *> &x)
Definition: onestep/localassembler.hh:178
Definition: onestep/localassembler.hh:139
GO::Traits::RangeField RangeField
The field type of the range (residual).
Definition: assemblerutilities.hh:51
virtual R d(int r) const =0
Return entries of the d Vector.
void setMethod(const OneStepParameters &method_)
Set the one step method parameters.
Definition: onestep/localassembler.hh:130
void setLocalPreStageEngine(PreStageEngine &prestage_engine_)
Definition: jacobianresidualengine.hh:53
Definition: assemblerutilities.hh:22
LA1::LocalPatternAssemblerEngine LocalExplicitPatternAssemblerEngine
Definition: onestep/localassembler.hh:55
void setDTAssemblingMode(DTAssemblingMode dt_mode_)
Definition: onestep/localassembler.hh:144
LocalPatternAssemblerEngine & localPatternAssemblerEngine(typename Traits::MatrixPattern &p)
Definition: onestep/localassembler.hh:169
Dune::PDELab::LocalAssemblerBase< typename GO::Traits::MatrixBackend, typename GO::Traits::TrialGridFunctionSpaceConstraints, typename GO::Traits::TestGridFunctionSpaceConstraints > Base
The base class.
Definition: onestep/localassembler.hh:46
Base parameter class for time stepping scheme parameters.
Definition: onestepparameter.hh:23
void static_checks()
Definition: onestep/localassembler.hh:66
void setSolution(const Solution &solution_)
Definition: onestep/residualengine.hh:81
LocalJacobianAssemblerEngine & localJacobianAssemblerEngine(typename Traits::Jacobian &a, const typename Traits::Solution &x)
Definition: onestep/localassembler.hh:199
For backward compatibility – Do not use this!
Definition: adaptivity.hh:27
OneStepLocalPreStageAssemblerEngine< OneStepLocalAssembler > LocalPreStageAssemblerEngine
Definition: onestep/localassembler.hh:51
void setStage(int stage_)
Set the current stage of the one step scheme.
Definition: onestep/localassembler.hh:135
GO::Traits::Jacobian Jacobian
The type of the jacobian.
Definition: assemblerutilities.hh:61
OneStepLocalPatternAssemblerEngine< OneStepLocalAssembler > LocalPatternAssemblerEngine
Definition: onestep/localassembler.hh:50
void setPattern(Pattern &pattern_)
Definition: onestep/patternengine.hh:61
LocalExplicitPatternAssemblerEngine & localExplicitPatternAssemblerEngine(typename Traits::MatrixPattern &p)
Definition: onestep/localassembler.hh:209
static const unsigned int value
Definition: gridfunctionspace/tags.hh:139
Definition: onestep/localassembler.hh:139
void setConstResidual(const Residual &const_residual_)
Definition: onestep/residualengine.hh:87
LA0 LocalAssemblerDT0
The types of the local assemblers of order one and zero.
Definition: onestep/localassembler.hh:37
void setLocalJacobianEngine(JacobianEngine &jacobian_engine_)
Definition: jacobianresidualengine.hh:58
LA1 LocalAssemblerDT1
Definition: onestep/localassembler.hh:38
void setResidual(Residual &residual_)
Definition: onestep/residualengine.hh:94
Dune::PDELab::LocalAssemblerTraits< GO > Traits
Definition: onestep/localassembler.hh:40
Real timeAtStage()
Access time at given stage.
Definition: onestep/localassembler.hh:154
LocalExplicitJacobianResidualAssemblerEngine & localExplicitJacobianResidualAssemblerEngine(typename Traits::Jacobian &a, typename Traits::Residual &r0, typename Traits::Residual &r1, const std::vector< typename Traits::Solution *> &x)
Definition: onestep/localassembler.hh:217
GO::Traits::Range Residual
The type of the range (residual).
Definition: assemblerutilities.hh:54
GO::Traits::Domain Solution
The type of the domain (solution).
Definition: assemblerutilities.hh:47
Dune::PDELab::TimeSteppingParameterInterface< Real > OneStepParameters
The type of the one step parameter object.
Definition: onestep/localassembler.hh:89
MatrixBackend::template Pattern< Jacobian, TestGridFunctionSpace, TrialGridFunctionSpace > MatrixPattern
The matrix pattern.
Definition: assemblerutilities.hh:68
void setConstResiduals(Residual &const_residual_0_, Residual &const_residual_1_)
Definition: prestageengine.hh:94
const GO::Traits::TrialGridFunctionSpaceConstraints & trialConstraints() const
get the constraints on the trial grid function space
Definition: assemblerutilities.hh:198
void setSolutions(const Solutions &solutions_)
Definition: prestageengine.hh:88
void preStep(Real time_, Real dt_, int stages_)
Definition: onestep/localassembler.hh:104
OneStepLocalAssembler(LA0 &la0_, LA1 &la1_, typename Traits::Residual &const_residual_)
Constructor with empty constraints.
Definition: onestep/localassembler.hh:92
Real timeAtStage(int stage_)
Access time at given stage.
Definition: onestep/localassembler.hh:149
const P & p
Definition: constraints.hh:147
Base class for local assembler.
Definition: assemblerutilities.hh:182
void setSolution(const Solution &solution_)
Definition: onestep/jacobianengine.hh:71
DTAssemblingMode
Definition: onestep/localassembler.hh:139
The local assembler for one step methods.
Definition: onestep/localassembler.hh:28
OneStepLocalJacobianAssemblerEngine< OneStepLocalAssembler > LocalJacobianAssemblerEngine
Definition: onestep/localassembler.hh:53