dune-pdelab  2.7-git
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  static constexpr bool isLinear = LOP::isLinear;
40 
42  typedef typename LA::Traits::Range Range;
43  typedef typename Range::ElementType RangeElement;
44 
46  typedef typename LA::Traits::Domain Domain;
47  typedef typename Domain::ElementType DomainElement;
48 
50  typedef typename LA::LFSU LFSU;
51  typedef typename LA::LFSUCache LFSUCache;
52  typedef typename LFSU::Traits::GridFunctionSpace GFSU;
53  typedef typename LA::LFSV LFSV;
54  typedef typename LA::LFSVCache LFSVCache;
55  typedef typename LFSV::Traits::GridFunctionSpace GFSV;
56 
57  typedef typename Domain::template ConstAliasedLocalView<LFSUCache> DomainView;
58  typedef typename Range::template AliasedLocalView<LFSVCache> RangeView;
59 
67  : local_assembler(local_assembler_),
68  lop(local_assembler_.localOperator())
69  {}
70 
73  bool requireSkeleton() const
74  { return local_assembler.doAlphaSkeleton(); }
76  { return local_assembler.doSkeletonTwoSided(); }
77  bool requireUVVolume() const
78  { return local_assembler.doAlphaVolume(); }
79  bool requireUVSkeleton() const
80  { return local_assembler.doAlphaSkeleton(); }
81  bool requireUVBoundary() const
82  { return local_assembler.doAlphaBoundary(); }
84  { return local_assembler.doAlphaVolumePostSkeleton(); }
86 
89  {
90  return local_assembler;
91  }
92 
94  const typename LocalAssembler::Traits::TrialGridFunctionSpaceConstraints& trialConstraints() const
95  {
96  return localAssembler().trialConstraints();
97  }
98 
100  const typename LocalAssembler::Traits::TestGridFunctionSpaceConstraints& testConstraints() const
101  {
102  return localAssembler().testConstraints();
103  }
104 
107  void setSolution(const Domain & solution_)
108  {
109  if (isLinear)
110  DUNE_THROW(Dune::Exception, "In the linear case the jacobian apply does not depend on the current solution and this method should never be called.");
111  global_solution_view_inside.attach(solution_);
112  global_solution_view_outside.attach(solution_);
113  }
114 
117  void setUpdate(const Domain & update_)
118  {
119  global_update_view_inside.attach(update_);
120  global_update_view_outside.attach(update_);
121  }
122 
125  void setResult(Range & result_)
126  {
127  global_result_view_inside.attach(result_);
128  global_result_view_outside.attach(result_);
129  }
130 
134  template<typename EG, typename LFSUC, typename LFSVC>
135  void onBindLFSUV(const EG & eg, const LFSUC & lfsu_cache, const LFSVC & lfsv_cache)
136  {
137  if (not isLinear)
138  global_solution_view_inside.bind(lfsu_cache);
139  global_update_view_inside.bind(lfsu_cache);
140  }
141 
142  template<typename EG, typename LFSVC>
143  void onBindLFSV(const EG & eg, const LFSVC & lfsv_cache)
144  {
145  global_result_view_inside.bind(lfsv_cache);
146  }
147 
148  template<typename IG, typename LFSUC, typename LFSVC>
149  void onBindLFSUVInside(const IG & ig, const LFSUC & lfsu_cache, const LFSVC & lfsv_cache)
150  {
151  if (not isLinear)
152  global_solution_view_inside.bind(lfsu_cache);
153  global_update_view_inside.bind(lfsu_cache);
154  }
155 
156  template<typename IG, typename LFSUC, typename LFSVC>
157  void onBindLFSUVOutside(const IG & ig,
158  const LFSUC & lfsu_s_cache, const LFSVC & lfsv_s_cache,
159  const LFSUC & lfsu_n_cache, const LFSVC & lfsv_n_cache)
160  {
161  if (not isLinear)
162  global_solution_view_outside.bind(lfsu_n_cache);
163  global_update_view_outside.bind(lfsu_n_cache);
164  }
165 
166  template<typename IG, typename LFSVC>
167  void onBindLFSVInside(const IG & ig, const LFSVC & lfsv_cache)
168  {
169  global_result_view_inside.bind(lfsv_cache);
170  }
171 
172  template<typename IG, typename LFSVC>
173  void onBindLFSVOutside(const IG & ig,
174  const LFSVC & lfsv_s_cache,
175  const LFSVC & lfsv_n_cache)
176  {
177  global_result_view_outside.bind(lfsv_n_cache);
178  }
179 
181 
185  template<typename EG, typename LFSVC>
186  void onUnbindLFSV(const EG & eg, const LFSVC & lfsv_cache)
187  {
188  global_result_view_inside.commit();
189  global_result_view_inside.unbind();
190  }
191 
192  template<typename IG, typename LFSVC>
193  void onUnbindLFSVInside(const IG & ig, const LFSVC & lfsv_cache)
194  {
195  global_result_view_inside.commit();
196  global_result_view_inside.unbind();
197  }
198 
199  template<typename IG, typename LFSVC>
200  void onUnbindLFSVOutside(const IG & ig,
201  const LFSVC & lfsv_s_cache,
202  const LFSVC & lfsv_n_cache)
203  {
204  global_result_view_outside.commit();
205  global_result_view_outside.unbind();
206  }
208 
211  template<typename LFSUC>
212  void loadCoefficientsLFSUInside(const LFSUC & lfsu_s_cache)
213  {}
214 
215  template<typename LFSUC>
216  void loadCoefficientsLFSUOutside(const LFSUC & lfsu_n_cache)
217  {}
218 
219  template<typename LFSUC>
220  void loadCoefficientsLFSUCoupling(const LFSUC & lfsu_c_cache)
221  {
222  DUNE_THROW(Dune::NotImplemented,"No coupling lfsu available for ");
223  }
225 
228 
229  void postAssembly(const GFSU& gfsu, const GFSV& gfsv)
230  {
231  if(local_assembler.doPostProcessing())
232  Dune::PDELab::constrain_residual(local_assembler.testConstraints(),
233  global_result_view_inside.container());
234  }
235 
237 
240 
247  template<typename EG>
248  bool assembleCell(const EG & eg)
249  {
250  return LocalAssembler::isNonOverlapping && eg.entity().partitionType() != Dune::InteriorEntity;
251  }
252 
253  template<typename EG, typename LFSUC, typename LFSVC>
254  void assembleUVVolume(const EG & eg, const LFSUC & lfsu_cache, const LFSVC & lfsv_cache)
255  {
256  global_result_view_inside.setWeight(local_assembler.weight());
258  jacobian_apply_volume(lop,eg,lfsu_cache.localFunctionSpace(),global_solution_view_inside,global_update_view_inside,lfsv_cache.localFunctionSpace(),global_result_view_inside);
259  }
260 
261  template<typename IG, typename LFSUC, typename LFSVC>
262  void assembleUVSkeleton(const IG & ig, const LFSUC & lfsu_s_cache, const LFSVC & lfsv_s_cache,
263  const LFSUC & lfsu_n_cache, const LFSVC & lfsv_n_cache)
264  {
265  global_result_view_inside.setWeight(local_assembler.weight());
266  global_result_view_outside.setWeight(local_assembler.weight());
269  lfsu_s_cache.localFunctionSpace(),global_solution_view_inside,global_update_view_inside,lfsv_s_cache.localFunctionSpace(),
270  lfsu_n_cache.localFunctionSpace(),global_solution_view_outside,global_update_view_outside,lfsv_n_cache.localFunctionSpace(),
271  global_result_view_inside,global_result_view_outside);
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  global_result_view_inside.setWeight(local_assembler.weight());
279  jacobian_apply_boundary(lop,ig,lfsu_s_cache.localFunctionSpace(),global_solution_view_inside,global_update_view_inside,lfsv_s_cache.localFunctionSpace(),global_result_view_inside);
280  }
281 
282  template<typename IG, typename LFSUC, typename LFSVC>
283  static void assembleUVEnrichedCoupling(const IG & ig,
284  const LFSUC & lfsu_s_cache, const LFSVC & lfsv_s_cache,
285  const LFSUC & lfsu_n_cache, const LFSVC & lfsv_n_cache,
286  const LFSUC & lfsu_coupling_cache, const LFSVC & lfsv_coupling_cache)
287  {
288  DUNE_THROW(Dune::NotImplemented,"Assembling of coupling spaces is not implemented for ");
289  }
290 
291  template<typename EG, typename LFSUC, typename LFSVC>
292  void assembleUVVolumePostSkeleton(const EG & eg, const LFSUC & lfsu_cache, const LFSVC & lfsv_cache)
293  {
294  global_result_view_inside.setWeight(local_assembler.weight());
296  jacobian_apply_volume_post_skeleton(lop,eg,lfsu_cache.localFunctionSpace(),global_solution_view_inside,global_update_view_inside,lfsv_cache.localFunctionSpace(),global_result_view_inside);
297  }
298 
300 
301  private:
304  const LocalAssembler & local_assembler;
305 
307  const LOP & lop;
308 
310  DomainView global_solution_view_inside;
311  DomainView global_solution_view_outside;
312 
314  DomainView global_update_view_inside;
315  DomainView global_update_view_outside;
316 
318  RangeView global_result_view_inside;
319  RangeView global_result_view_outside;
320 
323  typedef Dune::PDELab::TrialSpaceTag LocalTrialSpaceTag;
324  typedef Dune::PDELab::TestSpaceTag LocalTestSpaceTag;
325 
328 
330 
331  }; // End of class FastDGLocalJacobianAssemblerEngine
332 
333  }
334 }
335 #endif
const IG & ig
Definition: constraints.hh:149
void constrain_residual(const CG &cg, XG &xg)
transform residual into transformed basis: r -> r~
Definition: constraints.hh:904
For backward compatibility – Do not use this!
Definition: adaptivity.hh:28
Impl::LocalAssemblerCallSwitchHelper< LOP, doIt > LocalAssemblerCallSwitch
Definition: callswitch.hh:344
Definition: localfunctionspacetags.hh:48
Definition: localfunctionspacetags.hh:54
Base class for LocalAssemblerEngine implementations to avoid boilerplate code.
Definition: localassemblerenginebase.hh:22
The fast DG local assembler engine for DUNE grids which assembles the local application of the Jacobi...
Definition: fastdg/jacobianapplyengine.hh:23
LA LocalAssembler
The type of the wrapping local assembler.
Definition: fastdg/jacobianapplyengine.hh:33
void loadCoefficientsLFSUOutside(const LFSUC &lfsu_n_cache)
Definition: fastdg/jacobianapplyengine.hh:216
bool requireSkeletonTwoSided() const
Definition: fastdg/jacobianapplyengine.hh:75
void onBindLFSUV(const EG &eg, const LFSUC &lfsu_cache, const LFSVC &lfsv_cache)
Definition: fastdg/jacobianapplyengine.hh:135
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:50
void loadCoefficientsLFSUInside(const LFSUC &lfsu_s_cache)
Definition: fastdg/jacobianapplyengine.hh:212
void setResult(Range &result_)
Definition: fastdg/jacobianapplyengine.hh:125
LA::LFSV LFSV
Definition: fastdg/jacobianapplyengine.hh:53
void loadCoefficientsLFSUCoupling(const LFSUC &lfsu_c_cache)
Definition: fastdg/jacobianapplyengine.hh:220
Domain::template ConstAliasedLocalView< LFSUCache > DomainView
Definition: fastdg/jacobianapplyengine.hh:57
LFSV::Traits::GridFunctionSpace GFSV
Definition: fastdg/jacobianapplyengine.hh:55
void postAssembly(const GFSU &gfsu, const GFSV &gfsv)
Definition: fastdg/jacobianapplyengine.hh:229
Range::ElementType RangeElement
Definition: fastdg/jacobianapplyengine.hh:43
LA::Traits::Range Range
The type of the result vector.
Definition: fastdg/jacobianapplyengine.hh:42
void assembleUVVolume(const EG &eg, const LFSUC &lfsu_cache, const LFSVC &lfsv_cache)
Definition: fastdg/jacobianapplyengine.hh:254
static constexpr bool isLinear
Wheter the local operator is linear.
Definition: fastdg/jacobianapplyengine.hh:39
bool requireUVVolumePostSkeleton() const
Definition: fastdg/jacobianapplyengine.hh:83
bool requireUVBoundary() const
Definition: fastdg/jacobianapplyengine.hh:81
LA::Traits::Domain Domain
The type of the solution vector.
Definition: fastdg/jacobianapplyengine.hh:46
void onBindLFSVOutside(const IG &ig, const LFSVC &lfsv_s_cache, const LFSVC &lfsv_n_cache)
Definition: fastdg/jacobianapplyengine.hh:173
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:157
void assembleUVBoundary(const IG &ig, const LFSUC &lfsu_s_cache, const LFSVC &lfsv_s_cache)
Definition: fastdg/jacobianapplyengine.hh:275
bool requireUVVolume() const
Definition: fastdg/jacobianapplyengine.hh:77
void setSolution(const Domain &solution_)
Definition: fastdg/jacobianapplyengine.hh:107
bool needsConstraintsCaching(const TrialConstraintsContainer &cu, const TestConstraintsContainer &cv) const
Definition: fastdg/jacobianapplyengine.hh:27
Range::template AliasedLocalView< LFSVCache > RangeView
Definition: fastdg/jacobianapplyengine.hh:58
const LocalAssembler & localAssembler() const
Public access to the wrapping local assembler.
Definition: fastdg/jacobianapplyengine.hh:88
bool assembleCell(const EG &eg)
Definition: fastdg/jacobianapplyengine.hh:248
FastDGLocalJacobianApplyAssemblerEngine(const LocalAssembler &local_assembler_)
Constructor.
Definition: fastdg/jacobianapplyengine.hh:66
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:262
const LocalAssembler::Traits::TrialGridFunctionSpaceConstraints & trialConstraints() const
Trial space constraints.
Definition: fastdg/jacobianapplyengine.hh:94
void onBindLFSVInside(const IG &ig, const LFSVC &lfsv_cache)
Definition: fastdg/jacobianapplyengine.hh:167
void onBindLFSUVInside(const IG &ig, const LFSUC &lfsu_cache, const LFSVC &lfsv_cache)
Definition: fastdg/jacobianapplyengine.hh:149
LA::LFSVCache LFSVCache
Definition: fastdg/jacobianapplyengine.hh:54
void onUnbindLFSVOutside(const IG &ig, const LFSVC &lfsv_s_cache, const LFSVC &lfsv_n_cache)
Definition: fastdg/jacobianapplyengine.hh:200
const LocalAssembler::Traits::TestGridFunctionSpaceConstraints & testConstraints() const
Test space constraints.
Definition: fastdg/jacobianapplyengine.hh:100
LA::LFSUCache LFSUCache
Definition: fastdg/jacobianapplyengine.hh:51
void onBindLFSV(const EG &eg, const LFSVC &lfsv_cache)
Definition: fastdg/jacobianapplyengine.hh:143
void assembleUVVolumePostSkeleton(const EG &eg, const LFSUC &lfsu_cache, const LFSVC &lfsv_cache)
Definition: fastdg/jacobianapplyengine.hh:292
LFSU::Traits::GridFunctionSpace GFSU
Definition: fastdg/jacobianapplyengine.hh:52
Domain::ElementType DomainElement
Definition: fastdg/jacobianapplyengine.hh:47
bool requireSkeleton() const
Definition: fastdg/jacobianapplyengine.hh:73
void setUpdate(const Domain &update_)
Definition: fastdg/jacobianapplyengine.hh:117
void onUnbindLFSV(const EG &eg, const LFSVC &lfsv_cache)
Definition: fastdg/jacobianapplyengine.hh:186
void onUnbindLFSVInside(const IG &ig, const LFSVC &lfsv_cache)
Definition: fastdg/jacobianapplyengine.hh:193
bool requireUVSkeleton() const
Definition: fastdg/jacobianapplyengine.hh:79
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:283