dune-pdelab  2.5-dev
fastdg/nonlinearjacobianapplyengine.hh
Go to the documentation of this file.
1 #ifndef DUNE_PDELAB_GRIDOPERATOR_FASTDG_NONLINEARJACOBIANAPPLYENGINE_HH
2 #define DUNE_PDELAB_GRIDOPERATOR_FASTDG_NONLINEARJACOBIANAPPLYENGINE_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  FastDGLocalNonlinearJacobianApplyAssemblerEngine(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
86  {
87  return local_assembler;
88  }
89 
91  const typename LocalAssembler::Traits::TrialGridFunctionSpaceConstraints& trialConstraints() const
92  {
93  return localAssembler().trialConstraints();
94  }
95 
97  const typename LocalAssembler::Traits::TestGridFunctionSpaceConstraints& testConstraints() const
98  {
99  return localAssembler().testConstraints();
100  }
101 
104  void setResidual(Residual & residual_)
105  {
106  global_rl_view.attach(residual_);
107  global_rn_view.attach(residual_);
108  }
109 
112  void setSolution(const Solution & solution_)
113  {
114  global_sl_view.attach(solution_);
115  global_sn_view.attach(solution_);
116  }
117 
120  void setUpdate(const Solution & update_)
121  {
122  global_zl_view.attach(update_);
123  global_zn_view.attach(update_);
124  }
125 
129  template<typename EG, typename LFSUC, typename LFSVC>
130  void onBindLFSUV(const EG & eg, const LFSUC & lfsu_cache, const LFSVC & lfsv_cache)
131  {
132  global_sl_view.bind(lfsu_cache);
133  global_zl_view.bind(lfsu_cache);
134  }
135 
136  template<typename EG, typename LFSVC>
137  void onBindLFSV(const EG & eg, const LFSVC & lfsv_cache)
138  {
139  global_rl_view.bind(lfsv_cache);
140  }
141 
142  template<typename IG, typename LFSUC, typename LFSVC>
143  void onBindLFSUVInside(const IG & ig, const LFSUC & lfsu_cache, const LFSVC & lfsv_cache)
144  {
145  global_sl_view.bind(lfsu_cache);
146  global_zl_view.bind(lfsu_cache);
147  }
148 
149  template<typename IG, typename LFSUC, typename LFSVC>
150  void onBindLFSUVOutside(const IG & ig,
151  const LFSUC & lfsu_s_cache, const LFSVC & lfsv_s_cache,
152  const LFSUC & lfsu_n_cache, const LFSVC & lfsv_n_cache)
153  {
154  global_sn_view.bind(lfsu_n_cache);
155  global_zn_view.bind(lfsu_n_cache);
156  }
157 
158  template<typename IG, typename LFSVC>
159  void onBindLFSVInside(const IG & ig, const LFSVC & lfsv_cache)
160  {
161  global_rl_view.bind(lfsv_cache);
162  }
163 
164  template<typename IG, typename LFSVC>
165  void onBindLFSVOutside(const IG & ig,
166  const LFSVC & lfsv_s_cache,
167  const LFSVC & lfsv_n_cache)
168  {
169  global_rn_view.bind(lfsv_n_cache);
170  }
171 
173 
177  template<typename EG, typename LFSVC>
178  void onUnbindLFSV(const EG & eg, const LFSVC & lfsv_cache)
179  {
180  global_rl_view.commit();
181  global_rl_view.unbind();
182  }
183 
184  template<typename IG, typename LFSVC>
185  void onUnbindLFSVInside(const IG & ig, const LFSVC & lfsv_cache)
186  {
187  global_rl_view.commit();
188  global_rl_view.unbind();
189  }
190 
191  template<typename IG, typename LFSVC>
192  void onUnbindLFSVOutside(const IG & ig,
193  const LFSVC & lfsv_s_cache,
194  const LFSVC & lfsv_n_cache)
195  {
196  global_rn_view.commit();
197  global_rn_view.unbind();
198  }
200 
203  template<typename LFSUC>
204  void loadCoefficientsLFSUInside(const LFSUC & lfsu_s_cache)
205  {}
206 
207  template<typename LFSUC>
208  void loadCoefficientsLFSUOutside(const LFSUC & lfsu_n_cache)
209  {}
210 
211  template<typename LFSUC>
212  void loadCoefficientsLFSUCoupling(const LFSUC & lfsu_c_cache)
213  {
214  DUNE_THROW(Dune::NotImplemented,"No coupling lfsu available for ");
215  }
217 
220 
221  void postAssembly(const GFSU& gfsu, const GFSV& gfsv)
222  {
223  if(local_assembler.doPostProcessing())
224  Dune::PDELab::constrain_residual(local_assembler.testConstraints(),
225  global_rl_view.container());
226 
227  }
228 
230 
233 
240  template<typename EG>
241  bool assembleCell(const EG & eg)
242  {
243  return LocalAssembler::isNonOverlapping && eg.entity().partitionType() != Dune::InteriorEntity;
244  }
245 
246  template<typename EG, typename LFSUC, typename LFSVC>
247  void assembleUVVolume(const EG & eg, const LFSUC & lfsu_cache, const LFSVC & lfsv_cache)
248  {
249  global_rl_view.setWeight(local_assembler.weight());
251  nonlinear_jacobian_apply_volume(lop,eg,lfsu_cache.localFunctionSpace(),global_sl_view,global_zl_view,lfsv_cache.localFunctionSpace(),global_rl_view);
252  }
253 
254  template<typename IG, typename LFSUC, typename LFSVC>
255  void assembleUVSkeleton(const IG & ig, const LFSUC & lfsu_s_cache, const LFSVC & lfsv_s_cache,
256  const LFSUC & lfsu_n_cache, const LFSVC & lfsv_n_cache)
257  {
258  global_rl_view.setWeight(local_assembler.weight());
259  global_rn_view.setWeight(local_assembler.weight());
262  lfsu_s_cache.localFunctionSpace(),global_sl_view,global_zl_view,lfsv_s_cache.localFunctionSpace(),
263  lfsu_n_cache.localFunctionSpace(),global_sn_view,global_zn_view,lfsv_n_cache.localFunctionSpace(),
264  global_rl_view,global_rn_view);
265  }
266 
267  template<typename IG, typename LFSUC, typename LFSVC>
268  void assembleUVBoundary(const IG & ig, const LFSUC & lfsu_s_cache, const LFSVC & lfsv_s_cache)
269  {
270  global_rl_view.setWeight(local_assembler.weight());
272  nonlinear_jacobian_apply_boundary(lop,ig,lfsu_s_cache.localFunctionSpace(),global_sl_view,global_zl_view,lfsv_s_cache.localFunctionSpace(),global_rl_view);
273  }
274 
275  template<typename IG, typename LFSUC, typename LFSVC>
276  static void assembleUVEnrichedCoupling(const IG & ig,
277  const LFSUC & lfsu_s_cache, const LFSVC & lfsv_s_cache,
278  const LFSUC & lfsu_n_cache, const LFSVC & lfsv_n_cache,
279  const LFSUC & lfsu_coupling_cache, const LFSVC & lfsv_coupling_cache)
280  {DUNE_THROW(Dune::NotImplemented,"Assembling of coupling spaces is not implemented for ");}
281 
282  template<typename EG, typename LFSUC, typename LFSVC>
283  void assembleUVVolumePostSkeleton(const EG & eg, const LFSUC & lfsu_cache, const LFSVC & lfsv_cache)
284  {
285  global_rl_view.setWeight(local_assembler.weight());
287  nonlinear_jacobian_apply_volume_post_skeleton(lop,eg,lfsu_cache.localFunctionSpace(),global_sl_view,global_zl_view,lfsv_cache.localFunctionSpace(),global_rl_view);
288  }
289 
291 
292  private:
295  const LocalAssembler & local_assembler;
296 
298  const LOP & lop;
299 
301  ResidualView global_rl_view;
302  ResidualView global_rn_view;
303 
305  SolutionView global_sl_view;
306  SolutionView global_sn_view;
307 
309  SolutionView global_zl_view;
310  SolutionView global_zn_view;
311 
316 
319 
321 
322  }; // End of class FastDGLocalNonlinearJacobianApplyAssemblerEngine
323  }
324 }
325 #endif
const LocalAssembler & localAssembler() const
Public access to the wrapping local assembler.
Definition: fastdg/nonlinearjacobianapplyengine.hh:85
void assembleUVVolumePostSkeleton(const EG &eg, const LFSUC &lfsu_cache, const LFSVC &lfsv_cache)
Definition: fastdg/nonlinearjacobianapplyengine.hh:283
static void nonlinear_jacobian_apply_volume_post_skeleton(const LA &la, const EG &eg, const LFSU &lfsu, const X &x, const X &z, const LFSV &lfsv, Y &y)
Definition: callswitch.hh:159
void onBindLFSV(const EG &eg, const LFSVC &lfsv_cache)
Definition: fastdg/nonlinearjacobianapplyengine.hh:137
const LocalAssembler::Traits::TestGridFunctionSpaceConstraints & testConstraints() const
Test space constraints.
Definition: fastdg/nonlinearjacobianapplyengine.hh:97
Residual::template AliasedLocalView< LFSVCache > ResidualView
Definition: fastdg/nonlinearjacobianapplyengine.hh:55
void loadCoefficientsLFSUInside(const LFSUC &lfsu_s_cache)
Definition: fastdg/nonlinearjacobianapplyengine.hh:204
bool requireSkeletonTwoSided() const
Definition: fastdg/nonlinearjacobianapplyengine.hh:72
void setUpdate(const Solution &update_)
Definition: fastdg/nonlinearjacobianapplyengine.hh:120
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/nonlinearjacobianapplyengine.hh:36
Solution::ElementType SolutionElement
Definition: fastdg/nonlinearjacobianapplyengine.hh:44
void onBindLFSUV(const EG &eg, const LFSUC &lfsu_cache, const LFSVC &lfsv_cache)
Definition: fastdg/nonlinearjacobianapplyengine.hh:130
LA::LFSVCache LFSVCache
Definition: fastdg/nonlinearjacobianapplyengine.hh:51
Definition: localfunctionspacetags.hh:54
void onBindLFSVInside(const IG &ig, const LFSVC &lfsv_cache)
Definition: fastdg/nonlinearjacobianapplyengine.hh:159
bool requireUVSkeleton() const
Definition: fastdg/nonlinearjacobianapplyengine.hh:76
void setResidual(Residual &residual_)
Definition: fastdg/nonlinearjacobianapplyengine.hh:104
void postAssembly(const GFSU &gfsu, const GFSV &gfsv)
Definition: fastdg/nonlinearjacobianapplyengine.hh:221
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/nonlinearjacobianapplyengine.hh:276
bool requireUVVolume() const
Definition: fastdg/nonlinearjacobianapplyengine.hh:74
bool requireUVBoundary() const
Definition: fastdg/nonlinearjacobianapplyengine.hh:78
Base class for LocalAssemblerEngine implementations to avoid boilerplate code.
Definition: localassemblerenginebase.hh:21
static void nonlinear_jacobian_apply_boundary(const LA &la, const IG &ig, const LFSU &lfsu_s, const X &x_s, const X &z_s, const LFSV &lfsv_s, Y &y_s)
Definition: callswitch.hh:170
void onUnbindLFSVInside(const IG &ig, const LFSVC &lfsv_cache)
Definition: fastdg/nonlinearjacobianapplyengine.hh:185
LA::Traits::Solution Solution
The type of the solution vector.
Definition: fastdg/nonlinearjacobianapplyengine.hh:43
void onBindLFSUVInside(const IG &ig, const LFSUC &lfsu_cache, const LFSVC &lfsv_cache)
Definition: fastdg/nonlinearjacobianapplyengine.hh:143
void setSolution(const Solution &solution_)
Definition: fastdg/nonlinearjacobianapplyengine.hh:112
static void nonlinear_jacobian_apply_volume(const LA &la, const EG &eg, const LFSU &lfsu, const X &x, const X &z, const LFSV &lfsv, Y &y)
Definition: callswitch.hh:155
void loadCoefficientsLFSUOutside(const LFSUC &lfsu_n_cache)
Definition: fastdg/nonlinearjacobianapplyengine.hh:208
LFSU::Traits::GridFunctionSpace GFSU
Definition: fastdg/nonlinearjacobianapplyengine.hh:49
For backward compatibility – Do not use this!
Definition: adaptivity.hh:28
Solution::template ConstAliasedLocalView< LFSUCache > SolutionView
Definition: fastdg/nonlinearjacobianapplyengine.hh:54
void assembleUVBoundary(const IG &ig, const LFSUC &lfsu_s_cache, const LFSVC &lfsv_s_cache)
Definition: fastdg/nonlinearjacobianapplyengine.hh:268
LFSV::Traits::GridFunctionSpace GFSV
Definition: fastdg/nonlinearjacobianapplyengine.hh:52
const LocalAssembler::Traits::TrialGridFunctionSpaceConstraints & trialConstraints() const
Trial space constraints.
Definition: fastdg/nonlinearjacobianapplyengine.hh:91
void onUnbindLFSV(const EG &eg, const LFSVC &lfsv_cache)
Definition: fastdg/nonlinearjacobianapplyengine.hh:178
const IG & ig
Definition: constraints.hh:148
bool assembleCell(const EG &eg)
Definition: fastdg/nonlinearjacobianapplyengine.hh:241
LA::Traits::Residual Residual
The type of the residual vector.
Definition: fastdg/nonlinearjacobianapplyengine.hh:39
bool requireUVVolumePostSkeleton() const
Definition: fastdg/nonlinearjacobianapplyengine.hh:80
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/nonlinearjacobianapplyengine.hh:255
bool needsConstraintsCaching(const TrialConstraintsContainer &cu, const TestConstraintsContainer &cv) const
Definition: fastdg/nonlinearjacobianapplyengine.hh:27
LA::LFSUCache LFSUCache
Definition: fastdg/nonlinearjacobianapplyengine.hh:48
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/nonlinearjacobianapplyengine.hh:150
LA LocalAssembler
The type of the wrapping local assembler.
Definition: fastdg/nonlinearjacobianapplyengine.hh:33
Residual::ElementType ResidualElement
Definition: fastdg/nonlinearjacobianapplyengine.hh:40
void loadCoefficientsLFSUCoupling(const LFSUC &lfsu_c_cache)
Definition: fastdg/nonlinearjacobianapplyengine.hh:212
bool requireSkeleton() const
Definition: fastdg/nonlinearjacobianapplyengine.hh:70
LA::LFSU LFSU
The local function spaces.
Definition: fastdg/nonlinearjacobianapplyengine.hh:47
void onUnbindLFSVOutside(const IG &ig, const LFSVC &lfsv_s_cache, const LFSVC &lfsv_n_cache)
Definition: fastdg/nonlinearjacobianapplyengine.hh:192
void assembleUVVolume(const EG &eg, const LFSUC &lfsu_cache, const LFSVC &lfsv_cache)
Definition: fastdg/nonlinearjacobianapplyengine.hh:247
void onBindLFSVOutside(const IG &ig, const LFSVC &lfsv_s_cache, const LFSVC &lfsv_n_cache)
Definition: fastdg/nonlinearjacobianapplyengine.hh:165
FastDGLocalNonlinearJacobianApplyAssemblerEngine(const LocalAssembler &local_assembler_)
Constructor.
Definition: fastdg/nonlinearjacobianapplyengine.hh:63
Definition: localfunctionspacetags.hh:48
LA::LFSV LFSV
Definition: fastdg/nonlinearjacobianapplyengine.hh:50
The fast DG local assembler engine for DUNE grids which assembles the local application of the Jacobi...
Definition: fastdg/nonlinearjacobianapplyengine.hh:21
static void nonlinear_jacobian_apply_skeleton(const LA &la, const IG &ig, const LFSU &lfsu_s, const X &x_s, const X &z_s, const LFSV &lfsv_s, const LFSU &lfsu_n, const X &x_n, const X &z_n, const LFSV &lfsv_n, Y &y_s, Y &y_n)
Definition: callswitch.hh:163