dune-pdelab  2.7-git
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 
12 
14 
16 
18 
19 namespace Dune{
20  namespace PDELab{
21 
28  template<typename GO, typename LA0, typename LA1>
31  typename GO::Traits::MatrixBackend,
32  typename GO::Traits::TrialGridFunctionSpaceConstraints,
33  typename GO::Traits::TestGridFunctionSpaceConstraints>
34  {
35  public:
36 
38  typedef LA0 LocalAssemblerDT0;
39  typedef LA1 LocalAssemblerDT1;
40 
42 
45  typename GO::Traits::MatrixBackend,
46  typename GO::Traits::TrialGridFunctionSpaceConstraints,
47  typename GO::Traits::TestGridFunctionSpaceConstraints> Base;
48 
56 
57  typedef typename LA1::LocalPatternAssemblerEngine LocalExplicitPatternAssemblerEngine;
60 
68 
70  {
71  static_assert((std::is_same<typename LA0::Traits::Jacobian::Pattern,
72  typename LA1::Traits::Jacobian::Pattern>::value),
73  "Received two local assemblers which are non-compatible "
74  "due to different matrix pattern types");
75  static_assert((std::is_same<typename LA0::Traits::Jacobian,
76  typename LA1::Traits::Jacobian>::value),
77  "Received two local assemblers which are non-compatible "
78  "due to different jacobian types");
79  static_assert((std::is_same<typename LA0::Traits::Solution,
80  typename LA1::Traits::Solution>::value),
81  "Received two local assemblers which are non-compatible "
82  "due to different solution vector types");
83  static_assert((std::is_same<typename LA0::Traits::Residual,
84  typename LA1::Traits::Residual>::value),
85  "Received two local assemblers which are non-compatible "
86  "due to different residual vector types");
87  }
88 
90  typedef typename Traits::RangeField Real;
91 
94 
96  OneStepLocalAssembler (LA0 & la0_, LA1 & la1_, typename Traits::Residual & const_residual_)
97  : Base(la0_.trialConstraints(),la0_.testConstraints()),
98  la0(la0_), la1(la1_),
99  const_residual(const_residual_),
100  time(0.0), dt_mode(MultiplyOperator0ByDT), stage(0),
101  pattern_engine(*this), prestage_engine(*this), residual_engine(*this), jacobian_engine(*this),
102  explicit_jacobian_residual_engine(*this),
103  jacobian_apply_engine(*this)
104  { static_checks(); }
105 
109  void preStep(Real time_, Real dt_, int stages_)
110  {
111  time = time_;
112  dt = dt_;
113 
114  // This switch decides which term will be multiplied with dt
115  if(dt_mode == DivideOperator1ByDT)
116  {
117  dt_factor0 = 1.0;
118  dt_factor1 = 1.0 / dt;
119  }
120  else if(dt_mode == MultiplyOperator0ByDT)
121  {
122  dt_factor0 = dt;
123  dt_factor1 = 1.0;
124  }
125  else if(dt_mode == DoNotAssembleDT)
126  {
127  dt_factor0 = 1.0;
128  dt_factor1 = 1.0;
129  }
130  else
131  {
132  DUNE_THROW(Dune::Exception,"Unknown mode for assembling of time step size!");
133  }
134 
135  la0.preStep(time_,dt_, stages_);
136  la1.preStep(time_,dt_, stages_);
137  }
138 
140  void setMethod(const OneStepParameters & method_)
141  {
142  osp_method = & method_;
143  }
144 
146  void setStage(int stage_)
147  {
148  stage = stage_;
149  }
150 
152 
157  {
158  dt_mode = dt_mode_;
159  }
160 
162  Real timeAtStage(int stage_) const
163  {
164  return time+osp_method->d(stage_)*dt;
165  }
166 
169  {
170  return time+osp_method->d(stage)*dt;
171  }
172 
173  void setWeight(const Real weight)
174  {
175  la0.setWeight(weight);
176  la1.setWeight(weight);
177  }
178 
181 
185  (typename Traits::MatrixPattern & p)
186  {
187  pattern_engine.setPattern(p);
188  return pattern_engine;
189  }
190 
194  (const std::vector<typename Traits::Solution*> & x)
195  {
196  prestage_engine.setSolutions(x);
197  prestage_engine.setConstResidual(const_residual);
198  return prestage_engine;
199  }
200 
204  (typename Traits::Residual & r, const typename Traits::Solution & x)
205  {
206  residual_engine.setSolution(x);
207  residual_engine.setConstResidual(const_residual);
208  residual_engine.setResidual(r);
209  return residual_engine;
210  }
211 
215  (typename Traits::Jacobian & a, const typename Traits::Solution & x)
216  {
217  jacobian_engine.setSolution(x);
218  jacobian_engine.setJacobian(a);
219  return jacobian_engine;
220  }
221 
225  (typename Traits::MatrixPattern & p)
226  {
227  return la1.localPatternAssemblerEngine(p);
228  }
229 
233  (typename Traits::Jacobian & a,
234  typename Traits::Residual & r0, typename Traits::Residual & r1,
235  const std::vector<typename Traits::Solution*> & x)
236  {
237  // Init pre stage engine
238  prestage_engine.setSolutions( x );
239  prestage_engine.setConstResiduals(r0,r1);
240  explicit_jacobian_residual_engine.setLocalPreStageEngine(prestage_engine);
241 
242  // Init jacobian engine
243  explicit_jacobian_residual_engine.setLocalJacobianEngine
244  (la1.localJacobianAssemblerEngine(a,*(x[stage])));
245 
246  return explicit_jacobian_residual_engine;
247  }
248 
252  (const typename Traits::Domain & update, typename Traits::Range & result)
253  {
254  jacobian_apply_engine.setUpdate(update);
255  jacobian_apply_engine.setResult(result);
256  return jacobian_apply_engine;
257  }
258 
262  (const typename Traits::Domain & solution, const typename Traits::Domain & update, typename Traits::Range & result)
263  {
264  jacobian_apply_engine.setSolution(solution);
265  jacobian_apply_engine.setUpdate(update);
266  jacobian_apply_engine.setResult(result);
267  return jacobian_apply_engine;
268  }
269 
271 
272  private:
273 
277  LA0 & la0;
278  LA1 & la1;
280 
283  const OneStepParameters * osp_method;
284 
286  typename Traits::Residual & const_residual;
287 
289  Real time;
290 
292  Real dt;
293 
308  Real dt_factor0, dt_factor1;
309 
312  DTAssemblingMode dt_mode;
313 
315  int stage;
316 
319  LocalPatternAssemblerEngine pattern_engine;
320  LocalPreStageAssemblerEngine prestage_engine;
321  LocalResidualAssemblerEngine residual_engine;
322  LocalJacobianAssemblerEngine jacobian_engine;
323  LocalExplicitJacobianResidualAssemblerEngine explicit_jacobian_residual_engine;
324  LocalJacobianApplyAssemblerEngine jacobian_apply_engine;
326  };
327 
328  }
329 }
330 #endif // DUNE_PDELAB_GRIDOPERATOR_ONESTEP_LOCALASSEMBLER_HH
const P & p
Definition: constraints.hh:148
virtual R d(int r) const =0
Return entries of the d Vector.
For backward compatibility – Do not use this!
Definition: adaptivity.hh:28
Definition: assemblerutilities.hh:23
GO::Traits::Range Residual
The type of the range (residual).
Definition: assemblerutilities.hh:60
MatrixBackend::template Pattern< Jacobian, TestGridFunctionSpace, TrialGridFunctionSpace > MatrixPattern
The matrix pattern.
Definition: assemblerutilities.hh:74
GO::Traits::Jacobian Jacobian
The type of the jacobian.
Definition: assemblerutilities.hh:67
GO::Traits::Domain Solution
The type of the domain (solution).
Definition: assemblerutilities.hh:50
GO::Traits::RangeField RangeField
The field type of the range (residual).
Definition: assemblerutilities.hh:54
GO::Traits::Range Range
The type of the range (residual).
Definition: assemblerutilities.hh:57
GO::Traits::Domain Domain
The type of the domain (solution).
Definition: assemblerutilities.hh:47
Base class for local assembler.
Definition: assemblerutilities.hh:189
const GO::Traits::TrialGridFunctionSpaceConstraints & trialConstraints() const
get the constraints on the trial grid function space
Definition: assemblerutilities.hh:205
const GO::Traits::TestGridFunctionSpaceConstraints & testConstraints() const
get the constraints on the test grid function space
Definition: assemblerutilities.hh:211
void setSolution(const Domain &solution_)
Definition: onestep/jacobianapplyengine.hh:73
void setUpdate(const Domain &update_)
Definition: onestep/jacobianapplyengine.hh:80
void setResult(Range &result_)
Definition: onestep/jacobianapplyengine.hh:87
void setJacobian(Jacobian &jacobian_)
Definition: onestep/jacobianengine.hh:78
void setSolution(const Solution &solution_)
Definition: onestep/jacobianengine.hh:71
void setLocalPreStageEngine(PreStageEngine &prestage_engine_)
Definition: jacobianresidualengine.hh:53
void setLocalJacobianEngine(JacobianEngine &jacobian_engine_)
Definition: jacobianresidualengine.hh:58
The local assembler for one step methods.
Definition: onestep/localassembler.hh:34
Dune::PDELab::LocalAssemblerBase< typename GO::Traits::MatrixBackend, typename GO::Traits::TrialGridFunctionSpaceConstraints, typename GO::Traits::TestGridFunctionSpaceConstraints > Base
The base class.
Definition: onestep/localassembler.hh:47
Dune::PDELab::LocalAssemblerTraits< GO > Traits
Definition: onestep/localassembler.hh:41
void setStage(int stage_)
Set the current stage of the one step scheme.
Definition: onestep/localassembler.hh:146
LA1 LocalAssemblerDT1
Definition: onestep/localassembler.hh:39
void setWeight(const Real weight)
Definition: onestep/localassembler.hh:173
Real timeAtStage(int stage_) const
Access time at given stage.
Definition: onestep/localassembler.hh:162
OneStepExplicitLocalJacobianResidualAssemblerEngine< OneStepLocalAssembler > LocalExplicitJacobianResidualAssemblerEngine
Definition: onestep/localassembler.hh:59
Real timeAtStage() const
Access time at given stage.
Definition: onestep/localassembler.hh:168
void setMethod(const OneStepParameters &method_)
Set the one step method parameters.
Definition: onestep/localassembler.hh:140
void setDTAssemblingMode(DTAssemblingMode dt_mode_)
Definition: onestep/localassembler.hh:156
LA0 LocalAssemblerDT0
The types of the local assemblers of order one and zero.
Definition: onestep/localassembler.hh:38
OneStepLocalJacobianAssemblerEngine< OneStepLocalAssembler > LocalJacobianAssemblerEngine
Definition: onestep/localassembler.hh:54
LocalJacobianApplyAssemblerEngine & localJacobianApplyAssemblerEngine(const typename Traits::Domain &update, typename Traits::Range &result)
Definition: onestep/localassembler.hh:252
void static_checks()
Definition: onestep/localassembler.hh:69
LocalJacobianApplyAssemblerEngine & localJacobianApplyAssemblerEngine(const typename Traits::Domain &solution, const typename Traits::Domain &update, typename Traits::Range &result)
Definition: onestep/localassembler.hh:262
OneStepLocalPreStageAssemblerEngine< OneStepLocalAssembler > LocalPreStageAssemblerEngine
Definition: onestep/localassembler.hh:52
DTAssemblingMode
Definition: onestep/localassembler.hh:151
@ DoNotAssembleDT
Definition: onestep/localassembler.hh:151
@ DivideOperator1ByDT
Definition: onestep/localassembler.hh:151
@ MultiplyOperator0ByDT
Definition: onestep/localassembler.hh:151
OneStepLocalPatternAssemblerEngine< OneStepLocalAssembler > LocalPatternAssemblerEngine
Definition: onestep/localassembler.hh:51
OneStepLocalAssembler(LA0 &la0_, LA1 &la1_, typename Traits::Residual &const_residual_)
Constructor with empty constraints.
Definition: onestep/localassembler.hh:96
void preStep(Real time_, Real dt_, int stages_)
Definition: onestep/localassembler.hh:109
LA1::LocalPatternAssemblerEngine LocalExplicitPatternAssemblerEngine
Definition: onestep/localassembler.hh:57
OneStepLocalResidualAssemblerEngine< OneStepLocalAssembler > LocalResidualAssemblerEngine
Definition: onestep/localassembler.hh:53
LocalExplicitPatternAssemblerEngine & localExplicitPatternAssemblerEngine(typename Traits::MatrixPattern &p)
Definition: onestep/localassembler.hh:225
LocalPatternAssemblerEngine & localPatternAssemblerEngine(typename Traits::MatrixPattern &p)
Definition: onestep/localassembler.hh:185
LocalPreStageAssemblerEngine & localPreStageAssemblerEngine(const std::vector< typename Traits::Solution * > &x)
Definition: onestep/localassembler.hh:194
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:233
LocalResidualAssemblerEngine & localResidualAssemblerEngine(typename Traits::Residual &r, const typename Traits::Solution &x)
Definition: onestep/localassembler.hh:204
OneStepLocalJacobianApplyAssemblerEngine< OneStepLocalAssembler > LocalJacobianApplyAssemblerEngine
Definition: onestep/localassembler.hh:55
Dune::PDELab::TimeSteppingParameterInterface< Real > OneStepParameters
The type of the one step parameter object.
Definition: onestep/localassembler.hh:93
Traits::RangeField Real
The local operators type for real numbers e.g. time.
Definition: onestep/localassembler.hh:90
LocalJacobianAssemblerEngine & localJacobianAssemblerEngine(typename Traits::Jacobian &a, const typename Traits::Solution &x)
Definition: onestep/localassembler.hh:215
void setPattern(Pattern &pattern_)
Definition: onestep/patternengine.hh:62
void setSolutions(const Solutions &solutions_)
Definition: prestageengine.hh:88
void setConstResiduals(Residual &const_residual_0_, Residual &const_residual_1_)
Definition: prestageengine.hh:95
void setConstResidual(Residual &const_residual_)
Definition: prestageengine.hh:108
void setConstResidual(const Residual &const_residual_)
Definition: onestep/residualengine.hh:88
void setSolution(const Solution &solution_)
Definition: onestep/residualengine.hh:81
void setResidual(Residual &residual_)
Definition: onestep/residualengine.hh:96
Base parameter class for time stepping scheme parameters.
Definition: onestepparameter.hh:44
static const unsigned int value
Definition: gridfunctionspace/tags.hh:139