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 
67  {
68  static_assert((std::is_same<typename LA0::Traits::Jacobian::Pattern,
69  typename LA1::Traits::Jacobian::Pattern>::value),
70  "Received two local assemblers which are non-compatible "
71  "due to different matrix pattern types");
72  static_assert((std::is_same<typename LA0::Traits::Jacobian,
73  typename LA1::Traits::Jacobian>::value),
74  "Received two local assemblers which are non-compatible "
75  "due to different jacobian types");
76  static_assert((std::is_same<typename LA0::Traits::Solution,
77  typename LA1::Traits::Solution>::value),
78  "Received two local assemblers which are non-compatible "
79  "due to different solution vector types");
80  static_assert((std::is_same<typename LA0::Traits::Residual,
81  typename LA1::Traits::Residual>::value),
82  "Received two local assemblers which are non-compatible "
83  "due to different residual vector types");
84  }
85 
87  typedef typename Traits::RangeField Real;
88 
91 
93  OneStepLocalAssembler (LA0 & la0_, LA1 & la1_, typename Traits::Residual & const_residual_)
94  : Base(la0_.trialConstraints(),la0_.testConstraints()),
95  la0(la0_), la1(la1_),
96  const_residual(const_residual_),
97  time(0.0), dt_mode(MultiplyOperator0ByDT), stage(0),
98  pattern_engine(*this), prestage_engine(*this), residual_engine(*this), jacobian_engine(*this),
99  explicit_jacobian_residual_engine(*this)
100  { static_checks(); }
101 
105  void preStep(Real time_, Real dt_, int stages_)
106  {
107  time = time_;
108  dt = dt_;
109 
110  // This switch decides which term will be multiplied with dt
111  if(dt_mode == DivideOperator1ByDT)
112  {
113  dt_factor0 = 1.0;
114  dt_factor1 = 1.0 / dt;
115  }
116  else if(dt_mode == MultiplyOperator0ByDT)
117  {
118  dt_factor0 = dt;
119  dt_factor1 = 1.0;
120  }
121  else if(dt_mode == DoNotAssembleDT)
122  {
123  dt_factor0 = 1.0;
124  dt_factor1 = 1.0;
125  }
126  else
127  {
128  DUNE_THROW(Dune::Exception,"Unknown mode for assembling of time step size!");
129  }
130 
131  la0.preStep(time_,dt_, stages_);
132  la1.preStep(time_,dt_, stages_);
133  }
134 
136  void setMethod(const OneStepParameters & method_)
137  {
138  osp_method = & method_;
139  }
140 
142  void setStage(int stage_)
143  {
144  stage = stage_;
145  }
146 
148 
153  {
154  dt_mode = dt_mode_;
155  }
156 
158  Real timeAtStage(int stage_) const
159  {
160  return time+osp_method->d(stage_)*dt;
161  }
162 
164  Real timeAtStage() const
165  {
166  return time+osp_method->d(stage)*dt;
167  }
168 
169  void setWeight(const Real weight)
170  {
171  la0.setWeight(weight);
172  la1.setWeight(weight);
173  }
174 
177 
180  LocalPatternAssemblerEngine & localPatternAssemblerEngine
182  {
183  pattern_engine.setPattern(p);
184  return pattern_engine;
185  }
186 
189  LocalPreStageAssemblerEngine & localPreStageAssemblerEngine
190  (const std::vector<typename Traits::Solution*> & x)
191  {
192  prestage_engine.setSolutions(x);
193  prestage_engine.setConstResidual(const_residual);
194  return prestage_engine;
195  }
196 
199  LocalResidualAssemblerEngine & localResidualAssemblerEngine
200  (typename Traits::Residual & r, const typename Traits::Solution & x)
201  {
202  residual_engine.setSolution(x);
203  residual_engine.setConstResidual(const_residual);
204  residual_engine.setResidual(r);
205  return residual_engine;
206  }
207 
210  LocalJacobianAssemblerEngine & localJacobianAssemblerEngine
211  (typename Traits::Jacobian & a, const typename Traits::Solution & x)
212  {
213  jacobian_engine.setSolution(x);
214  jacobian_engine.setJacobian(a);
215  return jacobian_engine;
216  }
217 
220  LocalExplicitPatternAssemblerEngine & localExplicitPatternAssemblerEngine
221  (typename Traits::MatrixPattern & p)
222  {
223  return la1.localPatternAssemblerEngine(p);
224  }
225 
229  (typename Traits::Jacobian & a,
230  typename Traits::Residual & r0, typename Traits::Residual & r1,
231  const std::vector<typename Traits::Solution*> & x)
232  {
233  // Init pre stage engine
234  prestage_engine.setSolutions( x );
235  prestage_engine.setConstResiduals(r0,r1);
236  explicit_jacobian_residual_engine.setLocalPreStageEngine(prestage_engine);
237 
238  // Init jacobian engine
239  explicit_jacobian_residual_engine.setLocalJacobianEngine
240  (la1.localJacobianAssemblerEngine(a,*(x[stage])));
241 
242  return explicit_jacobian_residual_engine;
243  }
244 
246 
247  private:
248 
252  LA0 & la0;
253  LA1 & la1;
255 
258  const OneStepParameters * osp_method;
259 
261  typename Traits::Residual & const_residual;
262 
264  Real time;
265 
267  Real dt;
268 
283  Real dt_factor0, dt_factor1;
284 
287  DTAssemblingMode dt_mode;
288 
290  int stage;
291 
294  LocalPatternAssemblerEngine pattern_engine;
295  LocalPreStageAssemblerEngine prestage_engine;
296  LocalResidualAssemblerEngine residual_engine;
297  LocalJacobianAssemblerEngine jacobian_engine;
298  LocalExplicitJacobianResidualAssemblerEngine explicit_jacobian_residual_engine;
300  };
301 
302  }
303 }
304 #endif // DUNE_PDELAB_GRIDOPERATOR_ONESTEP_LOCALASSEMBLER_HH
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:229
LocalExplicitPatternAssemblerEngine & localExplicitPatternAssemblerEngine(typename Traits::MatrixPattern &p)
Definition: onestep/localassembler.hh:221
void setLocalPreStageEngine(PreStageEngine &prestage_engine_)
Definition: jacobianresidualengine.hh:53
void setResidual(Residual &residual_)
Definition: onestep/residualengine.hh:96
OneStepExplicitLocalJacobianResidualAssemblerEngine< OneStepLocalAssembler > LocalExplicitJacobianResidualAssemblerEngine
Definition: onestep/localassembler.hh:57
LocalPatternAssemblerEngine & localPatternAssemblerEngine(typename Traits::MatrixPattern &p)
Definition: onestep/localassembler.hh:181
GO::Traits::Range Residual
The type of the range (residual).
Definition: assemblerutilities.hh:54
Base class for local assembler.
Definition: assemblerutilities.hh:182
LocalJacobianAssemblerEngine & localJacobianAssemblerEngine(typename Traits::Jacobian &a, const typename Traits::Solution &x)
Definition: onestep/localassembler.hh:211
Definition: assemblerutilities.hh:22
GO::Traits::Jacobian Jacobian
The type of the jacobian.
Definition: assemblerutilities.hh:61
GO::Traits::RangeField RangeField
The field type of the range (residual).
Definition: assemblerutilities.hh:51
MatrixBackend::template Pattern< Jacobian, TestGridFunctionSpace, TrialGridFunctionSpace > MatrixPattern
The matrix pattern.
Definition: assemblerutilities.hh:68
Definition: onestep/localassembler.hh:147
Base parameter class for time stepping scheme parameters.
Definition: onestepparameter.hh:23
static const unsigned int value
Definition: gridfunctionspace/tags.hh:139
void setStage(int stage_)
Set the current stage of the one step scheme.
Definition: onestep/localassembler.hh:142
const GO::Traits::TestGridFunctionSpaceConstraints & testConstraints() const
get the constraints on the test grid function space
Definition: assemblerutilities.hh:205
Dune::PDELab::TimeSteppingParameterInterface< Real > OneStepParameters
The type of the one step parameter object.
Definition: onestep/localassembler.hh:90
void setConstResidual(Residual &const_residual_)
Definition: prestageengine.hh:108
void setSolution(const Solution &solution_)
Definition: onestep/residualengine.hh:81
LA1 LocalAssemblerDT1
Definition: onestep/localassembler.hh:38
void setWeight(const Real weight)
Definition: onestep/localassembler.hh:169
void setMethod(const OneStepParameters &method_)
Set the one step method parameters.
Definition: onestep/localassembler.hh:136
void setSolutions(const Solutions &solutions_)
Definition: prestageengine.hh:88
const GO::Traits::TrialGridFunctionSpaceConstraints & trialConstraints() const
get the constraints on the trial grid function space
Definition: assemblerutilities.hh:199
For backward compatibility – Do not use this!
Definition: adaptivity.hh:28
Definition: onestep/localassembler.hh:147
LocalPreStageAssemblerEngine & localPreStageAssemblerEngine(const std::vector< typename Traits::Solution *> &x)
Definition: onestep/localassembler.hh:190
LA1::LocalPatternAssemblerEngine LocalExplicitPatternAssemblerEngine
Definition: onestep/localassembler.hh:55
OneStepLocalPreStageAssemblerEngine< OneStepLocalAssembler > LocalPreStageAssemblerEngine
Definition: onestep/localassembler.hh:51
LA0 LocalAssemblerDT0
The types of the local assemblers of order one and zero.
Definition: onestep/localassembler.hh:37
void setDTAssemblingMode(DTAssemblingMode dt_mode_)
Definition: onestep/localassembler.hh:152
Dune::PDELab::LocalAssemblerTraits< GO > Traits
Definition: onestep/localassembler.hh:40
OneStepLocalJacobianAssemblerEngine< OneStepLocalAssembler > LocalJacobianAssemblerEngine
Definition: onestep/localassembler.hh:53
void setPattern(Pattern &pattern_)
Definition: onestep/patternengine.hh:62
GO::Traits::Domain Solution
The type of the domain (solution).
Definition: assemblerutilities.hh:47
void setSolution(const Solution &solution_)
Definition: onestep/jacobianengine.hh:71
void static_checks()
Definition: onestep/localassembler.hh:66
void setLocalJacobianEngine(JacobianEngine &jacobian_engine_)
Definition: jacobianresidualengine.hh:58
Definition: onestep/localassembler.hh:147
DTAssemblingMode
Definition: onestep/localassembler.hh:147
Real timeAtStage(int stage_) const
Access time at given stage.
Definition: onestep/localassembler.hh:158
void setConstResidual(const Residual &const_residual_)
Definition: onestep/residualengine.hh:88
void preStep(Real time_, Real dt_, int stages_)
Definition: onestep/localassembler.hh:105
OneStepLocalPatternAssemblerEngine< OneStepLocalAssembler > LocalPatternAssemblerEngine
Definition: onestep/localassembler.hh:50
void setJacobian(Jacobian &jacobian_)
Definition: onestep/jacobianengine.hh:78
OneStepLocalAssembler(LA0 &la0_, LA1 &la1_, typename Traits::Residual &const_residual_)
Constructor with empty constraints.
Definition: onestep/localassembler.hh:93
Dune::PDELab::LocalAssemblerBase< typename GO::Traits::MatrixBackend, typename GO::Traits::TrialGridFunctionSpaceConstraints, typename GO::Traits::TestGridFunctionSpaceConstraints > Base
The base class.
Definition: onestep/localassembler.hh:46
OneStepLocalResidualAssemblerEngine< OneStepLocalAssembler > LocalResidualAssemblerEngine
Definition: onestep/localassembler.hh:52
Real timeAtStage() const
Access time at given stage.
Definition: onestep/localassembler.hh:164
void setConstResiduals(Residual &const_residual_0_, Residual &const_residual_1_)
Definition: prestageengine.hh:95
The local assembler for one step methods.
Definition: onestep/localassembler.hh:28
const P & p
Definition: constraints.hh:147
LocalResidualAssemblerEngine & localResidualAssemblerEngine(typename Traits::Residual &r, const typename Traits::Solution &x)
Definition: onestep/localassembler.hh:200
Traits::RangeField Real
The local operators type for real numbers e.g. time.
Definition: onestep/localassembler.hh:87
virtual R d(int r) const =0
Return entries of the d Vector.