dune-pdelab  2.5-dev
fastdg/localassembler.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_FASTDG_LOCALASSEMBLER_HH
4 #define DUNE_PDELAB_GRIDOPERATOR_FASTDG_LOCALASSEMBLER_HH
5 
6 #include <memory>
7 
8 #include <dune/common/shared_ptr.hh>
9 
10 #include <dune/typetree/typetree.hh>
11 
19 
20 namespace Dune{
21  namespace PDELab{
22 
37  template<typename GO, typename LOP, bool nonoverlapping_mode = false>
39  public Dune::PDELab::LocalAssemblerBase<typename GO::Traits::MatrixBackend,
40  typename GO::Traits::TrialGridFunctionSpaceConstraints,
41  typename GO::Traits::TestGridFunctionSpaceConstraints>
42  {
43  public:
44 
47 
49  typedef typename Traits::Residual::ElementType RangeField;
50  typedef RangeField Real;
51 
54 
57 
60 
62  typedef LOP LocalOperator;
63 
64  static const bool isNonOverlapping = nonoverlapping_mode;
65 
68  // Types of local function spaces
73 
75 
83 
84  // friend declarations such that engines are able to call scatter_jacobian() and add_entry() from base class
87 
89  FastDGLocalAssembler (LOP & lop_, shared_ptr<typename GO::BorderDOFExchanger> border_dof_exchanger)
90  : lop(stackobject_to_shared_ptr(lop_)),
91  weight_(1.0),
92  doPreProcessing_(true),
93  doPostProcessing_(true),
94  pattern_engine(*this,border_dof_exchanger), residual_engine(*this), jacobian_engine(*this), jacobian_apply_engine(*this), nonlinear_jacobian_apply_engine(*this)
95  , _reconstruct_border_entries(isNonOverlapping)
96  {}
97 
99  FastDGLocalAssembler (LOP & lop_, const CU& cu_, const CV& cv_,
100  shared_ptr<typename GO::BorderDOFExchanger> border_dof_exchanger)
101  : Base(cu_, cv_),
102  lop(stackobject_to_shared_ptr(lop_)),
103  weight_(1.0),
104  doPreProcessing_(true),
105  doPostProcessing_(true),
106  pattern_engine(*this,border_dof_exchanger), residual_engine(*this), jacobian_engine(*this), jacobian_apply_engine(*this), nonlinear_jacobian_apply_engine(*this)
107  , _reconstruct_border_entries(isNonOverlapping)
108  {}
109 
112  {
113  return *lop;
114  }
116  const LOP &localOperator() const
117  {
118  return *lop;
119  }
120 
124  void setTime(Real time_)
125  {
126  lop->setTime(time_);
127  }
128 
130  RangeField weight() const { return weight_; }
131 
133  void setWeight(RangeField weight){
134  weight_ = weight;
135  }
136 
139  void preStage (Real time_, int r_) { lop->preStage(time_,r_); }
140  void preStep (Real time_, Real dt_, std::size_t stages_){ lop->preStep(time_,dt_,stages_); }
141  void postStep (){ lop->postStep(); }
142  void postStage (){ lop->postStage(); }
143  Real suggestTimestep (Real dt) const{return lop->suggestTimestep(dt); }
145 
147  {
148  return _reconstruct_border_entries;
149  }
150 
153 
156  LocalPatternAssemblerEngine & localPatternAssemblerEngine
158  {
159  pattern_engine.setPattern(p);
160  return pattern_engine;
161  }
162 
165  LocalResidualAssemblerEngine & localResidualAssemblerEngine
166  (typename Traits::Residual & r, const typename Traits::Solution & x)
167  {
168  residual_engine.setResidual(r);
169  residual_engine.setSolution(x);
170  return residual_engine;
171  }
172 
175  LocalJacobianAssemblerEngine & localJacobianAssemblerEngine
176  (typename Traits::Jacobian & a, const typename Traits::Solution & x)
177  {
178  jacobian_engine.setJacobian(a);
179  jacobian_engine.setSolution(x);
180  return jacobian_engine;
181  }
182 
185  LocalJacobianApplyAssemblerEngine & localJacobianApplyAssemblerEngine
186  (typename Traits::Residual & r, const typename Traits::Solution & x)
187  {
188  jacobian_apply_engine.setResidual(r);
189  jacobian_apply_engine.setSolution(x);
190  return jacobian_apply_engine;
191  }
192 
195  LocalNonlinearJacobianApplyAssemblerEngine & localNonlinearJacobianApplyAssemblerEngine
196  (typename Traits::Residual & r, const typename Traits::Solution & x, const typename Traits::Solution & z)
197  {
198  nonlinear_jacobian_apply_engine.setResidual(r);
199  nonlinear_jacobian_apply_engine.setSolution(x);
200  nonlinear_jacobian_apply_engine.setUpdate(z);
201  return nonlinear_jacobian_apply_engine;
202  }
203 
205 
210  static bool doAlphaVolume() { return LOP::doAlphaVolume; }
211  static bool doLambdaVolume() { return LOP::doLambdaVolume; }
212  static bool doAlphaSkeleton() { return LOP::doAlphaSkeleton; }
213  static bool doLambdaSkeleton() { return LOP::doLambdaSkeleton; }
214  static bool doAlphaBoundary() { return LOP::doAlphaBoundary; }
215  static bool doLambdaBoundary() { return LOP::doLambdaBoundary; }
216  static bool doAlphaVolumePostSkeleton() { return LOP::doAlphaVolumePostSkeleton; }
217  static bool doLambdaVolumePostSkeleton() { return LOP::doLambdaVolumePostSkeleton; }
218  static bool doSkeletonTwoSided() { return LOP::doSkeletonTwoSided; }
219  static bool doPatternVolume() { return LOP::doPatternVolume; }
220  static bool doPatternSkeleton() { return LOP::doPatternSkeleton; }
221  static bool doPatternBoundary() { return LOP::doPatternBoundary; }
222  static bool doPatternVolumePostSkeleton() { return LOP::doPatternVolumePostSkeleton; }
224 
226 
229  bool doPreProcessing() const { return doPreProcessing_; }
230 
235  void preProcessing(bool v)
236  {
237  doPreProcessing_ = v;
238  }
239 
241 
244  bool doPostProcessing() const { return doPostProcessing_; }
245 
250  void postProcessing(bool v)
251  {
252  doPostProcessing_ = v;
253  }
254 
255  template<typename GFSV, typename GC, typename C>
256  void set_trivial_rows(const GFSV& gfsv, GC& globalcontainer, const C& c) const
257  {
258  typedef typename C::const_iterator global_row_iterator;
259  for (global_row_iterator cit = c.begin(); cit != c.end(); ++cit)
260  globalcontainer.clear_row_block(cit->first,1);
261  }
262 
263  template<typename GFSV, typename GC>
264  void set_trivial_rows(const GFSV& gfsv, GC& globalcontainer, const EmptyTransformation& c) const
265  {
266  }
267 
268  template<typename GFSV, typename GC>
269  void handle_dirichlet_constraints(const GFSV& gfsv, GC& globalcontainer) const
270  {
271  globalcontainer.flush();
272  set_trivial_rows(gfsv,globalcontainer,*this->pconstraintsv);
273  globalcontainer.finalize();
274  }
275 
276  private:
277 
279  std::shared_ptr<LOP> lop;
280 
282  RangeField weight_;
283 
286  bool doPreProcessing_;
287 
290  bool doPostProcessing_;
291 
294  LocalPatternAssemblerEngine pattern_engine;
295  LocalResidualAssemblerEngine residual_engine;
296  LocalJacobianAssemblerEngine jacobian_engine;
297  LocalJacobianApplyAssemblerEngine jacobian_apply_engine;
298  LocalNonlinearJacobianApplyAssemblerEngine nonlinear_jacobian_apply_engine;
300 
301  bool _reconstruct_border_entries;
302 
303  }; // end class FastDGLocalAssembler
304  } // end namespace PDELab
305 } // end namespace Dune
306 #endif
static bool doAlphaSkeleton()
Definition: fastdg/localassembler.hh:212
Real suggestTimestep(Real dt) const
Definition: fastdg/localassembler.hh:143
Traits::TestGridFunctionSpace GFSV
Definition: fastdg/localassembler.hh:53
bool doPreProcessing() const
Query whether to do preprocessing in the engines.
Definition: fastdg/localassembler.hh:229
static const bool isNonOverlapping
Definition: fastdg/localassembler.hh:64
Traits::TestGridFunctionSpaceConstraints CV
Definition: fastdg/localassembler.hh:56
LocalJacobianApplyAssemblerEngine & localJacobianApplyAssemblerEngine(typename Traits::Residual &r, const typename Traits::Solution &x)
Definition: fastdg/localassembler.hh:186
void setSolution(const Solution &solution_)
Definition: fastdg/jacobianapplyengine.hh:109
GO::Traits::Range Residual
The type of the range (residual).
Definition: assemblerutilities.hh:54
GO::Traits::TrialGridFunctionSpaceConstraints TrialGridFunctionSpaceConstraints
The type of the trial grid function space constraints.
Definition: assemblerutilities.hh:33
Definition: constraintstransformation.hh:111
void preStage(Real time_, int r_)
Definition: fastdg/localassembler.hh:139
Base class for local assembler.
Definition: assemblerutilities.hh:182
void setUpdate(const Solution &update_)
Definition: fastdg/nonlinearjacobianapplyengine.hh:120
The local assembler for DUNE grids.
Definition: fastdg/localassembler.hh:38
void postProcessing(bool v)
Definition: fastdg/localassembler.hh:250
Definition: assemblerutilities.hh:22
static bool doPatternVolumePostSkeleton()
Definition: fastdg/localassembler.hh:222
static bool doLambdaVolumePostSkeleton()
Definition: fastdg/localassembler.hh:217
FastDGLocalResidualAssemblerEngine< FastDGLocalAssembler > LocalResidualAssemblerEngine
Definition: fastdg/localassembler.hh:79
GO::Traits::Jacobian Jacobian
The type of the jacobian.
Definition: assemblerutilities.hh:61
Definition: lfsindexcache.hh:977
void setSolution(const Solution &solution_)
Definition: fastdg/residualengine.hh:148
FastDGLocalPatternAssemblerEngine< FastDGLocalAssembler > LocalPatternAssemblerEngine
Definition: fastdg/localassembler.hh:78
static bool doPatternSkeleton()
Definition: fastdg/localassembler.hh:220
GO::Traits::TrialGridFunctionSpace TrialGridFunctionSpace
The trial grid function space.
Definition: assemblerutilities.hh:26
MatrixBackend::template Pattern< Jacobian, TestGridFunctionSpace, TrialGridFunctionSpace > MatrixPattern
The matrix pattern.
Definition: assemblerutilities.hh:68
static bool doAlphaVolume()
Query methods for the assembler engines. Theses methods do not belong to the assembler interface...
Definition: fastdg/localassembler.hh:210
LFSIndexCache< LFSU, CU, true > LFSUCache
Definition: fastdg/localassembler.hh:71
bool doPostProcessing() const
Query whether to do postprocessing in the engines.
Definition: fastdg/localassembler.hh:244
static bool doAlphaVolumePostSkeleton()
Definition: fastdg/localassembler.hh:216
void setResidual(Residual &residual_)
Definition: fastdg/nonlinearjacobianapplyengine.hh:104
void setTime(Real time_)
Definition: fastdg/localassembler.hh:124
void setResidual(Residual &residual_)
Definition: fastdg/jacobianapplyengine.hh:101
FastDGLocalNonlinearJacobianApplyAssemblerEngine< FastDGLocalAssembler > LocalNonlinearJacobianApplyAssemblerEngine
Definition: fastdg/localassembler.hh:82
bool reconstructBorderEntries() const
Definition: fastdg/localassembler.hh:146
FastDGLocalAssembler(LOP &lop_, const CU &cu_, const CV &cv_, shared_ptr< typename GO::BorderDOFExchanger > border_dof_exchanger)
Constructor for non trivial constraints.
Definition: fastdg/localassembler.hh:99
Traits::TrialGridFunctionSpace GFSU
Definition: fastdg/localassembler.hh:52
static bool doPatternBoundary()
Definition: fastdg/localassembler.hh:221
GO::Traits::TestGridFunctionSpace TestGridFunctionSpace
The test grid function space.
Definition: assemblerutilities.hh:29
RangeField weight() const
Obtain the weight that was set last.
Definition: fastdg/localassembler.hh:130
LocalNonlinearJacobianApplyAssemblerEngine & localNonlinearJacobianApplyAssemblerEngine(typename Traits::Residual &r, const typename Traits::Solution &x, const typename Traits::Solution &z)
Definition: fastdg/localassembler.hh:196
void setSolution(const Solution &solution_)
Definition: fastdg/nonlinearjacobianapplyengine.hh:112
Dune::PDELab::LocalFunctionSpace< GFSV, Dune::PDELab::TestSpaceTag > LFSV
Definition: fastdg/localassembler.hh:70
RangeField Real
Definition: fastdg/localassembler.hh:50
LOP LocalOperator
The local operator.
Definition: fastdg/localassembler.hh:62
GO::Traits::TestGridFunctionSpaceConstraints TestGridFunctionSpaceConstraints
The type of the test grid function space constraints.
Definition: assemblerutilities.hh:36
static bool doLambdaVolume()
Definition: fastdg/localassembler.hh:211
Traits::TrialGridFunctionSpaceConstraints CU
Definition: fastdg/localassembler.hh:55
For backward compatibility – Do not use this!
Definition: adaptivity.hh:28
Traits::Residual::ElementType RangeField
The local operators type for real numbers e.g. time.
Definition: fastdg/localassembler.hh:49
void set_trivial_rows(const GFSV &gfsv, GC &globalcontainer, const EmptyTransformation &c) const
Definition: fastdg/localassembler.hh:264
Dune::PDELab::LocalAssemblerTraits< GO > Traits
The traits class.
Definition: fastdg/localassembler.hh:46
void postStage()
Definition: fastdg/localassembler.hh:142
void set_trivial_rows(const GFSV &gfsv, GC &globalcontainer, const C &c) const
Definition: fastdg/localassembler.hh:256
static bool doAlphaBoundary()
Definition: fastdg/localassembler.hh:214
void setWeight(RangeField weight)
Notifies the assembler about the current weight of assembling.
Definition: fastdg/localassembler.hh:133
void setSolution(const Solution &solution_)
Definition: fastdg/jacobianengine.hh:142
Dune::PDELab::LocalFunctionSpace< GFSU, Dune::PDELab::TrialSpaceTag > LFSU
Definition: fastdg/localassembler.hh:69
static bool doPatternVolume()
Definition: fastdg/localassembler.hh:219
static bool doLambdaBoundary()
Definition: fastdg/localassembler.hh:215
GO::Traits::Domain Solution
The type of the domain (solution).
Definition: assemblerutilities.hh:47
void preProcessing(bool v)
Definition: fastdg/localassembler.hh:235
void preStep(Real time_, Real dt_, std::size_t stages_)
Definition: fastdg/localassembler.hh:140
FastDGLocalJacobianApplyAssemblerEngine< FastDGLocalAssembler > LocalJacobianApplyAssemblerEngine
Definition: fastdg/localassembler.hh:81
LocalPatternAssemblerEngine & localPatternAssemblerEngine(typename Traits::MatrixPattern &p)
Definition: fastdg/localassembler.hh:157
LocalJacobianAssemblerEngine & localJacobianAssemblerEngine(typename Traits::Jacobian &a, const typename Traits::Solution &x)
Definition: fastdg/localassembler.hh:176
FastDGLocalAssembler(LOP &lop_, shared_ptr< typename GO::BorderDOFExchanger > border_dof_exchanger)
Constructor with empty constraints.
Definition: fastdg/localassembler.hh:89
void setJacobian(Jacobian &jacobian_)
Definition: fastdg/jacobianengine.hh:132
LOP & localOperator()
get a reference to the local operator
Definition: fastdg/localassembler.hh:111
const LOP & localOperator() const
get a reference to the local operator
Definition: fastdg/localassembler.hh:116
void setResidual(Residual &residual_)
Definition: fastdg/residualengine.hh:140
static bool doSkeletonTwoSided()
Definition: fastdg/localassembler.hh:218
FastDGLocalJacobianAssemblerEngine< FastDGLocalAssembler > LocalJacobianAssemblerEngine
Definition: fastdg/localassembler.hh:80
static bool doLambdaSkeleton()
Definition: fastdg/localassembler.hh:213
Create a local function space from a global function space.
Definition: localfunctionspace.hh:698
void postStep()
Definition: fastdg/localassembler.hh:141
LFSIndexCache< LFSV, CV, true > LFSVCache
Definition: fastdg/localassembler.hh:72
Dune::PDELab::LocalAssemblerBase< typename Traits::MatrixBackend, CU, CV > Base
The base class of this local assembler.
Definition: fastdg/localassembler.hh:59
void handle_dirichlet_constraints(const GFSV &gfsv, GC &globalcontainer) const
Definition: fastdg/localassembler.hh:269
LocalResidualAssemblerEngine & localResidualAssemblerEngine(typename Traits::Residual &r, const typename Traits::Solution &x)
Definition: fastdg/localassembler.hh:166
const P & p
Definition: constraints.hh:147
void setPattern(Pattern &pattern_)
Definition: fastdg/patternengine.hh:92