dune-fem  2.8-git
diffusionmodel.hh
Go to the documentation of this file.
1 #ifndef ELLIPTC_MODEL_HH
2 #define ELLIPTC_MODEL_HH
3 
4 #include <cassert>
5 #include <cmath>
6 
7 #include <dune/common/visibility.hh>
8 
10 #include <dune/fem/io/parameter.hh>
14 #include <dune/fempy/quadrature/fempyquadratures.hh>
15 
16 #define VirtualDiffusionModelMethods(POINT) \
17  virtual void source ( const POINT &x,\
18  const DRangeType &value,\
19  const DJacobianRangeType &gradient,\
20  RRangeType &flux ) const = 0;\
21  virtual void linSource ( const DRangeType& uBar,\
22  const DJacobianRangeType &gradientBar,\
23  const POINT &x,\
24  const DRangeType &value,\
25  const DJacobianRangeType &gradient,\
26  RRangeType &flux ) const = 0;\
27  virtual void diffusiveFlux ( const POINT &x,\
28  const DRangeType &value,\
29  const DJacobianRangeType &gradient,\
30  RJacobianRangeType &flux ) const = 0;\
31  virtual void linDiffusiveFlux ( const DRangeType& uBar,\
32  const DJacobianRangeType& gradientBar,\
33  const POINT &x,\
34  const DRangeType &value,\
35  const DJacobianRangeType &gradient,\
36  RJacobianRangeType &flux ) const = 0;\
37  virtual void fluxDivergence( const POINT &x,\
38  const DRangeType &value,\
39  const DJacobianRangeType &jacobian,\
40  const DHessianRangeType &hessian,\
41  RRangeType &flux) const = 0;\
42  virtual void alpha(const POINT &x,\
43  const DRangeType &value,\
44  RRangeType &val) const = 0;\
45  virtual void linAlpha(const DRangeType &uBar,\
46  const POINT &x,\
47  const DRangeType &value,\
48  RRangeType &val) const = 0;\
49  virtual void dirichlet( int bndId, const POINT &x,\
50  RRangeType &value) const = 0;
51 
52 #define WrapperDiffusionModelMethods(POINT) \
53  virtual void source ( const POINT &x,\
54  const DRangeType &value,\
55  const DJacobianRangeType &gradient,\
56  RRangeType &flux ) const \
57  { impl().source(x, value, gradient, flux); } \
58  virtual void linSource ( const DRangeType& uBar,\
59  const DJacobianRangeType &gradientBar,\
60  const POINT &x,\
61  const DRangeType &value,\
62  const DJacobianRangeType &gradient,\
63  RRangeType &flux ) const \
64  { impl().linSource(uBar, gradientBar, x, value, gradient, flux); } \
65  virtual void diffusiveFlux ( const POINT &x,\
66  const DRangeType &value,\
67  const DJacobianRangeType &gradient,\
68  RJacobianRangeType &flux ) const \
69  { impl().diffusiveFlux(x, value, gradient, flux); } \
70  virtual void linDiffusiveFlux ( const DRangeType& uBar,\
71  const DJacobianRangeType& gradientBar,\
72  const POINT &x,\
73  const DRangeType &value,\
74  const DJacobianRangeType &gradient,\
75  RJacobianRangeType &flux ) const \
76  { impl().linDiffusiveFlux(uBar, gradientBar, x, value, gradient, flux); } \
77  virtual void fluxDivergence( const POINT &x,\
78  const DRangeType &value,\
79  const DJacobianRangeType &jacobian,\
80  const DHessianRangeType &hessian,\
81  RRangeType &flux) const \
82  { impl().fluxDivergence(x, value, jacobian, hessian, flux); } \
83  virtual void alpha(const POINT &x,\
84  const DRangeType &value,\
85  RRangeType &val) const \
86  { impl().alpha(x, value, val); } \
87  virtual void linAlpha(const DRangeType &uBar,\
88  const POINT &x,\
89  const DRangeType &value,\
90  RRangeType &val) const \
91  { impl().linAlpha(uBar, x, value, val); } \
92  virtual void dirichlet( int bndId, const POINT &x,\
93  RRangeType &value) const \
94  { impl().dirichlet(bndId, x, value); }
95 
96 template< class GridPart, int dimDomain, int dimRange=dimDomain, class RangeField = double >
98 {
99  typedef GridPart GridPartType;
100  static const int dimD = dimDomain;
101  static const int dimR = dimRange;
103  typedef RangeField RangeFieldType;
104 
105  typedef Dune::Fem::FunctionSpace< double, RangeFieldType,
106  GridPart::dimensionworld, dimD > DFunctionSpaceType;
107  typedef Dune::Fem::FunctionSpace< double, RangeFieldType,
108  GridPart::dimensionworld, dimR > RFunctionSpaceType;
118 
119  typedef typename GridPartType::template Codim<0>::EntityType EntityType;
120  typedef typename GridPartType::IntersectionType IntersectionType;
121  typedef typename EntityType::Geometry::LocalCoordinate LocalDomainType;
122 
123  // quadrature points from dune-fempy quadratures
124  template <class F,int d>
125  using Traits = Dune::FemPy::FempyQuadratureTraits<F,d>;
127  QuadraturePointWrapperType Point;
129  QuadraturePointWrapperType IntersectionPoint;
131  QuadraturePointWrapperType ElementPoint;
133  QuadraturePointWrapperType ElementIntersectionPoint;
134 
135  // quadrature points from dune-fem quadratures
137  QuadraturePointWrapperType OriginalPoint;
139  QuadraturePointWrapperType OriginalIntersectionPoint;
141  QuadraturePointWrapperType OriginalElementPoint;
143  QuadraturePointWrapperType OriginalElementIntersectionPoint;
144 
145  /*
146  static const bool isLinear;
147  static const bool isSymmetric;
148  */
149 
150 public:
152  { }
153  virtual ~DiffusionModel() {}
154 
155  virtual std::string name() const = 0;
156 
157  virtual bool init( const EntityType &entity) const = 0;
158 
159  // if the ModelImpl provides a time method (constant) then
160  // it maybe be set using this method
161  virtual void setTime( const double t ) const = 0;
162 
163  // if the ModelImpl provides a 'double& time()' method then
164  // it maybe retrieved by this method
165  virtual double time() const = 0;
166 
167  // virtual methods for fempy quadratures
172 
173  // virtual methods for fem quadratures
178 
180 
181  typedef std::array<int, dimRange> DirichletComponentType;
182  virtual bool hasDirichletBoundary () const = 0;
183  virtual bool hasNeumanBoundary () const = 0;
184  virtual bool isDirichletIntersection( const IntersectionType& inter, DirichletComponentType &dirichletComponent ) const = 0;
185 };
186 
187 namespace detail {
188 
189  template <class M>
190  class CheckTimeMethod
191  {
192  template <class T>
193  static std::true_type testSignature(double& (T::*)());
194 
195  template <class T>
196  static decltype(testSignature(&T::time)) test(std::nullptr_t);
197 
198  template <class T>
199  static std::false_type test(...);
200 
201  using type = decltype(test<M>(nullptr));
202  public:
203  static const bool value = type::value;
204  };
205 
206 
207  template <class Model, bool>
208  struct CallSetTime
209  {
210  static void setTime( Model&, const double ) {}
211  static double time( const Model& ) { return -1.0; }
212  };
213 
214  template <class Model>
215  struct CallSetTime< Model, true >
216  {
217  static void setTime( Model& model, const double t ) { model.time() = t; }
218  static double time ( const Model& model ) { return model.time(); }
219  };
220 } // end namespace detail
221 
222 
223 template < class ModelImpl >
224 struct DiffusionModelWrapper : public DiffusionModel<typename ModelImpl::GridPartType, ModelImpl::dimD, ModelImpl::dimR,
225  typename ModelImpl::RRangeFieldType>
226 {
227  typedef ModelImpl Impl;
228  typedef typename ModelImpl::GridPartType GridPartType;
229  static const int dimD = ModelImpl::dimD;
230  static const int dimR = ModelImpl::dimR;
232 
233  typedef typename Base::Point Point;
242  typedef typename Base::DomainType DomainType;
243  typedef typename Base::DRangeType DRangeType;
246  typedef typename Base::RRangeType RRangeType;
249  typedef typename Base::EntityType EntityType;
251 
252  template< class... Args, std::enable_if_t< std::is_constructible< ModelImpl, Args &&... >::value, int > = 0 >
253  explicit DiffusionModelWrapper ( Args &&... args )
254  : impl_( std::forward< Args >( args )... )
255  {}
256 
258  {
259  }
260 
261  // virtual methods for fempy quadratures
266 
267  // virtual methods for fem quadratures
272 
274 
275  // other virtual functions
276  virtual std::string name() const
277  {
278  return impl().name();
279  }
280 
281  // if the ModelImpl provides a 'double& time()' method then
282  // it maybe be set using this method
283  virtual void setTime( const double t ) const
284  {
285  detail::CallSetTime< ModelImpl, detail::CheckTimeMethod< ModelImpl >::value >
286  ::setTime( const_cast< ModelImpl& > (impl()), t );
287  }
288 
289  // if the ModelImpl provides a 'double& time()' method then
290  // it maybe be set using this method
291  virtual double time() const
292  {
293  return detail::CallSetTime< ModelImpl, detail::CheckTimeMethod< ModelImpl >::value >::time( impl() );
294  }
295 
296  typedef std::array<int, dimR> DirichletComponentType;
297  virtual bool hasDirichletBoundary () const
298  {
299  return impl().hasDirichletBoundary();
300  }
301  virtual bool hasNeumanBoundary () const
302  {
303  return impl().hasNeumanBoundary();
304  }
305  virtual bool isDirichletIntersection( const IntersectionType& inter, DirichletComponentType &dirichletComponent ) const
306  {
307  return impl().isDirichletIntersection(inter, dirichletComponent);
308  }
309  virtual bool init( const EntityType &entity) const
310  {
311  return impl().init(entity);
312  }
313  const ModelImpl &impl() const
314  {
315  return impl_;
316  }
317  ModelImpl &impl()
318  {
319  return impl_;
320  }
321  private:
322  ModelImpl impl_;
323 };
324 
325 #define VirtualPenaltyMethods(POINT) \
326  virtual RRangeType penalty ( const POINT &x,\
327  const DRangeType &value ) const = 0;\
328  virtual RRangeType linPenalty ( const POINT &x,\
329  const DRangeType &value ) const = 0;
330 #define WrapperPenaltyMethods(POINT) \
331  virtual RRangeType penalty ( const POINT &x,\
332  const DRangeType &value ) const \
333  { return impl().penalty(x, value); } \
334  virtual RRangeType linPenalty ( const POINT &x,\
335  const DRangeType &value ) const \
336  { return impl().linPenalty( x, value); }
337 
338 template< class GridPart, int dimDomain, int dimRange=dimDomain, class RangeField = double >
340 {
341  typedef GridPart GridPartType;
342  static const int dimD = dimDomain;
343  static const int dimR = dimRange;
345  typedef RangeField RangeFieldType;
346 
347  typedef Dune::Fem::FunctionSpace< double, RangeFieldType,
348  GridPart::dimensionworld, dimD > DFunctionSpaceType;
349  typedef Dune::Fem::FunctionSpace< double, RangeFieldType,
350  GridPart::dimensionworld, dimR > RFunctionSpaceType;
360 
361  typedef typename GridPartType::template Codim<0>::EntityType EntityType;
362  typedef typename GridPartType::IntersectionType IntersectionType;
363  typedef typename EntityType::Geometry::LocalCoordinate LocalDomainType;
364 
365  // quadrature points from dune-fempy quadratures
366  template <class F,int d>
367  using Traits = Dune::FemPy::FempyQuadratureTraits<F,d>;
369  QuadraturePointWrapperType Point;
371  QuadraturePointWrapperType IntersectionPoint;
373  QuadraturePointWrapperType ElementPoint;
375  QuadraturePointWrapperType ElementIntersectionPoint;
376 
377  // quadrature points from dune-fem quadratures
379  QuadraturePointWrapperType OriginalPoint;
381  QuadraturePointWrapperType OriginalIntersectionPoint;
383  QuadraturePointWrapperType OriginalElementPoint;
385  QuadraturePointWrapperType OriginalElementIntersectionPoint;
386 
387  /*
388  static const bool isLinear;
389  static const bool isSymmetric;
390  */
391 
392 public:
394  { }
395  virtual ~DGDiffusionModel() {}
396 
397  virtual std::string name() const = 0;
398 
399  virtual bool init( const EntityType &entity) const = 0;
400 
401  // if the ModelImpl provides a 'double& time()' method then
402  // it maybe be set using this method
403  virtual void setTime( const double t ) const = 0;
404 
405  // if the ModelImpl provides a 'double& time()' method then
406  // it maybe retrieved by this method
407  virtual double time() const = 0;
408 
409  // virtual methods for fempy quadratures
414 
415  // virtual methods for fem quadratures
420 
422 
423  typedef std::array<int, dimRange> DirichletComponentType;
424  virtual bool hasDirichletBoundary () const = 0;
425  virtual bool hasNeumanBoundary () const = 0;
426  virtual bool isDirichletIntersection( const IntersectionType& inter, DirichletComponentType &dirichletComponent ) const = 0;
427 
432 
433  // virtual methods for fem quadratures
438 
440 
441 };
442 
443 template < class ModelImpl >
444 struct DGDiffusionModelWrapper : public DGDiffusionModel<typename ModelImpl::GridPartType, ModelImpl::dimD, ModelImpl::dimR,
445  typename ModelImpl::RRangeFieldType>
446 {
447  typedef ModelImpl Impl;
448  typedef typename ModelImpl::GridPartType GridPartType;
449  static const int dimD = ModelImpl::dimD;
450  static const int dimR = ModelImpl::dimR;
452 
453  typedef typename Base::Point Point;
462  typedef typename Base::DomainType DomainType;
463  typedef typename Base::DRangeType DRangeType;
466  typedef typename Base::RRangeType RRangeType;
469  typedef typename Base::EntityType EntityType;
471 
472  template< class... Args, std::enable_if_t< std::is_constructible< ModelImpl, Args &&... >::value, int > = 0 >
473  explicit DGDiffusionModelWrapper ( Args &&... args )
474  : impl_( std::forward< Args >( args )... )
475  {}
476 
478  {
479  }
480 
481  // virtual methods for fempy quadratures
486 
487  // virtual methods for fem quadratures
492 
494 
499 
500  // virtual methods for fem quadratures
505 
507  // other virtual functions
508  virtual std::string name() const
509  {
510  return impl().name();
511  }
512 
513  // if the ModelImpl provides a 'double& time()' method then
514  // it maybe be set using this method
515  virtual void setTime( const double t ) const
516  {
517  detail::CallSetTime< ModelImpl, detail::CheckTimeMethod< ModelImpl >::value >
518  ::setTime( const_cast< ModelImpl& > (impl()), t );
519  }
520 
521  // if the ModelImpl provides a 'double& time()' method then
522  // it maybe be set using this method
523  virtual double time() const
524  {
525  return detail::CallSetTime< ModelImpl, detail::CheckTimeMethod< ModelImpl >::value >::time( impl() );
526  }
527 
528  typedef std::array<int, dimR> DirichletComponentType;
529  virtual bool hasDirichletBoundary () const
530  {
531  return impl().hasDirichletBoundary();
532  }
533  virtual bool hasNeumanBoundary () const
534  {
535  return impl().hasNeumanBoundary();
536  }
537  virtual bool isDirichletIntersection( const IntersectionType& inter, DirichletComponentType &dirichletComponent ) const
538  {
539  return impl().isDirichletIntersection(inter, dirichletComponent);
540  }
541  virtual bool init( const EntityType &entity) const
542  {
543  return impl().init(entity);
544  }
545  const ModelImpl &impl() const
546  {
547  return impl_;
548  }
549  ModelImpl &impl()
550  {
551  return impl_;
552  }
553  private:
554  ModelImpl impl_;
555 };
556 
557 #if 0 // doesn't seem to be used anymore - needs checking
558 template < class ModelTraits >
559 struct DiffusionModelEngine : public DiffusionModel<typename ModelTraits::GridPartType, ModelTraits::dimD, ModelTraits::dimR, typename ModelTraits::RRangeFieldType>
560 {
561  typedef typename ModelTraits::GridPartType GridPartType;
562  static const int dimD = ModelTraits::dimD;
563  static const int dimR = ModelTraits::dimR;
565 
566  typedef typename Base::Point Point;
567  typedef typename Base::IntersectionPoint IntersectionPoint;
568  typedef typename Base::ElementPoint ElementPoint;
569  typedef typename Base::ElementIntersectionPoint ElementIntersectionPoint;
570  typedef typename Base::LocalDomainType LocalDomainType;
571  typedef typename Base::DomainType DomainType;
572  typedef typename Base::DRangeType DRangeType;
573  typedef typename Base::DJacobianRangeType DJacobianRangeType;
574  typedef typename Base::DHessianRangeType DHessianRangeType;
575  typedef typename Base::RRangeType RRangeType;
576  typedef typename Base::RJacobianRangeType RJacobianRangeType;
577  typedef typename Base::RHessianRangeType RHessianRangeType;
578  typedef typename Base::EntityType EntityType;
579  typedef typename Base::IntersectionType IntersectionType;
580 
581  typedef typename ModelTraits::ModelType ModelImpl;
582 
583  DiffusionModelEngine(ModelImpl &model) : impl_(model) {}
584  ~DiffusionModelEngine()
585  {
586  }
587 
589  WrapperDiffusionModelMethods(ElementPoint);
590  WrapperDiffusionModelMethods(IntersectionPoint);
591  WrapperDiffusionModelMethods(ElementIntersectionPoint);
592  WrapperDiffusionModelMethods(LocalDomainType);
593 
594  // other virtual functions
595  virtual std::string name() const
596  {
597  return impl().name();
598  }
599  virtual bool hasDirichletBoundary () const
600  {
601  return impl().hasDirichletBoundary();
602  }
603  virtual bool hasNeumanBoundary () const
604  {
605  return impl().hasNeumanBoundary();
606  }
607  virtual bool isDirichletIntersection( const IntersectionType& inter, std::array<int, dimR> &dirichletComponent ) const
608  {
609  return impl().isDirichletIntersection(inter, dirichletComponent);
610  }
611  virtual bool init( const EntityType &entity) const
612  {
613  return impl().init(entity);
614  }
615  const ModelImpl &impl() const
616  {
617  return impl_;
618  }
619  ModelImpl &impl()
620  {
621  return impl_;
622  }
623  private:
624  ModelImpl &impl_;
625 };
626 #endif
627 
628 #endif // #ifndef ELLIPTC_MODEL_HH
#define WrapperDiffusionModelMethods(POINT)
Definition: diffusionmodel.hh:52
#define VirtualDiffusionModelMethods(POINT)
Definition: diffusionmodel.hh:16
#define VirtualPenaltyMethods(POINT)
Definition: diffusionmodel.hh:325
Definition: explicitfieldvector.hh:75
quadrature class supporting base function caching
Definition: cachingquadrature.hh:41
quadrature on the codim-0 reference element
Definition: elementquadrature.hh:18
Definition: diffusionmodel.hh:98
DiffusionModel()
Definition: diffusionmodel.hh:151
virtual VirtualDiffusionModelMethods(Point) VirtualDiffusionModelMethods(ElementPoint) VirtualDiffusionModelMethods(IntersectionPoint) VirtualDiffusionModelMethods(ElementIntersectionPoint) VirtualDiffusionModelMethods(OriginalPoint) VirtualDiffusionModelMethods(OriginalElementPoint) VirtualDiffusionModelMethods(OriginalIntersectionPoint) VirtualDiffusionModelMethods(OriginalElementIntersectionPoint) VirtualDiffusionModelMethods(LocalDomainType) typedef std bool hasDirichletBoundary() const =0
Dune::Fem::ElementQuadrature< GridPartType, 1 >::QuadraturePointWrapperType OriginalElementIntersectionPoint
Definition: diffusionmodel.hh:143
DiffusionModel< GridPartType, dimD, dimR, RangeField > ModelType
Definition: diffusionmodel.hh:102
Dune::Fem::ElementQuadrature< GridPartType, 1, Traits >::QuadraturePointWrapperType ElementIntersectionPoint
Definition: diffusionmodel.hh:133
DFunctionSpaceType::DomainFieldType DDomainFieldType
Definition: diffusionmodel.hh:113
RFunctionSpaceType::DomainFieldType rDomainFieldType
Definition: diffusionmodel.hh:117
Dune::Fem::CachingQuadrature< GridPartType, 1 >::QuadraturePointWrapperType OriginalIntersectionPoint
Definition: diffusionmodel.hh:139
virtual std::string name() const =0
virtual bool hasNeumanBoundary() const =0
virtual bool init(const EntityType &entity) const =0
Dune::Fem::ElementQuadrature< GridPartType, 0 >::QuadraturePointWrapperType OriginalElementPoint
Definition: diffusionmodel.hh:141
RFunctionSpaceType::JacobianRangeType RJacobianRangeType
Definition: diffusionmodel.hh:115
Dune::Fem::CachingQuadrature< GridPartType, 1, Traits >::QuadraturePointWrapperType IntersectionPoint
Definition: diffusionmodel.hh:129
Dune::Fem::CachingQuadrature< GridPartType, 0 >::QuadraturePointWrapperType OriginalPoint
Definition: diffusionmodel.hh:137
virtual double time() const =0
virtual ~DiffusionModel()
Definition: diffusionmodel.hh:153
Dune::Fem::FunctionSpace< double, RangeFieldType, GridPart::dimensionworld, dimR > RFunctionSpaceType
Definition: diffusionmodel.hh:108
GridPartType::IntersectionType IntersectionType
Definition: diffusionmodel.hh:120
DFunctionSpaceType::JacobianRangeType DJacobianRangeType
Definition: diffusionmodel.hh:111
virtual bool isDirichletIntersection(const IntersectionType &inter, DirichletComponentType &dirichletComponent) const =0
Dune::Fem::ElementQuadrature< GridPartType, 0, Traits >::QuadraturePointWrapperType ElementPoint
Definition: diffusionmodel.hh:131
GridPartType::template Codim< 0 >::EntityType EntityType
Definition: diffusionmodel.hh:119
Dune::Fem::CachingQuadrature< GridPartType, 0, Traits >::QuadraturePointWrapperType Point
Definition: diffusionmodel.hh:127
RFunctionSpaceType::HessianRangeType RHessianRangeType
Definition: diffusionmodel.hh:116
static const int dimD
Definition: diffusionmodel.hh:100
RangeField RangeFieldType
Definition: diffusionmodel.hh:103
EntityType::Geometry::LocalCoordinate LocalDomainType
Definition: diffusionmodel.hh:121
Dune::FemPy::FempyQuadratureTraits< F, d > Traits
Definition: diffusionmodel.hh:125
Dune::Fem::FunctionSpace< double, RangeFieldType, GridPart::dimensionworld, dimD > DFunctionSpaceType
Definition: diffusionmodel.hh:106
GridPart GridPartType
Definition: diffusionmodel.hh:99
DFunctionSpaceType::HessianRangeType DHessianRangeType
Definition: diffusionmodel.hh:112
DFunctionSpaceType::DomainType DomainType
Definition: diffusionmodel.hh:109
static const int dimR
Definition: diffusionmodel.hh:101
RFunctionSpaceType::RangeType RRangeType
Definition: diffusionmodel.hh:114
virtual void setTime(const double t) const =0
DFunctionSpaceType::RangeType DRangeType
Definition: diffusionmodel.hh:110
Definition: diffusionmodel.hh:226
Base::LocalDomainType LocalDomainType
Definition: diffusionmodel.hh:241
Base::DRangeType DRangeType
Definition: diffusionmodel.hh:243
DiffusionModelWrapper(Args &&... args)
Definition: diffusionmodel.hh:253
virtual bool hasNeumanBoundary() const
Definition: diffusionmodel.hh:301
WrapperDiffusionModelMethods(OriginalPoint)
virtual double time() const
Definition: diffusionmodel.hh:291
static const int dimR
Definition: diffusionmodel.hh:230
Base::RRangeType RRangeType
Definition: diffusionmodel.hh:246
const ModelImpl & impl() const
Definition: diffusionmodel.hh:313
ModelImpl & impl()
Definition: diffusionmodel.hh:317
virtual std::string name() const
Definition: diffusionmodel.hh:276
Base::DJacobianRangeType DJacobianRangeType
Definition: diffusionmodel.hh:244
~DiffusionModelWrapper()
Definition: diffusionmodel.hh:257
std::array< int, dimR > DirichletComponentType
Definition: diffusionmodel.hh:296
WrapperDiffusionModelMethods(OriginalIntersectionPoint)
Base::ElementPoint ElementPoint
Definition: diffusionmodel.hh:235
Base::RJacobianRangeType RJacobianRangeType
Definition: diffusionmodel.hh:247
Base::ElementIntersectionPoint ElementIntersectionPoint
Definition: diffusionmodel.hh:236
WrapperDiffusionModelMethods(IntersectionPoint)
Base::OriginalElementIntersectionPoint OriginalElementIntersectionPoint
Definition: diffusionmodel.hh:240
virtual bool hasDirichletBoundary() const
Definition: diffusionmodel.hh:297
WrapperDiffusionModelMethods(OriginalElementIntersectionPoint)
virtual bool init(const EntityType &entity) const
Definition: diffusionmodel.hh:309
Base::DomainType DomainType
Definition: diffusionmodel.hh:242
Base::EntityType EntityType
Definition: diffusionmodel.hh:249
virtual void setTime(const double t) const
Definition: diffusionmodel.hh:283
Base::OriginalElementPoint OriginalElementPoint
Definition: diffusionmodel.hh:239
Base::IntersectionType IntersectionType
Definition: diffusionmodel.hh:250
Base::Point Point
Definition: diffusionmodel.hh:233
WrapperDiffusionModelMethods(ElementPoint)
Base::IntersectionPoint IntersectionPoint
Definition: diffusionmodel.hh:234
WrapperDiffusionModelMethods(OriginalElementPoint)
WrapperDiffusionModelMethods(ElementIntersectionPoint)
Base::OriginalIntersectionPoint OriginalIntersectionPoint
Definition: diffusionmodel.hh:238
DiffusionModel< GridPartType, dimD, dimR, typename ModelImpl::RRangeFieldType > Base
Definition: diffusionmodel.hh:231
static const int dimD
Definition: diffusionmodel.hh:229
Base::OriginalPoint OriginalPoint
Definition: diffusionmodel.hh:237
WrapperDiffusionModelMethods(LocalDomainType)
virtual bool isDirichletIntersection(const IntersectionType &inter, DirichletComponentType &dirichletComponent) const
Definition: diffusionmodel.hh:305
Base::RHessianRangeType RHessianRangeType
Definition: diffusionmodel.hh:248
ModelImpl::GridPartType GridPartType
Definition: diffusionmodel.hh:228
Base::DHessianRangeType DHessianRangeType
Definition: diffusionmodel.hh:245
WrapperDiffusionModelMethods(Point)
Definition: diffusionmodel.hh:340
Dune::Fem::CachingQuadrature< GridPartType, 0, Traits >::QuadraturePointWrapperType Point
Definition: diffusionmodel.hh:369
GridPartType::IntersectionType IntersectionType
Definition: diffusionmodel.hh:362
DFunctionSpaceType::RangeType DRangeType
Definition: diffusionmodel.hh:352
DFunctionSpaceType::DomainFieldType DDomainFieldType
Definition: diffusionmodel.hh:355
virtual double time() const =0
Dune::Fem::FunctionSpace< double, RangeFieldType, GridPart::dimensionworld, dimD > DFunctionSpaceType
Definition: diffusionmodel.hh:348
GridPart GridPartType
Definition: diffusionmodel.hh:341
Dune::Fem::CachingQuadrature< GridPartType, 1, Traits >::QuadraturePointWrapperType IntersectionPoint
Definition: diffusionmodel.hh:371
GridPartType::template Codim< 0 >::EntityType EntityType
Definition: diffusionmodel.hh:361
static const int dimR
Definition: diffusionmodel.hh:343
Dune::Fem::ElementQuadrature< GridPartType, 0, Traits >::QuadraturePointWrapperType ElementPoint
Definition: diffusionmodel.hh:373
Dune::Fem::CachingQuadrature< GridPartType, 0 >::QuadraturePointWrapperType OriginalPoint
Definition: diffusionmodel.hh:379
virtual std::string name() const =0
RFunctionSpaceType::HessianRangeType RHessianRangeType
Definition: diffusionmodel.hh:358
RangeField RangeFieldType
Definition: diffusionmodel.hh:345
EntityType::Geometry::LocalCoordinate LocalDomainType
Definition: diffusionmodel.hh:363
Dune::Fem::ElementQuadrature< GridPartType, 1 >::QuadraturePointWrapperType OriginalElementIntersectionPoint
Definition: diffusionmodel.hh:385
virtual bool init(const EntityType &entity) const =0
DFunctionSpaceType::DomainType DomainType
Definition: diffusionmodel.hh:351
virtual ~DGDiffusionModel()
Definition: diffusionmodel.hh:395
Dune::FemPy::FempyQuadratureTraits< F, d > Traits
Definition: diffusionmodel.hh:367
Dune::Fem::ElementQuadrature< GridPartType, 0 >::QuadraturePointWrapperType OriginalElementPoint
Definition: diffusionmodel.hh:383
DFunctionSpaceType::JacobianRangeType DJacobianRangeType
Definition: diffusionmodel.hh:353
Dune::Fem::FunctionSpace< double, RangeFieldType, GridPart::dimensionworld, dimR > RFunctionSpaceType
Definition: diffusionmodel.hh:350
Dune::Fem::ElementQuadrature< GridPartType, 1, Traits >::QuadraturePointWrapperType ElementIntersectionPoint
Definition: diffusionmodel.hh:375
virtual bool hasNeumanBoundary() const =0
DGDiffusionModel< GridPartType, dimD, dimR, RangeField > ModelType
Definition: diffusionmodel.hh:344
virtual VirtualDiffusionModelMethods(Point) VirtualDiffusionModelMethods(ElementPoint) VirtualDiffusionModelMethods(IntersectionPoint) VirtualDiffusionModelMethods(ElementIntersectionPoint) VirtualDiffusionModelMethods(OriginalPoint) VirtualDiffusionModelMethods(OriginalElementPoint) VirtualDiffusionModelMethods(OriginalIntersectionPoint) VirtualDiffusionModelMethods(OriginalElementIntersectionPoint) VirtualDiffusionModelMethods(LocalDomainType) typedef std bool hasDirichletBoundary() const =0
Dune::Fem::CachingQuadrature< GridPartType, 1 >::QuadraturePointWrapperType OriginalIntersectionPoint
Definition: diffusionmodel.hh:381
RFunctionSpaceType::DomainFieldType rDomainFieldType
Definition: diffusionmodel.hh:359
static const int dimD
Definition: diffusionmodel.hh:342
DFunctionSpaceType::HessianRangeType DHessianRangeType
Definition: diffusionmodel.hh:354
virtual void setTime(const double t) const =0
RFunctionSpaceType::JacobianRangeType RJacobianRangeType
Definition: diffusionmodel.hh:357
virtual bool isDirichletIntersection(const IntersectionType &inter, DirichletComponentType &dirichletComponent) const =0
DGDiffusionModel()
Definition: diffusionmodel.hh:393
RFunctionSpaceType::RangeType RRangeType
Definition: diffusionmodel.hh:356
Definition: diffusionmodel.hh:446
static const int dimD
Definition: diffusionmodel.hh:449
static const int dimR
Definition: diffusionmodel.hh:450
virtual bool hasNeumanBoundary() const
Definition: diffusionmodel.hh:533
WrapperDiffusionModelMethods(OriginalElementIntersectionPoint)
WrapperDiffusionModelMethods(OriginalIntersectionPoint)
Base::LocalDomainType LocalDomainType
Definition: diffusionmodel.hh:461
Base::DRangeType DRangeType
Definition: diffusionmodel.hh:463
Base::EntityType EntityType
Definition: diffusionmodel.hh:469
Base::Point Point
Definition: diffusionmodel.hh:453
Base::ElementIntersectionPoint ElementIntersectionPoint
Definition: diffusionmodel.hh:456
WrapperPenaltyMethods(LocalDomainType)
Base::RRangeType RRangeType
Definition: diffusionmodel.hh:466
WrapperPenaltyMethods(OriginalIntersectionPoint)
WrapperDiffusionModelMethods(LocalDomainType)
Base::OriginalPoint OriginalPoint
Definition: diffusionmodel.hh:457
WrapperDiffusionModelMethods(IntersectionPoint)
DGDiffusionModelWrapper(Args &&... args)
Definition: diffusionmodel.hh:473
Base::OriginalElementIntersectionPoint OriginalElementIntersectionPoint
Definition: diffusionmodel.hh:460
Base::OriginalElementPoint OriginalElementPoint
Definition: diffusionmodel.hh:459
Base::DomainType DomainType
Definition: diffusionmodel.hh:462
Base::IntersectionPoint IntersectionPoint
Definition: diffusionmodel.hh:454
~DGDiffusionModelWrapper()
Definition: diffusionmodel.hh:477
Base::IntersectionType IntersectionType
Definition: diffusionmodel.hh:470
std::array< int, dimR > DirichletComponentType
Definition: diffusionmodel.hh:528
WrapperDiffusionModelMethods(OriginalElementPoint)
Base::OriginalIntersectionPoint OriginalIntersectionPoint
Definition: diffusionmodel.hh:458
virtual bool hasDirichletBoundary() const
Definition: diffusionmodel.hh:529
virtual void setTime(const double t) const
Definition: diffusionmodel.hh:515
Base::RJacobianRangeType RJacobianRangeType
Definition: diffusionmodel.hh:467
DGDiffusionModel< GridPartType, dimD, dimR, typename ModelImpl::RRangeFieldType > Base
Definition: diffusionmodel.hh:451
WrapperPenaltyMethods(OriginalElementPoint)
ModelImpl & impl()
Definition: diffusionmodel.hh:549
ModelImpl::GridPartType GridPartType
Definition: diffusionmodel.hh:448
Base::DHessianRangeType DHessianRangeType
Definition: diffusionmodel.hh:465
virtual double time() const
Definition: diffusionmodel.hh:523
WrapperDiffusionModelMethods(ElementIntersectionPoint)
const ModelImpl & impl() const
Definition: diffusionmodel.hh:545
WrapperDiffusionModelMethods(OriginalPoint)
Base::ElementPoint ElementPoint
Definition: diffusionmodel.hh:455
virtual bool init(const EntityType &entity) const
Definition: diffusionmodel.hh:541
WrapperPenaltyMethods(OriginalElementIntersectionPoint)
Base::RHessianRangeType RHessianRangeType
Definition: diffusionmodel.hh:468
WrapperDiffusionModelMethods(ElementPoint)
virtual std::string name() const
Definition: diffusionmodel.hh:508
Base::DJacobianRangeType DJacobianRangeType
Definition: diffusionmodel.hh:464
WrapperPenaltyMethods(Point) WrapperPenaltyMethods(ElementPoint) WrapperPenaltyMethods(IntersectionPoint) WrapperPenaltyMethods(ElementIntersectionPoint) WrapperPenaltyMethods(OriginalPoint)
virtual bool isDirichletIntersection(const IntersectionType &inter, DirichletComponentType &dirichletComponent) const
Definition: diffusionmodel.hh:537
A vector valued function space.
Definition: functionspace.hh:60
FunctionSpaceTraits::DomainFieldType DomainFieldType
Intrinsic type used for values in the domain field (usually a double)
Definition: functionspaceinterface.hh:60
FunctionSpaceTraits::RangeType RangeType
Type of range vector (using type of range field) has a Dune::FieldVector type interface.
Definition: functionspaceinterface.hh:71
FunctionSpaceTraits::LinearMappingType JacobianRangeType
Intrinsic type used for the jacobian values has a Dune::FieldMatrix type interface.
Definition: functionspaceinterface.hh:75
FunctionSpaceTraits::DomainType DomainType
Type of domain vector (using type of domain field) has a Dune::FieldVector type interface.
Definition: functionspaceinterface.hh:67