dune-pdelab  2.5-dev
default/jacobianapplyengine.hh
Go to the documentation of this file.
1 #ifndef DUNE_PDELAB_GRIDOPERATOR_DEFAULT_JACOBIANAPPLYENGINE_HH
2 #define DUNE_PDELAB_GRIDOPERATOR_DEFAULT_JACOBIANAPPLYENGINE_HH
3 
9 
10 namespace 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  typedef typename LA::Traits::Residual Residual;
40  typedef typename Residual::ElementType ResidualElement;
41 
43  typedef typename LA::Traits::Solution Solution;
44  typedef typename Solution::ElementType SolutionElement;
45 
47  typedef typename LA::LFSU LFSU;
48  typedef typename LA::LFSUCache LFSUCache;
49  typedef typename LFSU::Traits::GridFunctionSpace GFSU;
50  typedef typename LA::LFSV LFSV;
51  typedef typename LA::LFSVCache LFSVCache;
52  typedef typename LFSV::Traits::GridFunctionSpace GFSV;
53 
54  typedef typename Solution::template ConstLocalView<LFSUCache> SolutionView;
55  typedef typename Residual::template LocalView<LFSVCache> ResidualView;
56 
64  : local_assembler(local_assembler_),
65  lop(local_assembler_.localOperator()),
66  rl_view(rl,1.0),
67  rn_view(rn,1.0)
68  {}
69 
72  bool requireSkeleton() const
73  { return local_assembler.doAlphaSkeleton(); }
75  { return local_assembler.doSkeletonTwoSided(); }
76  bool requireUVVolume() const
77  { return local_assembler.doAlphaVolume(); }
78  bool requireUVSkeleton() const
79  { return local_assembler.doAlphaSkeleton(); }
80  bool requireUVBoundary() const
81  { return local_assembler.doAlphaBoundary(); }
83  { return local_assembler.doAlphaVolumePostSkeleton(); }
85 
88  {
89  return local_assembler;
90  }
91 
93  const typename LocalAssembler::Traits::TrialGridFunctionSpaceConstraints& trialConstraints() const
94  {
95  return localAssembler().trialConstraints();
96  }
97 
99  const typename LocalAssembler::Traits::TestGridFunctionSpaceConstraints& testConstraints() const
100  {
101  return localAssembler().testConstraints();
102  }
103 
106  void setResidual(Residual & residual_)
107  {
108  global_rl_view.attach(residual_);
109  global_rn_view.attach(residual_);
110  }
111 
114  void setSolution(const Solution & solution_)
115  {
116  global_sl_view.attach(solution_);
117  global_sn_view.attach(solution_);
118  }
119 
123  template<typename EG, typename LFSUC, typename LFSVC>
124  void onBindLFSUV(const EG & eg, const LFSUC & lfsu_cache, const LFSVC & lfsv_cache)
125  {
126  global_sl_view.bind(lfsu_cache);
127  xl.resize(lfsu_cache.size());
128  }
129 
130  template<typename EG, typename LFSVC>
131  void onBindLFSV(const EG & eg, const LFSVC & lfsv_cache)
132  {
133  global_rl_view.bind(lfsv_cache);
134  rl.assign(lfsv_cache.size(),0.0);
135  }
136 
137  template<typename IG, typename LFSUC, typename LFSVC>
138  void onBindLFSUVInside(const IG & ig, const LFSUC & lfsu_cache, const LFSVC & lfsv_cache)
139  {
140  global_sl_view.bind(lfsu_cache);
141  xl.resize(lfsu_cache.size());
142  }
143 
144  template<typename IG, typename LFSUC, typename LFSVC>
145  void onBindLFSUVOutside(const IG & ig,
146  const LFSUC & lfsu_s_cache, const LFSVC & lfsv_s_cache,
147  const LFSUC & lfsu_n_cache, const LFSVC & lfsv_n_cache)
148  {
149  global_sn_view.bind(lfsu_n_cache);
150  xn.resize(lfsu_n_cache.size());
151  }
152 
153  template<typename IG, typename LFSVC>
154  void onBindLFSVInside(const IG & ig, const LFSVC & lfsv_cache)
155  {
156  global_rl_view.bind(lfsv_cache);
157  rl.assign(lfsv_cache.size(),0.0);
158  }
159 
160  template<typename IG, typename LFSVC>
161  void onBindLFSVOutside(const IG & ig,
162  const LFSVC & lfsv_s_cache,
163  const LFSVC & lfsv_n_cache)
164  {
165  global_rn_view.bind(lfsv_n_cache);
166  rn.assign(lfsv_n_cache.size(),0.0);
167  }
168 
170 
174  template<typename EG, typename LFSVC>
175  void onUnbindLFSV(const EG & eg, const LFSVC & lfsv_cache)
176  {
177  global_rl_view.add(rl);
178  global_rl_view.commit();
179  }
180 
181  template<typename IG, typename LFSVC>
182  void onUnbindLFSVInside(const IG & ig, const LFSVC & lfsv_cache)
183  {
184  global_rl_view.add(rl);
185  global_rl_view.commit();
186  }
187 
188  template<typename IG, typename LFSVC>
189  void onUnbindLFSVOutside(const IG & ig,
190  const LFSVC & lfsv_s_cache,
191  const LFSVC & lfsv_n_cache)
192  {
193  global_rn_view.add(rn);
194  global_rn_view.commit();
195  }
197 
200  template<typename LFSUC>
201  void loadCoefficientsLFSUInside(const LFSUC & lfsu_s_cache)
202  {
203  global_sl_view.read(xl);
204  }
205  template<typename LFSUC>
206  void loadCoefficientsLFSUOutside(const LFSUC & lfsu_n_cache)
207  {
208  global_sn_view.read(xn);
209  }
210  template<typename LFSUC>
211  void loadCoefficientsLFSUCoupling(const LFSUC & lfsu_c_cache)
212  {
213  DUNE_THROW(Dune::NotImplemented,"No coupling lfsu available for ");
214  }
216 
219 
220  void postAssembly(const GFSU& gfsu, const GFSV& gfsv)
221  {
222  if(local_assembler.doPostProcessing())
223  Dune::PDELab::constrain_residual(local_assembler.testConstraints(),
224  global_rl_view.container());
225 
226  }
227 
229 
232 
239  template<typename EG>
240  bool assembleCell(const EG & eg)
241  {
242  return LocalAssembler::isNonOverlapping && eg.entity().partitionType() != Dune::InteriorEntity;
243  }
244 
245  template<typename EG, typename LFSUC, typename LFSVC>
246  void assembleUVVolume(const EG & eg, const LFSUC & lfsu_cache, const LFSVC & lfsv_cache)
247  {
248  rl_view.setWeight(local_assembler.weight());
250  jacobian_apply_volume(lop,eg,lfsu_cache.localFunctionSpace(),xl,lfsv_cache.localFunctionSpace(),rl_view);
251  }
252 
253  template<typename IG, typename LFSUC, typename LFSVC>
254  void assembleUVSkeleton(const IG & ig, const LFSUC & lfsu_s_cache, const LFSVC & lfsv_s_cache,
255  const LFSUC & lfsu_n_cache, const LFSVC & lfsv_n_cache)
256  {
257  rl_view.setWeight(local_assembler.weight());
258  rn_view.setWeight(local_assembler.weight());
261  lfsu_s_cache.localFunctionSpace(),xl,lfsv_s_cache.localFunctionSpace(),
262  lfsu_n_cache.localFunctionSpace(),xn,lfsv_n_cache.localFunctionSpace(),
263  rl_view,rn_view);
264  }
265 
266  template<typename IG, typename LFSUC, typename LFSVC>
267  void assembleUVBoundary(const IG & ig, const LFSUC & lfsu_s_cache, const LFSVC & lfsv_s_cache)
268  {
269  rl_view.setWeight(local_assembler.weight());
271  jacobian_apply_boundary(lop,ig,lfsu_s_cache.localFunctionSpace(),xl,lfsv_s_cache.localFunctionSpace(),rl_view);
272  }
273 
274  template<typename IG, typename LFSUC, typename LFSVC>
275  static void assembleUVEnrichedCoupling(const IG & ig,
276  const LFSUC & lfsu_s_cache, const LFSVC & lfsv_s_cache,
277  const LFSUC & lfsu_n_cache, const LFSVC & lfsv_n_cache,
278  const LFSUC & lfsu_coupling_cache, const LFSVC & lfsv_coupling_cache)
279  {
280  DUNE_THROW(Dune::NotImplemented,"Assembling of coupling spaces is not implemented for ");
281  }
282 
283  template<typename EG, typename LFSUC, typename LFSVC>
284  void assembleUVVolumePostSkeleton(const EG & eg, const LFSUC & lfsu_cache, const LFSVC & lfsv_cache)
285  {
286  rl_view.setWeight(local_assembler.weight());
288  jacobian_apply_volume_post_skeleton(lop,eg,lfsu_cache.localFunctionSpace(),xl,lfsv_cache.localFunctionSpace(),rl_view);
289  }
290 
292 
293  private:
296  const LocalAssembler & local_assembler;
297 
299  const LOP & lop;
300 
302  ResidualView global_rl_view;
303  ResidualView global_rn_view;
304 
306  SolutionView global_sl_view;
307  SolutionView global_sn_view;
308 
311  typedef Dune::PDELab::TrialSpaceTag LocalTrialSpaceTag;
312  typedef Dune::PDELab::TestSpaceTag LocalTestSpaceTag;
313 
316 
318  SolutionVector xl;
320  SolutionVector xn;
322  ResidualVector rl;
324  ResidualVector rn;
330 
331  }; // End of class DefaultLocalJacobianAssemblerEngine
332 
333  }
334 }
335 #endif // DUNE_PDELAB_GRIDOPERATOR_DEFAULT_JACOBIANAPPLYENGINE_HH
Dune::PDELab::TestSpaceTag
Definition: localfunctionspacetags.hh:54
assemblerutilities.hh
Dune::PDELab::DefaultLocalJacobianApplyAssemblerEngine::assembleUVBoundary
void assembleUVBoundary(const IG &ig, const LFSUC &lfsu_s_cache, const LFSVC &lfsv_s_cache)
Definition: default/jacobianapplyengine.hh:267
Dune::PDELab::DefaultLocalJacobianApplyAssemblerEngine::SolutionView
Solution::template ConstLocalView< LFSUCache > SolutionView
Definition: default/jacobianapplyengine.hh:54
Dune::PDELab::DefaultLocalJacobianApplyAssemblerEngine::onBindLFSVInside
void onBindLFSVInside(const IG &ig, const LFSVC &lfsv_cache)
Definition: default/jacobianapplyengine.hh:154
Dune::PDELab::DefaultLocalJacobianApplyAssemblerEngine::onBindLFSV
void onBindLFSV(const EG &eg, const LFSVC &lfsv_cache)
Definition: default/jacobianapplyengine.hh:131
Dune::PDELab::DefaultLocalJacobianApplyAssemblerEngine::LOP
LA::LocalOperator LOP
The type of the local operator.
Definition: default/jacobianapplyengine.hh:36
Dune::PDELab::DefaultLocalJacobianApplyAssemblerEngine::onBindLFSUVOutside
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: default/jacobianapplyengine.hh:145
Dune::PDELab::WeightedVectorAccumulationView::setWeight
void setWeight(weight_type weight)
Resets the weighting coefficient of the view.
Definition: localvector.hh:72
Dune::PDELab::LocalVector< SolutionElement, LocalTrialSpaceTag >
callswitch.hh
Dune::PDELab::DefaultLocalJacobianApplyAssemblerEngine::requireUVVolumePostSkeleton
bool requireUVVolumePostSkeleton() const
Definition: default/jacobianapplyengine.hh:82
Dune
For backward compatibility – Do not use this!
Definition: adaptivity.hh:28
Dune::PDELab::DefaultLocalJacobianApplyAssemblerEngine::GFSU
LFSU::Traits::GridFunctionSpace GFSU
Definition: default/jacobianapplyengine.hh:49
Dune::PDELab::DefaultLocalJacobianApplyAssemblerEngine::requireUVVolume
bool requireUVVolume() const
Definition: default/jacobianapplyengine.hh:76
Dune::PDELab::DefaultLocalJacobianApplyAssemblerEngine::testConstraints
const LocalAssembler::Traits::TestGridFunctionSpaceConstraints & testConstraints() const
Test space constraints.
Definition: default/jacobianapplyengine.hh:99
Dune::PDELab::DefaultLocalJacobianApplyAssemblerEngine::postAssembly
void postAssembly(const GFSU &gfsu, const GFSV &gfsv)
Definition: default/jacobianapplyengine.hh:220
Dune::PDELab::DefaultLocalJacobianApplyAssemblerEngine::setResidual
void setResidual(Residual &residual_)
Definition: default/jacobianapplyengine.hh:106
Dune::PDELab::DefaultLocalJacobianApplyAssemblerEngine::onBindLFSUV
void onBindLFSUV(const EG &eg, const LFSUC &lfsu_cache, const LFSVC &lfsv_cache)
Definition: default/jacobianapplyengine.hh:124
Dune::PDELab::DefaultLocalJacobianApplyAssemblerEngine
The local assembler engine for DUNE grids which assembles the local application of the Jacobian.
Definition: default/jacobianapplyengine.hh:21
Dune::PDELab::DefaultLocalJacobianApplyAssemblerEngine::LFSVCache
LA::LFSVCache LFSVCache
Definition: default/jacobianapplyengine.hh:51
Dune::PDELab::DefaultLocalJacobianApplyAssemblerEngine::onUnbindLFSVOutside
void onUnbindLFSVOutside(const IG &ig, const LFSVC &lfsv_s_cache, const LFSVC &lfsv_n_cache)
Definition: default/jacobianapplyengine.hh:189
Dune::PDELab::DefaultLocalJacobianApplyAssemblerEngine::trialConstraints
const LocalAssembler::Traits::TrialGridFunctionSpaceConstraints & trialConstraints() const
Trial space constraints.
Definition: default/jacobianapplyengine.hh:93
Dune::PDELab::DefaultLocalJacobianApplyAssemblerEngine::loadCoefficientsLFSUCoupling
void loadCoefficientsLFSUCoupling(const LFSUC &lfsu_c_cache)
Definition: default/jacobianapplyengine.hh:211
Dune::PDELab::LocalAssemblerCallSwitch::jacobian_apply_volume_post_skeleton
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:136
Dune::PDELab::DefaultLocalJacobianApplyAssemblerEngine::requireSkeleton
bool requireSkeleton() const
Definition: default/jacobianapplyengine.hh:72
Dune::PDELab::DefaultLocalJacobianApplyAssemblerEngine::assembleUVVolumePostSkeleton
void assembleUVVolumePostSkeleton(const EG &eg, const LFSUC &lfsu_cache, const LFSVC &lfsv_cache)
Definition: default/jacobianapplyengine.hh:284
Dune::PDELab::LocalVector::resize
void resize(size_type size)
Resize the container.
Definition: localvector.hh:257
Dune::PDELab::DefaultLocalJacobianApplyAssemblerEngine::localAssembler
const LocalAssembler & localAssembler() const
Public access to the wrapping local assembler.
Definition: default/jacobianapplyengine.hh:87
Dune::PDELab::LocalAssemblerCallSwitch::jacobian_apply_boundary
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:147
Dune::PDELab::DefaultLocalJacobianApplyAssemblerEngine::Solution
LA::Traits::Solution Solution
The type of the solution vector.
Definition: default/jacobianapplyengine.hh:43
Dune::PDELab::DefaultLocalJacobianApplyAssemblerEngine::GFSV
LFSV::Traits::GridFunctionSpace GFSV
Definition: default/jacobianapplyengine.hh:52
localvector.hh
ig
const IG & ig
Definition: constraints.hh:148
Dune::PDELab::DefaultLocalJacobianApplyAssemblerEngine::assembleUVEnrichedCoupling
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: default/jacobianapplyengine.hh:275
Dune::PDELab::DefaultLocalJacobianApplyAssemblerEngine::LFSUCache
LA::LFSUCache LFSUCache
Definition: default/jacobianapplyengine.hh:48
Dune::PDELab::DefaultLocalJacobianApplyAssemblerEngine::SolutionElement
Solution::ElementType SolutionElement
Definition: default/jacobianapplyengine.hh:44
Dune::PDELab::LocalAssemblerCallSwitch::jacobian_apply_volume
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:132
Dune::PDELab::constrain_residual
void constrain_residual(const CG &cg, XG &xg)
transform residual into transformed basis: r -> r~
Definition: constraints.hh:906
Dune::PDELab::DefaultLocalJacobianApplyAssemblerEngine::requireUVSkeleton
bool requireUVSkeleton() const
Definition: default/jacobianapplyengine.hh:78
Dune::PDELab::TrialSpaceTag
Definition: localfunctionspacetags.hh:48
Dune::PDELab::DefaultLocalJacobianApplyAssemblerEngine::onBindLFSVOutside
void onBindLFSVOutside(const IG &ig, const LFSVC &lfsv_s_cache, const LFSVC &lfsv_n_cache)
Definition: default/jacobianapplyengine.hh:161
Dune::PDELab::DefaultLocalJacobianApplyAssemblerEngine::Residual
LA::Traits::Residual Residual
The type of the residual vector.
Definition: default/jacobianapplyengine.hh:39
Dune::PDELab::DefaultLocalJacobianApplyAssemblerEngine::onBindLFSUVInside
void onBindLFSUVInside(const IG &ig, const LFSUC &lfsu_cache, const LFSVC &lfsv_cache)
Definition: default/jacobianapplyengine.hh:138
Dune::PDELab::DefaultLocalJacobianApplyAssemblerEngine::assembleUVVolume
void assembleUVVolume(const EG &eg, const LFSUC &lfsu_cache, const LFSVC &lfsv_cache)
Definition: default/jacobianapplyengine.hh:246
Dune::PDELab::DefaultLocalJacobianApplyAssemblerEngine::requireSkeletonTwoSided
bool requireSkeletonTwoSided() const
Definition: default/jacobianapplyengine.hh:74
localassemblerenginebase.hh
Dune::PDELab::DefaultLocalJacobianApplyAssemblerEngine::setSolution
void setSolution(const Solution &solution_)
Definition: default/jacobianapplyengine.hh:114
Dune::PDELab::DefaultLocalJacobianApplyAssemblerEngine::ResidualView
Residual::template LocalView< LFSVCache > ResidualView
Definition: default/jacobianapplyengine.hh:55
constraints.hh
Dune::PDELab::LocalAssemblerEngineBase
Base class for LocalAssemblerEngine implementations to avoid boilerplate code.
Definition: localassemblerenginebase.hh:21
Dune::PDELab::DefaultLocalJacobianApplyAssemblerEngine::loadCoefficientsLFSUOutside
void loadCoefficientsLFSUOutside(const LFSUC &lfsu_n_cache)
Definition: default/jacobianapplyengine.hh:206
Dune::PDELab::DefaultLocalJacobianApplyAssemblerEngine::LFSU
LA::LFSU LFSU
The local function spaces.
Definition: default/jacobianapplyengine.hh:47
Dune::PDELab::DefaultLocalJacobianApplyAssemblerEngine::LFSV
LA::LFSV LFSV
Definition: default/jacobianapplyengine.hh:50
Dune::PDELab::DefaultLocalJacobianApplyAssemblerEngine::assembleCell
bool assembleCell(const EG &eg)
Definition: default/jacobianapplyengine.hh:240
Dune::PDELab::DefaultLocalJacobianApplyAssemblerEngine::loadCoefficientsLFSUInside
void loadCoefficientsLFSUInside(const LFSUC &lfsu_s_cache)
Definition: default/jacobianapplyengine.hh:201
Dune::PDELab::LocalVector::assign
void assign(size_type size, const T &value)
Resize the container to size and assign the passed value to all entries.
Definition: localvector.hh:263
Dune::PDELab::DefaultLocalJacobianApplyAssemblerEngine::requireUVBoundary
bool requireUVBoundary() const
Definition: default/jacobianapplyengine.hh:80
Dune::PDELab::DefaultLocalJacobianApplyAssemblerEngine::DefaultLocalJacobianApplyAssemblerEngine
DefaultLocalJacobianApplyAssemblerEngine(const LocalAssembler &local_assembler_)
Constructor.
Definition: default/jacobianapplyengine.hh:63
Dune::PDELab::DefaultLocalJacobianApplyAssemblerEngine::LocalAssembler
LA LocalAssembler
The type of the wrapping local assembler.
Definition: default/jacobianapplyengine.hh:33
Dune::PDELab::DefaultLocalJacobianApplyAssemblerEngine::needsConstraintsCaching
bool needsConstraintsCaching(const TrialConstraintsContainer &cu, const TestConstraintsContainer &cv) const
Definition: default/jacobianapplyengine.hh:27
Dune::PDELab::LocalVector< ResidualElement, LocalTestSpaceTag >::WeightedAccumulationView
WeightedVectorAccumulationView< LocalVector > WeightedAccumulationView
An accumulate-only view of this container that automatically applies a weight to all contributions.
Definition: localvector.hh:198
Dune::PDELab::DefaultLocalJacobianApplyAssemblerEngine::assembleUVSkeleton
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: default/jacobianapplyengine.hh:254
Dune::PDELab::DefaultLocalJacobianApplyAssemblerEngine::onUnbindLFSVInside
void onUnbindLFSVInside(const IG &ig, const LFSVC &lfsv_cache)
Definition: default/jacobianapplyengine.hh:182
Dune::PDELab::DefaultLocalJacobianApplyAssemblerEngine::ResidualElement
Residual::ElementType ResidualElement
Definition: default/jacobianapplyengine.hh:40
Dune::PDELab::LocalAssemblerCallSwitch::jacobian_apply_skeleton
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:140
Dune::PDELab::DefaultLocalJacobianApplyAssemblerEngine::onUnbindLFSV
void onUnbindLFSV(const EG &eg, const LFSVC &lfsv_cache)
Definition: default/jacobianapplyengine.hh:175