dune-fem  2.4.1-rc
modelcaller.hh
Go to the documentation of this file.
1 #ifndef DUNE_FEM_PASS_LOCALDG_MODELCALLER_HH
2 #define DUNE_FEM_PASS_LOCALDG_MODELCALLER_HH
3 
4 #include <cstddef>
5 #include <vector>
6 
13 
14 namespace Dune
15 {
16 
17  namespace Fem
18  {
19 
20  // DGDiscreteModelCaller
21  // ---------------------
22 
29  template< class DiscreteModel, class Argument, class PassIds >
31  {
33 
34  public:
35  // discrete model type
36  typedef DiscreteModel DiscreteModelType;
37  // argument type
38  typedef Argument ArgumentType;
39 
41  typedef typename DiscreteModelType::Selector Selector;
42 
43  // entity type
44  typedef typename DiscreteModelType::EntityType EntityType;
45  // intersection type
46  typedef typename DiscreteModelType::IntersectionType IntersectionType;
47 
48  // type of volume quadrature
49  typedef typename DiscreteModelType::Traits::VolumeQuadratureType VolumeQuadratureType;
50  // type of face quadrature
51  typedef typename DiscreteModelType::Traits::FaceQuadratureType FaceQuadratureType;
52 
53  // type of mass factor
54  typedef typename DiscreteModelType::MassFactorType MassFactorType;
55 
56  // function space type
57  typedef typename DiscreteModelType::FunctionSpaceType FunctionSpaceType;
58  // range type
59  typedef typename FunctionSpaceType::RangeType RangeType;
60  // jacobian range type
61  typedef typename FunctionSpaceType::JacobianRangeType JacobianRangeType;
62 
63  private:
65  typedef typename FilterType::type DiscreteFunctionPointerTupleType;
66 
67  typedef typename TupleTypeTraits< DiscreteFunctionPointerTupleType >::PointeeTupleType DiscreteFunctionTupleType;
68 
69  protected:
71 
74 
75  public:
76  DGDiscreteModelCaller ( ArgumentType &argument, DiscreteModelType &discreteModel )
77  : discreteModel_( discreteModel ),
78  time_( 0. ),
79  discreteFunctions_( FilterType::apply( argument ) ),
80  localFunctionsInside_( Dune::DereferenceTuple< DiscreteFunctionPointerTupleType >::apply( discreteFunctions_ ) ),
81  localFunctionsOutside_( Dune::DereferenceTuple< DiscreteFunctionPointerTupleType >::apply( discreteFunctions_ ) )
82  {}
83 
84  // return true, if discrete model has flux
85  bool hasFlux () const { return discreteModel().hasFlux(); }
86  // return true, if discrete model has mass
87  bool hasMass () const { return discreteModel().hasMass(); }
88  // return true, if discrete model has source
89  bool hasSource () const { return discreteModel().hasSource(); }
90 
91  // set time
92  void setTime ( double time ) { time_ = time; }
93  // return time
94  double time () const { return time_; }
95 
96  // set inside local funtions to entity
97  void setEntity ( const EntityType &entity )
98  {
99  localFunctionsInside_.init( entity );
100  discreteModel().setEntity( entity );
101  }
102 
103  // set outside local functions to entity
104  void setNeighbor ( const EntityType &entity )
105  {
106  localFunctionsOutside_.init( entity );
107  discreteModel().setNeighbor( entity );
108  }
109 
110  // evaluate inside local functions in all quadrature points
111  void setEntity ( const EntityType &entity, const VolumeQuadratureType &quadrature )
112  {
113  setEntity( entity );
114  values_.resize( quadrature.nop() );
115  jacobians_.resize( quadrature.nop() );
118  }
119 
120  // evaluate outside local functions in all quadrature points
121  template< class QuadratureType >
122  void setNeighbor ( const EntityType &neighbor,
123  const QuadratureType &inside,
124  const QuadratureType &outside )
125  {
126  // we assume that setEntity() was called in advance!
127  valuesInside_.resize( inside.nop() );
129 
130  setNeighbor( neighbor );
131 
132  valuesOutside_.resize( outside.nop() );
134  }
135 
136  // please doc me
137  template< class QuadratureType >
138  void setBoundary ( const EntityType &entity, const QuadratureType &quadrature )
139  {
140  valuesInside_.resize( quadrature.nop() );
142  }
143 
144  // evaluate analytical flux
145  void analyticalFlux ( const EntityType &entity,
146  const VolumeQuadratureType &quadrature,
147  const int qp,
148  JacobianRangeType &flux )
149  {
150  assert( hasFlux() );
151  discreteModel().analyticalFlux( entity, time(), quadrature.point( qp ), values_[ qp ], flux );
152  }
153 
154  // evaluate analytical flux
155  double source ( const EntityType &entity,
156  const VolumeQuadratureType &quadrature,
157  const int qp,
158  RangeType &source )
159  {
160  assert( hasSource() );
161  return discreteModel().source( entity, time(), quadrature.point( qp ), values_[ qp ], jacobians_[ qp ], source );
162  }
163 
164 
165  // evaluate both analytical flux and source
166  double analyticalFluxAndSource ( const EntityType &entity,
167  const VolumeQuadratureType &quadrature,
168  const int qp,
169  JacobianRangeType &flux,
170  RangeType &source )
171  {
172  // we may only assume that hasSource() == true, cf. pass.hh
173  if( hasFlux() )
174  analyticalFlux( entity, quadrature, qp, flux );
175  return ThisType::source( entity, quadrature, qp, source );
176  }
177 
178  // evaluate numerical flux
179  template< class QuadratureType >
180  double numericalFlux ( const IntersectionType &intersection,
181  const QuadratureType &inside,
182  const QuadratureType &outside,
183  const int qp,
184  RangeType &gLeft, RangeType &gRight )
185  {
186  return discreteModel().numericalFlux( intersection, time(), inside.localPoint( qp ), valuesInside_[ qp ], valuesOutside_[ qp ], gLeft, gRight );
187  }
188 
189  // evaluate boundary flux
190  double boundaryFlux ( const IntersectionType &intersection,
191  const FaceQuadratureType &quadrature,
192  const int qp,
193  RangeType &gLeft )
194  {
195  return discreteModel().boundaryFlux( intersection, time(), quadrature.localPoint( qp ), valuesInside_[ qp ], gLeft );
196  }
197 
198  // evaluate mass
199  void mass ( const EntityType &entity,
200  const VolumeQuadratureType &quadrature,
201  const int qp,
202  MassFactorType &massFactor )
203  {
204  discreteModel().mass( entity, time(), quadrature.point( qp ), values_[ qp ], massFactor );
205  }
206 
207  protected:
208  DiscreteModelType &discreteModel () { return discreteModel_; }
209  const DiscreteModelType &discreteModel () const { return discreteModel_; }
210 
211  private:
212  // forbid copying and assignment
213  DGDiscreteModelCaller ( const ThisType & );
214  ThisType operator= ( const ThisType & );
215 
216  DiscreteModelType &discreteModel_;
217  double time_;
218 
219  DiscreteFunctionPointerTupleType discreteFunctions_;
220 
221  protected:
223  std::vector< Dune::TypeIndexedTuple< RangeTupleType, Selector > > values_, valuesInside_, valuesOutside_;
224  std::vector< Dune::TypeIndexedTuple< JacobianRangeTupleType, Selector > > jacobians_;
225  };
226 
227  } // namespace Fem
228 
229 } // namespace Dune
230 
231 #endif // #ifndef DUNE_FEM_PASS_LOCALDG_MODELCALLER_HH
std::vector< Dune::TypeIndexedTuple< JacobianRangeTupleType, Selector > > jacobians_
Definition: modelcaller.hh:224
bool hasSource() const
Definition: modelcaller.hh:89
DiscreteModelType::MassFactorType MassFactorType
Definition: modelcaller.hh:54
DiscreteModelType::Selector Selector
selector (tuple of integral constants)
Definition: modelcaller.hh:41
bool hasMass() const
Definition: modelcaller.hh:87
std::vector< Dune::TypeIndexedTuple< RangeTupleType, Selector > > valuesOutside_
Definition: modelcaller.hh:223
Dune::conditional< isPointerTuple, typename Dune::ForEachType< PointeeTypeEvaluator, Tuple >::Type, EmptyTuple< std::tuple_size< Tuple >::value > >::type PointeeTupleType
Definition: tupletypetraits.hh:120
DiscreteModelType & discreteModel()
Definition: modelcaller.hh:208
std::vector< Dune::TypeIndexedTuple< RangeTupleType, Selector > > valuesInside_
Definition: modelcaller.hh:223
void analyticalFlux(const EntityType &entity, const VolumeQuadratureType &quadrature, const int qp, JacobianRangeType &flux)
Definition: modelcaller.hh:145
DiscreteModelType::IntersectionType IntersectionType
Definition: modelcaller.hh:46
void setTime(double time)
Definition: modelcaller.hh:92
std::vector< Dune::TypeIndexedTuple< RangeTupleType, Selector > > values_
Definition: modelcaller.hh:223
void setEntity(const EntityType &entity, const VolumeQuadratureType &quadrature)
Definition: modelcaller.hh:111
DiscreteModel DiscreteModelType
Definition: modelcaller.hh:36
Dune::ForEachType< RangeTypeEvaluator, LocalFunctionTupleType >::Type RangeTupleType
Definition: localfunctiontuple.hh:120
LocalFunctionTupleType::RangeTupleType RangeTupleType
Definition: modelcaller.hh:72
void init(const EntityType &entity)
set local functions to given entity
Definition: localfunctiontuple.hh:141
double boundaryFlux(const IntersectionType &intersection, const FaceQuadratureType &quadrature, const int qp, RangeType &gLeft)
Definition: modelcaller.hh:190
double numericalFlux(const IntersectionType &intersection, const QuadratureType &inside, const QuadratureType &outside, const int qp, RangeType &gLeft, RangeType &gRight)
Definition: modelcaller.hh:180
void setBoundary(const EntityType &entity, const QuadratureType &quadrature)
Definition: modelcaller.hh:138
DiscreteModelType::Traits::VolumeQuadratureType VolumeQuadratureType
Definition: modelcaller.hh:49
Argument ArgumentType
Definition: modelcaller.hh:38
Dune::ForEachType< JacobianRangeTypeEvaluator, LocalFunctionTupleType >::Type JacobianRangeTupleType
Definition: localfunctiontuple.hh:122
Definition: coordinate.hh:4
FunctionSpaceType::JacobianRangeType JacobianRangeType
Definition: modelcaller.hh:61
void evaluateQuadrature(const QuadratureType &quadrature, TupleVectorType &vector) const
evaluate local functions for quadrature
Definition: localfunctiontuple.hh:191
DiscreteModelType::FunctionSpaceType FunctionSpaceType
Definition: modelcaller.hh:57
const DiscreteModelType & discreteModel() const
Definition: modelcaller.hh:209
double time() const
Definition: modelcaller.hh:94
Definition: tupleutility.hh:239
model caller for local DG pass
Definition: modelcaller.hh:30
double source(const EntityType &entity, const VolumeQuadratureType &quadrature, const int qp, RangeType &source)
Definition: modelcaller.hh:155
Dereference pointer tuple.
Definition: tupletypetraits.hh:286
LocalFunctionTuple< DiscreteFunctionTupleType, EntityType > LocalFunctionTupleType
Definition: modelcaller.hh:70
void setNeighbor(const EntityType &neighbor, const QuadratureType &inside, const QuadratureType &outside)
Definition: modelcaller.hh:122
LocalFunctionTupleType::JacobianRangeTupleType JacobianRangeTupleType
Definition: modelcaller.hh:73
FunctionSpaceType::RangeType RangeType
Definition: modelcaller.hh:59
LocalFunctionTupleType localFunctionsOutside_
Definition: modelcaller.hh:222
DiscreteModelType::EntityType EntityType
Definition: modelcaller.hh:44
DiscreteModelType::Traits::FaceQuadratureType FaceQuadratureType
Definition: modelcaller.hh:51
LocalFunctionTupleType localFunctionsInside_
Definition: modelcaller.hh:222
void mass(const EntityType &entity, const VolumeQuadratureType &quadrature, const int qp, MassFactorType &massFactor)
Definition: modelcaller.hh:199
NextType::type type
Definition: tupleutility.hh:261
double analyticalFluxAndSource(const EntityType &entity, const VolumeQuadratureType &quadrature, const int qp, JacobianRangeType &flux, RangeType &source)
Definition: modelcaller.hh:166
void setEntity(const EntityType &entity)
Definition: modelcaller.hh:97
DGDiscreteModelCaller(ArgumentType &argument, DiscreteModelType &discreteModel)
Definition: modelcaller.hh:76
void setNeighbor(const EntityType &entity)
Definition: modelcaller.hh:104
bool hasFlux() const
Definition: modelcaller.hh:85