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