dune-pdelab  2.4-dev
transportccfv.hh
Go to the documentation of this file.
1 // -*- tab-width: 4; indent-tabs-mode: nil -*-
2 #ifndef DUNE_PDELAB_TRANSPORTCCFV_HH
3 #define DUNE_PDELAB_TRANSPORTCCFV_HH
4 
5 #include<dune/common/fvector.hh>
6 #include<dune/geometry/referenceelements.hh>
7 
9 
10 #include"../common/function.hh"
11 #include"pattern.hh"
12 #include"flags.hh"
13 #include"idefault.hh"
14 
15 namespace Dune {
16  namespace PDELab {
17 
19  template<typename GV, typename RF>
21  {
23  typedef GV GridViewType;
24 
26  enum {
28  dimDomain = GV::dimension
29  };
30 
32  typedef typename GV::Grid::ctype DomainFieldType;
33 
35  typedef Dune::FieldVector<DomainFieldType,dimDomain> DomainType;
36 
38  typedef Dune::FieldVector<DomainFieldType,dimDomain-1> IntersectionDomainType;
39 
41  typedef RF RangeFieldType;
42 
44  typedef Dune::FieldVector<RF,GV::dimensionworld> RangeType;
45 
47  typedef typename GV::Traits::template Codim<0>::Entity ElementType;
48  typedef typename GV::Intersection IntersectionType;
49  };
50 
52  template<class T, class Imp>
54  {
55  public:
56  typedef T Traits;
57 
59  typename Traits::RangeType
60  v (const typename Traits::ElementType& e, const typename Traits::DomainType& x) const
61  {
62  return asImp().v(e,x);
63  }
64 
66  typename Traits::RangeFieldType
67  D (const typename Traits::ElementType& e, const typename Traits::DomainType& x) const
68  {
69  return asImp().D(e,x);
70  }
71 
73  typename Traits::RangeFieldType
74  q (const typename Traits::ElementType& e, const typename Traits::DomainType& x) const
75  {
76  return asImp().q(e,x);
77  }
78 
80 
85  int
86  bc (const typename Traits::IntersectionType& is, const typename Traits::IntersectionDomainType& x) const
87  {
88  return asImp().bc(is,x);
89  }
90 
92  typename Traits::RangeFieldType
93  g (const typename Traits::ElementType& e, const typename Traits::DomainType& x) const
94  {
95  return asImp().g(e,x);
96  }
97 
99  typename Traits::RangeFieldType
100  j (const typename Traits::ElementType& e, const typename Traits::DomainType& x) const
101  {
102  return asImp().j(e,x);
103  }
104 
105  private:
106  Imp& asImp () {return static_cast<Imp &> (*this);}
107  const Imp& asImp () const {return static_cast<const Imp &>(*this);}
108  };
109 
110 
115  template<typename T>
117  : public Dune::PDELab::BoundaryGridFunctionBase<Dune::PDELab::
118  BoundaryGridFunctionTraits<typename T::Traits::GridViewType,int,1,
119  Dune::FieldVector<int,1> >,
120  BoundaryConditionType_Transport<T> >
121  {
122  typename T::Traits::GridViewType gv;
123  const T& t;
124 
125  public:
126  typedef Dune::PDELab::BoundaryGridFunctionTraits<typename T::Traits::GridViewType,int,1,
127  Dune::FieldVector<int,1> > Traits;
129 
130  BoundaryConditionType_Transport (const typename T::Traits::GridViewType& gv_, const T& t_) : gv(gv_), t(t_) {}
131 
132  template<typename I>
134  const typename Traits::DomainType& x,
135  typename Traits::RangeType& y) const
136  {
137  y = t.bc(ig.intersection(),x);
138  }
139 
141  inline const typename T::Traits::GridViewType& getGridView ()
142  {
143  return gv;
144  }
145  };
146 
147 
152  template<typename T>
154  : public GridFunctionBase<Dune::PDELab::GridFunctionTraits<typename T::Traits::GridViewType,
155  typename T::Traits::RangeFieldType,
156  1,Dune::FieldVector<typename T::Traits::RangeFieldType,1> >
157  ,DirichletBoundaryCondition_Transport<T> >
158  {
159  public:
160  typedef Dune::PDELab::GridFunctionTraits<typename T::Traits::GridViewType,
161  typename T::Traits::RangeFieldType,
162  1,Dune::FieldVector<typename T::Traits::RangeFieldType,1> > Traits;
163 
165  DirichletBoundaryCondition_Transport (const typename Traits::GridViewType& g_, const T& t_) : g(g_), t(t_) {}
166 
168  inline void evaluate (const typename Traits::ElementType& e,
169  const typename Traits::DomainType& x,
170  typename Traits::RangeType& y) const
171  {
172  y = t.g(e,x);
173  }
174 
175  inline const typename Traits::GridViewType& getGridView () const
176  {
177  return g;
178  }
179 
180  private:
181  typename Traits::GridViewType g;
182  const T& t;
183  };
184 
198  template<typename TP>
200  public NumericalJacobianApplySkeleton<CCFVSpatialTransportOperator<TP> >,
201  public NumericalJacobianSkeleton<CCFVSpatialTransportOperator<TP> >,
202  public NumericalJacobianApplyBoundary<CCFVSpatialTransportOperator<TP> >,
203  public NumericalJacobianBoundary<CCFVSpatialTransportOperator<TP> >,
204  public NumericalJacobianApplyVolumePostSkeleton<CCFVSpatialTransportOperator<TP> >,
205  public NumericalJacobianVolumePostSkeleton<CCFVSpatialTransportOperator<TP> >,
206  public FullSkeletonPattern,
207  public FullVolumePattern,
209  public InstationaryLocalOperatorDefaultMethods<typename TP::Traits::RangeFieldType>
210  {
211  enum { dim = TP::Traits::GridViewType::dimension };
212 
213  public:
214  // pattern assembly flags
215  enum { doPatternVolume = true };
216  enum { doPatternSkeleton = true };
217 
218  // residual assembly flags
219  enum { doAlphaVolume = true };
220  enum { doAlphaSkeleton = true };
221  enum { doAlphaVolumePostSkeleton = true };
222  enum { doAlphaBoundary = true };
223  enum { doLambdaVolume = true };
224 
225  enum { doSkeletonTwoSided = true }; // need to see face from both sides for CFL calculation
226 
228  : tp(tp_)
229  {
230  }
231 
232  // volume integral depending on test and ansatz functions
233  template<typename EG, typename LFSU, typename X, typename LFSV, typename R>
234  void alpha_volume (const EG& eg, const LFSU& lfsu, const X& x, const LFSV& lfsv, R& r) const
235  {
236  cellinflux = 0.0; // prepare dt computation
237  }
238 
239  // jacobian of volume term
240  template<typename EG, typename LFSU, typename X, typename LFSV, typename M>
241  void jacobian_volume (const EG& eg, const LFSU& lfsu, const X& x, const LFSV& lfsv,
242  M& mat) const
243  {
244  // do nothing; alpha_volume only needed for dt computations
245  }
246 
247  // skeleton integral depending on test and ansatz functions
248  // each face is only visited ONCE!
249  template<typename IG, typename LFSU, typename X, typename LFSV, typename R>
250  void alpha_skeleton (const IG& ig,
251  const LFSU& lfsu_s, const X& x_s, const LFSV& lfsv_s,
252  const LFSU& lfsu_n, const X& x_n, const LFSV& lfsv_n,
253  R& r_s, R& r_n) const
254  {
255  // domain and range field type
256  typedef typename LFSU::Traits::FiniteElementType::
257  Traits::LocalBasisType::Traits::DomainFieldType DF;
258  typedef typename LFSU::Traits::FiniteElementType::
259  Traits::LocalBasisType::Traits::RangeFieldType RF;
260 
261  // face geometry
262  const Dune::FieldVector<DF,IG::dimension-1>&
263  face_local = Dune::ReferenceElements<DF,IG::dimension-1>::general(ig.geometry().type()).position(0,0);
264  RF face_volume = ig.geometry().volume();
265 
266  // face center in element coordinates
267  Dune::FieldVector<DF,IG::dimension> face_center_in_element = ig.geometryInInside().global(face_local);
268 
269  // evaluate velocity
270  typename TP::Traits::RangeType v(tp.v(*(ig.inside()),face_center_in_element));
271 
272  // the normal velocity
273  RF vn = v*ig.centerUnitOuterNormal();
274 
275  // convective flux
276  RF u_upwind=0;
277  if (vn>=0) u_upwind = x_s(lfsu_s,0); else u_upwind = x_n(lfsu_n,0);
278  r_s.accumulate(lfsu_s,0,(u_upwind*vn)*face_volume);
279  if (vn>=0)
280  cellinflux += vn*face_volume; // dt computation
281 
282  // evaluate diffusion coefficients
283  const Dune::FieldVector<DF,IG::dimension>&
284  inside_local = Dune::ReferenceElements<DF,IG::dimension>::general(ig.inside()->type()).position(0,0);
285  const Dune::FieldVector<DF,IG::dimension>&
286  outside_local = Dune::ReferenceElements<DF,IG::dimension>::general(ig.outside()->type()).position(0,0);
287  typename TP::Traits::RangeFieldType D_inside = tp.D(*(ig.inside()),inside_local);
288  typename TP::Traits::RangeFieldType D_outside = tp.D(*(ig.outside()),outside_local);
289  typename TP::Traits::RangeFieldType D_avg = 2.0/(1.0/(D_inside+1E-30) + 1.0/(D_outside+1E-30));
290 
291  // distance between cell centers in global coordinates
292  Dune::FieldVector<DF,IG::dimension>
293  inside_global = ig.inside()->geometry().center();
294  Dune::FieldVector<DF,IG::dimension>
295  outside_global = ig.outside()->geometry().center();
296  inside_global -= outside_global;
297  RF distance = inside_global.two_norm();
298 
299  // diffusive flux
300  // note: we do only one-sided evaluation here
301  r_s.accumulate(lfsu_s,0,-(D_avg*(x_n(lfsu_n,0)-x_s(lfsu_s,0))/distance)*face_volume);
302  }
303 
304  // post skeleton: compute time step allowable for cell; to be done later
305  template<typename EG, typename LFSU, typename X, typename LFSV, typename R>
306  void alpha_volume_post_skeleton(const EG& eg, const LFSU& lfsu, const X& x,
307  const LFSV& lfsv, R& r) const
308  {
309  // domain and range field type
310  typedef typename LFSU::Traits::FiniteElementType::
311  Traits::LocalBasisType::Traits::DomainFieldType DF;
312  const int dim = EG::Geometry::dimension;
313 
314  if (!first_stage) return; // time step calculation is only done in first stage
315 
316  // cell center
317  const Dune::FieldVector<DF,dim>&
318  inside_local = Dune::ReferenceElements<DF,dim>::general(eg.entity().type()).position(0,0);
319 
320  // compute optimal dt for this cell
321  typename TP::Traits::RangeFieldType cellcapacity = tp.c(eg.entity(),inside_local)*eg.geometry().volume();
322  typename TP::Traits::RangeFieldType celldt = cellcapacity/(cellinflux+1E-30);
323  dtmin = std::min(dtmin,celldt);
324  }
325 
326  // skeleton integral depending on test and ansatz functions
327  // We put the Dirchlet evaluation also in the alpha term to save some geometry evaluations
328  template<typename IG, typename LFSU, typename X, typename LFSV, typename R>
329  void alpha_boundary (const IG& ig,
330  const LFSU& lfsu_s, const X& x_s, const LFSV& lfsv_s,
331  R& r_s) const
332  {
333  // domain and range field type
334  typedef typename LFSU::Traits::FiniteElementType::
335  Traits::LocalBasisType::Traits::DomainFieldType DF;
336  typedef typename LFSU::Traits::FiniteElementType::
337  Traits::LocalBasisType::Traits::RangeFieldType RF;
338 
339  // face geometry
340  const Dune::FieldVector<DF,IG::dimension-1>&
341  face_local = Dune::ReferenceElements<DF,IG::dimension-1>::general(ig.geometry().type()).position(0,0);
342  RF face_volume = ig.geometry().volume();
343  Dune::FieldVector<DF,dim> face_center_in_element = ig.geometryInInside().global(face_local);
344 
345  // evaluate boundary condition type
346  int bc = tp.bc(ig.intersection(),face_local);
347 
348  // do things depending on boundary condition type
349  if (bc==0) // Neumann boundary
350  {
351  typename TP::Traits::RangeFieldType j = tp.j(*(ig.inside()),face_center_in_element);
352  r_s.accumulate(lfsu_s,0,j*face_volume);
353  return;
354  }
355 
356  // evaluate velocity
357  typename TP::Traits::RangeType v(tp.v(*(ig.inside()),face_center_in_element));
358 
359  // the normal velocity
360  RF vn = v*ig.centerUnitOuterNormal();
361  if (vn>=0)
362  cellinflux += vn*face_volume; // dt computation
363 
364  if (bc==2) // Outflow boundary
365  {
366  r_s.accumulate(lfsu_s,0,vn*x_s(lfsu_s,0)*face_volume);
367  return;
368  }
369 
370  if (bc==1) // Dirichlet boundary
371  {
372  typename TP::Traits::RangeFieldType g;
373  if (vn>=0) g=x_s(lfsu_s,0); else g=tp.g(*(ig.inside()),face_center_in_element);
374  const Dune::FieldVector<DF,IG::dimension>&
375  inside_local = Dune::ReferenceElements<DF,IG::dimension>::general(ig.inside()->type()).position(0,0);
376  typename TP::Traits::RangeFieldType D_inside = tp.D(*(ig.inside()),inside_local);
377  Dune::FieldVector<DF,IG::dimension>
378  inside_global = ig.inside()->geometry().center();
379  Dune::FieldVector<DF,IG::dimension>
380  outside_global = ig.geometry().center();
381  inside_global -= outside_global;
382  RF distance = inside_global.two_norm();
383  r_s.accumulate(lfsu_s,0,(g*vn - D_inside*(g-x_s(lfsu_s,0))/distance)*face_volume);
384  return;
385  }
386  }
387 
388  // volume integral depending only on test functions
389  template<typename EG, typename LFSV, typename R>
390  void lambda_volume (const EG& eg, const LFSV& lfsv, R& r) const
391  {
392  // domain and range field type
393  typedef typename LFSV::Traits::FiniteElementType::
394  Traits::LocalBasisType::Traits::DomainFieldType DF;
395  const int dim = EG::Geometry::dimension;
396 
397  // cell center
398  const Dune::FieldVector<DF,dim>&
399  inside_local = Dune::ReferenceElements<DF,dim>::general(eg.entity().type()).position(0,0);
400 
401  // evaluate source term
402  typename TP::Traits::RangeFieldType q = tp.q(eg.entity(),inside_local);
403 
404  r.accumulate(lfsv,0,-q*eg.geometry().volume());
405  }
406 
408  void setTime (typename TP::Traits::RangeFieldType t)
409  {
410  }
411 
413  void preStep (typename TP::Traits::RangeFieldType time, typename TP::Traits::RangeFieldType dt,
414  int stages)
415  {
416  }
417 
419  void preStage (typename TP::Traits::RangeFieldType time, int r)
420  {
421  if (r==1)
422  {
423  first_stage = true;
424  dtmin = 1E100;
425  }
426  else first_stage = false;
427  }
428 
430  void postStage ()
431  {
432  }
433 
435  typename TP::Traits::RangeFieldType suggestTimestep (typename TP::Traits::RangeFieldType dt) const
436  {
437  return dtmin;
438  }
439 
440  private:
441  TP& tp;
442  bool first_stage;
443  mutable typename TP::Traits::RangeFieldType dtmin; // accumulate minimum dt here
444  mutable typename TP::Traits::RangeFieldType cellinflux;
445  };
446 
447 
449  template<class T, class Imp>
451  {
452  public:
453  typedef T Traits;
454 
456  typename Traits::RangeFieldType
457  c (const typename Traits::ElementType& e, const typename Traits::DomainType& x) const
458  {
459  return asImp().c(e,x);
460  }
461 
462  private:
463  Imp& asImp () {return static_cast<Imp &> (*this);}
464  const Imp& asImp () const {return static_cast<const Imp &>(*this);}
465  };
466 
473  template<class TP>
474  class CCFVTemporalOperator : public NumericalJacobianApplyVolume<CCFVTemporalOperator<TP> >,
475  public FullVolumePattern,
477  public InstationaryLocalOperatorDefaultMethods<typename TP::Traits::RangeFieldType>
478  {
479  public:
480  // pattern assembly flags
481  enum { doPatternVolume = true };
482 
483  // residual assembly flags
484  enum { doAlphaVolume = true };
485 
487  : tp(tp_)
488  {
489  }
490 
491  // volume integral depending on test and ansatz functions
492  template<typename EG, typename LFSU, typename X, typename LFSV, typename R>
493  void alpha_volume (const EG& eg, const LFSU& lfsu, const X& x, const LFSV& lfsv, R& r) const
494  {
495  // domain and range field type
496  typedef typename LFSU::Traits::FiniteElementType::
497  Traits::LocalBasisType::Traits::DomainFieldType DF;
498 
499  // dimensions
500  const int dim = EG::Geometry::dimension;
501 
502  // cell center
503  const Dune::FieldVector<DF,dim>&
504  inside_local = Dune::ReferenceElements<DF,dim>::general(eg.entity().type()).position(0,0);
505 
506  // evaluate capacity
507  typename TP::Traits::RangeFieldType c = tp.c(eg.entity(),inside_local);
508 
509  // residual contribution
510  r.accumulate(lfsu,0,c*x(lfsu,0)*eg.geometry().volume());
511  }
512 
513  // jacobian of volume term
514  template<typename EG, typename LFSU, typename X, typename LFSV, typename M>
515  void jacobian_volume (const EG& eg, const LFSU& lfsu, const X& x, const LFSV& lfsv,
516  M& mat) const
517  {
518  // domain and range field type
519  typedef typename LFSU::Traits::FiniteElementType::
520  Traits::LocalBasisType::Traits::DomainFieldType DF;
521 
522  // dimensions
523  const int dim = EG::Geometry::dimension;
524 
525  // cell center
526  const Dune::FieldVector<DF,dim>&
527  inside_local = Dune::ReferenceElements<DF,dim>::general(eg.entity().type()).position(0,0);
528 
529  // evaluate capacity
530  typename TP::Traits::RangeFieldType c = tp.c(eg.entity(),inside_local);
531 
532  // residual contribution
533  mat.accumulate(lfsu,0,lfsu,0,c*eg.geometry().volume());
534  }
535 
536  private:
537  TP& tp;
538  };
539 
540  }
541 }
542 
543 #endif
Dune::FieldVector< DomainFieldType, dimDomain > DomainType
domain type
Definition: transportccfv.hh:35
GV::Grid::ctype DomainFieldType
Export type for domain field.
Definition: transportccfv.hh:32
Definition: transportccfv.hh:474
void lambda_volume(const EG &eg, const LFSV &lfsv, R &r) const
Definition: transportccfv.hh:390
Dune::PDELab::BoundaryGridFunctionTraits< typename T::Traits::GridViewType, int, 1, Dune::FieldVector< int, 1 > > Traits
Definition: transportccfv.hh:127
Traits::RangeFieldType c(const typename Traits::ElementType &e, const typename Traits::DomainType &x) const
source term
Definition: transportccfv.hh:457
const E & e
Definition: interpolate.hh:172
Implement jacobian_apply_volume() based on alpha_volume()
Definition: defaultimp.hh:321
const Traits::GridViewType & getGridView() const
Definition: transportccfv.hh:175
void alpha_boundary(const IG &ig, const LFSU &lfsu_s, const X &x_s, const LFSV &lfsv_s, R &r_s) const
Definition: transportccfv.hh:329
Traits::RangeFieldType q(const typename Traits::ElementType &e, const typename Traits::DomainType &x) const
source term
Definition: transportccfv.hh:74
traits class holding function signature, same as in local function
Definition: function.hh:230
BoundaryConditionType_Transport(const typename T::Traits::GridViewType &gv_, const T &t_)
Definition: transportccfv.hh:130
traits class for two phase parameter class
Definition: transportccfv.hh:20
void preStep(typename TP::Traits::RangeFieldType time, typename TP::Traits::RangeFieldType dt, int stages)
to be called once before each time step
Definition: transportccfv.hh:413
GV GridViewType
The type of the grid view the function lives on.
Definition: function.hh:114
void alpha_volume(const EG &eg, const LFSU &lfsu, const X &x, const LFSV &lfsv, R &r) const
Definition: transportccfv.hh:493
TP::Traits::RangeFieldType suggestTimestep(typename TP::Traits::RangeFieldType dt) const
to be called once before each stage
Definition: transportccfv.hh:435
const T::Traits::GridViewType & getGridView()
get a reference to the GridView
Definition: transportccfv.hh:141
Dune::PDELab::GridFunctionTraits< typename T::Traits::GridViewType, typename T::Traits::RangeFieldType, 1, Dune::FieldVector< typename T::Traits::RangeFieldType, 1 > > Traits
Definition: transportccfv.hh:162
const IG & ig
Definition: constraints.hh:147
GV::Traits::template Codim< 0 >::Entity ElementType
grid types
Definition: transportccfv.hh:47
leaf of a function tree
Definition: function.hh:576
Traits::RangeType v(const typename Traits::ElementType &e, const typename Traits::DomainType &x) const
velocity field
Definition: transportccfv.hh:60
Dune::FieldVector< DomainFieldType, dimDomain-1 > IntersectionDomainType
domain type
Definition: transportccfv.hh:38
void evaluate(const Dune::PDELab::IntersectionGeometry< I > &ig, const typename Traits::DomainType &x, typename Traits::RangeType &y) const
Definition: transportccfv.hh:133
Default flags for all local operators.
Definition: flags.hh:18
int bc(const typename Traits::IntersectionType &is, const typename Traits::IntersectionDomainType &x) const
boundary condition type function
Definition: transportccfv.hh:86
RF RangeFieldType
Export type for range field.
Definition: transportccfv.hh:41
sparsity pattern generator
Definition: pattern.hh:13
void jacobian_volume(const EG &eg, const LFSU &lfsu, const X &x, const LFSV &lfsv, M &mat) const
Definition: transportccfv.hh:515
Dune::FieldVector< RF, GV::dimensionworld > RangeType
range type
Definition: transportccfv.hh:44
T Traits
Definition: transportccfv.hh:56
T Traits
Definition: transportccfv.hh:453
GV GridViewType
the grid view
Definition: transportccfv.hh:23
A local operator for a cell-centered finite volume scheme for the transport equation.
Definition: transportccfv.hh:199
Dune::FieldVector< GV::Grid::ctype, GV::dimension-1 > DomainType
domain type in dim-size coordinates
Definition: function.hh:48
void preStage(typename TP::Traits::RangeFieldType time, int r)
to be called once before each stage
Definition: transportccfv.hh:419
const I & intersection() const
Definition: geometrywrapper.hh:268
CCFVTemporalOperator(TP &tp_)
Definition: transportccfv.hh:486
traits class holding the function signature, same as in local function
Definition: function.hh:175
Dune::PDELab::BoundaryGridFunctionBase< Traits, BoundaryConditionType_Transport< T > > BaseT
Definition: transportccfv.hh:128
dimension of the domain
Definition: transportccfv.hh:28
Implement jacobian_apply_skeleton() based on alpha_skeleton()
Definition: defaultimp.hh:445
Definition: adaptivity.hh:27
Wrap intersection.
Definition: geometrywrapper.hh:56
Traits::RangeFieldType D(const typename Traits::ElementType &e, const typename Traits::DomainType &x) const
scalar diffusion coefficient
Definition: transportccfv.hh:67
leaf of a function tree
Definition: function.hh:596
Implement jacobian_boundary() based on alpha_boundary()
Definition: defaultimp.hh:251
void alpha_volume_post_skeleton(const EG &eg, const LFSU &lfsu, const X &x, const LFSV &lfsv, R &r) const
Definition: transportccfv.hh:306
static const int dim
Definition: adaptivity.hh:83
void setTime(typename TP::Traits::RangeFieldType t)
set time in parameter class
Definition: transportccfv.hh:408
DirichletBoundaryCondition_Transport(const typename Traits::GridViewType &g_, const T &t_)
constructor
Definition: transportccfv.hh:165
base class for parameter class
Definition: transportccfv.hh:450
GV::Traits::template Codim< 0 >::Entity ElementType
codim 0 entity
Definition: function.hh:117
GV::Intersection IntersectionType
Definition: transportccfv.hh:48
void evaluate(const typename Traits::ElementType &e, const typename Traits::DomainType &x, typename Traits::RangeType &y) const
Evaluate the GridFunction at given position.
Definition: transportccfv.hh:168
sparsity pattern generator
Definition: pattern.hh:29
Traits::RangeFieldType j(const typename Traits::ElementType &e, const typename Traits::DomainType &x) const
Neumann boundary condition.
Definition: transportccfv.hh:100
CCFVSpatialTransportOperator(TP &tp_)
Definition: transportccfv.hh:227
void alpha_volume(const EG &eg, const LFSU &lfsu, const X &x, const LFSV &lfsv, R &r) const
Definition: transportccfv.hh:234
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: transportccfv.hh:250
Traits::RangeFieldType g(const typename Traits::ElementType &e, const typename Traits::DomainType &x) const
Dirichlet boundary condition on inflow.
Definition: transportccfv.hh:93
Implement jacobian_skeleton() based on alpha_skeleton()
Definition: defaultimp.hh:157
Definition: transportccfv.hh:116
void postStage()
to be called once at the end of each stage
Definition: transportccfv.hh:430
base class for parameter class
Definition: transportccfv.hh:53
const EG & eg
Definition: constraints.hh:280
Default class for additional methods in instationary local operators.
Definition: idefault.hh:89
Implement jacobian_apply_boundary() based on alpha_boundary()
Definition: defaultimp.hh:538
void jacobian_volume(const EG &eg, const LFSU &lfsu, const X &x, const LFSV &lfsv, M &mat) const
Definition: transportccfv.hh:241