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