dune-pdelab  2.4-dev
default/residualengine.hh
Go to the documentation of this file.
1 #ifndef DUNE_PDELAB_DEFAULT_RESIDUALENGINE_HH
2 #define DUNE_PDELAB_DEFAULT_RESIDUALENGINE_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 
63  DefaultLocalResidualAssemblerEngine(const LocalAssembler & local_assembler_)
64  : local_assembler(local_assembler_), lop(local_assembler_.lop),
65  rl_view(rl,1.0),
66  rn_view(rn,1.0)
67  {}
68 
71  bool requireSkeleton() const
72  { return ( local_assembler.doAlphaSkeleton() || local_assembler.doLambdaSkeleton() ); }
74  { return local_assembler.doSkeletonTwoSided(); }
75  bool requireUVVolume() const
76  { return local_assembler.doAlphaVolume(); }
77  bool requireVVolume() const
78  { return local_assembler.doLambdaVolume(); }
79  bool requireUVSkeleton() const
80  { return local_assembler.doAlphaSkeleton(); }
81  bool requireVSkeleton() const
82  { return local_assembler.doLambdaSkeleton(); }
83  bool requireUVBoundary() const
84  { return local_assembler.doAlphaBoundary(); }
85  bool requireVBoundary() const
86  { return local_assembler.doLambdaBoundary(); }
88  { return local_assembler.doAlphaVolumePostSkeleton(); }
90  { return local_assembler.doLambdaVolumePostSkeleton(); }
92 
94  const LocalAssembler & localAssembler() const { return local_assembler; }
95 
97  const typename LocalAssembler::Traits::TrialGridFunctionSpaceConstraints& trialConstraints() const
98  {
99  return localAssembler().trialConstraints();
100  }
101 
103  const typename LocalAssembler::Traits::TestGridFunctionSpaceConstraints& testConstraints() const
104  {
105  return localAssembler().testConstraints();
106  }
107 
110  void setResidual(Residual & residual_){
111  global_rl_view.attach(residual_);
112  global_rn_view.attach(residual_);
113  }
114 
117  void setSolution(const Solution & solution_){
118  global_sl_view.attach(solution_);
119  global_sn_view.attach(solution_);
120  }
121 
125  template<typename EG, typename LFSUC, typename LFSVC>
126  void onBindLFSUV(const EG & eg, const LFSUC & lfsu_cache, const LFSVC & lfsv_cache){
127  global_sl_view.bind(lfsu_cache);
128  xl.resize(lfsu_cache.size());
129  }
130 
131  template<typename EG, typename LFSVC>
132  void onBindLFSV(const EG & eg, const LFSVC & lfsv_cache){
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  global_sl_view.bind(lfsu_cache);
140  xl.resize(lfsu_cache.size());
141  }
142 
143  template<typename IG, typename LFSUC, typename LFSVC>
144  void onBindLFSUVOutside(const IG & ig,
145  const LFSUC & lfsu_s_cache, const LFSVC & lfsv_s_cache,
146  const LFSUC & lfsu_n_cache, const LFSVC & lfsv_n_cache)
147  {
148  global_sn_view.bind(lfsu_n_cache);
149  xn.resize(lfsu_n_cache.size());
150  }
151 
152  template<typename IG, typename LFSVC>
153  void onBindLFSVInside(const IG & ig, const LFSVC & lfsv_cache){
154  global_rl_view.bind(lfsv_cache);
155  rl.assign(lfsv_cache.size(),0.0);
156  }
157 
158  template<typename IG, typename LFSVC>
159  void onBindLFSVOutside(const IG & ig,
160  const LFSVC & lfsv_s_cache,
161  const LFSVC & lfsv_n_cache)
162  {
163  global_rn_view.bind(lfsv_n_cache);
164  rn.assign(lfsv_n_cache.size(),0.0);
165  }
166 
168 
172  template<typename EG, typename LFSVC>
173  void onUnbindLFSV(const EG & eg, const LFSVC & lfsv_cache){
174  global_rl_view.add(rl);
175  global_rl_view.commit();
176  }
177 
178  template<typename IG, typename LFSVC>
179  void onUnbindLFSVInside(const IG & ig, const LFSVC & lfsv_cache){
180  global_rl_view.add(rl);
181  global_rl_view.commit();
182  }
183 
184  template<typename IG, typename LFSVC>
185  void onUnbindLFSVOutside(const IG & ig,
186  const LFSVC & lfsv_s_cache,
187  const LFSVC & lfsv_n_cache)
188  {
189  global_rn_view.add(rn);
190  global_rn_view.commit();
191  }
193 
196  template<typename LFSUC>
197  void loadCoefficientsLFSUInside(const LFSUC & lfsu_s_cache){
198  global_sl_view.read(xl);
199  }
200  template<typename LFSUC>
201  void loadCoefficientsLFSUOutside(const LFSUC & lfsu_n_cache){
202  global_sn_view.read(xn);
203  }
204  template<typename LFSUC>
205  void loadCoefficientsLFSUCoupling(const LFSUC & lfsu_c_cache)
206  {DUNE_THROW(Dune::NotImplemented,"No coupling lfsu available for ");}
208 
211 
212  void postAssembly(const GFSU& gfsu, const GFSV& gfsv){
213  if(local_assembler.doPostProcessing)
214  {
215  Dune::PDELab::constrain_residual(*(local_assembler.pconstraintsv),global_rl_view.container());
216  }
217  }
218 
220 
223 
230  template<typename EG>
231  bool assembleCell(const EG & eg)
232  {
233  return LocalAssembler::isNonOverlapping && eg.entity().partitionType() != Dune::InteriorEntity;
234  }
235 
236  template<typename EG, typename LFSUC, typename LFSVC>
237  void assembleUVVolume(const EG & eg, const LFSUC & lfsu_cache, const LFSVC & lfsv_cache)
238  {
239  rl_view.setWeight(local_assembler.weight);
241  alpha_volume(lop,eg,lfsu_cache.localFunctionSpace(),xl,lfsv_cache.localFunctionSpace(),rl_view);
242  }
243 
244  template<typename EG, typename LFSVC>
245  void assembleVVolume(const EG & eg, const LFSVC & lfsv_cache)
246  {
247  rl_view.setWeight(local_assembler.weight);
249  lambda_volume(lop,eg,lfsv_cache.localFunctionSpace(),rl_view);
250  }
251 
252  template<typename IG, typename LFSUC, typename LFSVC>
253  void assembleUVSkeleton(const IG & ig, const LFSUC & lfsu_s_cache, const LFSVC & lfsv_s_cache,
254  const LFSUC & lfsu_n_cache, const LFSVC & lfsv_n_cache)
255  {
256  rl_view.setWeight(local_assembler.weight);
257  rn_view.setWeight(local_assembler.weight);
259  alpha_skeleton(lop,ig,
260  lfsu_s_cache.localFunctionSpace(),xl,lfsv_s_cache.localFunctionSpace(),
261  lfsu_n_cache.localFunctionSpace(),xn,lfsv_n_cache.localFunctionSpace(),
262  rl_view,rn_view);
263  }
264 
265  template<typename IG, typename LFSVC>
266  void assembleVSkeleton(const IG & ig, const LFSVC & lfsv_s_cache, const LFSVC & lfsv_n_cache)
267  {
268  rl_view.setWeight(local_assembler.weight);
269  rn_view.setWeight(local_assembler.weight);
271  lambda_skeleton(lop, ig, lfsv_s_cache.localFunctionSpace(), lfsv_n_cache.localFunctionSpace(), rl_view, rn_view);
272  }
273 
274  template<typename IG, typename LFSUC, typename LFSVC>
275  void assembleUVBoundary(const IG & ig, const LFSUC & lfsu_s_cache, const LFSVC & lfsv_s_cache)
276  {
277  rl_view.setWeight(local_assembler.weight);
279  alpha_boundary(lop,ig,lfsu_s_cache.localFunctionSpace(),xl,lfsv_s_cache.localFunctionSpace(),rl_view);
280  }
281 
282  template<typename IG, typename LFSVC>
283  void assembleVBoundary(const IG & ig, const LFSVC & lfsv_s_cache)
284  {
285  rl_view.setWeight(local_assembler.weight);
287  lambda_boundary(lop,ig,lfsv_s_cache.localFunctionSpace(),rl_view);
288  }
289 
290  template<typename IG, typename LFSUC, typename LFSVC>
291  static void assembleUVEnrichedCoupling(const IG & ig,
292  const LFSUC & lfsu_s_cache, const LFSVC & lfsv_s_cache,
293  const LFSUC & lfsu_n_cache, const LFSVC & lfsv_n_cache,
294  const LFSUC & lfsu_coupling_cache, const LFSVC & lfsv_coupling_cache)
295  {DUNE_THROW(Dune::NotImplemented,"Assembling of coupling spaces is not implemented for ");}
296 
297  template<typename IG, typename LFSVC>
298  static void assembleVEnrichedCoupling(const IG & ig,
299  const LFSVC & lfsv_s_cache,
300  const LFSVC & lfsv_n_cache,
301  const LFSVC & lfsv_coupling_cache)
302  {DUNE_THROW(Dune::NotImplemented,"Assembling of coupling spaces is not implemented for ");}
303 
304  template<typename EG, typename LFSUC, typename LFSVC>
305  void assembleUVVolumePostSkeleton(const EG & eg, const LFSUC & lfsu_cache, const LFSVC & lfsv_cache)
306  {
307  rl_view.setWeight(local_assembler.weight);
309  alpha_volume_post_skeleton(lop,eg,lfsu_cache.localFunctionSpace(),xl,lfsv_cache.localFunctionSpace(),rl_view);
310  }
311 
312  template<typename EG, typename LFSVC>
313  void assembleVVolumePostSkeleton(const EG & eg, const LFSVC & lfsv_cache)
314  {
315  rl_view.setWeight(local_assembler.weight);
317  lambda_volume_post_skeleton(lop,eg,lfsv_cache.localFunctionSpace(),rl_view);
318  }
319 
321 
322  private:
325  const LocalAssembler & local_assembler;
326 
328  const LOP & lop;
329 
331  ResidualView global_rl_view;
332  ResidualView global_rn_view;
333 
335  SolutionView global_sl_view;
336  SolutionView global_sn_view;
337 
340  typedef Dune::PDELab::TrialSpaceTag LocalTrialSpaceTag;
341  typedef Dune::PDELab::TestSpaceTag LocalTestSpaceTag;
342 
345 
347  SolutionVector xl;
349  SolutionVector xn;
351  ResidualVector rl;
353  ResidualVector rn;
359 
360  }; // End of class DefaultLocalResidualAssemblerEngine
361 
362  }
363 }
364 #endif
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/residualengine.hh:291
LA::LocalOperator LOP
The type of the local operator.
Definition: default/residualengine.hh:36
void assembleVBoundary(const IG &ig, const LFSVC &lfsv_s_cache)
Definition: default/residualengine.hh:283
void assembleVVolumePostSkeleton(const EG &eg, const LFSVC &lfsv_cache)
Definition: default/residualengine.hh:313
void onBindLFSUVInside(const IG &ig, const LFSUC &lfsu_cache, const LFSVC &lfsv_cache)
Definition: default/residualengine.hh:138
void resize(size_type size)
Resize the container.
Definition: localvector.hh:251
static void alpha_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, R &r_s, R &r_n)
Definition: callswitch.hh:47
void loadCoefficientsLFSUCoupling(const LFSUC &lfsu_c_cache)
Definition: default/residualengine.hh:205
void assembleUVVolumePostSkeleton(const EG &eg, const LFSUC &lfsu_cache, const LFSVC &lfsv_cache)
Definition: default/residualengine.hh:305
static void lambda_volume(const LA &la, const EG &eg, const LFSV &lfsv, R &r)
Definition: callswitch.hh:61
void onUnbindLFSVInside(const IG &ig, const LFSVC &lfsv_cache)
Definition: default/residualengine.hh:179
static void alpha_boundary(const LA &la, const IG &ig, const LFSU &lfsu_s, const X &x_s, const LFSV &lfsv_s, R &r_s)
Definition: callswitch.hh:54
Residual::template LocalView< LFSVCache > ResidualView
Definition: default/residualengine.hh:55
const LocalAssembler::Traits::TrialGridFunctionSpaceConstraints & trialConstraints() const
Trial space constraints.
Definition: default/residualengine.hh:97
void onBindLFSUV(const EG &eg, const LFSUC &lfsu_cache, const LFSVC &lfsv_cache)
Definition: default/residualengine.hh:126
Base class for LocalAssemblerEngine implementations to avoid boilerplate code.
Definition: localassemblerenginebase.hh:21
LA::LFSUCache LFSUCache
Definition: default/residualengine.hh:48
static void lambda_volume_post_skeleton(const LA &la, const EG &eg, const LFSV &lfsv, R &r)
Definition: callswitch.hh:65
bool requireUVBoundary() const
Definition: default/residualengine.hh:83
LA::Traits::Residual Residual
The type of the residual vector.
Definition: default/residualengine.hh:39
static void alpha_volume(const LA &la, const EG &eg, const LFSU &lfsu, const X &x, const LFSV &lfsv, R &r)
Definition: callswitch.hh:39
const IG & ig
Definition: constraints.hh:147
LA::Traits::Solution Solution
The type of the solution vector.
Definition: default/residualengine.hh:43
void onBindLFSV(const EG &eg, const LFSVC &lfsv_cache)
Definition: default/residualengine.hh:132
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
Definition: localfunctionspacetags.hh:48
void onBindLFSVOutside(const IG &ig, const LFSVC &lfsv_s_cache, const LFSVC &lfsv_n_cache)
Definition: default/residualengine.hh:159
Definition: localfunctionspacetags.hh:54
Solution::ElementType SolutionElement
Definition: default/residualengine.hh:44
const LocalAssembler::Traits::TestGridFunctionSpaceConstraints & testConstraints() const
Test space constraints.
Definition: default/residualengine.hh:103
void onUnbindLFSVOutside(const IG &ig, const LFSVC &lfsv_s_cache, const LFSVC &lfsv_n_cache)
Definition: default/residualengine.hh:185
bool needsConstraintsCaching(const TrialConstraintsContainer &cu, const TestConstraintsContainer &cv) const
Definition: default/residualengine.hh:27
bool assembleCell(const EG &eg)
Definition: default/residualengine.hh:231
Residual::ElementType ResidualElement
Definition: default/residualengine.hh:40
void assembleUVBoundary(const IG &ig, const LFSUC &lfsu_s_cache, const LFSVC &lfsv_s_cache)
Definition: default/residualengine.hh:275
void loadCoefficientsLFSUOutside(const LFSUC &lfsu_n_cache)
Definition: default/residualengine.hh:201
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/residualengine.hh:253
bool requireUVVolume() const
Definition: default/residualengine.hh:75
void assembleVVolume(const EG &eg, const LFSVC &lfsv_cache)
Definition: default/residualengine.hh:245
LA::LFSVCache LFSVCache
Definition: default/residualengine.hh:51
void onUnbindLFSV(const EG &eg, const LFSVC &lfsv_cache)
Definition: default/residualengine.hh:173
void loadCoefficientsLFSUInside(const LFSUC &lfsu_s_cache)
Definition: default/residualengine.hh:197
Definition: adaptivity.hh:27
The local assembler engine for DUNE grids which assembles the residual vector.
Definition: default/residualengine.hh:21
void setWeight(weight_type weight)
Resets the weighting coefficient of the view.
Definition: localvector.hh:72
LFSV::Traits::GridFunctionSpace GFSV
Definition: default/residualengine.hh:52
LA LocalAssembler
The type of the wrapping local assembler.
Definition: default/residualengine.hh:33
bool requireSkeleton() const
Definition: default/residualengine.hh:71
static void lambda_boundary(const LA &la, const IG &ig, const LFSV &lfsv, R &r)
Definition: callswitch.hh:75
static void lambda_skeleton(const LA &la, const IG &ig, const LFSV &lfsv_s, const LFSV &lfsv_n, R &r_s, R &r_n)
Definition: callswitch.hh:69
const LocalAssembler & localAssembler() const
Public access to the wrapping local assembler.
Definition: default/residualengine.hh:94
bool requireSkeletonTwoSided() const
Definition: default/residualengine.hh:73
bool requireVVolume() const
Definition: default/residualengine.hh:77
bool requireUVVolumePostSkeleton() const
Definition: default/residualengine.hh:87
bool requireVVolumePostSkeleton() const
Definition: default/residualengine.hh:89
bool requireVSkeleton() const
Definition: default/residualengine.hh:81
void onBindLFSVInside(const IG &ig, const LFSVC &lfsv_cache)
Definition: default/residualengine.hh:153
WeightedVectorAccumulationView< LocalVector > WeightedAccumulationView
An accumulate-only view of this container that automatically applies a weight to all contributions...
Definition: localvector.hh:198
static void alpha_volume_post_skeleton(const LA &la, const EG &eg, const LFSU &lfsu, const X &x, const LFSV &lfsv, R &r)
Definition: callswitch.hh:43
void constrain_residual(const CG &cg, XG &xg)
transform residual into transformed basis: r -> r~
Definition: constraints.hh:911
LA::LFSU LFSU
The local function spaces.
Definition: default/residualengine.hh:47
bool requireUVSkeleton() const
Definition: default/residualengine.hh:79
bool requireVBoundary() const
Definition: default/residualengine.hh:85
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/residualengine.hh:144
void postAssembly(const GFSU &gfsu, const GFSV &gfsv)
Definition: default/residualengine.hh:212
LFSU::Traits::GridFunctionSpace GFSU
Definition: default/residualengine.hh:49
DefaultLocalResidualAssemblerEngine(const LocalAssembler &local_assembler_)
Constructor.
Definition: default/residualengine.hh:63
void assembleUVVolume(const EG &eg, const LFSUC &lfsu_cache, const LFSVC &lfsv_cache)
Definition: default/residualengine.hh:237
void setResidual(Residual &residual_)
Definition: default/residualengine.hh:110
void setSolution(const Solution &solution_)
Definition: default/residualengine.hh:117
Solution::template ConstLocalView< LFSUCache > SolutionView
Definition: default/residualengine.hh:54
const EG & eg
Definition: constraints.hh:280
LA::LFSV LFSV
Definition: default/residualengine.hh:50
static void assembleVEnrichedCoupling(const IG &ig, const LFSVC &lfsv_s_cache, const LFSVC &lfsv_n_cache, const LFSVC &lfsv_coupling_cache)
Definition: default/residualengine.hh:298
void assembleVSkeleton(const IG &ig, const LFSVC &lfsv_s_cache, const LFSVC &lfsv_n_cache)
Definition: default/residualengine.hh:266