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