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