1 #ifndef DUNE_PDELAB_LOCALOPERATOR_STOKESPARAMETER_HH
2 #define DUNE_PDELAB_LOCALOPERATOR_STOKESPARAMETER_HH
4 #include <dune/common/parametertree.hh>
44 template<
typename GV,
typename RF>
60 typedef Dune::FieldVector<DomainField,dimDomain>
Domain;
78 typedef typename GV::Traits::template Codim<0>::Entity
Element;
89 template<
typename GF,
typename Entity,
typename Domain>
90 typename GF::Traits::RangeType
91 evaluateVelocityGridFunction(
const GF& gf,
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;
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,
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)
116 gf.child(d).evaluate(e,x,cy);
142 template <
typename GV,
typename RF,
typename F,
typename B,
typename V,
typename J,
bool navier = false,
bool tensor = false>
159 : _rho(config.get<double>(
"rho"))
160 , _mu(config.get<double>(
"mu"))
183 template<
typename EG>
187 typename F::Traits::RangeType fvalue;
188 return evaluateVelocityGridFunction(_f,e.entity(),x);
192 template<
typename IG>
197 typename B::Traits::RangeType y;
203 template<
typename EG>
211 template<
typename IG>
219 template<
typename EG>
227 template<
typename IG>
235 template<
typename EG>
239 typename V::Traits::RangeType y;
240 _v.evaluate(e.entity(),x,y);
245 template<
typename IG>
249 typename IG::EntityPointer ep = ig.inside();
250 typename V::Traits::RangeType y;
251 _v.evaluate(*ep,ig.geometryInInside().global(x),y);
256 template<
typename EG>
266 template<
typename IG>
275 template<
typename IG>
277 J::Traits::dimRange == 1 &&
278 (GV::dimension > 1) &&
286 typename J::Traits::RangeType r;
287 typename IG::EntityPointer ep = ig.inside();
288 _j.evaluate(*ep,ig.geometryInInside().global(x),r);
294 template<
typename IG>
296 J::Traits::dimRange == GV::dimension &&
304 typename IG::EntityPointer ep = ig.inside();
305 typename J::Traits::RangeType y;
306 _j.evaluate(*ep,ig.geometryInInside().global(x),y);
333 template<
typename PRM>
349 const Dune::FieldVector<typename I::ctype, I::dimension-1> & coord)
const
360 template<
typename PRM>
376 const Dune::FieldVector<typename I::ctype, I::dimension-1> & coord)
const
387 template<
typename PRM,
int rangeDim>
388 class NavierStokesFunctionAdapterBase
390 Dune::PDELab::GridFunctionTraits<
391 typename PRM::Traits::GridView,
392 typename PRM::Traits::RangeField,
394 Dune::FieldVector<typename PRM::Traits::RangeField,rangeDim>
396 NavierStokesFunctionAdapterBase<PRM,rangeDim>
402 typename PRM::Traits::GridView,
403 typename PRM::Traits::RangeField,
405 Dune::FieldVector<typename PRM::Traits::RangeField,rangeDim>
409 NavierStokesFunctionAdapterBase(PRM& prm)
413 void setTime(
const double time)
418 const PRM& parameters()
const
424 const typename Traits::GridViewType& getGridView ()
const
426 return _prm.gridView();
442 template<
typename PRM>
444 :
public NavierStokesFunctionAdapterBase<PRM,PRM::Traits::dimDomain>
448 typedef NavierStokesFunctionAdapterBase<PRM,PRM::Traits::dimDomain> Base;
450 using Base::parameters;
462 void evaluate (
const typename Traits::ElementType& e,
463 const typename Traits::DomainType& x,
464 typename Traits::RangeType& y)
const
466 y = parameters().g(e,x);
475 template <
typename PRM >
476 class NavierStokesDirichletFunctionAdapterFactory
480 NavierStokesDirichletFunctionAdapter<PRM>,
481 NavierStokesPressureDirichletFunctionAdapter<PRM> >
482 BoundaryDirichletFunction;
484 NavierStokesDirichletFunctionAdapterFactory(PRM & prm)
485 : v(prm),
p(prm), df(v,
p)
488 BoundaryDirichletFunction & dirichletFunction()
494 NavierStokesVelocityDirichletFunctionAdapter<PRM> v;
495 NavierStokesPressureDirichletFunctionAdapter<PRM>
p;
496 BoundaryDirichletFunction df;
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
Definition: stokesparameter.hh:36
const EG & eg
Definition: common/constraints.hh:277
Definition: stokesparameter.hh:37
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:35
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