1 #ifndef DUNE_PDELAB_DEFAULT_APPLYJACOBIANENGINE_HH
2 #define DUNE_PDELAB_DEFAULT_APPLYJACOBIANENGINE_HH
4 #include <dune/common/nullptr.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 ConstLocalView<LFSUCache>
SolutionView;
56 typedef typename Residual::template LocalView<LFSVCache>
ResidualView;
65 : local_assembler(local_assembler_), lop(local_assembler_.lop),
73 {
return local_assembler.doAlphaSkeleton(); }
75 {
return local_assembler.doSkeletonTwoSided(); }
77 {
return local_assembler.doAlphaVolume(); }
79 {
return local_assembler.doAlphaSkeleton(); }
81 {
return local_assembler.doAlphaBoundary(); }
83 {
return local_assembler.doAlphaVolumePostSkeleton(); }
90 const typename LocalAssembler::Traits::TrialGridFunctionSpaceConstraints&
trialConstraints()
const
96 const typename LocalAssembler::Traits::TestGridFunctionSpaceConstraints&
testConstraints()
const
104 global_rl_view.attach(residual_);
105 global_rn_view.attach(residual_);
111 global_sl_view.attach(solution_);
112 global_sn_view.attach(solution_);
118 template<
typename EG,
typename LFSUC,
typename LFSVC>
119 void onBindLFSUV(
const EG &
eg,
const LFSUC & lfsu_cache,
const LFSVC & lfsv_cache){
120 global_sl_view.bind(lfsu_cache);
121 xl.
resize(lfsu_cache.size());
124 template<
typename EG,
typename LFSVC>
126 global_rl_view.bind(lfsv_cache);
127 rl.
assign(lfsv_cache.size(),0.0);
130 template<
typename IG,
typename LFSUC,
typename LFSVC>
132 global_sl_view.bind(lfsu_cache);
133 xl.
resize(lfsu_cache.size());
136 template<
typename IG,
typename LFSUC,
typename LFSVC>
138 const LFSUC & lfsu_s_cache,
const LFSVC & lfsv_s_cache,
139 const LFSUC & lfsu_n_cache,
const LFSVC & lfsv_n_cache)
141 global_sn_view.bind(lfsu_n_cache);
142 xn.
resize(lfsu_n_cache.size());
145 template<
typename IG,
typename LFSVC>
147 global_rl_view.bind(lfsv_cache);
148 rl.
assign(lfsv_cache.size(),0.0);
151 template<
typename IG,
typename LFSVC>
153 const LFSVC & lfsv_s_cache,
154 const LFSVC & lfsv_n_cache)
156 global_rn_view.bind(lfsv_n_cache);
157 rn.
assign(lfsv_n_cache.size(),0.0);
165 template<
typename EG,
typename LFSVC>
167 global_rl_view.add(rl);
168 global_rl_view.commit();
171 template<
typename IG,
typename LFSVC>
173 global_rl_view.add(rl);
174 global_rl_view.commit();
177 template<
typename IG,
typename LFSVC>
179 const LFSVC & lfsv_s_cache,
180 const LFSVC & lfsv_n_cache)
182 global_rn_view.add(rn);
183 global_rn_view.commit();
189 template<
typename LFSUC>
191 global_sl_view.read(xl);
193 template<
typename LFSUC>
195 global_sn_view.read(xn);
197 template<
typename LFSUC>
199 {DUNE_THROW(Dune::NotImplemented,
"No coupling lfsu available for ");}
206 if(local_assembler.doPostProcessing){
222 template<
typename EG>
225 return LocalAssembler::isNonOverlapping && eg.entity().partitionType() != Dune::InteriorEntity;
228 template<
typename EG,
typename LFSUC,
typename LFSVC>
231 rl_view.
setWeight(local_assembler.weight);
233 jacobian_apply_volume(lop,eg,lfsu_cache.localFunctionSpace(),xl,lfsv_cache.localFunctionSpace(),rl_view);
236 template<
typename IG,
typename LFSUC,
typename LFSVC>
238 const LFSUC & lfsu_n_cache,
const LFSVC & lfsv_n_cache)
240 rl_view.
setWeight(local_assembler.weight);
241 rn_view.
setWeight(local_assembler.weight);
244 lfsu_s_cache.localFunctionSpace(),xl,lfsv_s_cache.localFunctionSpace(),
245 lfsu_n_cache.localFunctionSpace(),xn,lfsv_n_cache.localFunctionSpace(),
249 template<
typename IG,
typename LFSUC,
typename LFSVC>
252 rl_view.
setWeight(local_assembler.weight);
254 jacobian_apply_boundary(lop,ig,lfsu_s_cache.localFunctionSpace(),xl,lfsv_s_cache.localFunctionSpace(),rl_view);
257 template<
typename IG,
typename LFSUC,
typename LFSVC>
259 const LFSUC & lfsu_s_cache,
const LFSVC & lfsv_s_cache,
260 const LFSUC & lfsu_n_cache,
const LFSVC & lfsv_n_cache,
261 const LFSUC & lfsu_coupling_cache,
const LFSVC & lfsv_coupling_cache)
262 {DUNE_THROW(Dune::NotImplemented,
"Assembling of coupling spaces is not implemented for ");}
264 template<
typename EG,
typename LFSUC,
typename LFSVC>
267 rl_view.
setWeight(local_assembler.weight);
277 const LocalAssembler & local_assembler;
283 ResidualView global_rl_view;
284 ResidualView global_rn_view;
287 SolutionView global_sl_view;
288 SolutionView global_sn_view;
bool assembleCell(const EG &eg)
Definition: jacobianapplyengine.hh:223
DefaultLocalJacobianApplyAssemblerEngine(const LocalAssembler &local_assembler_)
Constructor.
Definition: jacobianapplyengine.hh:64
bool requireUVBoundary() const
Definition: jacobianapplyengine.hh:80
LA::LocalOperator LOP
The type of the local operator.
Definition: jacobianapplyengine.hh:37
void resize(size_type size)
Resize the container.
Definition: localvector.hh:251
LA LocalAssembler
The type of the wrapping local assembler.
Definition: jacobianapplyengine.hh:34
void onBindLFSVOutside(const IG &ig, const LFSVC &lfsv_s_cache, const LFSVC &lfsv_n_cache)
Definition: jacobianapplyengine.hh:152
void assembleUVVolumePostSkeleton(const EG &eg, const LFSUC &lfsu_cache, const LFSVC &lfsv_cache)
Definition: jacobianapplyengine.hh:265
Residual::ElementType ResidualElement
Definition: jacobianapplyengine.hh:41
LA::LFSV LFSV
Definition: jacobianapplyengine.hh:51
Solution::template ConstLocalView< LFSUCache > SolutionView
Definition: jacobianapplyengine.hh:55
void onBindLFSUVInside(const IG &ig, const LFSUC &lfsu_cache, const LFSVC &lfsv_cache)
Definition: jacobianapplyengine.hh:131
void loadCoefficientsLFSUOutside(const LFSUC &lfsu_n_cache)
Definition: jacobianapplyengine.hh:194
bool requireUVVolume() const
Definition: jacobianapplyengine.hh:76
Base class for LocalAssemblerEngine implementations to avoid boilerplate code.
Definition: localassemblerenginebase.hh:21
Residual::template LocalView< LFSVCache > ResidualView
Definition: jacobianapplyengine.hh:56
LA::LFSU LFSU
The local function spaces.
Definition: jacobianapplyengine.hh:48
static void jacobian_apply_skeleton(const LA &la, const IG &ig, const LFSU &lfsu_s, const X &x_s, const LFSV &lfsv_s, const LFSU &lfsu_n, const X &x_n, const LFSV &lfsv_n, Y &y_s, Y &y_n)
Definition: callswitch.hh:88
LA::LFSUCache LFSUCache
Definition: jacobianapplyengine.hh:49
static void jacobian_apply_volume_post_skeleton(const LA &la, const EG &eg, const LFSU &lfsu, const X &x, const LFSV &lfsv, Y &y)
Definition: callswitch.hh:84
bool needsConstraintsCaching(const TrialConstraintsContainer &cu, const TestConstraintsContainer &cv) const
Definition: jacobianapplyengine.hh:28
const IG & ig
Definition: constraints.hh:147
void assembleUVBoundary(const IG &ig, const LFSUC &lfsu_s_cache, const LFSVC &lfsv_s_cache)
Definition: jacobianapplyengine.hh:250
void assign(size_type size, const T &value)
Resize the container to size and assign the passed value to all entries.
Definition: localvector.hh:257
void onUnbindLFSV(const EG &eg, const LFSVC &lfsv_cache)
Definition: jacobianapplyengine.hh:166
Definition: localfunctionspacetags.hh:48
Definition: localfunctionspacetags.hh:54
void onUnbindLFSVOutside(const IG &ig, const LFSVC &lfsv_s_cache, const LFSVC &lfsv_n_cache)
Definition: jacobianapplyengine.hh:178
The local assembler engine for DUNE grids which assembles the local application of the Jacobian...
Definition: jacobianapplyengine.hh:22
bool requireUVVolumePostSkeleton() const
Definition: jacobianapplyengine.hh:82
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: jacobianapplyengine.hh:137
const LocalAssembler & localAssembler() const
Public access to the wrapping local assembler.
Definition: jacobianapplyengine.hh:87
const LocalAssembler::Traits::TrialGridFunctionSpaceConstraints & trialConstraints() const
Trial space constraints.
Definition: jacobianapplyengine.hh:90
void loadCoefficientsLFSUCoupling(const LFSUC &lfsu_c_cache)
Definition: jacobianapplyengine.hh:198
LFSV::Traits::GridFunctionSpace GFSV
Definition: jacobianapplyengine.hh:53
void onBindLFSUV(const EG &eg, const LFSUC &lfsu_cache, const LFSVC &lfsv_cache)
Definition: jacobianapplyengine.hh:119
const LocalAssembler::Traits::TestGridFunctionSpaceConstraints & testConstraints() const
Test space constraints.
Definition: jacobianapplyengine.hh:96
Definition: adaptivity.hh:27
bool requireUVSkeleton() const
Definition: jacobianapplyengine.hh:78
bool requireSkeleton() const
Definition: jacobianapplyengine.hh:72
void setWeight(weight_type weight)
Resets the weighting coefficient of the view.
Definition: localvector.hh:72
LA::Traits::Residual Residual
The type of the residual vector.
Definition: jacobianapplyengine.hh:40
static void jacobian_apply_volume(const LA &la, const EG &eg, const LFSU &lfsu, const X &x, const LFSV &lfsv, Y &y)
Definition: callswitch.hh:80
LFSU::Traits::GridFunctionSpace GFSU
Definition: jacobianapplyengine.hh:50
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: jacobianapplyengine.hh:237
void setSolution(const Solution &solution_)
Definition: jacobianapplyengine.hh:110
bool requireSkeletonTwoSided() const
Definition: jacobianapplyengine.hh:74
void postAssembly(const GFSU &gfsu, const GFSV &gfsv)
Definition: jacobianapplyengine.hh:205
WeightedVectorAccumulationView< LocalVector > WeightedAccumulationView
An accumulate-only view of this container that automatically applies a weight to all contributions...
Definition: localvector.hh:198
static void jacobian_apply_boundary(const LA &la, const IG &ig, const LFSU &lfsu_s, const X &x_s, const LFSV &lfsv_s, Y &y_s)
Definition: callswitch.hh:95
Solution::ElementType SolutionElement
Definition: jacobianapplyengine.hh:45
void onUnbindLFSVInside(const IG &ig, const LFSVC &lfsv_cache)
Definition: jacobianapplyengine.hh:172
void constrain_residual(const CG &cg, XG &xg)
transform residual into transformed basis: r -> r~
Definition: constraints.hh:911
void onBindLFSV(const EG &eg, const LFSVC &lfsv_cache)
Definition: jacobianapplyengine.hh:125
void loadCoefficientsLFSUInside(const LFSUC &lfsu_s_cache)
Definition: jacobianapplyengine.hh:190
LA::Traits::Solution Solution
The type of the solution vector.
Definition: jacobianapplyengine.hh:44
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: jacobianapplyengine.hh:258
void setResidual(Residual &residual_)
Definition: jacobianapplyengine.hh:103
const EG & eg
Definition: constraints.hh:280
LA::LFSVCache LFSVCache
Definition: jacobianapplyengine.hh:52
void assembleUVVolume(const EG &eg, const LFSUC &lfsu_cache, const LFSVC &lfsv_cache)
Definition: jacobianapplyengine.hh:229
void onBindLFSVInside(const IG &ig, const LFSVC &lfsv_cache)
Definition: jacobianapplyengine.hh:146