dune-pdelab 2.7-git
Loading...
Searching...
No Matches
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
18
19namespace Dune{
20 namespace PDELab{
21
36 template<typename GO, typename LOP, bool nonoverlapping_mode = false>
38 public Dune::PDELab::LocalAssemblerBase<typename GO::Traits::MatrixBackend,
39 typename GO::Traits::TrialGridFunctionSpaceConstraints,
40 typename GO::Traits::TestGridFunctionSpaceConstraints>
41 {
42 public:
43
46
48 typedef typename Traits::Residual::ElementType RangeField;
50
53
56
59
61 typedef LOP LocalOperator;
62
63 static const bool isNonOverlapping = nonoverlapping_mode;
64
67 // Types of local function spaces
72
74
81
82 // friend declarations such that engines are able to call scatter_jacobian() and add_entry() from base class
85
87 FastDGLocalAssembler (LOP & lop_, std::shared_ptr<typename GO::BorderDOFExchanger> border_dof_exchanger)
88 : lop(stackobject_to_shared_ptr(lop_)),
89 weight_(1.0),
90 doPreProcessing_(true),
91 doPostProcessing_(true),
92 pattern_engine(*this,border_dof_exchanger), residual_engine(*this), jacobian_engine(*this)
93 , jacobian_apply_engine(*this)
94 , _reconstruct_border_entries(isNonOverlapping)
95 {}
96
98 FastDGLocalAssembler (LOP & lop_, const CU& cu_, const CV& cv_,
99 std::shared_ptr<typename GO::BorderDOFExchanger> border_dof_exchanger)
100 : Base(cu_, cv_),
101 lop(stackobject_to_shared_ptr(lop_)),
102 weight_(1.0),
103 doPreProcessing_(true),
104 doPostProcessing_(true),
105 pattern_engine(*this,border_dof_exchanger), residual_engine(*this), jacobian_engine(*this)
106 , jacobian_apply_engine(*this)
107 , _reconstruct_border_entries(isNonOverlapping)
108 {}
109
112 {
113 return *lop;
114 }
115
117 const LOP &localOperator() const
118 {
119 return *lop;
120 }
121
125 void setTime(Real time_)
126 {
127 lop->setTime(time_);
128 }
129
132 {
133 return weight_;
134 }
135
138 weight_ = weight;
139 }
140
143 void preStage (Real time_, int r_) { lop->preStage(time_,r_); }
144 void preStep (Real time_, Real dt_, std::size_t stages_){ lop->preStep(time_,dt_,stages_); }
145 void postStep (){ lop->postStep(); }
146 void postStage (){ lop->postStage(); }
147 Real suggestTimestep (Real dt) const{return lop->suggestTimestep(dt); }
149
151 {
152 return _reconstruct_border_entries;
153 }
154
157
161 (typename Traits::MatrixPattern & p)
162 {
163 pattern_engine.setPattern(p);
164 return pattern_engine;
165 }
166
170 (typename Traits::Residual & r, const typename Traits::Solution & x)
171 {
172 residual_engine.setResidual(r);
173 residual_engine.setSolution(x);
174 return residual_engine;
175 }
176
180 (typename Traits::Jacobian & a, const typename Traits::Solution & x)
181 {
182 jacobian_engine.setJacobian(a);
183 jacobian_engine.setSolution(x);
184 return jacobian_engine;
185 }
186
190 (const typename Traits::Domain & update, typename Traits::Range & result)
191 {
192 jacobian_apply_engine.setUpdate(update);
193 jacobian_apply_engine.setResult(result);
194 return jacobian_apply_engine;
195 }
196
200 (const typename Traits::Domain & solution, const typename Traits::Domain & update, typename Traits::Range & result)
201 {
202 jacobian_apply_engine.setSolution(solution);
203 jacobian_apply_engine.setUpdate(update);
204 jacobian_apply_engine.setResult(result);
205 return jacobian_apply_engine;
206 }
207
209
214 static constexpr bool doAlphaVolume() { return LOP::doAlphaVolume; }
215 static constexpr bool doLambdaVolume() { return LOP::doLambdaVolume; }
216 static constexpr bool doAlphaSkeleton() { return LOP::doAlphaSkeleton; }
217 static constexpr bool doLambdaSkeleton() { return LOP::doLambdaSkeleton; }
218 static constexpr bool doAlphaBoundary() { return LOP::doAlphaBoundary; }
219 static constexpr bool doLambdaBoundary() { return LOP::doLambdaBoundary; }
220 static constexpr bool doAlphaVolumePostSkeleton() { return LOP::doAlphaVolumePostSkeleton; }
221 static constexpr bool doLambdaVolumePostSkeleton() { return LOP::doLambdaVolumePostSkeleton; }
222 static constexpr bool doSkeletonTwoSided() { return LOP::doSkeletonTwoSided; }
223 static constexpr bool doPatternVolume() { return LOP::doPatternVolume; }
224 static constexpr bool doPatternSkeleton() { return LOP::doPatternSkeleton; }
225 static constexpr bool doPatternBoundary() { return LOP::doPatternBoundary; }
226 static constexpr bool doPatternVolumePostSkeleton() { return LOP::doPatternVolumePostSkeleton; }
227 static constexpr bool isLinear() { return LOP::isLinear;}
229
231
234 bool doPreProcessing() const { return doPreProcessing_; }
235
240 void preProcessing(bool v)
241 {
242 doPreProcessing_ = v;
243 }
244
246
249 bool doPostProcessing() const { return doPostProcessing_; }
250
255 void postProcessing(bool v)
256 {
257 doPostProcessing_ = v;
258 }
259
260 template<typename GFSV, typename GC, typename C>
261 void set_trivial_rows(const GFSV& gfsv, GC& globalcontainer, const C& c) const
262 {
263 typedef typename C::const_iterator global_row_iterator;
264 for (global_row_iterator cit = c.begin(); cit != c.end(); ++cit)
265 globalcontainer.clear_row_block(cit->first,1);
266 }
267
268 template<typename GFSV, typename GC>
269 void set_trivial_rows(const GFSV& gfsv, GC& globalcontainer, const EmptyTransformation& c) const
270 {
271 }
272
273 template<typename GFSV, typename GC>
274 void handle_dirichlet_constraints(const GFSV& gfsv, GC& globalcontainer) const
275 {
276 globalcontainer.flush();
277 set_trivial_rows(gfsv,globalcontainer,*this->pconstraintsv);
278 globalcontainer.finalize();
279 }
280
281 private:
282
284 std::shared_ptr<LOP> lop;
285
287 RangeField weight_;
288
291 bool doPreProcessing_;
292
295 bool doPostProcessing_;
296
299 LocalPatternAssemblerEngine pattern_engine;
300 LocalResidualAssemblerEngine residual_engine;
301 LocalJacobianAssemblerEngine jacobian_engine;
302 LocalJacobianApplyAssemblerEngine jacobian_apply_engine;
304
305 bool _reconstruct_border_entries;
306 };
307
308 } // end namespace PDELab
309} // end namespace Dune
310#endif
const P & p
Definition constraints.hh:148
For backward compatibility – Do not use this!
Definition adaptivity.hh:28
Definition constraintstransformation.hh:112
Definition lfsindexcache.hh:979
Create a local function space from a global function space.
Definition localfunctionspace.hh:717
Definition assemblerutilities.hh:23
GO::Traits::Range Residual
The type of the range (residual).
Definition assemblerutilities.hh:60
MatrixBackend::template Pattern< Jacobian, TestGridFunctionSpace, TrialGridFunctionSpace > MatrixPattern
The matrix pattern.
Definition assemblerutilities.hh:74
GO::Traits::TrialGridFunctionSpace TrialGridFunctionSpace
The trial grid function space.
Definition assemblerutilities.hh:26
GO::Traits::Jacobian Jacobian
The type of the jacobian.
Definition assemblerutilities.hh:67
GO::Traits::TestGridFunctionSpace TestGridFunctionSpace
The test grid function space.
Definition assemblerutilities.hh:29
GO::Traits::Domain Solution
The type of the domain (solution).
Definition assemblerutilities.hh:50
GO::Traits::Range Range
The type of the range (residual).
Definition assemblerutilities.hh:57
GO::Traits::Domain Domain
The type of the domain (solution).
Definition assemblerutilities.hh:47
GO::Traits::TrialGridFunctionSpaceConstraints TrialGridFunctionSpaceConstraints
The type of the trial grid function space constraints.
Definition assemblerutilities.hh:33
GO::Traits::TestGridFunctionSpaceConstraints TestGridFunctionSpaceConstraints
The type of the test grid function space constraints.
Definition assemblerutilities.hh:36
Base class for local assembler.
Definition assemblerutilities.hh:189
The fast DG local assembler engine for DUNE grids which assembles the local application of the Jacobi...
Definition fastdg/jacobianapplyengine.hh:23
void setResult(Range &result_)
Definition fastdg/jacobianapplyengine.hh:125
void setSolution(const Domain &solution_)
Definition fastdg/jacobianapplyengine.hh:107
void setUpdate(const Domain &update_)
Definition fastdg/jacobianapplyengine.hh:117
The fast DG local assembler engine for DUNE grids which assembles the jacobian matrix.
Definition fastdg/jacobianengine.hh:27
void setSolution(const Solution &solution_)
Definition fastdg/jacobianengine.hh:142
void setJacobian(Jacobian &jacobian_)
Definition fastdg/jacobianengine.hh:132
The local assembler for DUNE grids.
Definition fastdg/localassembler.hh:41
LOP & localOperator()
get a reference to the local operator
Definition fastdg/localassembler.hh:111
static constexpr bool doAlphaVolumePostSkeleton()
Definition fastdg/localassembler.hh:220
static constexpr bool doPatternVolume()
Definition fastdg/localassembler.hh:223
static constexpr bool doLambdaBoundary()
Definition fastdg/localassembler.hh:219
Traits::TrialGridFunctionSpaceConstraints CU
Definition fastdg/localassembler.hh:54
FastDGLocalAssembler(LOP &lop_, const CU &cu_, const CV &cv_, std::shared_ptr< typename GO::BorderDOFExchanger > border_dof_exchanger)
Constructor for non trivial constraints.
Definition fastdg/localassembler.hh:98
FastDGLocalResidualAssemblerEngine< FastDGLocalAssembler > LocalResidualAssemblerEngine
Definition fastdg/localassembler.hh:78
static const bool isNonOverlapping
Definition fastdg/localassembler.hh:63
void postStage()
Definition fastdg/localassembler.hh:146
LFSIndexCache< LFSV, CV, true > LFSVCache
Definition fastdg/localassembler.hh:71
bool reconstructBorderEntries() const
Definition fastdg/localassembler.hh:150
LOP LocalOperator
The local operator.
Definition fastdg/localassembler.hh:61
Traits::Residual::ElementType RangeField
The local operators type for real numbers e.g. time.
Definition fastdg/localassembler.hh:48
Dune::PDELab::LocalFunctionSpace< GFSU, Dune::PDELab::TrialSpaceTag > LFSU
Definition fastdg/localassembler.hh:68
static constexpr bool doLambdaVolume()
Definition fastdg/localassembler.hh:215
void postStep()
Definition fastdg/localassembler.hh:145
static constexpr bool doPatternBoundary()
Definition fastdg/localassembler.hh:225
Traits::TestGridFunctionSpaceConstraints CV
Definition fastdg/localassembler.hh:55
void postProcessing(bool v)
Definition fastdg/localassembler.hh:255
static constexpr bool doAlphaSkeleton()
Definition fastdg/localassembler.hh:216
Dune::PDELab::LocalFunctionSpace< GFSV, Dune::PDELab::TestSpaceTag > LFSV
Definition fastdg/localassembler.hh:69
LocalJacobianApplyAssemblerEngine & localJacobianApplyAssemblerEngine(const typename Traits::Domain &solution, const typename Traits::Domain &update, typename Traits::Range &result)
Definition fastdg/localassembler.hh:200
static constexpr bool doAlphaVolume()
Query methods for the assembler engines. Theses methods do not belong to the assembler interface,...
Definition fastdg/localassembler.hh:214
const LOP & localOperator() const
get a reference to the local operator
Definition fastdg/localassembler.hh:117
static constexpr bool doLambdaVolumePostSkeleton()
Definition fastdg/localassembler.hh:221
RangeField Real
Definition fastdg/localassembler.hh:49
void preStage(Real time_, int r_)
Definition fastdg/localassembler.hh:143
LocalJacobianAssemblerEngine & localJacobianAssemblerEngine(typename Traits::Jacobian &a, const typename Traits::Solution &x)
Definition fastdg/localassembler.hh:180
static constexpr bool doAlphaBoundary()
Definition fastdg/localassembler.hh:218
void preProcessing(bool v)
Definition fastdg/localassembler.hh:240
FastDGLocalAssembler(LOP &lop_, std::shared_ptr< typename GO::BorderDOFExchanger > border_dof_exchanger)
Constructor with empty constraints.
Definition fastdg/localassembler.hh:87
static constexpr bool doLambdaSkeleton()
Definition fastdg/localassembler.hh:217
void setTime(Real time_)
Definition fastdg/localassembler.hh:125
bool doPreProcessing() const
Query whether to do preprocessing in the engines.
Definition fastdg/localassembler.hh:234
LFSIndexCache< LFSU, CU, true > LFSUCache
Definition fastdg/localassembler.hh:70
bool doPostProcessing() const
Query whether to do postprocessing in the engines.
Definition fastdg/localassembler.hh:249
Traits::TestGridFunctionSpace GFSV
Definition fastdg/localassembler.hh:52
Traits::TrialGridFunctionSpace GFSU
Definition fastdg/localassembler.hh:51
LocalPatternAssemblerEngine & localPatternAssemblerEngine(typename Traits::MatrixPattern &p)
Definition fastdg/localassembler.hh:161
static constexpr bool doPatternSkeleton()
Definition fastdg/localassembler.hh:224
LocalJacobianApplyAssemblerEngine & localJacobianApplyAssemblerEngine(const typename Traits::Domain &update, typename Traits::Range &result)
Definition fastdg/localassembler.hh:190
RangeField weight() const
Obtain the weight that was set last.
Definition fastdg/localassembler.hh:131
void preStep(Real time_, Real dt_, std::size_t stages_)
Definition fastdg/localassembler.hh:144
Dune::PDELab::LocalAssemblerBase< typename Traits::MatrixBackend, CU, CV > Base
The base class of this local assembler.
Definition fastdg/localassembler.hh:58
void set_trivial_rows(const GFSV &gfsv, GC &globalcontainer, const C &c) const
Definition fastdg/localassembler.hh:261
FastDGLocalPatternAssemblerEngine< FastDGLocalAssembler > LocalPatternAssemblerEngine
Definition fastdg/localassembler.hh:77
FastDGLocalJacobianAssemblerEngine< FastDGLocalAssembler > LocalJacobianAssemblerEngine
Definition fastdg/localassembler.hh:79
Real suggestTimestep(Real dt) const
Definition fastdg/localassembler.hh:147
void setWeight(RangeField weight)
Notifies the assembler about the current weight of assembling.
Definition fastdg/localassembler.hh:137
static constexpr bool doPatternVolumePostSkeleton()
Definition fastdg/localassembler.hh:226
Dune::PDELab::LocalAssemblerTraits< GO > Traits
The traits class.
Definition fastdg/localassembler.hh:45
static constexpr bool doSkeletonTwoSided()
Definition fastdg/localassembler.hh:222
void set_trivial_rows(const GFSV &gfsv, GC &globalcontainer, const EmptyTransformation &c) const
Definition fastdg/localassembler.hh:269
FastDGLocalJacobianApplyAssemblerEngine< FastDGLocalAssembler > LocalJacobianApplyAssemblerEngine
Definition fastdg/localassembler.hh:80
LocalResidualAssemblerEngine & localResidualAssemblerEngine(typename Traits::Residual &r, const typename Traits::Solution &x)
Definition fastdg/localassembler.hh:170
static constexpr bool isLinear()
Definition fastdg/localassembler.hh:227
void handle_dirichlet_constraints(const GFSV &gfsv, GC &globalcontainer) const
Definition fastdg/localassembler.hh:274
The fast DG local assembler engine for DUNE grids which creates the matrix pattern.
Definition fastdg/patternengine.hh:24
void setPattern(Pattern &pattern_)
Definition fastdg/patternengine.hh:92
The fast DG local assembler engine for DUNE grids which assembles the residual vector.
Definition fastdg/residualengine.hh:24
void setSolution(const Solution &solution_)
Definition fastdg/residualengine.hh:148
void setResidual(Residual &residual_)
Definition fastdg/residualengine.hh:140