dune-pdelab  2.4.1
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  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  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>
193  typename Traits::BoundaryCondition::Type
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,x,y);
241  return y;
242  }
243 
245  template<typename IG>
246  typename Traits::VelocityRange
247  DUNE_DEPRECATED_MSG("Deprecated in DUNE-PDELab 2.4, use entity-based version instead!")
248  g(const IG& ig, const typename Traits::IntersectionDomain& x) const
249  {
250  auto e = ig.inside();
251  typename V::Traits::RangeType y;
252  _v.evaluate(e,ig.geometryInInside().global(x),y);
253  return y;
254  }
255 
257  template<typename EG>
258  typename Traits::RangeField
259  g2(const EG& e, const typename Traits::Domain& x) const
260  {
261  return 0;
262  }
263 
264 #ifdef DOXYGEN
265 
267  template<typename IG>
268  typename Traits::VelocityRange>
269  j(const IG& ig,
270  const typename Traits::IntersectionDomain& x,
271  const typename Traits::Domain& normal) const;
272 
273 #else // DOXYGEN
274 
276  template<typename IG>
277  typename enable_if<
278  J::Traits::dimRange == 1 &&
279  (GV::dimension > 1) &&
280  AlwaysTrue<IG>::value, // required to force lazy evaluation
281  typename Traits::VelocityRange
282  >::type
283  j(const IG& ig,
284  const typename Traits::IntersectionDomain& x,
285  typename Traits::Domain normal) const
286  {
287  typename J::Traits::RangeType r;
288  auto e = ig.inside();
289  _j.evaluate(e,ig.geometryInInside().global(x),r);
290  normal *= r;
291  return normal;
292  }
293 
295  template<typename IG>
296  typename enable_if<
297  J::Traits::dimRange == GV::dimension &&
298  AlwaysTrue<IG>::value, // required to force lazy evaluation
299  typename Traits::VelocityRange
300  >::type
301  j(const IG& ig,
302  const typename Traits::IntersectionDomain& x,
303  const typename Traits::Domain& normal) const
304  {
305  auto e = ig.inside();
306  typename J::Traits::RangeType y;
307  _j.evaluate(e,ig.geometryInInside().global(x),y);
308  return y;
309  }
310 
311 #endif // DOXYGEN
312 
313  void setTime(RF time)
314  {
315  _f.setTime(time);
316  _b.setTime(time);
317  _v.setTime(time);
318  _j.setTime(time);
319  }
320 
321  private:
322  const RF _rho;
323  const RF _mu;
324  const F& _f;
325  const B& _b;
326  const V& _v;
327  const J& _j;
328  };
329 
330 
334  template<typename PRM>
337  {
338  private:
339  const PRM & prm_;
340 
341  public:
342 
345  : prm_(_prm) { }
346 
348  template<typename I>
349  bool isDirichlet(const I & intersection,
350  const Dune::FieldVector<typename I::ctype, I::dimension-1> & coord) const
351  {
352  StokesBoundaryCondition::Type bctype = prm_.bctype(intersection,coord);
354  }
355 
356  };
357 
361  template<typename PRM>
364  {
365  private:
366  const PRM & prm_;
367 
368  public:
369 
372  : prm_(_prm) { }
373 
375  template<typename I>
376  bool isDirichlet(const I & intersection,
377  const Dune::FieldVector<typename I::ctype, I::dimension-1> & coord) const
378  { return false; }
379  };
380 
381 
382 
383 #ifndef DOXYGEN
384 
388  template<typename PRM, int rangeDim>
389  class NavierStokesFunctionAdapterBase
391  Dune::PDELab::GridFunctionTraits<
392  typename PRM::Traits::GridView,
393  typename PRM::Traits::RangeField,
394  rangeDim,
395  Dune::FieldVector<typename PRM::Traits::RangeField,rangeDim>
396  >,
397  NavierStokesFunctionAdapterBase<PRM,rangeDim>
398  >
399  {
400  public:
403  typename PRM::Traits::GridView,
404  typename PRM::Traits::RangeField,
405  rangeDim,
406  Dune::FieldVector<typename PRM::Traits::RangeField,rangeDim>
407  > Traits;
408 
410  NavierStokesFunctionAdapterBase(PRM& prm)
411  : _prm(prm)
412  {}
413 
414  void setTime(const double time)
415  {
416  _prm.setTime(time);
417  }
418 
419  const PRM& parameters() const
420  {
421  return _prm;
422  }
423 
425  const typename Traits::GridViewType& getGridView () const
426  {
427  return _prm.gridView();
428  }
429 
430  private:
431  PRM& _prm;
432  };
433 
434 
435 #endif // DOXYGEN
436 
443  template<typename PRM>
445  : public NavierStokesFunctionAdapterBase<PRM,PRM::Traits::dimDomain>
446  {
447 
449  typedef NavierStokesFunctionAdapterBase<PRM,PRM::Traits::dimDomain> Base;
450 
451  using Base::parameters;
452 
453  public:
454 
455  typedef typename Base::Traits Traits;
456 
459  : Base(prm)
460  {}
461 
463  void evaluate (const typename Traits::ElementType& e,
464  const typename Traits::DomainType& x,
465  typename Traits::RangeType& y) const
466  {
467  y = parameters().g(e,x);
468  }
469 };
470 
471 
472 
473 #if 0
474 
476 template < typename PRM >
477 class NavierStokesDirichletFunctionAdapterFactory
478 {
479 public:
481  NavierStokesDirichletFunctionAdapter<PRM>,
482  NavierStokesPressureDirichletFunctionAdapter<PRM> >
483  BoundaryDirichletFunction;
484 
485 NavierStokesDirichletFunctionAdapterFactory(PRM & prm)
486 : v(prm), p(prm), df(v,p)
487 {}
488 
489 BoundaryDirichletFunction & dirichletFunction()
490 {
491  return df;
492 }
493 
494 private:
495 NavierStokesVelocityDirichletFunctionAdapter<PRM> v;
496 NavierStokesPressureDirichletFunctionAdapter<PRM> p;
497 BoundaryDirichletFunction df;
498 };
499 #endif
500 
501 
502 } // end namespace PDELab
503 } // end namespace Dune
504 
505 #endif
Traits::BoundaryCondition::Type bctype(const IG &is, const typename Traits::IntersectionDomain &x) const
boundary condition type from local intersection coordinate
Definition: stokesparameter.hh:194
void evaluate(const typename Traits::ElementType &e, const typename Traits::DomainType &x, typename Traits::RangeType &y) const
Evaluate dirichlet function.
Definition: stokesparameter.hh:463
Dune::FieldVector< DomainField, dimDomain-1 > IntersectionDomain
domain type
Definition: stokesparameter.hh:63
static const unsigned int value
Definition: gridfunctionspace/tags.hh:177
composite functions
Definition: function.hh:799
Definition: stokesparameter.hh:444
StokesBoundaryCondition BoundaryCondition
boundary type value
Definition: stokesparameter.hh:75
GV::Intersection Intersection
grid intersection type
Definition: stokesparameter.hh:80
RF RangeField
Export type for range field.
Definition: stokesparameter.hh:66
const IG & ig
Definition: constraints.hh:148
StokesVelocityDirichletConstraints(const PRM &_prm)
Constructor.
Definition: stokesparameter.hh:344
GV GridView
the grid view
Definition: stokesparameter.hh:48
bool isDirichlet(const I &intersection, const Dune::FieldVector< typename I::ctype, I::dimension-1 > &coord) const
Definition: stokesparameter.hh:376
Definition: adaptivity.hh:27
Definition: stokesparameter.hh:45
Definition: stokesparameter.hh:32
NavierStokesDefaultParameters(const RF &mu, const RF &rho, F &f, B &b, V &v, J &j)
Definition: stokesparameter.hh:167
bool isDirichlet(const I &intersection, const Dune::FieldVector< typename I::ctype, I::dimension-1 > &coord) const
Definition: stokesparameter.hh:349
const Entity & e
Definition: localfunctionspace.hh:111
Base::Traits Traits
Definition: stokesparameter.hh:455
Traits::VelocityRange f(const EG &e, const typename Traits::Domain &x) const
source term
Definition: stokesparameter.hh:185
GV::Grid::ctype DomainField
Export type for domain field.
Definition: stokesparameter.hh:57
Definition: stokesparameter.hh:335
StokesPressureDirichletConstraints(const PRM &_prm)
Constructor.
Definition: stokesparameter.hh:371
NavierStokesParameterTraits< GV, RF > Traits
Type traits.
Definition: stokesparameter.hh:151
Traits::RangeField rho(const IG &ig, const typename Traits::IntersectionDomain &x) const
Density value from local intersection coordinate.
Definition: stokesparameter.hh:229
Definition: stokesparameter.hh:34
Definition: stokesparameter.hh:143
Definition: constraintsparameters.hh:24
Traits::RangeField rho(const EG &eg, const typename Traits::Domain &x) const
Density value from local cell coordinate.
Definition: stokesparameter.hh:221
Dune::FieldVector< RF, 1 > PressureRange
pressure range type
Definition: stokesparameter.hh:72
const P & p
Definition: constraints.hh:147
Traits::VelocityRange g(const EG &e, const typename Traits::Domain &x) const
Dirichlet boundary condition value from local cell coordinate.
Definition: stokesparameter.hh:237
Type
Definition: stokesparameter.hh:33
GV::Traits::template Codim< 0 >::Entity Element
grid element type
Definition: stokesparameter.hh:78
Traits::RangeField mu(const IG &ig, const typename Traits::IntersectionDomain &x) const
Dynamic viscosity value from local intersection coordinate.
Definition: stokesparameter.hh:213
void setTime(RF time)
Definition: stokesparameter.hh:313
NavierStokesDefaultParameters(const Dune::ParameterTree &config, F &f, B &b, V &v, J &j)
Constructor.
Definition: stokesparameter.hh:154
Traits::RangeField mu(const EG &e, const typename Traits::Domain &x) const
Dynamic viscosity value from local cell coordinate.
Definition: stokesparameter.hh:205
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:362
leaf of a function tree
Definition: function.hh:576
Traits::RangeField g2(const EG &e, const typename Traits::Domain &x) const
pressure source term
Definition: stokesparameter.hh:259
NavierStokesVelocityFunctionAdapter(PRM &prm)
Constructor.
Definition: stokesparameter.hh:458
traits class holding the function signature, same as in local function
Definition: function.hh:175