Go to the documentation of this file. 1 #ifndef DUNE_PDELAB_GRIDOPERATOR_FASTDG_RESIDUALENGINE_HH
2 #define DUNE_PDELAB_GRIDOPERATOR_FASTDG_RESIDUALENGINE_HH
27 template<
typename TrialConstra
intsContainer,
typename TestConstra
intsContainer>
37 typedef typename LA::LocalOperator
LOP;
40 typedef typename LA::Traits::Residual
Residual;
44 typedef typename LA::Traits::Solution
Solution;
48 typedef typename LA::LFSU
LFSU;
50 typedef typename LFSU::Traits::GridFunctionSpace
GFSU;
51 typedef typename LA::LFSV
LFSV;
53 typedef typename LFSV::Traits::GridFunctionSpace
GFSV;
55 typedef typename Solution::template ConstAliasedLocalView<LFSUCache>
SolutionView;
56 typedef typename Residual::template AliasedLocalView<LFSVCache>
ResidualView;
65 : local_assembler(local_assembler_),
66 lop(local_assembler_.localOperator())
87 local_assembler(other.local_assembler), lop(other.lop),
88 global_rl_view(other.global_rl_view),
89 global_rn_view(other.global_rn_view),
90 global_sl_view(other.global_sl_view),
91 global_sn_view(other.global_sn_view)
99 {
return ( local_assembler.doAlphaSkeleton() || local_assembler.doLambdaSkeleton() ); }
101 {
return local_assembler.doSkeletonTwoSided(); }
103 {
return local_assembler.doAlphaVolume(); }
105 {
return local_assembler.doLambdaVolume(); }
107 {
return local_assembler.doAlphaSkeleton(); }
109 {
return local_assembler.doLambdaSkeleton(); }
111 {
return local_assembler.doAlphaBoundary(); }
113 {
return local_assembler.doLambdaBoundary(); }
115 {
return local_assembler.doAlphaVolumePostSkeleton(); }
117 {
return local_assembler.doLambdaVolumePostSkeleton(); }
123 return local_assembler;
127 const typename LocalAssembler::Traits::TrialGridFunctionSpaceConstraints&
trialConstraints()
const
133 const typename LocalAssembler::Traits::TestGridFunctionSpaceConstraints&
testConstraints()
const
142 global_rl_view.attach(residual_);
143 global_rn_view.attach(residual_);
150 global_sl_view.attach(solution_);
151 global_sn_view.attach(solution_);
157 template<
typename EG,
typename LFSUC,
typename LFSVC>
158 void onBindLFSUV(
const EG & eg,
const LFSUC & lfsu_cache,
const LFSVC & lfsv_cache)
160 global_sl_view.bind(lfsu_cache);
164 template<
typename EG,
typename LFSVC>
167 global_rl_view.bind(lfsv_cache);
171 template<
typename IG,
typename LFSUC,
typename LFSVC>
174 global_sl_view.bind(lfsu_cache);
178 template<
typename IG,
typename LFSUC,
typename LFSVC>
180 const LFSUC & lfsu_s_cache,
const LFSVC & lfsv_s_cache,
181 const LFSUC & lfsu_n_cache,
const LFSVC & lfsv_n_cache)
183 global_sn_view.bind(lfsu_n_cache);
187 template<
typename IG,
typename LFSVC>
190 global_rl_view.bind(lfsv_cache);
194 template<
typename IG,
typename LFSVC>
196 const LFSVC & lfsv_s_cache,
197 const LFSVC & lfsv_n_cache)
199 global_rn_view.bind(lfsv_n_cache);
208 template<
typename EG,
typename LFSVC>
212 global_rl_view.unbind();
213 global_rl_view.commit();
216 template<
typename IG,
typename LFSVC>
220 global_rl_view.commit();
221 global_rl_view.unnbind();
224 template<
typename IG,
typename LFSVC>
226 const LFSVC & lfsv_s_cache,
227 const LFSVC & lfsv_n_cache)
230 global_rn_view.commit();
231 global_rn_view.unbind();
237 template<
typename LFSUC>
242 template<
typename LFSUC>
247 template<
typename LFSUC>
250 DUNE_THROW(Dune::NotImplemented,
"No coupling lfsu available for ");
259 if(local_assembler.doPostProcessing())
261 global_rl_view.container());
276 template<
typename EG>
279 return LocalAssembler::isNonOverlapping && eg.entity().partitionType() != Dune::InteriorEntity;
282 template<
typename EG,
typename LFSUC,
typename LFSVC>
286 global_rl_view.setWeight(local_assembler.weight());
288 alpha_volume(lop,eg,lfsu_cache.localFunctionSpace(),global_sl_view,lfsv_cache.localFunctionSpace(),global_rl_view);
291 template<
typename EG,
typename LFSVC>
295 global_rl_view.setWeight(local_assembler.weight());
297 lambda_volume(lop,eg,lfsv_cache.localFunctionSpace(),global_rl_view);
300 template<
typename IG,
typename LFSUC,
typename LFSVC>
302 const LFSUC & lfsu_n_cache,
const LFSVC & lfsv_n_cache)
306 global_rl_view.setWeight(local_assembler.weight());
307 global_rn_view.setWeight(local_assembler.weight());
310 lfsu_s_cache.localFunctionSpace(),global_sl_view,lfsv_s_cache.localFunctionSpace(),
311 lfsu_n_cache.localFunctionSpace(),global_sn_view,lfsv_n_cache.localFunctionSpace(),
312 global_rl_view,global_rn_view);
315 template<
typename IG,
typename LFSVC>
320 global_rl_view.setWeight(local_assembler.weight());
321 global_rn_view.setWeight(local_assembler.weight());
323 lambda_skeleton(lop,
ig, lfsv_s_cache.localFunctionSpace(), lfsv_n_cache.localFunctionSpace(), global_rl_view, global_rn_view);
326 template<
typename IG,
typename LFSUC,
typename LFSVC>
330 global_rl_view.setWeight(local_assembler.weight());
332 alpha_boundary(lop,
ig,lfsu_s_cache.localFunctionSpace(),global_sl_view,lfsv_s_cache.localFunctionSpace(),global_rl_view);
335 template<
typename IG,
typename LFSVC>
339 global_rl_view.setWeight(local_assembler.weight());
344 template<
typename IG,
typename LFSUC,
typename LFSVC>
346 const LFSUC & lfsu_s_cache,
const LFSVC & lfsv_s_cache,
347 const LFSUC & lfsu_n_cache,
const LFSVC & lfsv_n_cache,
348 const LFSUC & lfsu_coupling_cache,
const LFSVC & lfsv_coupling_cache)
350 DUNE_THROW(Dune::NotImplemented,
"Assembling of coupling spaces is not implemented for ");
353 template<
typename IG,
typename LFSVC>
355 const LFSVC & lfsv_s_cache,
356 const LFSVC & lfsv_n_cache,
357 const LFSVC & lfsv_coupling_cache)
359 DUNE_THROW(Dune::NotImplemented,
"Assembling of coupling spaces is not implemented for ");
362 template<
typename EG,
typename LFSUC,
typename LFSVC>
366 global_rl_view.setWeight(local_assembler.weight());
368 alpha_volume_post_skeleton(lop,eg,lfsu_cache.localFunctionSpace(),global_sl_view,lfsv_cache.localFunctionSpace(),global_rl_view);
371 template<
typename EG,
typename LFSVC>
375 global_rl_view.setWeight(local_assembler.weight());
bool requireSkeletonTwoSided() const
Definition: fastdg/residualengine.hh:100
void assembleUVSkeleton(const IG &ig, const LFSUC &lfsu_s_cache, const LFSVC &lfsv_s_cache, const LFSUC &lfsu_n_cache, const LFSVC &lfsv_n_cache)
Definition: fastdg/residualengine.hh:301
bool needsConstraintsCaching(const TrialConstraintsContainer &cu, const TestConstraintsContainer &cv) const
Definition: fastdg/residualengine.hh:28
Definition: localfunctionspacetags.hh:48
void onBindLFSVOutside(const IG &ig, const LFSVC &lfsv_s_cache, const LFSVC &lfsv_n_cache)
Definition: fastdg/residualengine.hh:195
Solution::ElementType SolutionElement
Definition: fastdg/residualengine.hh:45
LA::LFSU LFSU
The local function spaces.
Definition: fastdg/residualengine.hh:48
void assembleUVBoundary(const IG &ig, const LFSUC &lfsu_s_cache, const LFSVC &lfsv_s_cache)
Definition: fastdg/residualengine.hh:327
void setSolution(const Solution &solution_)
Definition: fastdg/residualengine.hh:148
void postAssembly(const GFSU &gfsu, const GFSV &gfsv)
Definition: fastdg/residualengine.hh:257
void assembleUVVolumePostSkeleton(const EG &eg, const LFSUC &lfsu_cache, const LFSVC &lfsv_cache)
Definition: fastdg/residualengine.hh:363
bool requireVVolumePostSkeleton() const
Definition: fastdg/residualengine.hh:116
LA::LFSUCache LFSUCache
Definition: fastdg/residualengine.hh:49
void setResidual(Residual &residual_)
Definition: fastdg/residualengine.hh:140
For backward compatibility – Do not use this!
Definition: adaptivity.hh:28
LA::Traits::Solution Solution
The type of the solution vector.
Definition: fastdg/residualengine.hh:44
void onBindLFSVInside(const IG &ig, const LFSVC &lfsv_cache)
Definition: fastdg/residualengine.hh:188
void loadCoefficientsLFSUInside(const LFSUC &lfsu_s_cache)
Definition: fastdg/residualengine.hh:238
bool requireUVVolume() const
Definition: fastdg/residualengine.hh:102
Residual::template AliasedLocalView< LFSVCache > ResidualView
Definition: fastdg/residualengine.hh:56
LA::LocalOperator LOP
The type of the local operator.
Definition: fastdg/residualengine.hh:37
bool requireVSkeleton() const
Definition: fastdg/residualengine.hh:108
void onBindLFSUVInside(const IG &ig, const LFSUC &lfsu_cache, const LFSVC &lfsv_cache)
Definition: fastdg/residualengine.hh:172
void constrain_residual(const CG &cg, XG &xg)
transform residual into transformed basis: r -> r~
Definition: constraints.hh:904
void assembleVVolumePostSkeleton(const EG &eg, const LFSVC &lfsv_cache)
Definition: fastdg/residualengine.hh:372
LA LocalAssembler
The type of the wrapping local assembler.
Definition: fastdg/residualengine.hh:34
FastDGLocalResidualAssemblerEngine(const LocalAssembler &local_assembler_)
Constructor.
Definition: fastdg/residualengine.hh:64
void assembleVVolume(const EG &eg, const LFSVC &lfsv_cache)
Definition: fastdg/residualengine.hh:292
bool requireUVSkeleton() const
Definition: fastdg/residualengine.hh:106
bool assembleCell(const EG &eg)
Definition: fastdg/residualengine.hh:277
static void assembleVEnrichedCoupling(const IG &ig, const LFSVC &lfsv_s_cache, const LFSVC &lfsv_n_cache, const LFSVC &lfsv_coupling_cache)
Definition: fastdg/residualengine.hh:354
void onBindLFSUV(const EG &eg, const LFSUC &lfsu_cache, const LFSVC &lfsv_cache)
Definition: fastdg/residualengine.hh:158
bool requireVBoundary() const
Definition: fastdg/residualengine.hh:112
LA::LFSVCache LFSVCache
Definition: fastdg/residualengine.hh:52
Solution::template ConstAliasedLocalView< LFSUCache > SolutionView
Definition: fastdg/residualengine.hh:55
Definition: localfunctionspacetags.hh:54
void onBindLFSV(const EG &eg, const LFSVC &lfsv_cache)
Definition: fastdg/residualengine.hh:165
bool requireUVBoundary() const
Definition: fastdg/residualengine.hh:110
Residual::ElementType ResidualElement
Definition: fastdg/residualengine.hh:41
void onUnbindLFSVInside(const IG &ig, const LFSVC &lfsv_cache)
Definition: fastdg/residualengine.hh:217
Impl::LocalAssemblerCallSwitchHelper< LOP, doIt > LocalAssemblerCallSwitch
Definition: callswitch.hh:344
Base class for LocalAssemblerEngine implementations to avoid boilerplate code.
Definition: localassemblerenginebase.hh:21
void onUnbindLFSV(const EG &eg, const LFSVC &lfsv_cache)
Definition: fastdg/residualengine.hh:209
The local assembler for DUNE grids.
Definition: fastdg/localassembler.hh:37
const LocalAssembler::Traits::TrialGridFunctionSpaceConstraints & trialConstraints() const
Trial space constraints.
Definition: fastdg/residualengine.hh:127
const LocalAssembler::Traits::TestGridFunctionSpaceConstraints & testConstraints() const
Test space constraints.
Definition: fastdg/residualengine.hh:133
void assembleVBoundary(const IG &ig, const LFSVC &lfsv_s_cache)
Definition: fastdg/residualengine.hh:336
LFSV::Traits::GridFunctionSpace GFSV
Definition: fastdg/residualengine.hh:53
void loadCoefficientsLFSUCoupling(const LFSUC &lfsu_c_cache)
Definition: fastdg/residualengine.hh:248
void assembleUVVolume(const EG &eg, const LFSUC &lfsu_cache, const LFSVC &lfsv_cache)
Definition: fastdg/residualengine.hh:283
void loadCoefficientsLFSUOutside(const LFSUC &lfsu_n_cache)
Definition: fastdg/residualengine.hh:243
bool requireVVolume() const
Definition: fastdg/residualengine.hh:104
const LocalAssembler & localAssembler() const
Public access to the wrapping local assembler.
Definition: fastdg/residualengine.hh:121
The fast DG local assembler engine for DUNE grids which assembles the residual vector.
Definition: fastdg/residualengine.hh:22
static void assembleUVEnrichedCoupling(const IG &ig, const LFSUC &lfsu_s_cache, const LFSVC &lfsv_s_cache, const LFSUC &lfsu_n_cache, const LFSVC &lfsv_n_cache, const LFSUC &lfsu_coupling_cache, const LFSVC &lfsv_coupling_cache)
Definition: fastdg/residualengine.hh:345
const IG & ig
Definition: constraints.hh:149
bool requireSkeleton() const
Definition: fastdg/residualengine.hh:98
LFSU::Traits::GridFunctionSpace GFSU
Definition: fastdg/residualengine.hh:50
void onBindLFSUVOutside(const IG &ig, const LFSUC &lfsu_s_cache, const LFSVC &lfsv_s_cache, const LFSUC &lfsu_n_cache, const LFSVC &lfsv_n_cache)
Definition: fastdg/residualengine.hh:179
LA::LFSV LFSV
Definition: fastdg/residualengine.hh:51
void assembleVSkeleton(const IG &ig, const LFSVC &lfsv_s_cache, const LFSVC &lfsv_n_cache)
Definition: fastdg/residualengine.hh:316
void onUnbindLFSVOutside(const IG &ig, const LFSVC &lfsv_s_cache, const LFSVC &lfsv_n_cache)
Definition: fastdg/residualengine.hh:225
bool requireUVVolumePostSkeleton() const
Definition: fastdg/residualengine.hh:114
LA::Traits::Residual Residual
The type of the residual vector.
Definition: fastdg/residualengine.hh:40