dune-pdelab  2.5-dev
default/localassembler.hh
Go to the documentation of this file.
1 #ifndef DUNE_PDELAB_GRIDOPERATOR_DEFAULT_LOCALASSEMBLER_HH
2 #define DUNE_PDELAB_GRIDOPERATOR_DEFAULT_LOCALASSEMBLER_HH
3 
4 #include <dune/typetree/typetree.hh>
5 
13 
14 namespace Dune{
15  namespace PDELab{
16 
31  template<typename GO, typename LOP, bool nonoverlapping_mode = false>
33  public Dune::PDELab::LocalAssemblerBase<typename GO::Traits::MatrixBackend,
34  typename GO::Traits::TrialGridFunctionSpaceConstraints,
35  typename GO::Traits::TestGridFunctionSpaceConstraints>
36  {
37 
38  // The GridOperator has to be a friend to modify the do{Pre,Post}Processing flags
39  template<typename, typename, typename,
40  typename, typename, typename, typename,
41  typename, typename, int>
42  friend class GridOperator;
43 
44  public:
45 
48 
50  typedef typename Traits::Residual::ElementType RangeField;
51  typedef RangeField Real;
52 
55 
58 
61 
63  typedef LOP LocalOperator;
64 
65  static const bool isNonOverlapping = nonoverlapping_mode;
66 
69  // Types of local function spaces
74 
76 
84 
91 
93  DefaultLocalAssembler (LOP & lop_, shared_ptr<typename GO::BorderDOFExchanger> border_dof_exchanger)
94  : lop(lop_), weight(1.0), doPreProcessing(true), doPostProcessing(true),
95  pattern_engine(*this,border_dof_exchanger), residual_engine(*this), jacobian_engine(*this)
96  , jacobian_apply_engine(*this)
97  , nonlinear_jacobian_apply_engine(*this)
98  , _reconstruct_border_entries(isNonOverlapping)
99  {}
100 
102  DefaultLocalAssembler (LOP & lop_, const CU& cu_, const CV& cv_,
103  shared_ptr<typename GO::BorderDOFExchanger> border_dof_exchanger)
104  : Base(cu_, cv_),
105  lop(lop_), weight(1.0), doPreProcessing(true), doPostProcessing(true),
106  pattern_engine(*this,border_dof_exchanger), residual_engine(*this), jacobian_engine(*this)
107  , jacobian_apply_engine(*this)
108  , nonlinear_jacobian_apply_engine(*this)
109  , _reconstruct_border_entries(isNonOverlapping)
110  {}
111 
115  void setTime(Real time_){
116  lop.setTime(time_);
117  }
118 
120  void setWeight(RangeField weight_){
121  weight = weight_;
122  }
123 
126  void preStage (Real time_, int r_) { lop.preStage(time_,r_); }
127  void preStep (Real time_, Real dt_, std::size_t stages_){ lop.preStep(time_,dt_,stages_); }
128  void postStep (){ lop.postStep(); }
129  void postStage (){ lop.postStage(); }
130  Real suggestTimestep (Real dt) const{return lop.suggestTimestep(dt); }
132 
134  {
135  return _reconstruct_border_entries;
136  }
137 
140 
143  LocalPatternAssemblerEngine & localPatternAssemblerEngine
145  {
146  pattern_engine.setPattern(p);
147  return pattern_engine;
148  }
149 
152  LocalResidualAssemblerEngine & localResidualAssemblerEngine
153  (typename Traits::Residual & r, const typename Traits::Solution & x)
154  {
155  residual_engine.setResidual(r);
156  residual_engine.setSolution(x);
157  return residual_engine;
158  }
159 
162  LocalJacobianAssemblerEngine & localJacobianAssemblerEngine
163  (typename Traits::Jacobian & a, const typename Traits::Solution & x)
164  {
165  jacobian_engine.setJacobian(a);
166  jacobian_engine.setSolution(x);
167  return jacobian_engine;
168  }
169 
172  LocalJacobianApplyAssemblerEngine & localJacobianApplyAssemblerEngine
173  (typename Traits::Residual & r, const typename Traits::Solution & x)
174  {
175  jacobian_apply_engine.setResidual(r);
176  jacobian_apply_engine.setSolution(x);
177  return jacobian_apply_engine;
178  }
179 
182  LocalNonlinearJacobianApplyAssemblerEngine & localNonlinearJacobianApplyAssemblerEngine
183  (typename Traits::Residual & r, const typename Traits::Solution & x, const typename Traits::Solution & z)
184  {
185  nonlinear_jacobian_apply_engine.setResidual(r);
186  nonlinear_jacobian_apply_engine.setSolution(x);
187  nonlinear_jacobian_apply_engine.setUpdate(z);
188  return nonlinear_jacobian_apply_engine;
189  }
190 
192 
197  static bool doAlphaVolume() { return LOP::doAlphaVolume; }
198  static bool doLambdaVolume() { return LOP::doLambdaVolume; }
199  static bool doAlphaSkeleton() { return LOP::doAlphaSkeleton; }
200  static bool doLambdaSkeleton() { return LOP::doLambdaSkeleton; }
201  static bool doAlphaBoundary() { return LOP::doAlphaBoundary; }
202  static bool doLambdaBoundary() { return LOP::doLambdaBoundary; }
203  static bool doAlphaVolumePostSkeleton() { return LOP::doAlphaVolumePostSkeleton; }
204  static bool doLambdaVolumePostSkeleton() { return LOP::doLambdaVolumePostSkeleton; }
205  static bool doSkeletonTwoSided() { return LOP::doSkeletonTwoSided; }
206  static bool doPatternVolume() { return LOP::doPatternVolume; }
207  static bool doPatternSkeleton() { return LOP::doPatternSkeleton; }
208  static bool doPatternBoundary() { return LOP::doPatternBoundary; }
209  static bool doPatternVolumePostSkeleton() { return LOP::doPatternVolumePostSkeleton; }
211 
216  void preProcessing(bool v)
217  {
218  doPreProcessing = v;
219  }
220 
225  void postProcessing(bool v)
226  {
227  doPostProcessing = v;
228  }
229 
230  private:
231 
233  LOP & lop;
234 
236  RangeField weight;
237 
240  bool doPreProcessing;
241 
244  bool doPostProcessing;
245 
248  LocalPatternAssemblerEngine pattern_engine;
249  LocalResidualAssemblerEngine residual_engine;
250  LocalJacobianAssemblerEngine jacobian_engine;
251  LocalJacobianApplyAssemblerEngine jacobian_apply_engine;
252  LocalNonlinearJacobianApplyAssemblerEngine nonlinear_jacobian_apply_engine;
254 
255  bool _reconstruct_border_entries;
256 
257  };
258 
259  }
260 }
261 #endif // DUNE_PDELAB_GRIDOPERATOR_DEFAULT_LOCALASSEMBLER_HH
DefaultLocalPatternAssemblerEngine< DefaultLocalAssembler > LocalPatternAssemblerEngine
Definition: default/localassembler.hh:79
GO::Traits::TrialGridFunctionSpace TrialGridFunctionSpace
The trial grid function space.
Definition: assemblerutilities.hh:26
void setResidual(Residual &residual_)
Definition: jacobianapplyengine.hh:102
Traits::TestGridFunctionSpace GFSV
Definition: default/localassembler.hh:54
Traits::Residual::ElementType RangeField
The local operators type for real numbers e.g. time.
Definition: default/localassembler.hh:50
void setSolution(const Solution &solution_)
Definition: nonlinearjacobianapplyengine.hh:109
static bool doPatternVolumePostSkeleton()
Definition: default/localassembler.hh:209
DefaultLocalAssembler(LOP &lop_, shared_ptr< typename GO::BorderDOFExchanger > border_dof_exchanger)
Constructor with empty constraints.
Definition: default/localassembler.hh:93
LFSIndexCache< LFSV, CV > LFSVCache
Definition: default/localassembler.hh:73
Traits::TestGridFunctionSpaceConstraints CV
Definition: default/localassembler.hh:57
void setPattern(Pattern &pattern_)
Definition: default/patternengine.hh:92
LOP LocalOperator
The local operator.
Definition: default/localassembler.hh:63
void setResidual(Residual &residual_)
Definition: default/residualengine.hh:110
Dune::PDELab::LocalFunctionSpace< GFSU, Dune::PDELab::TrialSpaceTag > LFSU
Definition: default/localassembler.hh:70
RangeField Real
Definition: default/localassembler.hh:51
Definition: lfsindexcache.hh:948
static bool doLambdaSkeleton()
Definition: default/localassembler.hh:200
static bool doLambdaVolume()
Definition: default/localassembler.hh:198
LFSIndexCache< LFSU, CU > LFSUCache
Definition: default/localassembler.hh:72
void preStage(Real time_, int r_)
Definition: default/localassembler.hh:126
Real suggestTimestep(Real dt) const
Definition: default/localassembler.hh:130
Definition: assemblerutilities.hh:22
void setSolution(const Solution &solution_)
Definition: default/residualengine.hh:117
static bool doLambdaVolumePostSkeleton()
Definition: default/localassembler.hh:204
DefaultLocalAssembler(LOP &lop_, const CU &cu_, const CV &cv_, shared_ptr< typename GO::BorderDOFExchanger > border_dof_exchanger)
Constructor for non trivial constraints.
Definition: default/localassembler.hh:102
LocalJacobianAssemblerEngine & localJacobianAssemblerEngine(typename Traits::Jacobian &a, const typename Traits::Solution &x)
Definition: default/localassembler.hh:163
void setJacobian(Jacobian &jacobian_)
Definition: default/jacobianengine.hh:106
For backward compatibility – Do not use this!
Definition: adaptivity.hh:27
LocalPatternAssemblerEngine & localPatternAssemblerEngine(typename Traits::MatrixPattern &p)
Definition: default/localassembler.hh:144
Create a local function space from a global function space.
Definition: localfunctionspace.hh:670
GO::Traits::TestGridFunctionSpaceConstraints TestGridFunctionSpaceConstraints
The type of the test grid function space constraints.
Definition: assemblerutilities.hh:36
Standard grid operator implementation.
Definition: gridoperator.hh:57
LocalResidualAssemblerEngine & localResidualAssemblerEngine(typename Traits::Residual &r, const typename Traits::Solution &x)
Definition: default/localassembler.hh:153
GO::Traits::Jacobian Jacobian
The type of the jacobian.
Definition: assemblerutilities.hh:61
void setResidual(Residual &residual_)
Definition: nonlinearjacobianapplyengine.hh:102
bool reconstructBorderEntries() const
Definition: default/localassembler.hh:133
static bool doPatternVolume()
Definition: default/localassembler.hh:206
static const bool isNonOverlapping
Definition: default/localassembler.hh:65
void postStep()
Definition: default/localassembler.hh:128
The local assembler for DUNE grids.
Definition: default/localassembler.hh:32
void setTime(Real time_)
Definition: default/localassembler.hh:115
static bool doAlphaVolumePostSkeleton()
Definition: default/localassembler.hh:203
Dune::PDELab::LocalAssemblerTraits< GO > Traits
The traits class.
Definition: default/localassembler.hh:47
void setWeight(RangeField weight_)
Notifies the assembler about the current weight of assembling.
Definition: default/localassembler.hh:120
GO::Traits::Range Residual
The type of the range (residual).
Definition: assemblerutilities.hh:54
DefaultLocalJacobianApplyAssemblerEngine< DefaultLocalAssembler > LocalJacobianApplyAssemblerEngine
Definition: default/localassembler.hh:82
DefaultLocalNonlinearJacobianApplyAssemblerEngine< DefaultLocalAssembler > LocalNonlinearJacobianApplyAssemblerEngine
Definition: default/localassembler.hh:83
GO::Traits::Domain Solution
The type of the domain (solution).
Definition: assemblerutilities.hh:47
DefaultLocalResidualAssemblerEngine< DefaultLocalAssembler > LocalResidualAssemblerEngine
Definition: default/localassembler.hh:80
void setSolution(const Solution &solution_)
Definition: jacobianapplyengine.hh:109
DefaultLocalJacobianAssemblerEngine< DefaultLocalAssembler > LocalJacobianAssemblerEngine
Definition: default/localassembler.hh:81
MatrixBackend::template Pattern< Jacobian, TestGridFunctionSpace, TrialGridFunctionSpace > MatrixPattern
The matrix pattern.
Definition: assemblerutilities.hh:68
static bool doLambdaBoundary()
Definition: default/localassembler.hh:202
static bool doSkeletonTwoSided()
Definition: default/localassembler.hh:205
GO::Traits::TrialGridFunctionSpaceConstraints TrialGridFunctionSpaceConstraints
The type of the trial grid function space constraints.
Definition: assemblerutilities.hh:33
void postProcessing(bool v)
Definition: default/localassembler.hh:225
Traits::TrialGridFunctionSpace GFSU
Definition: default/localassembler.hh:53
const P & p
Definition: constraints.hh:147
Base class for local assembler.
Definition: assemblerutilities.hh:182
void setUpdate(const Solution &update_)
Definition: nonlinearjacobianapplyengine.hh:116
static bool doPatternBoundary()
Definition: default/localassembler.hh:208
void setSolution(const Solution &solution_)
Definition: default/jacobianengine.hh:115
void preProcessing(bool v)
Definition: default/localassembler.hh:216
Traits::TrialGridFunctionSpaceConstraints CU
Definition: default/localassembler.hh:56
Dune::PDELab::LocalAssemblerBase< typename Traits::MatrixBackend, CU, CV > Base
The base class of this local assembler.
Definition: default/localassembler.hh:60
static bool doAlphaVolume()
Query methods for the assembler engines. Theses methods do not belong to the assembler interface...
Definition: default/localassembler.hh:197
void postStage()
Definition: default/localassembler.hh:129
static bool doAlphaBoundary()
Definition: default/localassembler.hh:201
void preStep(Real time_, Real dt_, std::size_t stages_)
Definition: default/localassembler.hh:127
static bool doPatternSkeleton()
Definition: default/localassembler.hh:207
LocalJacobianApplyAssemblerEngine & localJacobianApplyAssemblerEngine(typename Traits::Residual &r, const typename Traits::Solution &x)
Definition: default/localassembler.hh:173
static bool doAlphaSkeleton()
Definition: default/localassembler.hh:199
GO::Traits::TestGridFunctionSpace TestGridFunctionSpace
The test grid function space.
Definition: assemblerutilities.hh:29
LocalNonlinearJacobianApplyAssemblerEngine & localNonlinearJacobianApplyAssemblerEngine(typename Traits::Residual &r, const typename Traits::Solution &x, const typename Traits::Solution &z)
Definition: default/localassembler.hh:183
Dune::PDELab::LocalFunctionSpace< GFSV, Dune::PDELab::TestSpaceTag > LFSV
Definition: default/localassembler.hh:71