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