dune-pdelab  2.5-dev
gridoperator/onestep.hh
Go to the documentation of this file.
1 
2 #ifndef DUNE_PDELAB_GRIDOPERATOR_ONESTEP_HH
3 #define DUNE_PDELAB_GRIDOPERATOR_ONESTEP_HH
4 
5 #include <tuple>
6 
11 
12 namespace Dune{
13  namespace PDELab{
14 
15  template<typename GO0, typename GO1, bool implicit = true>
17  {
18  public:
19 
21  typedef typename GO0::Pattern Pattern;
22 
24  typedef typename GO0::Traits::Assembler Assembler;
25 
28  typedef typename GO0::Traits::LocalAssembler LocalAssemblerDT0;
29  typedef typename GO1::Traits::LocalAssembler LocalAssemblerDT1;
31 
34 
36  typedef typename GO0::BorderDOFExchanger BorderDOFExchanger;
37 
40  <typename GO0::Traits::TrialGridFunctionSpace,
41  typename GO0::Traits::TestGridFunctionSpace,
42  typename GO0::Traits::MatrixBackend,
43  typename GO0::Traits::DomainField,
44  typename GO0::Traits::RangeField,
45  typename GO0::Traits::JacobianField,
46  typename GO0::Traits::TrialGridFunctionSpaceConstraints,
47  typename GO0::Traits::TestGridFunctionSpaceConstraints,
48  Assembler,
49  LocalAssembler> Traits;
50 
53  typedef typename Traits::Domain Domain;
54  typedef typename Traits::Range Range;
55  typedef typename Traits::Jacobian Jacobian;
57 
58  template <typename MFT>
60  typedef Jacobian Type;
61  };
62 
64  typedef typename LocalAssembler::Real Real;
65 
68 
70  OneStepGridOperator(GO0 & go0_, GO1 & go1_)
71  : global_assembler(go0_.assembler()),
72  go0(go0_), go1(go1_),
73  la0(go0_.localAssembler()), la1(go1_.localAssembler()),
74  const_residual( go0_.testGridFunctionSpace() ),
75  local_assembler(la0,la1, const_residual)
76  {
77  GO0::setupGridOperators(std::tie(go0_,go1_));
78  if(!implicit)
79  local_assembler.setDTAssemblingMode(LocalAssembler::DoNotAssembleDT);
80  }
81 
86  {
87  if(!implicit)
88  DUNE_THROW(Dune::Exception,"This function should not be called in explicit mode");
89  local_assembler.setDTAssemblingMode(LocalAssembler::DivideOperator1ByDT);
90  }
92  {
93  if(!implicit)
94  DUNE_THROW(Dune::Exception,"This function should not be called in explicit mode");
95  local_assembler.setDTAssemblingMode(LocalAssembler::MultiplyOperator0ByDT);
96  }
97 
100  {
101  return global_assembler.trialGridFunctionSpace();
102  }
103 
106  {
107  return global_assembler.testGridFunctionSpace();
108  }
109 
111  typename Traits::TrialGridFunctionSpace::Traits::SizeType globalSizeU () const
112  {
113  return trialGridFunctionSpace().globalSize();
114  }
115 
117  typename Traits::TestGridFunctionSpace::Traits::SizeType globalSizeV () const
118  {
119  return testGridFunctionSpace().globalSize();
120  }
121 
122  Assembler & assembler() const { return global_assembler; }
123 
124  LocalAssembler & localAssembler() const { return local_assembler; }
125 
127  void fill_pattern(Pattern & p) const {
128  if(implicit){
129  typedef typename LocalAssembler::LocalPatternAssemblerEngine PatternEngine;
130  PatternEngine & pattern_engine = local_assembler.localPatternAssemblerEngine(p);
131  global_assembler.assemble(pattern_engine);
132  } else {
133  typedef typename LocalAssembler::LocalExplicitPatternAssemblerEngine PatternEngine;
134  PatternEngine & pattern_engine = local_assembler.localExplicitPatternAssemblerEngine(p);
135  global_assembler.assemble(pattern_engine);
136  }
137  }
138 
140  void preStage(unsigned int stage, const std::vector<Domain*> & x){
141  if(!implicit){DUNE_THROW(Dune::Exception,"This function should not be called in explicit mode");}
142 
143  typedef typename LocalAssembler::LocalPreStageAssemblerEngine PreStageEngine;
144  local_assembler.setStage(stage);
145  PreStageEngine & prestage_engine = local_assembler.localPreStageAssemblerEngine(x);
146  global_assembler.assemble(prestage_engine);
147  //Dune::printvector(std::cout,const_residual.base(),"const residual","row",4,9,1);
148  }
149 
151  void residual(const Domain & x, Range & r) const {
152  if(!implicit){DUNE_THROW(Dune::Exception,"This function should not be called in explicit mode");}
153 
154  typedef typename LocalAssembler::LocalResidualAssemblerEngine ResidualEngine;
155  ResidualEngine & residual_engine = local_assembler.localResidualAssemblerEngine(r,x);
156  global_assembler.assemble(residual_engine);
157  //Dune::printvector(std::cout,r.base(),"residual","row",4,9,1);
158  }
159 
161  void jacobian(const Domain & x, Jacobian & a) const {
162  if(!implicit){DUNE_THROW(Dune::Exception,"This function should not be called in explicit mode");}
163 
164  typedef typename LocalAssembler::LocalJacobianAssemblerEngine JacobianEngine;
165  JacobianEngine & jacobian_engine = local_assembler.localJacobianAssemblerEngine(a,x);
166  global_assembler.assemble(jacobian_engine);
167  //printmatrix(std::cout,a.base(),"global stiffness matrix","row",9,1);
168  }
169 
171  void explicit_jacobian_residual(unsigned int stage, const std::vector<Domain*> & x,
172  Jacobian & a, Range & r1, Range & r0)
173  {
174  if(implicit){DUNE_THROW(Dune::Exception,"This function should not be called in implicit mode");}
175 
176  local_assembler.setStage(stage);
177 
179  ExplicitJacobianResidualEngine;
180 
181  ExplicitJacobianResidualEngine & jacobian_residual_engine
182  = local_assembler.localExplicitJacobianResidualAssemblerEngine(a,r0,r1,x);
183 
184  global_assembler.assemble(jacobian_residual_engine);
185  }
186 
188  template<typename F, typename X>
189  void interpolate (unsigned stage, const X& xold, F& f, X& x) const
190  {
191  // Set time in boundary value function
192  f.setTime(local_assembler.timeAtStage(stage));
193 
194  go0.localAssembler().setTime(local_assembler.timeAtStage(stage));
195 
196  // Interpolate
197  go0.interpolate(xold,f,x);
198 
199  // Copy non-constrained dofs from old time step
200  Dune::PDELab::copy_nonconstrained_dofs(local_assembler.trialConstraints(),xold,x);
201  }
202 
205  {
206  local_assembler.setMethod(method_);
207  }
208 
210  void preStep (const TimeSteppingParameterInterface<Real>& method_, Real time_, Real dt_)
211  {
212  local_assembler.setMethod(method_);
213  local_assembler.preStep(time_,dt_,method_.s());
214  }
215 
217  void postStep ()
218  {
219  la0.postStep();
220  la1.postStep();
221  }
222 
224  void postStage ()
225  {
226  la0.postStage();
227  la1.postStage();
228  }
229 
231  Real suggestTimestep (Real dt) const
232  {
233  Real suggested_dt = std::min(la0.suggestTimestep(dt),la1.suggestTimestep(dt));
234  if (trialGridFunctionSpace().gridView().comm().size()>1)
235  suggested_dt = trialGridFunctionSpace().gridView().comm().min(suggested_dt);
236  return suggested_dt;
237  }
238 
239  void update()
240  {
241  go0.update();
242  go1.update();
243  const_residual = Range(go0.testGridFunctionSpace());
244  }
245 
246  const typename Traits::MatrixBackend& matrixBackend() const
247  {
248  return go0.matrixBackend();
249  }
250 
252  {
253  return local_assembler.trialConstraints();
254  }
255 
256  private:
257  Assembler & global_assembler;
258  GO0 & go0;
259  GO1 & go1;
260  LocalAssemblerDT0 & la0;
261  LocalAssemblerDT1 & la1;
262  Range const_residual;
263  mutable LocalAssembler local_assembler;
264  };
265 
266  }
267 }
268 #endif // DUNE_PDELAB_GRIDOPERATOR_ONESTEP_HH
Traits class for the grid operator.
Definition: gridoperatorutilities.hh:33
Definition: gridoperator/onestep.hh:16
OneStepLocalAssembler< OneStepGridOperator, LocalAssemblerDT0, LocalAssemblerDT1 > LocalAssembler
The local assembler type.
Definition: gridoperator/onestep.hh:33
Traits::TrialGridFunctionSpace::Traits::SizeType globalSizeU() const
Get dimension of space u.
Definition: gridoperator/onestep.hh:111
void jacobian(const Domain &x, Jacobian &a) const
Assemble jacobian.
Definition: gridoperator/onestep.hh:161
void setMethod(const TimeSteppingParameterInterface< Real > &method_)
set time stepping method
Definition: gridoperator/onestep.hh:204
void interpolate(unsigned stage, const X &xold, F &f, X &x) const
Interpolate constrained values from given function f.
Definition: gridoperator/onestep.hh:189
GO0::Pattern Pattern
The sparsity pattern container for the jacobian matrix.
Definition: gridoperator/onestep.hh:21
OneStepGridOperator(GO0 &go0_, GO1 &go1_)
Constructor for non trivial constraints.
Definition: gridoperator/onestep.hh:70
Jacobian Type
Definition: gridoperator/onestep.hh:60
void multiplySpatialTermByDeltaT()
Definition: gridoperator/onestep.hh:91
Base parameter class for time stepping scheme parameters.
Definition: onestepparameter.hh:23
void divideMassTermByDeltaT()
Definition: gridoperator/onestep.hh:85
const LocalAssembler::Traits::TrialGridFunctionSpaceConstraints trialConstraints() const
Definition: gridoperator/onestep.hh:251
void copy_nonconstrained_dofs(const CG &cg, const XG &xgin, XG &xgout)
Definition: constraints.hh:989
Dune::PDELab::Backend::Matrix< MB, Domain, Range, JF > Jacobian
The type of the jacobian.
Definition: gridoperatorutilities.hh:72
GO1::Traits::LocalAssembler LocalAssemblerDT1
Definition: gridoperator/onestep.hh:29
Dune::PDELab::GridOperatorTraits< typename GO0::Traits::TrialGridFunctionSpace, typename GO0::Traits::TestGridFunctionSpace, typename GO0::Traits::MatrixBackend, typename GO0::Traits::DomainField, typename GO0::Traits::RangeField, typename GO0::Traits::JacobianField, typename GO0::Traits::TrialGridFunctionSpaceConstraints, typename GO0::Traits::TestGridFunctionSpaceConstraints, Assembler, LocalAssembler > Traits
The grid operator traits.
Definition: gridoperator/onestep.hh:49
Traits::RangeField Real
The local operators type for real numbers e.g. time.
Definition: onestep/localassembler.hh:86
LocalAssembler::Real Real
The type for real number e.g. time.
Definition: gridoperator/onestep.hh:64
const Traits::TrialGridFunctionSpace & trialGridFunctionSpace() const
Get the trial grid function space.
Definition: gridoperator/onestep.hh:99
GO0::Traits::Assembler Assembler
The global UDG assembler type.
Definition: gridoperator/onestep.hh:24
Traits::TestGridFunctionSpace::Traits::SizeType globalSizeV() const
Get dimension of space v.
Definition: gridoperator/onestep.hh:117
Traits::Jacobian Jacobian
Definition: gridoperator/onestep.hh:55
Dune::PDELab::Backend::Vector< GFSU, DF > Domain
The type of the domain (solution).
Definition: gridoperatorutilities.hh:58
void explicit_jacobian_residual(unsigned int stage, const std::vector< Domain *> &x, Jacobian &a, Range &r1, Range &r0)
Assemble jacobian and residual simultaneously for explicit treatment.
Definition: gridoperator/onestep.hh:171
GO::Traits::TrialGridFunctionSpaceConstraints TrialGridFunctionSpaceConstraints
The type of the trial grid function space constraints.
Definition: assemblerutilities.hh:33
const P & p
Definition: constraints.hh:147
const Traits::TestGridFunctionSpace & testGridFunctionSpace() const
Get the test grid function space.
Definition: gridoperator/onestep.hh:105
GO0::BorderDOFExchanger BorderDOFExchanger
The BorderDOFExchanger.
Definition: gridoperator/onestep.hh:36
void fill_pattern(Pattern &p) const
Fill pattern of jacobian matrix.
Definition: gridoperator/onestep.hh:127
Real suggestTimestep(Real dt) const
to be called once before each stage
Definition: gridoperator/onestep.hh:231
void preStage(unsigned int stage, const std::vector< Domain *> &x)
Assemble constant part of residual.
Definition: gridoperator/onestep.hh:140
LocalAssembler & localAssembler() const
Definition: gridoperator/onestep.hh:124
Traits::Domain Domain
Definition: gridoperator/onestep.hh:53
MB MatrixBackend
The matrix backend of the grid operator.
Definition: gridoperatorutilities.hh:51
For backward compatibility – Do not use this!
Definition: adaptivity.hh:27
LocalAssemblerDT1 ::LocalPatternAssemblerEngine LocalExplicitPatternAssemblerEngine
Definition: onestep/localassembler.hh:55
Assembler & assembler() const
Definition: gridoperator/onestep.hh:122
GO0::Traits::LocalAssembler LocalAssemblerDT0
Definition: gridoperator/onestep.hh:28
const Traits::MatrixBackend & matrixBackend() const
Definition: gridoperator/onestep.hh:246
void postStage()
to be called after stage is completed
Definition: gridoperator/onestep.hh:224
LocalAssembler::OneStepParameters OneStepParameters
The type of the one step method parameters.
Definition: gridoperator/onestep.hh:67
GFSU TrialGridFunctionSpace
The trial grid function space.
Definition: gridoperatorutilities.hh:37
void residual(const Domain &x, Range &r) const
Assemble residual.
Definition: gridoperator/onestep.hh:151
virtual unsigned s() const =0
Return number of stages of the method.
void update()
Definition: gridoperator/onestep.hh:239
Definition: gridoperator/onestep.hh:59
GFSV TestGridFunctionSpace
The test grid function space.
Definition: gridoperatorutilities.hh:40
Traits::Range Range
Definition: gridoperator/onestep.hh:54
Dune::PDELab::Backend::Vector< GFSV, RF > Range
The type of the range (residual).
Definition: gridoperatorutilities.hh:65
void preStep(const TimeSteppingParameterInterface< Real > &method_, Real time_, Real dt_)
parametrize assembler with a time-stepping method
Definition: gridoperator/onestep.hh:210
void postStep()
to be called after step is completed
Definition: gridoperator/onestep.hh:217