DUNE PDELab (git)

jacobianapplyengine.hh
1#ifndef DUNE_PDELAB_GRIDOPERATOR_FASTDG_JACOBIANAPPLYENGINE_HH
2#define DUNE_PDELAB_GRIDOPERATOR_FASTDG_JACOBIANAPPLYENGINE_HH
3
5#include <dune/pdelab/gridoperator/common/assemblerutilities.hh>
6#include <dune/pdelab/gridoperator/common/localassemblerenginebase.hh>
7#include <dune/pdelab/constraints/common/constraints.hh>
8#include <dune/pdelab/localoperator/callswitch.hh>
9
10namespace Dune{
11 namespace PDELab{
12
20 template<typename LA>
23 {
24 public:
25
26 template<typename TrialConstraintsContainer, typename TestConstraintsContainer>
27 bool needsConstraintsCaching(const TrialConstraintsContainer& cu, const TestConstraintsContainer& cv) const
28 {
29 return false;
30 }
31
33 typedef LA LocalAssembler;
34
36 typedef typename LA::LocalOperator LOP;
37
39 static constexpr bool isLinear = LOP::isLinear;
40
42 typedef typename LA::Traits::Range Range;
43 typedef typename Range::ElementType RangeElement;
44
46 typedef typename LA::Traits::Domain Domain;
47 typedef typename Domain::ElementType DomainElement;
48
50 typedef typename LA::LFSU LFSU;
51 typedef typename LA::LFSUCache LFSUCache;
52 typedef typename LFSU::Traits::GridFunctionSpace GFSU;
53 typedef typename LA::LFSV LFSV;
54 typedef typename LA::LFSVCache LFSVCache;
55 typedef typename LFSV::Traits::GridFunctionSpace GFSV;
56
57 typedef typename Domain::template ConstAliasedLocalView<LFSUCache> DomainView;
58 typedef typename Range::template AliasedLocalView<LFSVCache> RangeView;
59
67 : local_assembler(local_assembler_),
68 lop(local_assembler_.localOperator())
69 {}
70
73 bool requireSkeleton() const
74 { return local_assembler.doAlphaSkeleton(); }
75 bool requireSkeletonTwoSided() const
76 { return local_assembler.doSkeletonTwoSided(); }
77 bool requireUVVolume() const
78 { return local_assembler.doAlphaVolume(); }
79 bool requireUVSkeleton() const
80 { return local_assembler.doAlphaSkeleton(); }
81 bool requireUVBoundary() const
82 { return local_assembler.doAlphaBoundary(); }
83 bool requireUVVolumePostSkeleton() const
84 { return local_assembler.doAlphaVolumePostSkeleton(); }
86
89 {
90 return local_assembler;
91 }
92
94 const typename LocalAssembler::Traits::TrialGridFunctionSpaceConstraints& trialConstraints() const
95 {
96 return localAssembler().trialConstraints();
97 }
98
100 const typename LocalAssembler::Traits::TestGridFunctionSpaceConstraints& testConstraints() const
101 {
102 return localAssembler().testConstraints();
103 }
104
107 void setSolution(const Domain & solution_)
108 {
109 if (isLinear)
110 DUNE_THROW(Dune::Exception, "In the linear case the jacobian apply does not depend on the current solution and this method should never be called.");
111 global_solution_view_inside.attach(solution_);
112 global_solution_view_outside.attach(solution_);
113 }
114
117 void setUpdate(const Domain & update_)
118 {
119 global_update_view_inside.attach(update_);
120 global_update_view_outside.attach(update_);
121 }
122
125 void setResult(Range & result_)
126 {
127 global_result_view_inside.attach(result_);
128 global_result_view_outside.attach(result_);
129 }
130
134 template<typename EG, typename LFSUC, typename LFSVC>
135 void onBindLFSUV(const EG & eg, const LFSUC & lfsu_cache, const LFSVC & lfsv_cache)
136 {
137 if (not isLinear)
138 global_solution_view_inside.bind(lfsu_cache);
139 global_update_view_inside.bind(lfsu_cache);
140 }
141
142 template<typename EG, typename LFSVC>
143 void onBindLFSV(const EG & eg, const LFSVC & lfsv_cache)
144 {
145 global_result_view_inside.bind(lfsv_cache);
146 }
147
148 template<typename IG, typename LFSUC, typename LFSVC>
149 void onBindLFSUVInside(const IG & ig, const LFSUC & lfsu_cache, const LFSVC & lfsv_cache)
150 {
151 if (not isLinear)
152 global_solution_view_inside.bind(lfsu_cache);
153 global_update_view_inside.bind(lfsu_cache);
154 }
155
156 template<typename IG, typename LFSUC, typename LFSVC>
157 void onBindLFSUVOutside(const IG & ig,
158 const LFSUC & lfsu_s_cache, const LFSVC & lfsv_s_cache,
159 const LFSUC & lfsu_n_cache, const LFSVC & lfsv_n_cache)
160 {
161 if (not isLinear)
162 global_solution_view_outside.bind(lfsu_n_cache);
163 global_update_view_outside.bind(lfsu_n_cache);
164 }
165
166 template<typename IG, typename LFSVC>
167 void onBindLFSVInside(const IG & ig, const LFSVC & lfsv_cache)
168 {
169 global_result_view_inside.bind(lfsv_cache);
170 }
171
172 template<typename IG, typename LFSVC>
173 void onBindLFSVOutside(const IG & ig,
174 const LFSVC & lfsv_s_cache,
175 const LFSVC & lfsv_n_cache)
176 {
177 global_result_view_outside.bind(lfsv_n_cache);
178 }
179
181
185 template<typename EG, typename LFSVC>
186 void onUnbindLFSV(const EG & eg, const LFSVC & lfsv_cache)
187 {
188 global_result_view_inside.commit();
189 global_result_view_inside.unbind();
190 }
191
192 template<typename IG, typename LFSVC>
193 void onUnbindLFSVInside(const IG & ig, const LFSVC & lfsv_cache)
194 {
195 global_result_view_inside.commit();
196 global_result_view_inside.unbind();
197 }
198
199 template<typename IG, typename LFSVC>
200 void onUnbindLFSVOutside(const IG & ig,
201 const LFSVC & lfsv_s_cache,
202 const LFSVC & lfsv_n_cache)
203 {
204 global_result_view_outside.commit();
205 global_result_view_outside.unbind();
206 }
208
211 template<typename LFSUC>
212 void loadCoefficientsLFSUInside(const LFSUC & lfsu_s_cache)
213 {}
214
215 template<typename LFSUC>
216 void loadCoefficientsLFSUOutside(const LFSUC & lfsu_n_cache)
217 {}
218
219 template<typename LFSUC>
220 void loadCoefficientsLFSUCoupling(const LFSUC & lfsu_c_cache)
221 {
222 DUNE_THROW(Dune::NotImplemented,"No coupling lfsu available for ");
223 }
225
228
229 void postAssembly(const GFSU& gfsu, const GFSV& gfsv)
230 {
231 if(local_assembler.doPostProcessing())
232 Dune::PDELab::constrain_residual(local_assembler.testConstraints(),
233 global_result_view_inside.container());
234 }
235
237
244 template <typename EG> bool skipEntity(const EG &eg)
245 {
246 return localAssembler().skipEntity(eg) ||
247 (LocalAssembler::isNonOverlapping &&
248 eg.entity().partitionType() != Dune::InteriorEntity);
249 }
250
257 template<typename IG>
258 bool skipIntersection(const IG & ig)
259 {
260 return localAssembler().skipIntersection(ig);
261 }
262
263 template<typename EG, typename LFSUC, typename LFSVC>
264 void assembleUVVolume(const EG & eg, const LFSUC & lfsu_cache, const LFSVC & lfsv_cache)
265 {
266 global_result_view_inside.setWeight(local_assembler.weight());
267 Dune::PDELab::LocalAssemblerCallSwitch<LOP,LOP::doAlphaVolume>::
268 jacobian_apply_volume(lop,eg,lfsu_cache.localFunctionSpace(),global_solution_view_inside,global_update_view_inside,lfsv_cache.localFunctionSpace(),global_result_view_inside);
269 }
270
271 template<typename IG, typename LFSUC, typename LFSVC>
272 void assembleUVSkeleton(const IG & ig, const LFSUC & lfsu_s_cache, const LFSVC & lfsv_s_cache,
273 const LFSUC & lfsu_n_cache, const LFSVC & lfsv_n_cache)
274 {
275 global_result_view_inside.setWeight(local_assembler.weight());
276 global_result_view_outside.setWeight(local_assembler.weight());
277 Dune::PDELab::LocalAssemblerCallSwitch<LOP,LOP::doAlphaSkeleton>::
278 jacobian_apply_skeleton(lop,ig,
279 lfsu_s_cache.localFunctionSpace(),global_solution_view_inside,global_update_view_inside,lfsv_s_cache.localFunctionSpace(),
280 lfsu_n_cache.localFunctionSpace(),global_solution_view_outside,global_update_view_outside,lfsv_n_cache.localFunctionSpace(),
281 global_result_view_inside,global_result_view_outside);
282 }
283
284 template<typename IG, typename LFSUC, typename LFSVC>
285 void assembleUVBoundary(const IG & ig, const LFSUC & lfsu_s_cache, const LFSVC & lfsv_s_cache)
286 {
287 global_result_view_inside.setWeight(local_assembler.weight());
288 Dune::PDELab::LocalAssemblerCallSwitch<LOP,LOP::doAlphaBoundary>::
289 jacobian_apply_boundary(lop,ig,lfsu_s_cache.localFunctionSpace(),global_solution_view_inside,global_update_view_inside,lfsv_s_cache.localFunctionSpace(),global_result_view_inside);
290 }
291
292 template<typename IG, typename LFSUC, typename LFSVC>
293 static void assembleUVEnrichedCoupling(const IG & ig,
294 const LFSUC & lfsu_s_cache, const LFSVC & lfsv_s_cache,
295 const LFSUC & lfsu_n_cache, const LFSVC & lfsv_n_cache,
296 const LFSUC & lfsu_coupling_cache, const LFSVC & lfsv_coupling_cache)
297 {
298 DUNE_THROW(Dune::NotImplemented,"Assembling of coupling spaces is not implemented for ");
299 }
300
301 template<typename EG, typename LFSUC, typename LFSVC>
302 void assembleUVVolumePostSkeleton(const EG & eg, const LFSUC & lfsu_cache, const LFSVC & lfsv_cache)
303 {
304 global_result_view_inside.setWeight(local_assembler.weight());
305 Dune::PDELab::LocalAssemblerCallSwitch<LOP,LOP::doAlphaVolumePostSkeleton>::
306 jacobian_apply_volume_post_skeleton(lop,eg,lfsu_cache.localFunctionSpace(),global_solution_view_inside,global_update_view_inside,lfsv_cache.localFunctionSpace(),global_result_view_inside);
307 }
308
310
311 private:
314 const LocalAssembler & local_assembler;
315
317 const LOP & lop;
318
320 DomainView global_solution_view_inside;
321 DomainView global_solution_view_outside;
322
324 DomainView global_update_view_inside;
325 DomainView global_update_view_outside;
326
328 RangeView global_result_view_inside;
329 RangeView global_result_view_outside;
330
333 typedef Dune::PDELab::TrialSpaceTag LocalTrialSpaceTag;
334 typedef Dune::PDELab::TestSpaceTag LocalTestSpaceTag;
335
338
340
341 }; // End of class FastDGLocalJacobianAssemblerEngine
342
343 }
344}
345#endif
Base class for Dune-Exceptions.
Definition: exceptions.hh:96
Default exception for dummy implementations.
Definition: exceptions.hh:263
The fast DG local assembler engine for DUNE grids which assembles the local application of the Jacobi...
Definition: jacobianapplyengine.hh:23
LA LocalAssembler
The type of the wrapping local assembler.
Definition: jacobianapplyengine.hh:33
bool skipEntity(const EG &eg)
Definition: jacobianapplyengine.hh:244
void onBindLFSUV(const EG &eg, const LFSUC &lfsu_cache, const LFSVC &lfsv_cache)
Definition: jacobianapplyengine.hh:135
LA::LocalOperator LOP
The type of the local operator.
Definition: jacobianapplyengine.hh:36
LA::LFSU LFSU
The local function spaces.
Definition: jacobianapplyengine.hh:50
void loadCoefficientsLFSUInside(const LFSUC &lfsu_s_cache)
Definition: jacobianapplyengine.hh:212
void setResult(Range &result_)
Definition: jacobianapplyengine.hh:125
bool skipIntersection(const IG &ig)
Definition: jacobianapplyengine.hh:258
void postAssembly(const GFSU &gfsu, const GFSV &gfsv)
Definition: jacobianapplyengine.hh:229
LA::Traits::Range Range
The type of the result vector.
Definition: jacobianapplyengine.hh:42
static constexpr bool isLinear
Wheter the local operator is linear.
Definition: jacobianapplyengine.hh:39
LA::Traits::Domain Domain
The type of the solution vector.
Definition: jacobianapplyengine.hh:46
void setSolution(const Domain &solution_)
Definition: jacobianapplyengine.hh:107
const LocalAssembler & localAssembler() const
Public access to the wrapping local assembler.
Definition: jacobianapplyengine.hh:88
FastDGLocalJacobianApplyAssemblerEngine(const LocalAssembler &local_assembler_)
Constructor.
Definition: jacobianapplyengine.hh:66
const LocalAssembler::Traits::TrialGridFunctionSpaceConstraints & trialConstraints() const
Trial space constraints.
Definition: jacobianapplyengine.hh:94
const LocalAssembler::Traits::TestGridFunctionSpaceConstraints & testConstraints() const
Test space constraints.
Definition: jacobianapplyengine.hh:100
bool requireSkeleton() const
Definition: jacobianapplyengine.hh:73
void setUpdate(const Domain &update_)
Definition: jacobianapplyengine.hh:117
void onUnbindLFSV(const EG &eg, const LFSVC &lfsv_cache)
Definition: jacobianapplyengine.hh:186
Base class for LocalAssemblerEngine implementations to avoid boilerplate code.
Definition: localassemblerenginebase.hh:22
#define DUNE_THROW(E, m)
Definition: exceptions.hh:218
@ InteriorEntity
all interior entities
Definition: gridenums.hh:31
void constrain_residual(const CG &cg, XG &xg)
transform residual into transformed basis: r -> r~
Definition: constraints.hh:904
Dune namespace.
Definition: alignedallocator.hh:13
Definition: localfunctionspacetags.hh:54
Definition: localfunctionspacetags.hh:48
Creative Commons License   |  Legal Statements / Impressum  |  Hosted by TU Dresden  |  generated with Hugo v0.111.3 (Nov 23, 23:29, 2024)