dune-pdelab  2.7-git
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,
37  SlipVelocity = 3
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,TypeTree::StaticDegree<GF>::value>
107  evaluateVelocityGridFunction(const GF& gf,
108  const Entity& e,
109  const Domain& x)
110  {
111  static_assert(Domain::dimension == TypeTree::StaticDegree<GF>::value,"dimension of function range does not match grid dimension");
112  FieldVector<typename GF::template Child<0>::Type::Traits::RangeFieldType,TypeTree::StaticDegree<GF>::value> 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<RF>("rho"))
160  , _mu(config.get<RF>("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,x);
189  }
190 
192  template<typename IG>
194  bctype(const IG& is,
195  const typename Traits::IntersectionDomain& x) const
196  {
197  typename B::Traits::BoundaryCondition::Type 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 EG>
246  typename Traits::RangeField
247  g2(const EG& e, const typename Traits::Domain& x) const
248  {
249  return 0;
250  }
251 
252 #ifdef DOXYGEN
253 
255  template<typename IG>
256  typename Traits::VelocityRange>
257  j(const IG& ig,
258  const typename Traits::IntersectionDomain& x,
259  const typename Traits::Domain& normal) const;
260 
261 #else // DOXYGEN
262 
264  template<typename IG>
265  typename std::enable_if<
266  J::Traits::dimRange == 1 &&
267  (GV::dimension > 1) &&
268  AlwaysTrue<IG>::value, // required to force lazy evaluation
269  typename Traits::VelocityRange
270  >::type
271  j(const IG& ig,
272  const typename Traits::IntersectionDomain& x,
273  typename Traits::Domain normal) const
274  {
275  typename J::Traits::RangeType r;
276  auto e = ig.inside();
277  _j.evaluate(e,ig.geometryInInside().global(x),r);
278  normal *= r;
279  return normal;
280  }
281 
283  template<typename IG>
284  typename std::enable_if<
285  J::Traits::dimRange == GV::dimension &&
286  AlwaysTrue<IG>::value, // required to force lazy evaluation
287  typename Traits::VelocityRange
288  >::type
289  j(const IG& ig,
290  const typename Traits::IntersectionDomain& x,
291  const typename Traits::Domain& normal) const
292  {
293  auto e = ig.inside();
294  typename J::Traits::RangeType y;
295  _j.evaluate(e,ig.geometryInInside().global(x),y);
296  return y;
297  }
298 
299 #endif // DOXYGEN
300 
301  void setTime(RF time)
302  {
303  _f.setTime(time);
304  _b.setTime(time);
305  _v.setTime(time);
306  _j.setTime(time);
307  }
308 
309  private:
310  const RF _rho;
311  const RF _mu;
312  F& _f;
313  B& _b;
314  V& _v;
315  J& _j;
316  };
317 
318 
322  template<typename PRM>
325  {
326  private:
327  const PRM & prm_;
328 
329  public:
330 
333  : prm_(_prm) { }
334 
336  template<typename I>
337  bool isDirichlet(const I & intersection,
338  const Dune::FieldVector<typename I::ctype, I::mydimension> & coord) const
339  {
340  StokesBoundaryCondition::Type bctype = prm_.bctype(intersection,coord);
342  }
343 
344  };
345 
349  template<typename PRM>
352  {
353  private:
354  const PRM & prm_;
355 
356  public:
357 
360  : prm_(_prm) { }
361 
363  template<typename I>
364  bool isDirichlet(const I & intersection,
365  const Dune::FieldVector<typename I::ctype, I::mydimension> & coord) const
366  { return false; }
367  };
368 
369 
370 
371 #ifndef DOXYGEN
372 
376  template<typename PRM, int rangeDim>
377  class NavierStokesFunctionAdapterBase
379  Dune::PDELab::GridFunctionTraits<
380  typename PRM::Traits::GridView,
381  typename PRM::Traits::RangeField,
382  rangeDim,
383  Dune::FieldVector<typename PRM::Traits::RangeField,rangeDim>
384  >,
385  NavierStokesFunctionAdapterBase<PRM,rangeDim>
386  >
387  {
388  public:
391  typename PRM::Traits::GridView,
392  typename PRM::Traits::RangeField,
393  rangeDim,
394  Dune::FieldVector<typename PRM::Traits::RangeField,rangeDim>
395  > Traits;
396 
398  NavierStokesFunctionAdapterBase(PRM& prm)
399  : _prm(prm)
400  {}
401 
402  void setTime(const double time)
403  {
404  _prm.setTime(time);
405  }
406 
407  const PRM& parameters() const
408  {
409  return _prm;
410  }
411 
413  const typename Traits::GridViewType& getGridView () const
414  {
415  return _prm.gridView();
416  }
417 
418  private:
419  PRM& _prm;
420  };
421 
422 
423 #endif // DOXYGEN
424 
431  template<typename PRM>
433  : public NavierStokesFunctionAdapterBase<PRM,PRM::Traits::dimDomain>
434  {
435 
437  typedef NavierStokesFunctionAdapterBase<PRM,PRM::Traits::dimDomain> Base;
438 
439  using Base::parameters;
440 
441  public:
442 
443  typedef typename Base::Traits Traits;
444 
447  : Base(prm)
448  {}
449 
451  void evaluate (const typename Traits::ElementType& e,
452  const typename Traits::DomainType& x,
453  typename Traits::RangeType& y) const
454  {
455  y = parameters().g(e,x);
456  }
457  };
458 
459 
460 
461 #if 0
464  template < typename PRM >
465  class NavierStokesDirichletFunctionAdapterFactory
466  {
467  public:
469  NavierStokesDirichletFunctionAdapter<PRM>,
470  NavierStokesPressureDirichletFunctionAdapter<PRM> >
471  BoundaryDirichletFunction;
472 
473  NavierStokesDirichletFunctionAdapterFactory(PRM & prm)
474  : v(prm), p(prm), df(v,p)
475  {}
476 
477  BoundaryDirichletFunction & dirichletFunction()
478  {
479  return df;
480  }
481 
482  private:
483  NavierStokesVelocityDirichletFunctionAdapter<PRM> v;
484  NavierStokesPressureDirichletFunctionAdapter<PRM> p;
485  BoundaryDirichletFunction df;
486  };
487 #endif
488 
489 
490  } // end namespace PDELab
491 } // end namespace Dune
492 
493 #endif // DUNE_PDELAB_LOCALOPERATOR_STOKESPARAMETER_HH
const P & p
Definition: constraints.hh:148
const IG & ig
Definition: constraints.hh:149
const Entity & e
Definition: localfunctionspace.hh:121
For backward compatibility – Do not use this!
Definition: adaptivity.hh:28
traits class holding the function signature, same as in local function
Definition: function.hh:183
leaf of a function tree
Definition: function.hh:302
composite functions
Definition: function.hh:539
Definition: constraintsparameters.hh:26
Definition: stokesparameter.hh:32
Type
Definition: stokesparameter.hh:33
@ SlipVelocity
Definition: stokesparameter.hh:37
@ StressNeumann
Definition: stokesparameter.hh:36
@ VelocityDirichlet
Definition: stokesparameter.hh:35
@ DoNothing
Definition: stokesparameter.hh:34
Definition: stokesparameter.hh:46
GV GridView
the grid view
Definition: stokesparameter.hh:48
GV::Grid::ctype DomainField
Export type for domain field.
Definition: stokesparameter.hh:57
Dune::FieldVector< DomainField, dimDomain > Domain
domain type
Definition: stokesparameter.hh:60
RF RangeField
Export type for range field.
Definition: stokesparameter.hh:66
@ dimDomain
dimension of the domain
Definition: stokesparameter.hh:53
GV::Traits::template Codim< 0 >::Entity Element
grid element type
Definition: stokesparameter.hh:78
Dune::FieldVector< RF, GV::dimensionworld > VelocityRange
deformation range type
Definition: stokesparameter.hh:69
StokesBoundaryCondition BoundaryCondition
boundary type value
Definition: stokesparameter.hh:75
Dune::FieldVector< RF, 1 > PressureRange
pressure range type
Definition: stokesparameter.hh:72
GV::Intersection Intersection
grid intersection type
Definition: stokesparameter.hh:80
Dune::FieldVector< DomainField, dimDomain-1 > IntersectionDomain
domain type
Definition: stokesparameter.hh:63
Definition: stokesparameter.hh:144
void setTime(RF time)
Definition: stokesparameter.hh:301
Traits::RangeField rho(const IG &ig, const typename Traits::IntersectionDomain &x) const
Density value from local intersection coordinate.
Definition: stokesparameter.hh:229
Traits::VelocityRange g(const EG &e, const typename Traits::Domain &x) const
Dirichlet boundary condition value from local cell coordinate.
Definition: stokesparameter.hh:237
Traits::RangeField rho(const EG &eg, const typename Traits::Domain &x) const
Density value from local cell coordinate.
Definition: stokesparameter.hh:221
NavierStokesDefaultParameters(const Dune::ParameterTree &config, F &f, B &b, V &v, J &j)
Constructor.
Definition: stokesparameter.hh:154
static const bool assemble_navier
Definition: stokesparameter.hh:147
Traits::VelocityRange j(const IG &ig, const typename Traits::IntersectionDomain &x, const typename Traits::Domain &normal) const
Neumann boundary condition (stress)
NavierStokesParameterTraits< GV, RF > Traits
Type traits.
Definition: stokesparameter.hh:151
static const bool assemble_full_tensor
Definition: stokesparameter.hh:148
Traits::VelocityRange f(const EG &e, const typename Traits::Domain &x) const
source term
Definition: stokesparameter.hh:185
NavierStokesDefaultParameters(const RF &mu, const RF &rho, F &f, B &b, V &v, J &j)
Definition: stokesparameter.hh:167
Traits::RangeField mu(const EG &e, const typename Traits::Domain &x) const
Dynamic viscosity value from local cell coordinate.
Definition: stokesparameter.hh:205
Traits::RangeField mu(const IG &ig, const typename Traits::IntersectionDomain &x) const
Dynamic viscosity value from local intersection coordinate.
Definition: stokesparameter.hh:213
Traits::RangeField g2(const EG &e, const typename Traits::Domain &x) const
pressure source term
Definition: stokesparameter.hh:247
Traits::BoundaryCondition::Type bctype(const IG &is, const typename Traits::IntersectionDomain &x) const
boundary condition type from local intersection coordinate
Definition: stokesparameter.hh:194
Definition: stokesparameter.hh:325
StokesVelocityDirichletConstraints(const PRM &_prm)
Constructor.
Definition: stokesparameter.hh:332
bool isDirichlet(const I &intersection, const Dune::FieldVector< typename I::ctype, I::mydimension > &coord) const
Definition: stokesparameter.hh:337
Definition: stokesparameter.hh:352
StokesPressureDirichletConstraints(const PRM &_prm)
Constructor.
Definition: stokesparameter.hh:359
bool isDirichlet(const I &intersection, const Dune::FieldVector< typename I::ctype, I::mydimension > &coord) const
Definition: stokesparameter.hh:364
Definition: stokesparameter.hh:434
void evaluate(const typename Traits::ElementType &e, const typename Traits::DomainType &x, typename Traits::RangeType &y) const
Evaluate dirichlet function.
Definition: stokesparameter.hh:451
NavierStokesVelocityFunctionAdapter(PRM &prm)
Constructor.
Definition: stokesparameter.hh:446
Base::Traits Traits
Definition: stokesparameter.hh:443
static const unsigned int value
Definition: gridfunctionspace/tags.hh:139