dune-fem  2.4.1-rc
discretemodel.hh
Go to the documentation of this file.
1 #ifndef DUNE_FEM_LOCALDG_DISCRETEMODEL_HH
2 #define DUNE_FEM_LOCALDG_DISCRETEMODEL_HH
3 
4 #include <cassert>
5 
6 #include <dune/common/bartonnackmanifcheck.hh>
7 
9 
10 namespace Dune
11 {
12  namespace Fem
13  {
14 
15  // DGDiscreteModelInterface
16  // ------------------------
17 
32  template <class DGDiscreteModelTraits>
34  {
36  typedef DGDiscreteModelTraits Traits;
38  typedef typename Traits::DGDiscreteModelType DGDiscreteModelType;
39 
41  typedef typename Traits::DiscreteFunctionSpaceType DiscreteFunctionSpaceType;
43  typedef typename DiscreteFunctionSpaceType::FunctionSpaceType FunctionSpaceType;
45  typedef typename FunctionSpaceType::DomainType DomainType;
47  typedef typename FunctionSpaceType::RangeType RangeType;
49  typedef typename FunctionSpaceType::RangeFieldType RangeFieldType;
51  typedef typename FunctionSpaceType::JacobianRangeType JacobianRangeType;
52 
54  typedef typename DiscreteFunctionSpaceType::GridPartType GridPartType;
56  typedef typename DiscreteFunctionSpaceType::EntityType EntityType;
58  typedef typename DiscreteFunctionSpaceType::IntersectionType IntersectionType;
59 
61  typedef typename EntityType::Geometry::LocalCoordinate LocalCoordinateType;
62 
64  enum { dimRange = DiscreteFunctionSpaceType::dimRange };
65 
67  template< class MatrixType > struct ExtractMatrix;
68 
69  template< template < class, int, int > class MatrixType,
70  class ctype, int dimD, int dimR
71  >
72  struct ExtractMatrix<MatrixType< ctype, dimD, dimR> >
73  {
74  typedef MatrixType< RangeFieldType, dimRange, dimRange > Type;
75  };
76 
79 
80  public:
85  bool hasFlux () const { return asImp().hasFlux(); }
86 
91  bool hasSource () const { return asImp().hasSource(); }
92 
94  bool hasMass () const { return asImp().hasMass(); }
95 
119  template <class ArgumentTuple, class FaceDomainType>
120  double numericalFlux ( const IntersectionType &intersection,
121  const double time,
122  const FaceDomainType &x,
123  const ArgumentTuple &uLeft,
124  const ArgumentTuple &uRight,
125  RangeType &gLeft,
126  RangeType &gRight)
127  {
128  CHECK_INTERFACE_IMPLEMENTATION( asImp().numericalFlux( intersection, time, x, uLeft, uRight, gLeft, gRight) );
129  return asImp().numericalFlux( intersection, time, x, uLeft, uRight, gLeft, gRight );
130  }
131 
132 
147  template <class ArgumentTuple, class FaceDomainType>
148  double boundaryFlux ( const IntersectionType &intersection,
149  const double time,
150  const FaceDomainType &x,
151  const ArgumentTuple &uLeft,
152  RangeType &gLeft )
153  {
154  CHECK_INTERFACE_IMPLEMENTATION( asImp().boundaryFlux( intersection, time, x, uLeft, gLeft ) );
155  return asImp().boundaryFlux( intersection, time, x, uLeft, gLeft );
156  }
157 
158 
169  template <class ArgumentTuple>
170  void analyticalFlux ( const EntityType& en,
171  const double time,
172  const LocalCoordinateType& x,
173  const ArgumentTuple& u,
174  JacobianRangeType& f )
175  {
176  CHECK_AND_CALL_INTERFACE_IMPLEMENTATION(
177  asImp().analyticalFlux(en, time, x, u, f) );
178  }
179 
193  template <class ArgumentTuple, class JacobianTuple>
194  double source ( const EntityType& en,
195  const double time,
196  const LocalCoordinateType& x,
197  const ArgumentTuple& u,
198  const JacobianTuple& jac,
199  RangeType& s )
200  {
201  CHECK_INTERFACE_IMPLEMENTATION(
202  asImp().source(en, time, x, u, jac, s) );
203  return asImp().source(en, time, x, u, jac, s);
204  }
205 
216  template <class ArgumentTuple>
217  void mass ( const EntityType& en,
218  const double time,
219  const LocalCoordinateType& x,
220  const ArgumentTuple& u,
221  MassFactorType& m )
222  {
223  CHECK_AND_CALL_INTERFACE_IMPLEMENTATION(
224  asImp().mass(en, time, x, u, m) );
225  }
226 
233  void setEntity ( EntityType& en )
234  {
235  CHECK_AND_CALL_INTERFACE_IMPLEMENTATION( asImp().setEntity(en) );
236  }
237 
244  void setNeighbor ( EntityType& nb )
245  {
246  CHECK_AND_CALL_INTERFACE_IMPLEMENTATION( asImp().setNeighbor(nb) );
247  }
248 
249  protected:
250  DGDiscreteModelType& asImp() { return static_cast<DGDiscreteModelType&>(*this); }
251  const DGDiscreteModelType& asImp() const { return static_cast<const DGDiscreteModelType&>(*this); }
252  };
253 
254 
255 
256  // DGDiscreteModelDefault
257  // ----------------------
258 
265  template< class DGDiscreteModelTraits
266  , int N1 = -1
267  , int N2 = -1
268  , int N3 = -1
269  , int N4 = -1
270  , int N5 = -1
271  , int N6 = -1
272  , int N7 = -1
273  , int N8 = -1
274  , int N9 = -1
275  >
277  public DGDiscreteModelInterface<DGDiscreteModelTraits>
278  {
280 
281  public:
306 
307  public:
311  typedef typename BaseType::RangeType RangeType;
314 
319  inline bool hasFlux () const { return false; }
320 
325  inline bool hasSource () const { return false; }
326 
331  inline bool hasMass () const { return false; }
332 
336  template <class ArgumentTuple, class FaceDomainType>
337  double numericalFlux ( const IntersectionType& it,
338  const double time,
339  const FaceDomainType& x,
340  const ArgumentTuple& uLeft,
341  const ArgumentTuple& uRight,
342  RangeType& gLeft,
343  RangeType& gRight )
344  {
345  assert(!this->asImp().hasFlux());
346  gLeft = 0.0;
347  gRight = 0.0;
348  return 0.0;
349  }
350 
354  template <class ArgumentTuple, class FaceDomainType>
355  double boundaryFlux ( const IntersectionType &intersection,
356  const double time,
357  const FaceDomainType &x,
358  const ArgumentTuple &uLeft,
359  RangeType &gLeft )
360  {
361  assert( !this->asImp().hasFlux() );
362  gLeft = 0.0;
363  return 0.0;
364  }
365 
369  template< class ArgumentTuple >
370  void analyticalFlux ( const EntityType &entity,
371  const double time,
372  const LocalCoordinateType &x,
373  const ArgumentTuple &u,
374  JacobianRangeType &f )
375  {
376  assert( !this->asImp().hasFlux() );
377  f = 0.0;
378  }
379 
383  template <class ArgumentTuple, class JacobianTuple>
384  double source ( const EntityType& en,
385  const double time,
386  const LocalCoordinateType& x,
387  const ArgumentTuple& u,
388  const JacobianTuple& jac,
389  RangeType& s )
390  {
391  assert(!this->asImp().hasSource());
392  s = 0.0;
393  return 0.0;
394  }
395 
399  template <class ArgumentTuple>
400  void mass (const EntityType& en,
401  const double time,
402  const LocalCoordinateType& x,
403  const ArgumentTuple& u,
404  MassFactorType& m )
405  {
406  enum { rows = MassFactorType :: rows };
407  enum { cols = MassFactorType :: cols };
408  // default implementation sets mass factor to identity
409  for(int i=0; i<rows; ++i)
410  {
411  // set diagonal to 1
412  m[i][i] = 1;
413  // and all other values to 0
414  for(int j=i+1; j<cols; ++j)
415  m[i][j] = m[j][i] = 0;
416  }
417  }
418 
420  void setEntity( const EntityType& en )
421  { }
422 
424  void setNeighbor( const EntityType& nb )
425  { }
426  };
427 
428  // DGDGAdaptiveDiscreteModel
429  // ---------------------------------------
430 
432  {
433  public:
435  template <class Adaptation>
436  void setAdaptation( Adaptation&,
437  const double weight = 1.0 ) {}
438 
440  template <class Adaptation, class ThreadFilter>
441  void setAdaptation( Adaptation&,
442  const ThreadFilter&,
443  const double weight = 1.0 ) {}
444 
447  };
448 
449 
450  // DGDiscreteModelDefaultWithInsideOutside
451  // ---------------------------------------
452 
456  template <class DGDiscreteModelTraits
457  , int N1 = -1
458  , int N2 = -1
459  , int N3 = -1
460  , int N4 = -1
461  , int N5 = -1
462  , int N6 = -1
463  , int N7 = -1
464  , int N8 = -1
465  , int N9 = -1
466  >
468  : public DGDiscreteModelDefault< DGDiscreteModelTraits, N1, N2, N3, N4, N5, N6, N7, N8, N9 >,
470  {
472 
473  public:
475 
477  : enVol_(-1.0) , nbVol_(-1.0) , en_(0) , nb_(0)
478  {}
479 
481  : enVol_(other.enVol_) , nbVol_(other.nbVol_),
482  en_(other.en_) , nb_(other.nb_)
483  {}
484 
489  void setEntity ( const EntityType& en )
490  {
491  en_ = &en;
492  enVol_ = en.geometry().volume();
493  }
494 
499  void setNeighbor ( const EntityType& nb )
500  {
501  nb_ = &nb;
502  nbVol_ = nb.geometry().volume();
503  }
504 
509  const EntityType &inside () const
510  {
511  assert( en_ );
512  return *en_;
513  }
514 
519  const EntityType &outside () const
520  {
521  assert( nb_ );
522  return *nb_;
523  }
524 
526  double enVolume () const
527  {
528  assert(enVol_ > 0.0);
529  return enVol_;
530  }
531 
533  double nbVolume () const
534  {
535  assert( nbVol_ > 0.0 );
536  return nbVol_;
537  }
538 
539  private:
540  double enVol_;
541  double nbVol_;
542 
543  const EntityType* en_;
544  const EntityType* nb_;
545  };
546 
547  } // namespace Fem
548 
549 } // namespace Dune
550 
551 #endif // #ifndef DUNE_FEM_LOCALDG_DISCRETEMODEL_HH
bool hasSource() const
Definition: discretemodel.hh:325
double source(const EntityType &en, const double time, const LocalCoordinateType &x, const ArgumentTuple &u, const JacobianTuple &jac, RangeType &s)
Empty implementation that fails if problem claims to have a source term.
Definition: discretemodel.hh:384
Dune::Fem::Selector< N1, N2, N3, N4, N5, N6, N7, N8, N9 >::Type Selector
Definition: discretemodel.hh:305
void setEntity(const EntityType &en)
empty implementation
Definition: discretemodel.hh:420
double enVolume() const
return volume of entity
Definition: discretemodel.hh:526
BaseType::EntityType EntityType
Definition: discretemodel.hh:474
void setNeighbor(EntityType &nb)
Passes the active neigbor entity to the model. This can be used, to set local functions required as d...
Definition: discretemodel.hh:244
Traits::DGDiscreteModelType DGDiscreteModelType
Implementation type for Barton-Nackman trick.
Definition: discretemodel.hh:38
void analyticalFlux(const EntityType &entity, const double time, const LocalCoordinateType &x, const ArgumentTuple &u, JacobianRangeType &f)
Empty implementation that fails if problem claims to have a flux contribution.
Definition: discretemodel.hh:370
void removeAdaptation()
remove pointer to adaptation handle
Definition: discretemodel.hh:446
void setEntity(const EntityType &en)
method setting pointer of inside entity and getting volume
Definition: discretemodel.hh:489
const DGDiscreteModelType & asImp() const
Definition: discretemodel.hh:251
ExtractMatrix< JacobianRangeType >::Type MassFactorType
Type of mass factor.
Definition: discretemodel.hh:78
void analyticalFlux(const EntityType &en, const double time, const LocalCoordinateType &x, const ArgumentTuple &u, JacobianRangeType &f)
Computes the analytical flux of the problem. Analytical flux of the problem.
Definition: discretemodel.hh:170
Mass factor type.
Definition: discretemodel.hh:67
BaseType::JacobianRangeType JacobianRangeType
Definition: discretemodel.hh:312
DiscreteFunctionSpaceType::EntityType EntityType
Element (codim 0 entity) of the grid.
Definition: discretemodel.hh:56
void setNeighbor(const EntityType &nb)
method seting pointer of outside entity and getting volume
Definition: discretemodel.hh:499
BaseType::EntityType EntityType
Definition: discretemodel.hh:308
double numericalFlux(const IntersectionType &it, const double time, const FaceDomainType &x, const ArgumentTuple &uLeft, const ArgumentTuple &uRight, RangeType &gLeft, RangeType &gRight)
Empty implementation that fails if problem claims to have a flux contribution.
Definition: discretemodel.hh:337
Interface class for problem definition in the LDG context.
Definition: discretemodel.hh:33
FunctionSpaceType::RangeType RangeType
Vector type of the function space&#39;s range.
Definition: discretemodel.hh:47
DGDiscreteModelType & asImp()
Definition: discretemodel.hh:250
Default implementation of the DGDiscreteModelInterface where methods for the fluxes and the source te...
Definition: discretemodel.hh:276
BaseType::IntersectionType IntersectionType
Definition: discretemodel.hh:309
bool hasMass() const
Returns true if problem has a mass matrix factor disctinct from the identity.
Definition: discretemodel.hh:94
bool hasSource() const
Returns true if problem has a source term. If you program this method to return true, make sure to implement method source well.
Definition: discretemodel.hh:91
BaseType::LocalCoordinateType LocalCoordinateType
Definition: discretemodel.hh:310
void setEntity(EntityType &en)
Passes the active entity to the model. This can be used, to set local functions required as data func...
Definition: discretemodel.hh:233
Traits::DiscreteFunctionSpaceType DiscreteFunctionSpaceType
Discrete function space.
Definition: discretemodel.hh:41
DiscreteFunctionSpaceType::GridPartType GridPartType
Type of GridPart.
Definition: discretemodel.hh:54
double numericalFlux(const IntersectionType &intersection, const double time, const FaceDomainType &x, const ArgumentTuple &uLeft, const ArgumentTuple &uRight, RangeType &gLeft, RangeType &gRight)
Computes the numerical flux at a cell interface.
Definition: discretemodel.hh:120
FunctionSpaceType::DomainType DomainType
Coordinate type (world coordinates)
Definition: discretemodel.hh:45
Definition: coordinate.hh:4
FunctionSpaceType::JacobianRangeType JacobianRangeType
Matrix type of the function space&#39;s jacobian.
Definition: discretemodel.hh:51
DGDiscreteModelTraits Traits
Traits class defined by the user.
Definition: discretemodel.hh:36
bool hasMass() const
Definition: discretemodel.hh:331
void setAdaptation(Adaptation &, const double weight=1.0)
default method for setting adaptation handle to discrete model
Definition: discretemodel.hh:436
Dune::Fem::ElementTuple< N1, N2, N3, N4, N5, N6, N7, N8, N9 >::Type Type
tuple consisting of Dune::integral_constant< int, N_i >
Definition: selector.hh:71
DGDiscreteModelDefaultWithInsideOutside(const DGDiscreteModelDefaultWithInsideOutside &other)
Definition: discretemodel.hh:480
double boundaryFlux(const IntersectionType &intersection, const double time, const FaceDomainType &x, const ArgumentTuple &uLeft, RangeType &gLeft)
Computes the flux at the boundary Special kind of numerical flux. The intersection iterator provides ...
Definition: discretemodel.hh:148
DiscreteFunctionSpaceType::FunctionSpaceType FunctionSpaceType
Function space type.
Definition: discretemodel.hh:43
Definition: threadfilter.hh:20
double boundaryFlux(const IntersectionType &intersection, const double time, const FaceDomainType &x, const ArgumentTuple &uLeft, RangeType &gLeft)
Empty implementation that fails if problem claims to have a flux contribution.
Definition: discretemodel.hh:355
double source(const EntityType &en, const double time, const LocalCoordinateType &x, const ArgumentTuple &u, const JacobianTuple &jac, RangeType &s)
Implements the source term of the problem. Source term contribution. The source term can contain part...
Definition: discretemodel.hh:194
BaseType::MassFactorType MassFactorType
Definition: discretemodel.hh:313
DGDiscreteModelDefaultWithInsideOutside()
Definition: discretemodel.hh:476
EntityType::Geometry::LocalCoordinate LocalCoordinateType
Type of local coordinate.
Definition: discretemodel.hh:61
const EntityType & inside() const
method returning reference to inside entity
Definition: discretemodel.hh:509
bool hasFlux() const
Definition: discretemodel.hh:319
void mass(const EntityType &en, const double time, const LocalCoordinateType &x, const ArgumentTuple &u, MassFactorType &m)
Implements the mass factor term of the problem. Mass term contribution.
Definition: discretemodel.hh:217
Definition: discretemodel.hh:64
void mass(const EntityType &en, const double time, const LocalCoordinateType &x, const ArgumentTuple &u, MassFactorType &m)
empty implementation for mass factor default implementation sets this factor to 1.0
Definition: discretemodel.hh:400
bool hasFlux() const
Returns true if problem has a flux contribution. If you program this method to return true...
Definition: discretemodel.hh:85
BaseType::RangeType RangeType
Definition: discretemodel.hh:311
const EntityType & outside() const
method returning reference to outside entity
Definition: discretemodel.hh:519
DiscreteFunctionSpaceType::IntersectionType IntersectionType
Intersection type.
Definition: discretemodel.hh:58
double nbVolume() const
return volume of neighbor
Definition: discretemodel.hh:533
void setAdaptation(Adaptation &, const ThreadFilter &, const double weight=1.0)
default method for setting adaptation handle and thead filter to discrete model
Definition: discretemodel.hh:441
FunctionSpaceType::RangeFieldType RangeFieldType
Vector type of the function space&#39;s range field type.
Definition: discretemodel.hh:49
void setNeighbor(const EntityType &nb)
empty implementation
Definition: discretemodel.hh:424
Definition: discretemodel.hh:431
MatrixType< RangeFieldType, dimRange, dimRange > Type
Definition: discretemodel.hh:74