dune-pdelab  2.5-dev
convectiondiffusionccfv.hh
Go to the documentation of this file.
1 // -*- tab-width: 2; indent-tabs-mode: nil -*-
2 #ifndef DUNE_PDELAB_LOCALOPERATOR_CONVECTIONDIFFUSIONCCFV_HH
3 #define DUNE_PDELAB_LOCALOPERATOR_CONVECTIONDIFFUSIONCCFV_HH
4 
5 #include<dune/common/exceptions.hh>
6 #include<dune/common/fvector.hh>
7 #include<dune/common/typetraits.hh>
8 #include<dune/geometry/referenceelements.hh>
9 
17 
18 namespace Dune {
19  namespace PDELab {
20 
36  template<typename TP>
38  public Dune::PDELab::NumericalJacobianApplyBoundary<ConvectionDiffusionCCFV<TP> >,
39  public FullSkeletonPattern,
40  public FullVolumePattern,
42  public InstationaryLocalOperatorDefaultMethods<typename TP::Traits::RangeFieldType>
43  {
44  using BCType = typename ConvectionDiffusionBoundaryConditions::Type;
45 
46  public:
47  // pattern assembly flags
48  enum { doPatternVolume = true };
49  enum { doPatternSkeleton = true };
50 
51  // residual assembly flags
52  enum { doAlphaVolume = true };
53  enum { doAlphaSkeleton = true };
54  enum { doAlphaBoundary = true };
55  enum { doLambdaVolume = true };
56  enum { doLambdaSkeleton = false };
57  enum { doLambdaBoundary = false };
58 
61  , param(param_)
62  {}
63 
64  // volume integral depending on test and ansatz functions
65  template<typename EG, typename LFSU, typename X, typename LFSV, typename R>
66  void alpha_volume (const EG& eg, const LFSU& lfsu, const X& x, const LFSV& lfsv, R& r) const
67  {
68  // get cell
69  const auto& cell = eg.entity();
70 
71  // cell center
72  auto geo = eg.geometry();
73  auto ref_el = referenceElement(geo);
74  auto local_inside = ref_el.position(0,0);
75 
76  // evaluate reaction term
77  auto c = param.c(cell,local_inside);
78 
79  // and accumulate
80  r.accumulate(lfsu,0,(c*x(lfsu,0))*geo.volume());
81  }
82 
83  // apply jacobian of volume term
84  template<typename EG, typename LFSU, typename X, typename LFSV, typename Y>
85  void jacobian_apply_volume (const EG& eg, const LFSU& lfsu, const X& x, const LFSV& lfsv, Y& y) const
86  {
87  alpha_volume(eg,lfsu,x,lfsv,y);
88  }
89 
90  // jacobian of volume term
91  template<typename EG, typename LFSU, typename X, typename LFSV, typename M>
92  void jacobian_volume (const EG& eg, const LFSU& lfsu, const X& x, const LFSV& lfsv,
93  M& mat) const
94  {
95  // get cell
96  const auto& cell = eg.entity();
97 
98  // cell center
99  auto geo = eg.geometry();
100  auto ref_el = referenceElement(geo);
101  auto local_inside = ref_el.position(0,0);
102 
103  // evaluate reaction term
104  auto c = param.c(cell,local_inside);
105 
106  // and accumulate
107  mat.accumulate(lfsu,0,lfsu,0,c*geo.volume());
108  }
109 
110  // skeleton integral depending on test and ansatz functions
111  // each face is only visited ONCE!
112  template<typename IG, typename LFSU, typename X, typename LFSV, typename R>
113  void alpha_skeleton (const IG& ig,
114  const LFSU& lfsu_s, const X& x_s, const LFSV& lfsv_s,
115  const LFSU& lfsu_n, const X& x_n, const LFSV& lfsv_n,
116  R& r_s, R& r_n) const
117  {
118  // define types
119  using RF = typename LFSU::Traits::FiniteElementType::
120  Traits::LocalBasisType::Traits::RangeFieldType;
121 
122  // dimensions
123  const auto dim = IG::Entity::dimension;
124 
125  // get cell entities from both sides of the intersection
126  auto cell_inside = ig.inside();
127  auto cell_outside = ig.outside();
128 
129  // get geometries
130  auto geo = ig.geometry();
131  auto geo_inside = cell_inside.geometry();
132  auto geo_outside = cell_outside.geometry();
133 
134  // get geometry of intersection in local coordinates of neighbor cells
135  auto geo_in_inside = ig.geometryInInside();
136 
137  // center in face's reference element
138  auto ref_el = referenceElement(geo);
139  auto face_local = ref_el.position(0,0);
140 
141  // face volume for integration
142  auto face_volume = geo.integrationElement(face_local) * ref_el.volume();
143 
144  // cell centers in references elements
145  auto ref_el_inside = referenceElement(geo_inside);
146  auto ref_el_outside = referenceElement(geo_outside);
147  auto local_inside = ref_el_inside.position(0,0);
148  auto local_outside = ref_el_outside.position(0,0);
149 
150  // evaluate diffusion coefficient from either side and take harmonic average
151  auto tensor_inside = param.A(cell_inside,local_inside);
152  auto tensor_outside = param.A(cell_outside,local_outside);
153  auto n_F = ig.centerUnitOuterNormal();
154  Dune::FieldVector<RF,dim> An_F;
155  tensor_inside.mv(n_F,An_F);
156  auto k_inside = n_F*An_F;
157  tensor_outside.mv(n_F,An_F);
158  auto k_outside = n_F*An_F;
159  auto k_avg = 2.0/(1.0/(k_inside+1E-30) + 1.0/(k_outside+1E-30));
160 
161  // evaluate convective term
162  auto iplocal_s = geo_in_inside.global(face_local);
163  auto b = param.b(cell_inside,iplocal_s);
164  auto vn = b*n_F;
165  RF u_upwind=0;
166  if (vn>=0) u_upwind = x_s(lfsu_s,0); else u_upwind = x_n(lfsu_n,0);
167 
168  // cell centers in global coordinates
169  auto global_inside = geo_inside.global(local_inside);
170  auto global_outside = geo_outside.global(local_outside);
171 
172  // distance between the two cell centers
173  global_inside -= global_outside;
174  auto distance = global_inside.two_norm();
175 
176  // contribution to residual on inside element, other residual is computed by symmetric call
177  r_s.accumulate(lfsu_s,0,(u_upwind*vn)*face_volume+k_avg*(x_s(lfsu_s,0)-x_n(lfsu_n,0))*face_volume/distance);
178  r_n.accumulate(lfsu_n,0,-(u_upwind*vn)*face_volume-k_avg*(x_s(lfsu_s,0)-x_n(lfsu_n,0))*face_volume/distance);
179  }
180 
181  // apply jacobian of skeleton term
182  template<typename IG, typename LFSU, typename X, typename LFSV, typename Y>
183  void jacobian_apply_skeleton (const IG& ig,
184  const LFSU& lfsu_s, const X& x_s, const LFSV& lfsv_s,
185  const LFSU& lfsu_n, const X& x_n, const LFSV& lfsv_n,
186  Y& y_s, Y& y_n) const
187  {
188  alpha_skeleton(ig,lfsu_s,x_s,lfsv_s,lfsu_n,x_n,lfsv_n,y_s,y_n);
189  }
190 
191  template<typename IG, typename LFSU, typename X, typename LFSV, typename M>
192  void jacobian_skeleton (const IG& ig,
193  const LFSU& lfsu_s, const X& x_s, const LFSV& lfsv_s,
194  const LFSU& lfsu_n, const X& x_n, const LFSV& lfsv_n,
195  M& mat_ss, M& mat_sn,
196  M& mat_ns, M& mat_nn) const
197  {
198  // define types
199  using RF = typename LFSU::Traits::FiniteElementType::
200  Traits::LocalBasisType::Traits::RangeFieldType;
201 
202  // dimensions
203  const auto dim = IG::Entity::dimension;
204 
205  // get cell entities from both sides of the intersection
206  auto cell_inside = ig.inside();
207  auto cell_outside = ig.outside();
208 
209  // get geometries
210  auto geo = ig.geometry();
211  auto geo_inside = cell_inside.geometry();
212  auto geo_outside = cell_outside.geometry();
213 
214  // get geometry of intersection in local coordinates of neighbor cells
215  auto geo_in_inside = ig.geometryInInside();
216 
217  // center in face's reference element
218  auto ref_el = referenceElement(geo);
219  auto face_local = ref_el.position(0,0);
220 
221  // face volume for integration
222  auto face_volume = geo.integrationElement(face_local) * ref_el.volume();
223 
224  // cell centers in references elements
225  auto ref_el_inside = referenceElement(geo_inside);
226  auto ref_el_outside = referenceElement(geo_outside);
227  auto local_inside = ref_el_inside.position(0,0);
228  auto local_outside = ref_el_outside.position(0,0);
229 
230  // evaluate diffusion coefficient from either side and take harmonic average
231  auto tensor_inside = param.A(cell_inside,local_inside);
232  auto tensor_outside = param.A(cell_outside,local_outside);
233  auto n_F = ig.centerUnitOuterNormal();
234  Dune::FieldVector<RF,dim> An_F;
235  tensor_inside.mv(n_F,An_F);
236  auto k_inside = n_F*An_F;
237  tensor_outside.mv(n_F,An_F);
238  auto k_outside = n_F*An_F;
239  auto k_avg = 2.0/(1.0/(k_inside+1E-30) + 1.0/(k_outside+1E-30));
240 
241  // evaluate convective term
242  auto iplocal_s = geo_in_inside.global(face_local);
243  auto b = param.b(cell_inside,iplocal_s);
244  auto vn = b*n_F;
245 
246  // cell centers in global coordinates
247  auto global_inside = geo_inside.global(local_inside);
248  auto global_outside = geo_outside.global(local_outside);
249 
250  // distance between the two cell centers
251  global_inside -= global_outside;
252  auto distance = global_inside.two_norm();
253 
254  // contribution to jacobians on inside element and outside element for test and trial function
255  mat_ss.accumulate(lfsu_s,0,lfsu_s,0, k_avg*face_volume/distance );
256  mat_ns.accumulate(lfsu_n,0,lfsu_s,0, -k_avg*face_volume/distance );
257  mat_sn.accumulate(lfsu_s,0,lfsu_n,0, -k_avg*face_volume/distance );
258  mat_nn.accumulate(lfsu_n,0,lfsu_n,0, k_avg*face_volume/distance );
259  if (vn>=0)
260  {
261  mat_ss.accumulate(lfsu_s,0,lfsu_s,0, vn*face_volume );
262  mat_ns.accumulate(lfsu_n,0,lfsu_s,0, -vn*face_volume );
263  }
264  else
265  {
266  mat_sn.accumulate(lfsu_s,0,lfsu_n,0, vn*face_volume );
267  mat_nn.accumulate(lfsu_n,0,lfsu_n,0, -vn*face_volume );
268  }
269  }
270 
271 
272  // post skeleton: compute time step allowable for cell; to be done later
273  template<typename EG, typename LFSU, typename X, typename LFSV, typename R>
274  void alpha_volume_post_skeleton(const EG& eg, const LFSU& lfsu, const X& x,
275  const LFSV& lfsv, R& r) const
276  {
277  if (not first_stage) return; // time step calculation is only done in first stage
278 
279  // get cell
280  const auto& cell = eg.entity();
281 
282  // cell center
283  auto geo = eg.geometry();
284  auto ref_el = referenceElement(geo);
285  auto local_inside = ref_el.position(0,0);
286 
287  // compute optimal dt for this cell
288  auto cellcapacity = param.d(cell,local_inside)*geo.volume();
289  auto celldt = cellcapacity/(cellinflux+1E-30);
290  dtmin = std::min(dtmin,celldt);
291  }
292 
293 
294  // skeleton integral depending on test and ansatz functions
295  // We put the Dirchlet evaluation also in the alpha term to save some geometry evaluations
296  template<typename IG, typename LFSU, typename X, typename LFSV, typename R>
297  void alpha_boundary (const IG& ig,
298  const LFSU& lfsu_s, const X& x_s, const LFSV& lfsv_s,
299  R& r_s) const
300  {
301  // define types
302  using RF = typename LFSU::Traits::FiniteElementType::
303  Traits::LocalBasisType::Traits::RangeFieldType;
304 
305  // dimensions
306  const auto dim = IG::Entity::dimension;
307 
308  // get cell entities from both sides of the intersection
309  auto cell_inside = ig.inside();
310 
311  // get geometries
312  auto geo = ig.geometry();
313  auto geo_inside = cell_inside.geometry();
314 
315  // get geometry of intersection in local coordinates of neighbor cells
316  auto geo_in_inside = ig.geometryInInside();
317 
318  // center in face's reference element
319  auto ref_el = referenceElement(geo);
320  auto face_local = ref_el.position(0,0);
321 
322  // face volume for integration
323  auto face_volume = geo.integrationElement(face_local) * ref_el.volume();
324 
325  // cell centers in references elements
326  auto ref_el_inside = referenceElement(geo_inside);
327  auto local_inside = ref_el_inside.position(0,0);
328 
329  // evaluate boundary condition type
330  auto bctype = param.bctype(ig.intersection(),face_local);
331 
333  {
334  // Dirichlet boundary
335  // distance between cell center and face center
336  auto global_inside = geo_inside.global(local_inside);
337  auto global_outside = geo.global(face_local);
338  global_inside -= global_outside;
339  auto distance = global_inside.two_norm();
340 
341  // evaluate diffusion coefficient
342  auto tensor_inside = param.A(cell_inside,local_inside);
343  auto n_F = ig.centerUnitOuterNormal();
344  Dune::FieldVector<RF,dim> An_F;
345  tensor_inside.mv(n_F,An_F);
346  auto k_inside = n_F*An_F;
347 
348  // evaluate boundary condition function
349  auto iplocal_s = geo_in_inside.global(face_local);
350  auto g = param.g(cell_inside,iplocal_s);
351 
352  // velocity needed for convection term
353  auto b = param.b(cell_inside,iplocal_s);
354  auto n = ig.centerUnitOuterNormal();
355 
356  // contribution to residual on inside element, assumes that Dirichlet boundary is inflow
357  r_s.accumulate(lfsu_s,0,(b*n)*g*face_volume + k_inside*(x_s(lfsu_s,0)-g)*face_volume/distance);
358 
359  return;
360  }
361 
363  {
364  // Neumann boundary
365  // evaluate flux boundary condition
366 
367  //evaluate boundary function
368  auto j = param.j(ig.intersection(),face_local);
369 
370  // contribution to residual on inside element
371  r_s.accumulate(lfsu_s,0,j*face_volume);
372 
373  return;
374  }
375 
377  {
378  // evaluate velocity field and outer unit normal
379  auto iplocal_s = geo_in_inside.global(face_local);
380  auto b = param.b(cell_inside,iplocal_s);
381  auto n = ig.centerUnitOuterNormal();
382 
383  // evaluate outflow boundary condition
384  auto o = param.o(ig.intersection(),face_local);
385 
386  // integrate o
387  r_s.accumulate(lfsu_s,0,((b*n)*x_s(lfsu_s,0) + o)*face_volume);
388 
389  return;
390  }
391  }
392 
393  template<typename IG, typename LFSU, typename X, typename LFSV, typename M>
394  void jacobian_boundary (const IG& ig,
395  const LFSU& lfsu_s, const X& x_s, const LFSV& lfsv_s,
396  M& mat_ss) const
397  {
398  // define types
399  using RF = typename LFSU::Traits::FiniteElementType::
400  Traits::LocalBasisType::Traits::RangeFieldType;
401 
402  // dimensions
403  const auto dim = IG::Entity::dimension;
404 
405  // get cell entities from both sides of the intersection
406  auto cell_inside = ig.inside();
407 
408  // get geometries
409  auto geo = ig.geometry();
410  auto geo_inside = cell_inside.geometry();
411 
412  // get geometry of intersection in local coordinates of neighbor cells
413  auto geo_in_inside = ig.geometryInInside();
414 
415  // center in face's reference element
416  auto ref_el = referenceElement(geo);
417  auto face_local = ref_el.position(0,0);
418 
419  // face volume for integration
420  auto face_volume = geo.integrationElement(face_local) * ref_el.volume();
421 
422  // cell centers in references elements
423  auto ref_el_inside = referenceElement(geo_inside);
424  auto local_inside = ref_el_inside.position(0,0);
425 
426  // evaluate boundary condition type
427  auto bctype = param.bctype(ig.intersection(),face_local);
428 
430  {
431  // Dirichlet boundary
432  // distance between cell center and face center
433  auto global_inside = geo_inside.global(local_inside);
434  auto global_outside = geo.global(face_local);
435  global_inside -= global_outside;
436  auto distance = global_inside.two_norm();
437 
438  // evaluate diffusion coefficient
439  auto tensor_inside = param.A(cell_inside,local_inside);
440  auto n_F = ig.centerUnitOuterNormal();
441  Dune::FieldVector<RF,dim> An_F;
442  tensor_inside.mv(n_F,An_F);
443  auto k_inside = n_F*An_F;
444 
445  // contribution to jacobian on inside element for test and trial function
446  mat_ss.accumulate(lfsu_s,0,lfsv_s,0, k_inside*face_volume/distance );
447 
448  return;
449  }
450 
452  {
453  // evaluate velocity field and outer unit normal
454  auto iplocal_s = geo_in_inside.global(face_local);
455  auto b = param.b(cell_inside,iplocal_s);
456  auto n = ig.centerUnitOuterNormal();
457 
458  // integrate o
459  mat_ss.accumulate(lfsu_s,0,lfsv_s,0, (b*n)*face_volume );
460 
461  return;
462  }
463  }
464 
465  // volume integral depending only on test functions
466  template<typename EG, typename LFSV, typename R>
467  void lambda_volume (const EG& eg, const LFSV& lfsv, R& r) const
468  {
469  // get cell
470  const auto& cell = eg.entity();
471 
472  // cell center
473  auto geo = eg.geometry();
474  auto ref_el = referenceElement(geo);
475  auto local_inside = ref_el.position(0,0);
476 
477  // evaluate source and sink term
478  auto f = param.f(cell,local_inside);
479 
480  r.accumulate(lfsv,0,-f*geo.volume());
481  }
482 
484  void setTime (typename TP::Traits::RangeFieldType t)
485  {
486  param.setTime(t);
487  }
488 
490  void preStep (typename TP::Traits::RangeFieldType time, typename TP::Traits::RangeFieldType dt,
491  int stages)
492  {
493  }
494 
496  void preStage (typename TP::Traits::RangeFieldType time, int r)
497  {
498  if (r==1)
499  {
500  first_stage = true;
501  dtmin = 1E100;
502  }
503  else first_stage = false;
504  }
505 
507  void postStage ()
508  {
509  }
510 
512  typename TP::Traits::RangeFieldType suggestTimestep (typename TP::Traits::RangeFieldType dt) const
513  {
514  return dtmin;
515  }
516 
517  private:
518  TP& param;
519  bool first_stage;
520  mutable typename TP::Traits::RangeFieldType dtmin; // accumulate minimum dt here
521  mutable typename TP::Traits::RangeFieldType cellinflux;
522  };
523 
524 
525 
526 
533  template<class TP>
535  public FullVolumePattern,
537  public InstationaryLocalOperatorDefaultMethods<typename TP::Traits::RangeFieldType>
538  {
539  public:
540  // pattern assembly flags
541  enum { doPatternVolume = true };
542 
543  // residual assembly flags
544  enum { doAlphaVolume = true };
545 
547  : param(param_)
548  {
549  }
550 
551  // volume integral depending on test and ansatz functions
552  template<typename EG, typename LFSU, typename X, typename LFSV, typename R>
553  void alpha_volume (const EG& eg, const LFSU& lfsu, const X& x, const LFSV& lfsv, R& r) const
554  {
555  // get cell
556  const auto& cell = eg.entity();
557 
558  // cell center
559  auto geo = eg.geometry();
560  auto ref_el = referenceElement(geo);
561  auto local_inside = ref_el.position(0,0);
562 
563  // capacity term
564  auto capacity = param.d(cell,local_inside);
565 
566  // residual contribution
567  r.accumulate(lfsu,0,capacity*x(lfsu,0)*geo.volume());
568  }
569 
570  // apply jacobian of volume term
571  template<typename EG, typename LFSU, typename X, typename LFSV, typename Y>
572  void jacobian_apply_volume (const EG& eg, const LFSU& lfsu, const X& x, const LFSV& lfsv, Y& y) const
573  {
574  alpha_volume(eg,lfsu,x,lfsv,y);
575  }
576 
577  // jacobian of volume term
578  template<typename EG, typename LFSU, typename X, typename LFSV, typename M>
579  void jacobian_volume (const EG& eg, const LFSU& lfsu, const X& x, const LFSV& lfsv,
580  M& mat) const
581  {
582  // get cell
583  const auto& cell = eg.entity();
584 
585  // cell center
586  auto geo = eg.geometry();
587  auto ref_el = referenceElement(geo);
588  auto local_inside = ref_el.position(0,0);
589 
590  // capacity term
591  auto capacity = param.d(cell,local_inside);
592 
593  // jacobian contribution
594  mat.accumulate(lfsu,0,lfsu,0,capacity*geo.volume());
595  }
596 
597  private:
598  TP& param;
599  };
600 
601 
603  } // namespace PDELab
604 } // namespace Dune
605 
606 #endif // DUNE_PDELAB_LOCALOPERATOR_CONVECTIONDIFFUSIONCCFV_HH
static const int dim
Definition: adaptivity.hh:84
const Entity & e
Definition: localfunctionspace.hh:120
Default flags for all local operators.
Definition: flags.hh:18
sparsity pattern generator
Definition: pattern.hh:29
void jacobian_volume(const EG &eg, const LFSU &lfsu, const X &x, const LFSV &lfsv, M &mat) const
Definition: convectiondiffusionccfv.hh:92
void jacobian_apply_volume(const EG &eg, const LFSU &lfsu, const X &x, const LFSV &lfsv, Y &y) const
Definition: convectiondiffusionccfv.hh:572
Definition: convectiondiffusionccfv.hh:54
TP::Traits::RangeFieldType suggestTimestep(typename TP::Traits::RangeFieldType dt) const
to be called once before each stage
Definition: convectiondiffusionccfv.hh:512
Default class for additional methods in instationary local operators.
Definition: idefault.hh:89
void preStep(typename TP::Traits::RangeFieldType time, typename TP::Traits::RangeFieldType dt, int stages)
to be called once before each time step
Definition: convectiondiffusionccfv.hh:490
sparsity pattern generator
Definition: pattern.hh:13
Definition: convectiondiffusionccfv.hh:57
Definition: convectiondiffusionccfv.hh:52
Definition: convectiondiffusionccfv.hh:55
void jacobian_boundary(const IG &ig, const LFSU &lfsu_s, const X &x_s, const LFSV &lfsv_s, M &mat_ss) const
Definition: convectiondiffusionccfv.hh:394
void setTime(typename TP::Traits::RangeFieldType t)
set time in parameter class
Definition: convectiondiffusionccfv.hh:484
Definition: convectiondiffusionparameter.hh:65
Definition: convectiondiffusionccfv.hh:49
Definition: convectiondiffusionparameter.hh:65
void alpha_volume(const EG &eg, const LFSU &lfsu, const X &x, const LFSV &lfsv, R &r) const
Definition: convectiondiffusionccfv.hh:66
void alpha_boundary(const IG &ig, const LFSU &lfsu_s, const X &x_s, const LFSV &lfsv_s, R &r_s) const
Definition: convectiondiffusionccfv.hh:297
For backward compatibility – Do not use this!
Definition: adaptivity.hh:28
void alpha_volume_post_skeleton(const EG &eg, const LFSU &lfsu, const X &x, const LFSV &lfsv, R &r) const
Definition: convectiondiffusionccfv.hh:274
ConvectionDiffusionCCFVTemporalOperator(TP &param_)
Definition: convectiondiffusionccfv.hh:546
Implements linear and nonlinear versions of jacobian_apply_boundary() based on alpha_boundary() ...
Definition: numericaljacobianapply.hh:285
const IG & ig
Definition: constraints.hh:148
ConvectionDiffusionCCFV(TP &param_)
Definition: convectiondiffusionccfv.hh:59
void preStage(typename TP::Traits::RangeFieldType time, int r)
to be called once before each stage
Definition: convectiondiffusionccfv.hh:496
Definition: convectiondiffusionccfv.hh:534
Definition: convectiondiffusionccfv.hh:53
void jacobian_volume(const EG &eg, const LFSU &lfsu, const X &x, const LFSV &lfsv, M &mat) const
Definition: convectiondiffusionccfv.hh:579
void alpha_volume(const EG &eg, const LFSU &lfsu, const X &x, const LFSV &lfsv, R &r) const
Definition: convectiondiffusionccfv.hh:553
void alpha_skeleton(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) const
Definition: convectiondiffusionccfv.hh:113
void jacobian_skeleton(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, M &mat_ss, M &mat_sn, M &mat_ns, M &mat_nn) const
Definition: convectiondiffusionccfv.hh:192
void jacobian_apply_volume(const EG &eg, const LFSU &lfsu, const X &x, const LFSV &lfsv, Y &y) const
Definition: convectiondiffusionccfv.hh:85
Type
Definition: convectiondiffusionparameter.hh:65
Definition: convectiondiffusionccfv.hh:56
void postStage()
to be called once at the end of each stage
Definition: convectiondiffusionccfv.hh:507
Definition: convectiondiffusionccfv.hh:48
Definition: convectiondiffusionparameter.hh:65
void jacobian_apply_skeleton(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) const
Definition: convectiondiffusionccfv.hh:183
void lambda_volume(const EG &eg, const LFSV &lfsv, R &r) const
Definition: convectiondiffusionccfv.hh:467
Definition: convectiondiffusionccfv.hh:37