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,
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)
82  }
83 
88  {
89  if(not implicit)
90  DUNE_THROW(Dune::Exception,"This function should not be called in explicit mode");
92  }
94  {
95  if(not implicit)
96  DUNE_THROW(Dune::Exception,"This function should not be called in explicit mode");
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  // using Dune::PDELab::Backend::native;
156  // Dune::printvector(std::cout,native(const_residual),"const residual","row",4,9,1);
157  }
158 
160  void residual(const Domain & x, Range & r) const
161  {
162  if(not implicit)
163  DUNE_THROW(Dune::Exception,"This function should not be called in explicit mode");
164 
165  typedef typename LocalAssembler::LocalResidualAssemblerEngine ResidualEngine;
166  ResidualEngine & residual_engine = local_assembler.localResidualAssemblerEngine(r,x);
167  global_assembler.assemble(residual_engine);
168  // using Dune::PDELab::Backend::native;
169  // Dune::printvector(std::cout,native(r),"residual","row",4,9,1);
170  }
171 
173  void jacobian(const Domain & x, Jacobian & a) const
174  {
175  if(not implicit)
176  DUNE_THROW(Dune::Exception,"This function should not be called in explicit mode");
177 
178  typedef typename LocalAssembler::LocalJacobianAssemblerEngine JacobianEngine;
179  JacobianEngine & jacobian_engine = local_assembler.localJacobianAssemblerEngine(a,x);
180  global_assembler.assemble(jacobian_engine);
181  // using Dune::PDELab::Backend::native;
182  // printmatrix(std::cout,native(a),"global stiffness matrix","row",9,1);
183  }
184 
186  void explicit_jacobian_residual(unsigned int stage, const std::vector<Domain*> & x,
187  Jacobian & a, Range & r1, Range & r0)
188  {
189  if(implicit)
190  DUNE_THROW(Dune::Exception,"This function should not be called in implicit mode");
191 
192  local_assembler.setStage(stage);
193 
195  ExplicitJacobianResidualEngine;
196 
197  ExplicitJacobianResidualEngine & jacobian_residual_engine
198  = local_assembler.localExplicitJacobianResidualAssemblerEngine(a,r0,r1,x);
199 
200  global_assembler.assemble(jacobian_residual_engine);
201  }
202 
204  template<typename F, typename X>
205  void interpolate (unsigned stage, const X& xold, F& f, X& x) const
206  {
207  // Set time in boundary value function
208  f.setTime(local_assembler.timeAtStage(stage));
209 
210  go0.localAssembler().setTime(local_assembler.timeAtStage(stage));
211 
212  // Interpolate
213  go0.interpolate(xold,f,x);
214 
215  // Copy non-constrained dofs from old time step
217  }
218 
221  {
222  local_assembler.setMethod(method_);
223  }
224 
226  void preStep (const TimeSteppingParameterInterface<Real>& method_, Real time_, Real dt_)
227  {
228  local_assembler.setMethod(method_);
229  local_assembler.preStep(time_,dt_,method_.s());
230  }
231 
233  void postStep ()
234  {
235  la0.postStep();
236  la1.postStep();
237  }
238 
240  void postStage ()
241  {
242  la0.postStage();
243  la1.postStage();
244  }
245 
248  {
249  Real suggested_dt = std::min(la0.suggestTimestep(dt),la1.suggestTimestep(dt));
250  if (trialGridFunctionSpace().gridView().comm().size()>1)
251  suggested_dt = trialGridFunctionSpace().gridView().comm().min(suggested_dt);
252  return suggested_dt;
253  }
254 
255  void update()
256  {
257  go0.update();
258  go1.update();
259  const_residual = Range(go0.testGridFunctionSpace());
260  }
261 
262  void make_consistent(Jacobian& a) const
263  {
264  go0.make_consistent(a);
265  }
266 
267  const typename Traits::MatrixBackend& matrixBackend() const
268  {
269  return go0.matrixBackend();
270  }
271 
273  {
274  return local_assembler.trialConstraints();
275  }
276 
277  private:
278  Assembler & global_assembler;
279  GO0 & go0;
280  GO1 & go1;
281  LocalAssemblerDT0 & la0;
282  LocalAssemblerDT1 & la1;
283  Range const_residual;
284  mutable LocalAssembler local_assembler;
285  };
286 
287  }
288 }
289 #endif // DUNE_PDELAB_GRIDOPERATOR_ONESTEP_HH
Dune::PDELab::OneStepGridOperator::divideMassTermByDeltaT
void divideMassTermByDeltaT()
Definition: gridoperator/onestep.hh:87
Dune::PDELab::OneStepLocalAssembler::localPreStageAssemblerEngine
LocalPreStageAssemblerEngine & localPreStageAssemblerEngine(const std::vector< typename Traits::Solution * > &x)
Definition: onestep/localassembler.hh:190
Dune::PDELab::TimeSteppingParameterInterface::s
virtual unsigned s() const =0
Return number of stages of the method.
onestep.hh
Dune::PDELab::GridOperatorTraits::TrialGridFunctionSpace
GFSU TrialGridFunctionSpace
The trial grid function space.
Definition: gridoperatorutilities.hh:37
Dune::PDELab::OneStepGridOperator::assembler
Assembler & assembler() const
Definition: gridoperator/onestep.hh:124
Dune::PDELab::copy_nonconstrained_dofs
void copy_nonconstrained_dofs(const CG &cg, const XG &xgin, XG &xgout)
Definition: constraints.hh:989
Dune::PDELab::OneStepLocalAssembler< OneStepGridOperator, LocalAssemblerDT0, LocalAssemblerDT1 >::DoNotAssembleDT
@ DoNotAssembleDT
Definition: onestep/localassembler.hh:147
Dune::PDELab::OneStepGridOperator::matrixBackend
const Traits::MatrixBackend & matrixBackend() const
Definition: gridoperator/onestep.hh:267
Dune::PDELab::OneStepGridOperator::BorderDOFExchanger
GO0::BorderDOFExchanger BorderDOFExchanger
The BorderDOFExchanger.
Definition: gridoperator/onestep.hh:37
Dune::PDELab::GridOperatorTraits::TestGridFunctionSpace
GFSV TestGridFunctionSpace
The test grid function space.
Definition: gridoperatorutilities.hh:40
Dune::PDELab::OneStepLocalAssembler< OneStepGridOperator, LocalAssemblerDT0, LocalAssemblerDT1 >::LocalExplicitPatternAssemblerEngine
LocalAssemblerDT1 ::LocalPatternAssemblerEngine LocalExplicitPatternAssemblerEngine
Definition: onestep/localassembler.hh:55
Dune::PDELab::GridOperatorTraits::Range
Dune::PDELab::Backend::Vector< GFSV, RF > Range
The type of the range (residual).
Definition: gridoperatorutilities.hh:65
Dune::PDELab::OneStepLocalPatternAssemblerEngine< OneStepLocalAssembler >
Dune::PDELab::OneStepGridOperator::localAssembler
LocalAssembler & localAssembler() const
Definition: gridoperator/onestep.hh:126
Dune::PDELab::OneStepGridOperator::multiplySpatialTermByDeltaT
void multiplySpatialTermByDeltaT()
Definition: gridoperator/onestep.hh:93
Dune::PDELab::OneStepLocalAssembler< OneStepGridOperator, LocalAssemblerDT0, LocalAssemblerDT1 >::Real
Traits::RangeField Real
The local operators type for real numbers e.g. time.
Definition: onestep/localassembler.hh:87
Dune
For backward compatibility – Do not use this!
Definition: adaptivity.hh:28
Dune::PDELab::TimeSteppingParameterInterface
Base parameter class for time stepping scheme parameters.
Definition: onestepparameter.hh:23
Dune::PDELab::GridOperatorTraits
Traits class for the grid operator.
Definition: gridoperatorutilities.hh:33
Dune::PDELab::OneStepLocalResidualAssemblerEngine< OneStepLocalAssembler >
Dune::PDELab::OneStepGridOperator::interpolate
void interpolate(unsigned stage, const X &xold, F &f, X &x) const
Interpolate constrained values from given function f.
Definition: gridoperator/onestep.hh:205
Dune::PDELab::OneStepGridOperator::Jacobian
Traits::Jacobian Jacobian
Definition: gridoperator/onestep.hh:56
Dune::PDELab::OneStepLocalAssembler::setStage
void setStage(int stage_)
Set the current stage of the one step scheme.
Definition: onestep/localassembler.hh:142
Dune::PDELab::OneStepGridOperator::Traits
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
gridoperatorutilities.hh
Dune::PDELab::OneStepGridOperator::preStage
void preStage(unsigned int stage, const std::vector< Domain * > &x)
Assemble constant part of residual.
Definition: gridoperator/onestep.hh:146
Dune::PDELab::LocalAssemblerTraits::TrialGridFunctionSpaceConstraints
GO::Traits::TrialGridFunctionSpaceConstraints TrialGridFunctionSpaceConstraints
The type of the trial grid function space constraints.
Definition: assemblerutilities.hh:33
Dune::PDELab::OneStepLocalAssembler::localPatternAssemblerEngine
LocalPatternAssemblerEngine & localPatternAssemblerEngine(typename Traits::MatrixPattern &p)
Definition: onestep/localassembler.hh:181
Dune::PDELab::OneStepGridOperator::LocalAssemblerDT1
GO1::Traits::LocalAssembler LocalAssemblerDT1
Definition: gridoperator/onestep.hh:30
Dune::PDELab::OneStepGridOperator::Range
Traits::Range Range
Definition: gridoperator/onestep.hh:55
Dune::PDELab::OneStepLocalJacobianAssemblerEngine< OneStepLocalAssembler >
Dune::PDELab::OneStepGridOperator::trialConstraints
const LocalAssembler::Traits::TrialGridFunctionSpaceConstraints trialConstraints() const
Definition: gridoperator/onestep.hh:272
Dune::PDELab::GridOperatorTraits::MatrixBackend
MB MatrixBackend
The matrix backend of the grid operator.
Definition: gridoperatorutilities.hh:51
Dune::PDELab::OneStepGridOperator::postStage
void postStage()
to be called after stage is completed
Definition: gridoperator/onestep.hh:240
Dune::PDELab::OneStepLocalAssembler< OneStepGridOperator, LocalAssemblerDT0, LocalAssemblerDT1 >::DivideOperator1ByDT
@ DivideOperator1ByDT
Definition: onestep/localassembler.hh:147
Dune::PDELab::OneStepGridOperator::update
void update()
Definition: gridoperator/onestep.hh:255
Dune::PDELab::OneStepGridOperator::postStep
void postStep()
to be called after step is completed
Definition: gridoperator/onestep.hh:233
Dune::PDELab::OneStepLocalAssembler::localExplicitJacobianResidualAssemblerEngine
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
Dune::PDELab::OneStepGridOperator::jacobian
void jacobian(const Domain &x, Jacobian &a) const
Assemble jacobian.
Definition: gridoperator/onestep.hh:173
Dune::PDELab::OneStepGridOperator::make_consistent
void make_consistent(Jacobian &a) const
Definition: gridoperator/onestep.hh:262
Dune::PDELab::OneStepGridOperator::setMethod
void setMethod(const TimeSteppingParameterInterface< Real > &method_)
set time stepping method
Definition: gridoperator/onestep.hh:220
Dune::PDELab::OneStepGridOperator::OneStepGridOperator
OneStepGridOperator(GO0 &go0_, GO1 &go1_)
Constructor for non trivial constraints.
Definition: gridoperator/onestep.hh:72
Dune::PDELab::OneStepGridOperator::LocalAssemblerDT0
GO0::Traits::LocalAssembler LocalAssemblerDT0
Definition: gridoperator/onestep.hh:29
Dune::PDELab::OneStepGridOperator::explicit_jacobian_residual
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:186
Dune::PDELab::OneStepGridOperator::Assembler
GO0::Traits::Assembler Assembler
The global UDG assembler type.
Definition: gridoperator/onestep.hh:25
Dune::PDELab::OneStepGridOperator::MatrixContainer
Definition: gridoperator/onestep.hh:60
Dune::PDELab::OneStepGridOperator::suggestTimestep
Real suggestTimestep(Real dt) const
to be called once before each stage
Definition: gridoperator/onestep.hh:247
Dune::PDELab::OneStepLocalPreStageAssemblerEngine< OneStepLocalAssembler >
Dune::PDELab::GridOperatorTraits::Domain
Dune::PDELab::Backend::Vector< GFSU, DF > Domain
The type of the domain (solution).
Definition: gridoperatorutilities.hh:58
Dune::PDELab::OneStepLocalAssembler::timeAtStage
Real timeAtStage(int stage_) const
Access time at given stage.
Definition: onestep/localassembler.hh:158
Dune::PDELab::OneStepGridOperator::testGridFunctionSpace
const Traits::TestGridFunctionSpace & testGridFunctionSpace() const
Get the test grid function space.
Definition: gridoperator/onestep.hh:107
Dune::PDELab::OneStepLocalAssembler::localExplicitPatternAssemblerEngine
LocalExplicitPatternAssemblerEngine & localExplicitPatternAssemblerEngine(typename Traits::MatrixPattern &p)
Definition: onestep/localassembler.hh:221
Dune::PDELab::OneStepGridOperator::residual
void residual(const Domain &x, Range &r) const
Assemble residual.
Definition: gridoperator/onestep.hh:160
Dune::PDELab::GridOperatorTraits::Jacobian
Dune::PDELab::Backend::Matrix< MB, Domain, Range, JF > Jacobian
The type of the jacobian.
Definition: gridoperatorutilities.hh:72
constraints.hh
Dune::PDELab::OneStepGridOperator::Real
LocalAssembler::Real Real
The type for real number e.g. time.
Definition: gridoperator/onestep.hh:66
Dune::PDELab::OneStepLocalAssembler< OneStepGridOperator, LocalAssemblerDT0, LocalAssemblerDT1 >
Dune::PDELab::OneStepGridOperator::globalSizeU
Traits::TrialGridFunctionSpace::Traits::SizeType globalSizeU() const
Get dimension of space u.
Definition: gridoperator/onestep.hh:113
Dune::PDELab::OneStepLocalAssembler::setMethod
void setMethod(const OneStepParameters &method_)
Set the one step method parameters.
Definition: onestep/localassembler.hh:136
Dune::PDELab::OneStepGridOperator
Definition: gridoperator/onestep.hh:17
Dune::PDELab::OneStepGridOperator::Pattern
GO0::Pattern Pattern
The sparsity pattern container for the jacobian matrix.
Definition: gridoperator/onestep.hh:22
Dune::PDELab::OneStepGridOperator::OneStepParameters
LocalAssembler::OneStepParameters OneStepParameters
The type of the one step method parameters.
Definition: gridoperator/onestep.hh:69
Dune::PDELab::OneStepGridOperator::globalSizeV
Traits::TestGridFunctionSpace::Traits::SizeType globalSizeV() const
Get dimension of space v.
Definition: gridoperator/onestep.hh:119
Dune::PDELab::OneStepGridOperator::trialGridFunctionSpace
const Traits::TrialGridFunctionSpace & trialGridFunctionSpace() const
Get the trial grid function space.
Definition: gridoperator/onestep.hh:101
Dune::PDELab::OneStepExplicitLocalJacobianResidualAssemblerEngine< OneStepLocalAssembler >
Dune::PDELab::LocalAssemblerBase::trialConstraints
const CU & trialConstraints() const
get the constraints on the trial grid function space
Definition: assemblerutilities.hh:199
Dune::PDELab::OneStepLocalAssembler::localJacobianAssemblerEngine
LocalJacobianAssemblerEngine & localJacobianAssemblerEngine(typename Traits::Jacobian &a, const typename Traits::Solution &x)
Definition: onestep/localassembler.hh:211
localassembler.hh
Dune::PDELab::OneStepGridOperator::LocalAssembler
OneStepLocalAssembler< OneStepGridOperator, LocalAssemblerDT0, LocalAssemblerDT1 > LocalAssembler
The local assembler type.
Definition: gridoperator/onestep.hh:34
Dune::PDELab::OneStepLocalAssembler< OneStepGridOperator, LocalAssemblerDT0, LocalAssemblerDT1 >::MultiplyOperator0ByDT
@ MultiplyOperator0ByDT
Definition: onestep/localassembler.hh:147
Dune::PDELab::OneStepGridOperator::MatrixContainer::Type
Jacobian Type
Definition: gridoperator/onestep.hh:62
Dune::PDELab::OneStepLocalAssembler::preStep
void preStep(Real time_, Real dt_, int stages_)
Definition: onestep/localassembler.hh:105
Dune::PDELab::OneStepGridOperator::Domain
Traits::Domain Domain
Definition: gridoperator/onestep.hh:54
Dune::PDELab::OneStepGridOperator::fill_pattern
void fill_pattern(Pattern &p) const
Fill pattern of jacobian matrix.
Definition: gridoperator/onestep.hh:129
Dune::PDELab::OneStepLocalAssembler::localResidualAssemblerEngine
LocalResidualAssemblerEngine & localResidualAssemblerEngine(typename Traits::Residual &r, const typename Traits::Solution &x)
Definition: onestep/localassembler.hh:200
Dune::PDELab::OneStepGridOperator::preStep
void preStep(const TimeSteppingParameterInterface< Real > &method_, Real time_, Real dt_)
parametrize assembler with a time-stepping method
Definition: gridoperator/onestep.hh:226
Dune::PDELab::OneStepLocalAssembler::setDTAssemblingMode
void setDTAssemblingMode(DTAssemblingMode dt_mode_)
Definition: onestep/localassembler.hh:152
p
const P & p
Definition: constraints.hh:147