dune-pdelab  2.5-dev
nonlinearjacobianapplyengine.hh
Go to the documentation of this file.
1 #ifndef DUNE_PDELAB_GRIDOPERATOR_DEFAULT_NONLINEARJACOBIANAPPLYENGINE_HH
2 #define DUNE_PDELAB_GRIDOPERATOR_DEFAULT_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 ConstLocalView<LFSUCache> SolutionView;
55  typedef typename Residual::template LocalView<LFSVCache> ResidualView;
56 
63  DefaultLocalNonlinearJacobianApplyAssemblerEngine(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(); }
74  { return local_assembler.doSkeletonTwoSided(); }
75  bool requireUVVolume() const
76  { return local_assembler.doAlphaVolume(); }
77  bool requireUVSkeleton() const
78  { return local_assembler.doAlphaSkeleton(); }
79  bool requireUVBoundary() const
80  { return local_assembler.doAlphaBoundary(); }
82  { return local_assembler.doAlphaVolumePostSkeleton(); }
84 
86  const LocalAssembler & localAssembler() const { return local_assembler; }
87 
89  const typename LocalAssembler::Traits::TrialGridFunctionSpaceConstraints& trialConstraints() const
90  {
91  return localAssembler().trialConstraints();
92  }
93 
95  const typename LocalAssembler::Traits::TestGridFunctionSpaceConstraints& testConstraints() const
96  {
97  return localAssembler().testConstraints();
98  }
99 
102  void setResidual(Residual & residual_){
103  global_rl_view.attach(residual_);
104  global_rn_view.attach(residual_);
105  }
106 
109  void setSolution(const Solution & solution_){
110  global_sl_view.attach(solution_);
111  global_sn_view.attach(solution_);
112  }
113 
116  void setUpdate(const Solution & update_){
117  global_zl_view.attach(update_);
118  global_zn_view.attach(update_);
119  }
120 
124  template<typename EG, typename LFSUC, typename LFSVC>
125  void onBindLFSUV(const EG & eg, const LFSUC & lfsu_cache, const LFSVC & lfsv_cache){
126  global_sl_view.bind(lfsu_cache);
127  xl.resize(lfsu_cache.size());
128  global_zl_view.bind(lfsu_cache);
129  zl.resize(lfsu_cache.size());
130  }
131 
132  template<typename EG, typename LFSVC>
133  void onBindLFSV(const EG & eg, const LFSVC & lfsv_cache){
134  global_rl_view.bind(lfsv_cache);
135  rl.assign(lfsv_cache.size(),0.0);
136  }
137 
138  template<typename IG, typename LFSUC, typename LFSVC>
139  void onBindLFSUVInside(const IG & ig, const LFSUC & lfsu_cache, const LFSVC & lfsv_cache){
140  global_sl_view.bind(lfsu_cache);
141  xl.resize(lfsu_cache.size());
142  global_zl_view.bind(lfsu_cache);
143  zl.resize(lfsu_cache.size());
144  }
145 
146  template<typename IG, typename LFSUC, typename LFSVC>
147  void onBindLFSUVOutside(const IG & ig,
148  const LFSUC & lfsu_s_cache, const LFSVC & lfsv_s_cache,
149  const LFSUC & lfsu_n_cache, const LFSVC & lfsv_n_cache)
150  {
151  global_sn_view.bind(lfsu_n_cache);
152  xn.resize(lfsu_n_cache.size());
153  global_zn_view.bind(lfsu_n_cache);
154  zn.resize(lfsu_n_cache.size());
155  }
156 
157  template<typename IG, typename LFSVC>
158  void onBindLFSVInside(const IG & ig, const LFSVC & lfsv_cache){
159  global_rl_view.bind(lfsv_cache);
160  rl.assign(lfsv_cache.size(),0.0);
161  }
162 
163  template<typename IG, typename LFSVC>
164  void onBindLFSVOutside(const IG & ig,
165  const LFSVC & lfsv_s_cache,
166  const LFSVC & lfsv_n_cache)
167  {
168  global_rn_view.bind(lfsv_n_cache);
169  rn.assign(lfsv_n_cache.size(),0.0);
170  }
171 
173 
177  template<typename EG, typename LFSVC>
178  void onUnbindLFSV(const EG & eg, const LFSVC & lfsv_cache){
179  global_rl_view.add(rl);
180  global_rl_view.commit();
181  }
182 
183  template<typename IG, typename LFSVC>
184  void onUnbindLFSVInside(const IG & ig, const LFSVC & lfsv_cache){
185  global_rl_view.add(rl);
186  global_rl_view.commit();
187  }
188 
189  template<typename IG, typename LFSVC>
190  void onUnbindLFSVOutside(const IG & ig,
191  const LFSVC & lfsv_s_cache,
192  const LFSVC & lfsv_n_cache)
193  {
194  global_rn_view.add(rn);
195  global_rn_view.commit();
196  }
198 
201  template<typename LFSUC>
202  void loadCoefficientsLFSUInside(const LFSUC & lfsu_s_cache){
203  global_sl_view.read(xl);
204  global_zl_view.read(zl);
205  }
206  template<typename LFSUC>
207  void loadCoefficientsLFSUOutside(const LFSUC & lfsu_n_cache){
208  global_sn_view.read(xn);
209  global_zn_view.read(zn);
210  }
211  template<typename LFSUC>
212  void loadCoefficientsLFSUCoupling(const LFSUC & lfsu_c_cache)
213  {DUNE_THROW(Dune::NotImplemented,"No coupling lfsu available for ");}
215 
218 
219  void postAssembly(const GFSU& gfsu, const GFSV& gfsv){
220  if(local_assembler.doPostProcessing){
221  Dune::PDELab::constrain_residual(*(local_assembler.pconstraintsv),global_rl_view.container());
222  }
223  }
224 
226 
229 
236  template<typename EG>
237  bool assembleCell(const EG & eg)
238  {
239  return LocalAssembler::isNonOverlapping && eg.entity().partitionType() != Dune::InteriorEntity;
240  }
241 
242  template<typename EG, typename LFSUC, typename LFSVC>
243  void assembleUVVolume(const EG & eg, const LFSUC & lfsu_cache, const LFSVC & lfsv_cache)
244  {
245  rl_view.setWeight(local_assembler.weight);
247  nonlinear_jacobian_apply_volume(lop,eg,lfsu_cache.localFunctionSpace(),xl,zl,lfsv_cache.localFunctionSpace(),rl_view);
248  }
249 
250  template<typename IG, typename LFSUC, typename LFSVC>
251  void assembleUVSkeleton(const IG & ig, const LFSUC & lfsu_s_cache, const LFSVC & lfsv_s_cache,
252  const LFSUC & lfsu_n_cache, const LFSVC & lfsv_n_cache)
253  {
254  rl_view.setWeight(local_assembler.weight);
255  rn_view.setWeight(local_assembler.weight);
258  lfsu_s_cache.localFunctionSpace(),xl,zl,lfsv_s_cache.localFunctionSpace(),
259  lfsu_n_cache.localFunctionSpace(),xn,zn,lfsv_n_cache.localFunctionSpace(),
260  rl_view,rn_view);
261  }
262 
263  template<typename IG, typename LFSUC, typename LFSVC>
264  void assembleUVBoundary(const IG & ig, const LFSUC & lfsu_s_cache, const LFSVC & lfsv_s_cache)
265  {
266  rl_view.setWeight(local_assembler.weight);
268  nonlinear_jacobian_apply_boundary(lop,ig,lfsu_s_cache.localFunctionSpace(),xl,zl,lfsv_s_cache.localFunctionSpace(),rl_view);
269  }
270 
271  template<typename IG, typename LFSUC, typename LFSVC>
272  static void assembleUVEnrichedCoupling(const IG & ig,
273  const LFSUC & lfsu_s_cache, const LFSVC & lfsv_s_cache,
274  const LFSUC & lfsu_n_cache, const LFSVC & lfsv_n_cache,
275  const LFSUC & lfsu_coupling_cache, const LFSVC & lfsv_coupling_cache)
276  {DUNE_THROW(Dune::NotImplemented,"Assembling of coupling spaces is not implemented for ");}
277 
278  template<typename EG, typename LFSUC, typename LFSVC>
279  void assembleUVVolumePostSkeleton(const EG & eg, const LFSUC & lfsu_cache, const LFSVC & lfsv_cache)
280  {
281  rl_view.setWeight(local_assembler.weight);
283  nonlinear_jacobian_apply_volume_post_skeleton(lop,eg,lfsu_cache.localFunctionSpace(),xl,zl,lfsv_cache.localFunctionSpace(),rl_view);
284  }
285 
287 
288  private:
291  const LocalAssembler & local_assembler;
292 
294  const LOP & lop;
295 
297  ResidualView global_rl_view;
298  ResidualView global_rn_view;
299 
301  SolutionView global_sl_view;
302  SolutionView global_sn_view;
303 
305  SolutionView global_zl_view;
306  SolutionView global_zn_view;
307 
312 
315 
317  SolutionVector xl;
319  SolutionVector xn;
321  SolutionVector zl;
323  SolutionVector zn;
325  ResidualVector rl;
327  ResidualVector rn;
333 
334  }; // End of class DefaultLocalJacobianAssemblerEngine
335 
336  }
337 }
338 
339 #endif // DUNE_PDELAB_GRIDOPERATOR_DEFAULT_NONLINEARJACOBIANAPPLYENGINE_HH
bool requireUVVolumePostSkeleton() const
Definition: nonlinearjacobianapplyengine.hh:81
const IG & ig
Definition: constraints.hh:148
bool assembleCell(const EG &eg)
Definition: nonlinearjacobianapplyengine.hh:237
LFSV::Traits::GridFunctionSpace GFSV
Definition: nonlinearjacobianapplyengine.hh:52
const LocalAssembler & localAssembler() const
Public access to the wrapping local assembler.
Definition: nonlinearjacobianapplyengine.hh:86
void assembleUVVolumePostSkeleton(const EG &eg, const LFSUC &lfsu_cache, const LFSVC &lfsv_cache)
Definition: nonlinearjacobianapplyengine.hh:279
void onUnbindLFSVOutside(const IG &ig, const LFSVC &lfsv_s_cache, const LFSVC &lfsv_n_cache)
Definition: nonlinearjacobianapplyengine.hh:190
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
LA LocalAssembler
The type of the wrapping local assembler.
Definition: nonlinearjacobianapplyengine.hh:33
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: nonlinearjacobianapplyengine.hh:251
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: nonlinearjacobianapplyengine.hh:147
void onBindLFSUVInside(const IG &ig, const LFSUC &lfsu_cache, const LFSVC &lfsv_cache)
Definition: nonlinearjacobianapplyengine.hh:139
void assembleUVBoundary(const IG &ig, const LFSUC &lfsu_s_cache, const LFSVC &lfsv_s_cache)
Definition: nonlinearjacobianapplyengine.hh:264
The local assembler engine for DUNE grids which assembles the local application of the Jacobian...
Definition: nonlinearjacobianapplyengine.hh:21
void onUnbindLFSV(const EG &eg, const LFSVC &lfsv_cache)
Definition: nonlinearjacobianapplyengine.hh:178
Definition: localfunctionspacetags.hh:54
void loadCoefficientsLFSUOutside(const LFSUC &lfsu_n_cache)
Definition: nonlinearjacobianapplyengine.hh:207
void onBindLFSV(const EG &eg, const LFSVC &lfsv_cache)
Definition: nonlinearjacobianapplyengine.hh:133
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
LA::LFSUCache LFSUCache
Definition: nonlinearjacobianapplyengine.hh:48
void onBindLFSUV(const EG &eg, const LFSUC &lfsu_cache, const LFSVC &lfsv_cache)
Definition: nonlinearjacobianapplyengine.hh:125
LA::LocalOperator LOP
The type of the local operator.
Definition: nonlinearjacobianapplyengine.hh:36
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
Definition: localfunctionspacetags.hh:48
const LocalAssembler::Traits::TestGridFunctionSpaceConstraints & testConstraints() const
Test space constraints.
Definition: nonlinearjacobianapplyengine.hh:95
void assembleUVVolume(const EG &eg, const LFSUC &lfsu_cache, const LFSVC &lfsv_cache)
Definition: nonlinearjacobianapplyengine.hh:243
Solution::ElementType SolutionElement
Definition: nonlinearjacobianapplyengine.hh:44
void loadCoefficientsLFSUCoupling(const LFSUC &lfsu_c_cache)
Definition: nonlinearjacobianapplyengine.hh:212
void setSolution(const Solution &solution_)
Definition: nonlinearjacobianapplyengine.hh:109
bool requireUVBoundary() const
Definition: nonlinearjacobianapplyengine.hh:79
void setWeight(weight_type weight)
Resets the weighting coefficient of the view.
Definition: localvector.hh:72
LA::Traits::Residual Residual
The type of the residual vector.
Definition: nonlinearjacobianapplyengine.hh:39
bool requireUVSkeleton() const
Definition: nonlinearjacobianapplyengine.hh:77
void loadCoefficientsLFSUInside(const LFSUC &lfsu_s_cache)
Definition: nonlinearjacobianapplyengine.hh:202
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
LA::Traits::Solution Solution
The type of the solution vector.
Definition: nonlinearjacobianapplyengine.hh:43
void postAssembly(const GFSU &gfsu, const GFSV &gfsv)
Definition: nonlinearjacobianapplyengine.hh:219
void setResidual(Residual &residual_)
Definition: nonlinearjacobianapplyengine.hh:102
DefaultLocalNonlinearJacobianApplyAssemblerEngine(const LocalAssembler &local_assembler_)
Constructor.
Definition: nonlinearjacobianapplyengine.hh:63
void onBindLFSVInside(const IG &ig, const LFSVC &lfsv_cache)
Definition: nonlinearjacobianapplyengine.hh:158
Residual::template LocalView< LFSVCache > ResidualView
Definition: nonlinearjacobianapplyengine.hh:55
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
For backward compatibility – Do not use this!
Definition: adaptivity.hh:27
LFSU::Traits::GridFunctionSpace GFSU
Definition: nonlinearjacobianapplyengine.hh:49
WeightedVectorAccumulationView< LocalVector > WeightedAccumulationView
An accumulate-only view of this container that automatically applies a weight to all contributions...
Definition: localvector.hh:198
void onBindLFSVOutside(const IG &ig, const LFSVC &lfsv_s_cache, const LFSVC &lfsv_n_cache)
Definition: nonlinearjacobianapplyengine.hh:164
bool needsConstraintsCaching(const TrialConstraintsContainer &cu, const TestConstraintsContainer &cv) const
Definition: nonlinearjacobianapplyengine.hh:27
void constrain_residual(const CG &cg, XG &xg)
transform residual into transformed basis: r -> r~
Definition: constraints.hh:906
void onUnbindLFSVInside(const IG &ig, const LFSVC &lfsv_cache)
Definition: nonlinearjacobianapplyengine.hh:184
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: nonlinearjacobianapplyengine.hh:272
void resize(size_type size)
Resize the container.
Definition: localvector.hh:251
Residual::ElementType ResidualElement
Definition: nonlinearjacobianapplyengine.hh:40
LA::LFSVCache LFSVCache
Definition: nonlinearjacobianapplyengine.hh:51
LA::LFSU LFSU
The local function spaces.
Definition: nonlinearjacobianapplyengine.hh:47
Base class for LocalAssemblerEngine implementations to avoid boilerplate code.
Definition: localassemblerenginebase.hh:21
bool requireUVVolume() const
Definition: nonlinearjacobianapplyengine.hh:75
bool requireSkeletonTwoSided() const
Definition: nonlinearjacobianapplyengine.hh:73
Solution::template ConstLocalView< LFSUCache > SolutionView
Definition: nonlinearjacobianapplyengine.hh:54
bool requireSkeleton() const
Definition: nonlinearjacobianapplyengine.hh:71
LA::LFSV LFSV
Definition: nonlinearjacobianapplyengine.hh:50
const LocalAssembler::Traits::TrialGridFunctionSpaceConstraints & trialConstraints() const
Trial space constraints.
Definition: nonlinearjacobianapplyengine.hh:89
void setUpdate(const Solution &update_)
Definition: nonlinearjacobianapplyengine.hh:116