2 #ifndef DUNE_PDELAB_ONESTEP_OPERATOR_HH 3 #define DUNE_PDELAB_ONESTEP_OPERATOR_HH 13 template<
typename GO0,
typename GO1,
bool implicit = true>
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,
56 template <
typename MFT>
73 local_assembler(la0,la1, const_residual)
75 GO0::setupGridOperators(Dune::tie(go0_,go1_));
86 DUNE_THROW(Dune::Exception,
"This function should not be called in explicit mode");
92 DUNE_THROW(Dune::Exception,
"This function should not be called in explicit mode");
99 return global_assembler.trialGridFunctionSpace();
105 return global_assembler.testGridFunctionSpace();
109 typename Traits::TrialGridFunctionSpace::Traits::SizeType
globalSizeU ()
const 115 typename Traits::TestGridFunctionSpace::Traits::SizeType
globalSizeV ()
const 120 Assembler &
assembler()
const {
return global_assembler; }
128 PatternEngine & pattern_engine = local_assembler.localPatternAssemblerEngine(p);
129 global_assembler.assemble(pattern_engine);
132 PatternEngine & pattern_engine = local_assembler.localExplicitPatternAssemblerEngine(p);
133 global_assembler.assemble(pattern_engine);
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");}
142 local_assembler.setStage(stage);
143 PreStageEngine & prestage_engine = local_assembler.localPreStageAssemblerEngine(x);
144 global_assembler.assemble(prestage_engine);
150 if(!implicit){DUNE_THROW(Dune::Exception,
"This function should not be called in explicit mode");}
153 ResidualEngine & residual_engine = local_assembler.localResidualAssemblerEngine(r,x);
154 global_assembler.assemble(residual_engine);
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");}
163 JacobianEngine & jacobian_engine = local_assembler.localJacobianAssemblerEngine(a,x);
164 global_assembler.assemble(jacobian_engine);
170 Jacobian & a, Range & r1, Range & r0)
172 if(implicit){DUNE_THROW(Dune::Exception,
"This function should not be called in implicit mode");}
174 local_assembler.setStage(stage);
177 ExplicitJacobianResidualEngine;
179 ExplicitJacobianResidualEngine & jacobian_residual_engine
180 = local_assembler.localExplicitJacobianResidualAssemblerEngine(a,r0,r1,x);
182 global_assembler.assemble(jacobian_residual_engine);
186 template<
typename F,
typename X>
187 void interpolate (
unsigned stage,
const X& xold, F& f, X& x)
const 190 f.setTime(local_assembler.timeAtStage(stage));
192 go0.localAssembler().setTime(local_assembler.timeAtStage(stage));
195 go0.interpolate(xold,f,x);
204 local_assembler.setMethod(method_);
210 local_assembler.setMethod(method_);
211 local_assembler.preStep(time_,dt_,method_.
s());
231 Real suggested_dt = std::min(la0.suggestTimestep(dt),la1.suggestTimestep(dt));
241 const_residual =
Range(go0.testGridFunctionSpace());
246 return go0.matrixBackend();
250 Assembler & global_assembler;
253 LocalAssemblerDT0 & la0;
254 LocalAssemblerDT1 & la1;
255 Range const_residual;
256 mutable LocalAssembler local_assembler;
void update()
Definition: gridoperator/onestep.hh:237
GFSV TestGridFunctionSpace
The test grid function space.
Definition: gridoperatorutilities.hh:40
Traits class for the grid operator.
Definition: gridoperatorutilities.hh:33
Traits::TestGridFunctionSpace::Traits::SizeType globalSizeV() const
Get dimension of space v.
Definition: gridoperator/onestep.hh:115
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
Traits::Jacobian Jacobian
Definition: gridoperator/onestep.hh:53
OneStepLocalAssembler< OneStepGridOperator, LocalAssemblerDT0, LocalAssemblerDT1 > LocalAssembler
The local assembler type.
Definition: gridoperator/onestep.hh:31
OneStepGridOperator(GO0 &go0_, GO1 &go1_)
Constructor for non trivial constraints.
Definition: gridoperator/onestep.hh:68
GO0::BorderDOFExchanger BorderDOFExchanger
The BorderDOFExchanger.
Definition: gridoperator/onestep.hh:34
Assembler & assembler() const
Definition: gridoperator/onestep.hh:120
MB MatrixBackend
The matrix backend of the grid operator.
Definition: gridoperatorutilities.hh:51
const Traits::TestGridFunctionSpace & testGridFunctionSpace() const
Get the test grid function space.
Definition: gridoperator/onestep.hh:103
void multiplySpatialTermByDeltaT()
Definition: gridoperator/onestep.hh:89
GO0::Pattern Pattern
The sparsity pattern container for the jacobian matrix.
Definition: gridoperator/onestep.hh:19
Definition: adaptivity.hh:27
Definition: onestep/localassembler.hh:139
virtual unsigned s() const =0
Return number of stages of the method.
GO0::Traits::LocalAssembler LocalAssemblerDT0
Definition: gridoperator/onestep.hh:26
LocalAssemblerDT1::LocalPatternAssemblerEngine LocalExplicitPatternAssemblerEngine
Definition: onestep/localassembler.hh:55
Traits::RangeField Real
The local operators type for real numbers e.g. time.
Definition: onestep/localassembler.hh:86
void residual(const Domain &x, Range &r) const
Assemble residual.
Definition: gridoperator/onestep.hh:149
Jacobian Type
Definition: gridoperator/onestep.hh:58
Dune::PDELab::Backend::Vector< GFSU, DF > Domain
The type of the domain (solution).
Definition: gridoperatorutilities.hh:58
void copy_nonconstrained_dofs(const CG &cg, const XG &xgin, XG &xgout)
Definition: constraints.hh:989
void setMethod(const TimeSteppingParameterInterface< Real > &method_)
set time stepping method
Definition: gridoperator/onestep.hh:202
Base parameter class for time stepping scheme parameters.
Definition: onestepparameter.hh:23
LocalAssembler::OneStepParameters OneStepParameters
The type of the one step method parameters.
Definition: gridoperator/onestep.hh:65
GFSU TrialGridFunctionSpace
The trial grid function space.
Definition: gridoperatorutilities.hh:37
Traits::Range Range
Definition: gridoperator/onestep.hh:52
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
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
const Traits::TrialGridFunctionSpace & trialGridFunctionSpace() const
Get the trial grid function space.
Definition: gridoperator/onestep.hh:97
void preStep(const TimeSteppingParameterInterface< Real > &method_, Real time_, Real dt_)
parametrize assembler with a time-stepping method
Definition: gridoperator/onestep.hh:208
Definition: onestep/localassembler.hh:139
void divideMassTermByDeltaT()
Definition: gridoperator/onestep.hh:83
LocalAssembler::Real Real
The type for real number e.g. time.
Definition: gridoperator/onestep.hh:62
Definition: onestep/localassembler.hh:139
const P & p
Definition: constraints.hh:147
void jacobian(const Domain &x, Jacobian &a) const
Assemble jacobian.
Definition: gridoperator/onestep.hh:159
GO1::Traits::LocalAssembler LocalAssemblerDT1
Definition: gridoperator/onestep.hh:27
const Traits::MatrixBackend & matrixBackend() const
Definition: gridoperator/onestep.hh:244
void fill_pattern(Pattern &p) const
Fill pattern of jacobian matrix.
Definition: gridoperator/onestep.hh:125
Definition: gridoperator/onestep.hh:14
void postStage()
to be called after stage is completed
Definition: gridoperator/onestep.hh:222
Traits::TrialGridFunctionSpace::Traits::SizeType globalSizeU() const
Get dimension of space u.
Definition: gridoperator/onestep.hh:109
void postStep()
to be called after step is completed
Definition: gridoperator/onestep.hh:215
LocalAssembler & localAssembler() const
Definition: gridoperator/onestep.hh:122
Dune::PDELab::Backend::Vector< GFSV, RF > Range
The type of the range (residual).
Definition: gridoperatorutilities.hh:65
Definition: gridoperator/onestep.hh:57
Real suggestTimestep(Real dt) const
to be called once before each stage
Definition: gridoperator/onestep.hh:229
void preStage(unsigned int stage, const std::vector< Domain * > &x)
Assemble constant part of residual.
Definition: gridoperator/onestep.hh:138
GO0::Traits::Assembler Assembler
The global UDG assembler type.
Definition: gridoperator/onestep.hh:22
void interpolate(unsigned stage, const X &xold, F &f, X &x) const
Interpolate constrained values from given function f.
Definition: gridoperator/onestep.hh:187