dune-pdelab  2.0.0
stokesparameter.hh
Go to the documentation of this file.
1 #ifndef DUNE_PDELAB_LOCALOPERATOR_STOKESPARAMETER_HH
2 #define DUNE_PDELAB_LOCALOPERATOR_STOKESPARAMETER_HH
3 
4 #include <dune/common/parametertree.hh>
7 
8 namespace Dune {
9  namespace PDELab {
10 
33  enum Type {
34  DoNothing = 0,
38  };
39  };
40 
44  template<typename GV, typename RF>
46  {
48  typedef GV GridView;
49 
51  enum {
53  dimDomain = GV::dimension
54  };
55 
57  typedef typename GV::Grid::ctype DomainField;
58 
60  typedef Dune::FieldVector<DomainField,dimDomain> Domain;
61 
63  typedef Dune::FieldVector<DomainField,dimDomain-1> IntersectionDomain;
64 
66  typedef RF RangeField;
67 
69  typedef Dune::FieldVector<RF,GV::dimensionworld> VelocityRange;
70 
72  typedef Dune::FieldVector<RF,1> PressureRange;
73 
76 
78  typedef typename GV::Traits::template Codim<0>::Entity Element;
80  typedef typename GV::Intersection Intersection;
81  };
82 
83  namespace {
84 
89  template<typename GF, typename Entity, typename Domain>
90  typename GF::Traits::RangeType
91  evaluateVelocityGridFunction(const GF& gf,
92  const Entity& e,
93  const Domain& x)
94  {
95  dune_static_assert(int(GF::Traits::dimRange) == int(Domain::dimension),"dimension of function range does not match grid dimension");
96  typename GF::Traits::RangeType y;
97  gf.evaluate(e,x,y);
98  return y;
99  }
100 
105  template<typename GF, typename Entity, typename Domain>
106  FieldVector<typename GF::template Child<0>::Type::Traits::RangeFieldType,GF::CHILDREN>
107  evaluateVelocityGridFunction(const GF& gf,
108  const Entity& e,
109  const Domain& x)
110  {
111  dune_static_assert(Domain::dimension == GF::CHILDREN,"dimension of function range does not match grid dimension");
112  FieldVector<typename GF::template Child<0>::Type::Traits::RangeFieldType,GF::CHILDREN> y;
113  typename GF::template Child<0>::Type::Traits::RangeType cy;
114  for (int d = 0; d < Domain::dimension; ++d)
115  {
116  gf.child(d).evaluate(e,x,cy);
117  y[d] = cy;
118  }
119  return y;
120  }
121 
122  }
123 
142  template <typename GV, typename RF, typename F, typename B, typename V, typename J, bool navier = false, bool tensor = false>
144  {
145  public:
146 
147  static const bool assemble_navier = navier;
148  static const bool assemble_full_tensor = tensor;
149 
152 
154  NavierStokesDefaultParameters(const Dune::ParameterTree& config,
155  F& f,
156  B& b,
157  V& v,
158  J& j)
159  : _rho(config.get<double>("rho"))
160  , _mu(config.get<double>("mu"))
161  , _f(f)
162  , _b(b)
163  , _v(v)
164  , _j(j)
165  {}
166 
168  const RF& rho,
169  F& f,
170  B& b,
171  V& v,
172  J& j)
173  : _rho(rho)
174  , _mu(mu)
175  , _f(f)
176  , _b(b)
177  , _v(v)
178  , _j(j)
179  {}
180 
181 
183  template<typename EG>
184  typename Traits::VelocityRange
185  f(const EG& e, const typename Traits::Domain& x) const
186  {
187  typename F::Traits::RangeType fvalue;
188  return evaluateVelocityGridFunction(_f,e.entity(),x);
189  }
190 
192  template<typename IG>
194  bctype(const IG& is,
195  const typename Traits::IntersectionDomain& x) const
196  {
197  typename B::Traits::RangeType y;
198  _b.evaluate(is,x,y);
199  return y;
200  }
201 
203  template<typename EG>
204  typename Traits::RangeField
205  mu(const EG& e, const typename Traits::Domain& x) const
206  {
207  return _mu;
208  }
209 
211  template<typename IG>
212  typename Traits::RangeField
213  mu(const IG& ig, const typename Traits::IntersectionDomain& x) const
214  {
215  return _mu;
216  }
217 
219  template<typename EG>
220  typename Traits::RangeField
221  rho(const EG& eg, const typename Traits::Domain& x) const
222  {
223  return _rho;
224  }
225 
227  template<typename IG>
228  typename Traits::RangeField
229  rho(const IG& ig, const typename Traits::IntersectionDomain& x) const
230  {
231  return _rho;
232  }
233 
235  template<typename EG>
236  typename Traits::VelocityRange
237  g(const EG& e, const typename Traits::Domain& x) const
238  {
239  typename V::Traits::RangeType y;
240  _v.evaluate(e.entity(),x,y);
241  return y;
242  }
243 
245  template<typename IG>
246  typename Traits::VelocityRange
247  g(const IG& ig, const typename Traits::IntersectionDomain& x) const
248  {
249  typename IG::EntityPointer ep = ig.inside();
250  typename V::Traits::RangeType y;
251  _v.evaluate(*ep,ig.geometryInInside().global(x),y);
252  return y;
253  }
254 
256  template<typename EG>
257  typename Traits::RangeField
258  g2(const EG& e, const typename Traits::Domain& x) const
259  {
260  return 0;
261  }
262 
263 #ifdef DOXYGEN
264 
266  template<typename IG>
267  typename Traits::VelocityRange>
268  j(const IG& ig,
269  const typename Traits::IntersectionDomain& x,
270  const typename Traits::Domain& normal) const;
271 
272 #else // DOXYGEN
273 
275  template<typename IG>
276  typename enable_if<
277  J::Traits::dimRange == 1 &&
278  (GV::dimension > 1) &&
279  AlwaysTrue<IG>::value, // required to force lazy evaluation
280  typename Traits::VelocityRange
281  >::type
282  j(const IG& ig,
283  const typename Traits::IntersectionDomain& x,
284  typename Traits::Domain normal) const
285  {
286  typename J::Traits::RangeType r;
287  typename IG::EntityPointer ep = ig.inside();
288  _j.evaluate(*ep,ig.geometryInInside().global(x),r);
289  normal *= r;
290  return normal;
291  }
292 
294  template<typename IG>
295  typename enable_if<
296  J::Traits::dimRange == GV::dimension &&
297  AlwaysTrue<IG>::value, // required to force lazy evaluation
298  typename Traits::VelocityRange
299  >::type
300  j(const IG& ig,
301  const typename Traits::IntersectionDomain& x,
302  const typename Traits::Domain& normal) const
303  {
304  typename IG::EntityPointer ep = ig.inside();
305  typename J::Traits::RangeType y;
306  _j.evaluate(*ep,ig.geometryInInside().global(x),y);
307  return y;
308  }
309 
310 #endif // DOXYGEN
311 
312  void setTime(RF time)
313  {
314  _f.setTime(time);
315  _b.setTime(time);
316  _v.setTime(time);
317  _j.setTime(time);
318  }
319 
320  private:
321  const RF _rho;
322  const RF _mu;
323  const F& _f;
324  const B& _b;
325  const V& _v;
326  const J& _j;
327  };
328 
329 
333  template<typename PRM>
336  {
337  private:
338  const PRM & prm_;
339 
340  public:
341 
344  : prm_(_prm) { }
345 
347  template<typename I>
348  bool isDirichlet(const I & intersection,
349  const Dune::FieldVector<typename I::ctype, I::dimension-1> & coord) const
350  {
351  StokesBoundaryCondition::Type bctype = prm_.bctype(intersection,coord);
353  }
354 
355  };
356 
360  template<typename PRM>
363  {
364  private:
365  const PRM & prm_;
366 
367  public:
368 
371  : prm_(_prm) { }
372 
374  template<typename I>
375  bool isDirichlet(const I & intersection,
376  const Dune::FieldVector<typename I::ctype, I::dimension-1> & coord) const
377  { return false; }
378  };
379 
380 
381 
382 #ifndef DOXYGEN
383 
387  template<typename PRM, int rangeDim>
388  class NavierStokesFunctionAdapterBase
390  Dune::PDELab::GridFunctionTraits<
391  typename PRM::Traits::GridView,
392  typename PRM::Traits::RangeField,
393  rangeDim,
394  Dune::FieldVector<typename PRM::Traits::RangeField,rangeDim>
395  >,
396  NavierStokesFunctionAdapterBase<PRM,rangeDim>
397  >
398  {
399  public:
402  typename PRM::Traits::GridView,
403  typename PRM::Traits::RangeField,
404  rangeDim,
405  Dune::FieldVector<typename PRM::Traits::RangeField,rangeDim>
406  > Traits;
407 
409  NavierStokesFunctionAdapterBase(PRM& prm)
410  : _prm(prm)
411  {}
412 
413  void setTime(const double time)
414  {
415  _prm.setTime(time);
416  }
417 
418  const PRM& parameters() const
419  {
420  return _prm;
421  }
422 
424  const typename Traits::GridViewType& getGridView () const
425  {
426  return _prm.gridView();
427  }
428 
429  private:
430  PRM& _prm;
431  };
432 
433 
434 #endif // DOXYGEN
435 
442  template<typename PRM>
444  : public NavierStokesFunctionAdapterBase<PRM,PRM::Traits::dimDomain>
445  {
446 
448  typedef NavierStokesFunctionAdapterBase<PRM,PRM::Traits::dimDomain> Base;
449 
450  using Base::parameters;
451 
452  public:
453 
454  typedef typename Base::Traits Traits;
455 
458  : Base(prm)
459  {}
460 
462  void evaluate (const typename Traits::ElementType& e,
463  const typename Traits::DomainType& x,
464  typename Traits::RangeType& y) const
465  {
466  y = parameters().g(e,x);
467  }
468 };
469 
470 
471 
472 #if 0
473 
475 template < typename PRM >
476 class NavierStokesDirichletFunctionAdapterFactory
477 {
478 public:
480  NavierStokesDirichletFunctionAdapter<PRM>,
481  NavierStokesPressureDirichletFunctionAdapter<PRM> >
482  BoundaryDirichletFunction;
483 
484 NavierStokesDirichletFunctionAdapterFactory(PRM & prm)
485 : v(prm), p(prm), df(v,p)
486 {}
487 
488 BoundaryDirichletFunction & dirichletFunction()
489 {
490  return df;
491 }
492 
493 private:
494 NavierStokesVelocityDirichletFunctionAdapter<PRM> v;
495 NavierStokesPressureDirichletFunctionAdapter<PRM> p;
496 BoundaryDirichletFunction df;
497 };
498 #endif
499 
500 
501 } // end namespace PDELab
502 } // end namespace Dune
503 
504 #endif
static const bool assemble_navier
Definition: stokesparameter.hh:147
Traits::VelocityRange f(const EG &e, const typename Traits::Domain &x) const
source term
Definition: stokesparameter.hh:185
Dune::FieldVector< DomainField, dimDomain-1 > IntersectionDomain
domain type
Definition: stokesparameter.hh:63
StokesBoundaryCondition BoundaryCondition
boundary type value
Definition: stokesparameter.hh:75
Definition: stokesparameter.hh:443
Traits::VelocityRange g(const IG &ig, const typename Traits::IntersectionDomain &x) const
Dirichlet boundary condition value from local intersection coordinate.
Definition: stokesparameter.hh:247
bool isDirichlet(const I &intersection, const Dune::FieldVector< typename I::ctype, I::dimension-1 > &coord) const
Definition: stokesparameter.hh:348
NavierStokesVelocityFunctionAdapter(PRM &prm)
Constructor.
Definition: stokesparameter.hh:457
Traits::VelocityRange j(const IG &ig, const typename Traits::IntersectionDomain &x, const typename Traits::Domain &normal) const
Neumann boundary condition (stress)
Type
Definition: stokesparameter.hh:33
RF RangeField
Export type for range field.
Definition: stokesparameter.hh:66
NavierStokesParameterTraits< GV, RF > Traits
Type traits.
Definition: stokesparameter.hh:151
Traits::RangeField mu(const IG &ig, const typename Traits::IntersectionDomain &x) const
Dynamic viscosity value from local intersection coordinate.
Definition: stokesparameter.hh:213
GV GridView
the grid view
Definition: stokesparameter.hh:48
NavierStokesDefaultParameters(const Dune::ParameterTree &config, F &f, B &b, V &v, J &j)
Constructor.
Definition: stokesparameter.hh:154
Definition: stokesparameter.hh:45
GV::Intersection Intersection
grid intersection type
Definition: stokesparameter.hh:80
composite functions
Definition: function.hh:800
static const bool assemble_full_tensor
Definition: stokesparameter.hh:148
Definition: common/constraintsparameters.hh:24
StokesVelocityDirichletConstraints(const PRM &_prm)
Constructor.
Definition: stokesparameter.hh:343
NavierStokesDefaultParameters(const RF &mu, const RF &rho, F &f, B &b, V &v, J &j)
Definition: stokesparameter.hh:167
Traits::RangeField rho(const IG &ig, const typename Traits::IntersectionDomain &x) const
Density value from local intersection coordinate.
Definition: stokesparameter.hh:229
void setTime(RF time)
Definition: stokesparameter.hh:312
Base::Traits Traits
Definition: stokesparameter.hh:454
Traits::RangeField rho(const EG &eg, const typename Traits::Domain &x) const
Density value from local cell coordinate.
Definition: stokesparameter.hh:221
Definition: stokesparameter.hh:32
Traits::RangeField g2(const EG &e, const typename Traits::Domain &x) const
pressure source term
Definition: stokesparameter.hh:258
GV::Grid::ctype DomainField
Export type for domain field.
Definition: stokesparameter.hh:57
Traits::VelocityRange g(const EG &e, const typename Traits::Domain &x) const
Dirichlet boundary condition value from local cell coordinate.
Definition: stokesparameter.hh:237
dimension of the domain
Definition: stokesparameter.hh:53
Traits::BoundaryCondition::Type bctype(const IG &is, const typename Traits::IntersectionDomain &x) const
boundary condition type from local intersection coordinate
Definition: stokesparameter.hh:194
const IG & ig
Definition: common/constraints.hh:146
leaf of a function tree
Definition: function.hh:577
Definition: stokesparameter.hh:334
Definition: stokesparameter.hh:143
StokesPressureDirichletConstraints(const PRM &_prm)
Constructor.
Definition: stokesparameter.hh:370
Dune::FieldVector< RF, 1 > PressureRange
pressure range type
Definition: stokesparameter.hh:72
Traits::RangeField mu(const EG &e, const typename Traits::Domain &x) const
Dynamic viscosity value from local cell coordinate.
Definition: stokesparameter.hh:205
static const unsigned int value
Definition: gridfunctionspace/tags.hh:175
traits class holding the function signature, same as in local function
Definition: function.hh:176
GV::Traits::template Codim< 0 >::Entity Element
grid element type
Definition: stokesparameter.hh:78
const EG & eg
Definition: common/constraints.hh:277
Dune::FieldVector< DomainField, dimDomain > Domain
domain type
Definition: stokesparameter.hh:60
Dune::FieldVector< RF, GV::dimensionworld > VelocityRange
deformation range type
Definition: stokesparameter.hh:69
Definition: stokesparameter.hh:34
R p(int i, D x)
Definition: qkdg.hh:62
Definition: stokesparameter.hh:361
const E & e
Definition: interpolate.hh:172
bool isDirichlet(const I &intersection, const Dune::FieldVector< typename I::ctype, I::dimension-1 > &coord) const
Definition: stokesparameter.hh:375
void evaluate(const typename Traits::ElementType &e, const typename Traits::DomainType &x, typename Traits::RangeType &y) const
Evaluate dirichlet function.
Definition: stokesparameter.hh:462