dune-pdelab  2.5-dev
fastdg/jacobianapplyengine.hh
Go to the documentation of this file.
1 #ifndef DUNE_PDELAB_GRIDOPERATOR_FASTDG_JACOBIANAPPLYENGINE_HH
2 #define DUNE_PDELAB_GRIDOPERATOR_FASTDG_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 ConstAliasedLocalView<LFSUCache> SolutionView;
55  typedef typename Residual::template AliasedLocalView<LFSVCache> ResidualView;
56 
63  FastDGLocalJacobianApplyAssemblerEngine(const LocalAssembler & local_assembler_)
64  : local_assembler(local_assembler_),
65  lop(local_assembler_.localOperator())
66  {}
67 
70  bool requireSkeleton() const
71  { return local_assembler.doAlphaSkeleton(); }
73  { return local_assembler.doSkeletonTwoSided(); }
74  bool requireUVVolume() const
75  { return local_assembler.doAlphaVolume(); }
76  bool requireUVSkeleton() const
77  { return local_assembler.doAlphaSkeleton(); }
78  bool requireUVBoundary() const
79  { return local_assembler.doAlphaBoundary(); }
81  { return local_assembler.doAlphaVolumePostSkeleton(); }
83 
85  const LocalAssembler & localAssembler() const { return local_assembler; }
86 
88  const typename LocalAssembler::Traits::TrialGridFunctionSpaceConstraints& trialConstraints() const
89  {
90  return localAssembler().trialConstraints();
91  }
92 
94  const typename LocalAssembler::Traits::TestGridFunctionSpaceConstraints& testConstraints() const
95  {
96  return localAssembler().testConstraints();
97  }
98 
101  void setResidual(Residual & residual_)
102  {
103  global_rl_view.attach(residual_);
104  global_rn_view.attach(residual_);
105  }
106 
109  void setSolution(const Solution & solution_)
110  {
111  global_sl_view.attach(solution_);
112  global_sn_view.attach(solution_);
113  }
114 
118  template<typename EG, typename LFSUC, typename LFSVC>
119  void onBindLFSUV(const EG & eg, const LFSUC & lfsu_cache, const LFSVC & lfsv_cache)
120  {
121  global_sl_view.bind(lfsu_cache);
122  }
123 
124  template<typename EG, typename LFSVC>
125  void onBindLFSV(const EG & eg, const LFSVC & lfsv_cache)
126  {
127  global_rl_view.bind(lfsv_cache);
128  }
129 
130  template<typename IG, typename LFSUC, typename LFSVC>
131  void onBindLFSUVInside(const IG & ig, const LFSUC & lfsu_cache, const LFSVC & lfsv_cache)
132  {
133  global_sl_view.bind(lfsu_cache);
134  }
135 
136  template<typename IG, typename LFSUC, typename LFSVC>
137  void onBindLFSUVOutside(const IG & ig,
138  const LFSUC & lfsu_s_cache, const LFSVC & lfsv_s_cache,
139  const LFSUC & lfsu_n_cache, const LFSVC & lfsv_n_cache)
140  {
141  global_sn_view.bind(lfsu_n_cache);
142  }
143 
144  template<typename IG, typename LFSVC>
145  void onBindLFSVInside(const IG & ig, const LFSVC & lfsv_cache)
146  {
147  global_rl_view.bind(lfsv_cache);
148  }
149 
150  template<typename IG, typename LFSVC>
151  void onBindLFSVOutside(const IG & ig,
152  const LFSVC & lfsv_s_cache,
153  const LFSVC & lfsv_n_cache)
154  {
155  global_rn_view.bind(lfsv_n_cache);
156  }
157 
159 
163  template<typename EG, typename LFSVC>
164  void onUnbindLFSV(const EG & eg, const LFSVC & lfsv_cache)
165  {
166  global_rl_view.commit();
167  global_rl_view.unbind();
168  }
169 
170  template<typename IG, typename LFSVC>
171  void onUnbindLFSVInside(const IG & ig, const LFSVC & lfsv_cache)
172  {
173  global_rl_view.commit();
174  global_rl_view.unbind();
175  }
176 
177  template<typename IG, typename LFSVC>
178  void onUnbindLFSVOutside(const IG & ig,
179  const LFSVC & lfsv_s_cache,
180  const LFSVC & lfsv_n_cache)
181  {
182  global_rn_view.commit();
183  global_rn_view.unbind();
184  }
186 
189  template<typename LFSUC>
190  void loadCoefficientsLFSUInside(const LFSUC & lfsu_s_cache)
191  {}
192 
193  template<typename LFSUC>
194  void loadCoefficientsLFSUOutside(const LFSUC & lfsu_n_cache)
195  {}
196 
197  template<typename LFSUC>
198  void loadCoefficientsLFSUCoupling(const LFSUC & lfsu_c_cache)
199  {
200  DUNE_THROW(Dune::NotImplemented,"No coupling lfsu available for ");
201  }
203 
206 
207  void postAssembly(const GFSU& gfsu, const GFSV& gfsv)
208  {
209  if(local_assembler.doPostProcessing())
210  Dune::PDELab::constrain_residual(local_assembler.testConstraints(),
211  global_rl_view.container());
212  }
213 
215 
218 
225  template<typename EG>
226  bool assembleCell(const EG & eg)
227  {
228  return LocalAssembler::isNonOverlapping && eg.entity().partitionType() != Dune::InteriorEntity;
229  }
230 
231  template<typename EG, typename LFSUC, typename LFSVC>
232  void assembleUVVolume(const EG & eg, const LFSUC & lfsu_cache, const LFSVC & lfsv_cache)
233  {
234  global_rl_view.setWeight(local_assembler.weight());
236  jacobian_apply_volume(lop,eg,lfsu_cache.localFunctionSpace(),global_sl_view,lfsv_cache.localFunctionSpace(),global_rl_view);
237  }
238 
239  template<typename IG, typename LFSUC, typename LFSVC>
240  void assembleUVSkeleton(const IG & ig, const LFSUC & lfsu_s_cache, const LFSVC & lfsv_s_cache,
241  const LFSUC & lfsu_n_cache, const LFSVC & lfsv_n_cache)
242  {
243  global_rl_view.setWeight(local_assembler.weight());
244  global_rn_view.setWeight(local_assembler.weight());
247  lfsu_s_cache.localFunctionSpace(),global_sl_view,lfsv_s_cache.localFunctionSpace(),
248  lfsu_n_cache.localFunctionSpace(),global_sn_view,lfsv_n_cache.localFunctionSpace(),
249  global_rl_view,global_rn_view);
250  }
251 
252  template<typename IG, typename LFSUC, typename LFSVC>
253  void assembleUVBoundary(const IG & ig, const LFSUC & lfsu_s_cache, const LFSVC & lfsv_s_cache)
254  {
255  global_rl_view.setWeight(local_assembler.weight());
257  jacobian_apply_boundary(lop,ig,lfsu_s_cache.localFunctionSpace(),global_sl_view,lfsv_s_cache.localFunctionSpace(),global_rl_view);
258  }
259 
260  template<typename IG, typename LFSUC, typename LFSVC>
261  static void assembleUVEnrichedCoupling(const IG & ig,
262  const LFSUC & lfsu_s_cache, const LFSVC & lfsv_s_cache,
263  const LFSUC & lfsu_n_cache, const LFSVC & lfsv_n_cache,
264  const LFSUC & lfsu_coupling_cache, const LFSVC & lfsv_coupling_cache)
265  {
266  DUNE_THROW(Dune::NotImplemented,"Assembling of coupling spaces is not implemented for ");
267  }
268 
269  template<typename EG, typename LFSUC, typename LFSVC>
270  void assembleUVVolumePostSkeleton(const EG & eg, const LFSUC & lfsu_cache, const LFSVC & lfsv_cache)
271  {
272  global_rl_view.setWeight(local_assembler.weight());
274  jacobian_apply_volume_post_skeleton(lop,eg,lfsu_cache.localFunctionSpace(),global_sl_view,lfsv_cache.localFunctionSpace(),global_rl_view);
275  }
276 
278 
279  private:
282  const LocalAssembler & local_assembler;
283 
285  const LOP & lop;
286 
288  ResidualView global_rl_view;
289  ResidualView global_rn_view;
290 
292  SolutionView global_sl_view;
293  SolutionView global_sn_view;
294 
299 
302 
304 
305  }; // End of class FastDGLocalJacobianAssemblerEngine
306 
307  }
308 }
309 #endif
const LocalAssembler & localAssembler() const
Public access to the wrapping local assembler.
Definition: fastdg/jacobianapplyengine.hh:85
The fast DG local assembler engine for DUNE grids which assembles the local application of the Jacobi...
Definition: fastdg/jacobianapplyengine.hh:21
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
void setSolution(const Solution &solution_)
Definition: fastdg/jacobianapplyengine.hh:109
void loadCoefficientsLFSUCoupling(const LFSUC &lfsu_c_cache)
Definition: fastdg/jacobianapplyengine.hh:198
void onUnbindLFSV(const EG &eg, const LFSVC &lfsv_cache)
Definition: fastdg/jacobianapplyengine.hh:164
void constrain_residual(const CG &cg, XG &xg)
transform residual into transformed basis: r -> r~
Definition: constraints.hh:906
LA::LocalOperator LOP
The type of the local operator.
Definition: fastdg/jacobianapplyengine.hh:36
LA::LFSU LFSU
The local function spaces.
Definition: fastdg/jacobianapplyengine.hh:47
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
bool requireUVVolume() const
Definition: fastdg/jacobianapplyengine.hh:74
Solution::template ConstAliasedLocalView< LFSUCache > SolutionView
Definition: fastdg/jacobianapplyengine.hh:54
void onBindLFSVInside(const IG &ig, const LFSVC &lfsv_cache)
Definition: fastdg/jacobianapplyengine.hh:145
Definition: localfunctionspacetags.hh:54
LA::Traits::Solution Solution
The type of the solution vector.
Definition: fastdg/jacobianapplyengine.hh:43
Solution::ElementType SolutionElement
Definition: fastdg/jacobianapplyengine.hh:44
bool assembleCell(const EG &eg)
Definition: fastdg/jacobianapplyengine.hh:226
void loadCoefficientsLFSUOutside(const LFSUC &lfsu_n_cache)
Definition: fastdg/jacobianapplyengine.hh:194
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
void setResidual(Residual &residual_)
Definition: fastdg/jacobianapplyengine.hh:101
Residual::ElementType ResidualElement
Definition: fastdg/jacobianapplyengine.hh:40
Base class for LocalAssemblerEngine implementations to avoid boilerplate code.
Definition: localassemblerenginebase.hh:21
FastDGLocalJacobianApplyAssemblerEngine(const LocalAssembler &local_assembler_)
Constructor.
Definition: fastdg/jacobianapplyengine.hh:63
bool needsConstraintsCaching(const TrialConstraintsContainer &cu, const TestConstraintsContainer &cv) const
Definition: fastdg/jacobianapplyengine.hh:27
void onBindLFSUVInside(const IG &ig, const LFSUC &lfsu_cache, const LFSVC &lfsv_cache)
Definition: fastdg/jacobianapplyengine.hh:131
bool requireUVVolumePostSkeleton() const
Definition: fastdg/jacobianapplyengine.hh:80
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/jacobianapplyengine.hh:261
LA::LFSV LFSV
Definition: fastdg/jacobianapplyengine.hh:50
bool requireUVBoundary() const
Definition: fastdg/jacobianapplyengine.hh:78
void loadCoefficientsLFSUInside(const LFSUC &lfsu_s_cache)
Definition: fastdg/jacobianapplyengine.hh:190
LA::LFSVCache LFSVCache
Definition: fastdg/jacobianapplyengine.hh:51
void assembleUVVolumePostSkeleton(const EG &eg, const LFSUC &lfsu_cache, const LFSVC &lfsv_cache)
Definition: fastdg/jacobianapplyengine.hh:270
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/jacobianapplyengine.hh:240
For backward compatibility – Do not use this!
Definition: adaptivity.hh:28
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
bool requireUVSkeleton() const
Definition: fastdg/jacobianapplyengine.hh:76
LA::LFSUCache LFSUCache
Definition: fastdg/jacobianapplyengine.hh:48
void onBindLFSVOutside(const IG &ig, const LFSVC &lfsv_s_cache, const LFSVC &lfsv_n_cache)
Definition: fastdg/jacobianapplyengine.hh:151
void onUnbindLFSVInside(const IG &ig, const LFSVC &lfsv_cache)
Definition: fastdg/jacobianapplyengine.hh:171
const IG & ig
Definition: constraints.hh:148
LFSV::Traits::GridFunctionSpace GFSV
Definition: fastdg/jacobianapplyengine.hh:52
void onBindLFSV(const EG &eg, const LFSVC &lfsv_cache)
Definition: fastdg/jacobianapplyengine.hh:125
void assembleUVBoundary(const IG &ig, const LFSUC &lfsu_s_cache, const LFSVC &lfsv_s_cache)
Definition: fastdg/jacobianapplyengine.hh:253
LA LocalAssembler
The type of the wrapping local assembler.
Definition: fastdg/jacobianapplyengine.hh:33
void onBindLFSUV(const EG &eg, const LFSUC &lfsu_cache, const LFSVC &lfsv_cache)
Definition: fastdg/jacobianapplyengine.hh:119
const LocalAssembler::Traits::TestGridFunctionSpaceConstraints & testConstraints() const
Test space constraints.
Definition: fastdg/jacobianapplyengine.hh:94
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/jacobianapplyengine.hh:137
const LocalAssembler::Traits::TrialGridFunctionSpaceConstraints & trialConstraints() const
Trial space constraints.
Definition: fastdg/jacobianapplyengine.hh:88
Definition: localfunctionspacetags.hh:48
LA::Traits::Residual Residual
The type of the residual vector.
Definition: fastdg/jacobianapplyengine.hh:39
LFSU::Traits::GridFunctionSpace GFSU
Definition: fastdg/jacobianapplyengine.hh:49
Residual::template AliasedLocalView< LFSVCache > ResidualView
Definition: fastdg/jacobianapplyengine.hh:55
bool requireSkeletonTwoSided() const
Definition: fastdg/jacobianapplyengine.hh:72
bool requireSkeleton() const
Definition: fastdg/jacobianapplyengine.hh:70
void postAssembly(const GFSU &gfsu, const GFSV &gfsv)
Definition: fastdg/jacobianapplyengine.hh:207
void onUnbindLFSVOutside(const IG &ig, const LFSVC &lfsv_s_cache, const LFSVC &lfsv_n_cache)
Definition: fastdg/jacobianapplyengine.hh:178
void assembleUVVolume(const EG &eg, const LFSUC &lfsu_cache, const LFSVC &lfsv_cache)
Definition: fastdg/jacobianapplyengine.hh:232