1 #ifndef DUNE_PDELAB_LOCALOPERATOR_STOKESDGPARAMETER_HH
2 #define DUNE_PDELAB_LOCALOPERATOR_STOKESDGPARAMETER_HH
4 #include <dune/common/parametertreeparser.hh>
25 template<
typename GV,
typename RF,
typename F,
typename B,
typename V,
typename J,
26 typename IP = DefaultInteriorPenalty<typename V::Traits::RangeFieldType> >
34 void initFromString(
const std::string & method)
36 std::string
s = method;
37 std::transform(s.begin(), s.end(), s.begin(), tolower);
40 if (s.find(
"nipg") != std::string::npos)
46 if (s.find(
"sipg") != std::string::npos)
60 if (3 == sscanf(s.c_str(),
"%d %lg %lg", &_epsilon, &sigma, &beta))
63 DUNE_THROW(Dune::Exception,
"Unknown DG type " << method);
68 F &
f, B & b, V & v, J &
j, IP & ip)
69 :
Base(mu,rho,f,b,v,j),
72 initFromString(method);
82 F &
f, B & b, V & v, J &
j, IP & ip)
83 :
Base(mu,1.0,f,b,v,j),
86 initFromString(method);
90 F &
f, B & b, V & v, J &
j, IP & ip)
91 :
Base(configuration,f,b,v,j),
93 _epsilon(configuration.get<int>(
"epsilon"))
113 return _ip.getFaceIP(ig);
124 return _ip.getFaceIP(ig);
153 template<
typename GV,
typename RF,
typename F,
typename B,
typename V,
typename J,
154 typename IP = DefaultInteriorPenalty<typename V::Traits::RangeField> >
165 F &
f, B & b, V & v, J &
j, IP & ip)
166 :
Base(method,mu,rho,f,b,v,j,ip)
170 F &
f, B & b, V & v, J &
j, IP & ip)
171 :
Base(configuration,f,b,v,j,ip)
177 namespace NavierStokesDGImp{
194 template<
typename PRM,
typename Dummy =
void>
197 template<
typename IntersectionGeometry>
198 static typename PRM::Traits::RangeField
202 const typename PRM::Traits::IntersectionDomain& )
209 template<
typename PRM>
211 <PRM,typename Dune::enable_if<PRM::enable_variable_slip>::type>
213 template<
typename IntersectionGeometry>
214 static typename PRM::Traits::RangeField
218 const typename PRM::Traits::IntersectionDomain& x)
220 return prm.boundarySlip(ig,x);
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
Compile-time switch for the boundary slip factor.
Definition: stokesdgparameter.hh:195
Traits::VelocityRange j(const IG &ig, const typename Traits::IntersectionDomain &x, const typename Traits::Domain &normal) const
Neumann boundary condition (stress)
RF RangeField
Export type for range field.
Definition: stokesparameter.hh:66
Definition: stokesparameter.hh:45
static PRM::Traits::RangeField boundarySlip(const PRM &, const IntersectionGeometry &, const typename PRM::Traits::IntersectionDomain &)
Definition: stokesdgparameter.hh:200
NavierStokesDGParameters(const std::string &method, const RF mu, const RF rho, F &f, B &b, V &v, J &j, IP &ip)
Definition: stokesdgparameter.hh:164
Traits::RangeField rho(const EG &eg, const typename Traits::Domain &x) const
Density value from local cell coordinate.
Definition: stokesparameter.hh:221
StokesDGParameters(const std::string &method, const RF mu, F &f, B &b, V &v, J &j, IP &ip)
Definition: stokesdgparameter.hh:81
Traits::RangeField getFaceIP(const I &ig) const
Interior penalty parameter.
Definition: stokesdgparameter.hh:111
const IG & ig
Definition: common/constraints.hh:146
Definition: stokesparameter.hh:143
Traits::RangeField incompressibilityScaling(typename Traits::RangeField dt) const
Rescaling factor for the incompressibility equation.
Definition: stokesdgparameter.hh:99
Traits::RangeField mu(const EG &e, const typename Traits::Domain &x) const
Dynamic viscosity value from local cell coordinate.
Definition: stokesparameter.hh:205
int epsilonIPSymmetryFactor()
Definition: stokesdgparameter.hh:130
StokesDGParameters(const std::string &method, const RF mu, const RF rho, F &f, B &b, V &v, J &j, IP &ip)
Definition: stokesdgparameter.hh:67
NavierStokesDGParameters(const Dune::ParameterTree &configuration, F &f, B &b, V &v, J &j, IP &ip)
Definition: stokesdgparameter.hh:169
Parameter class for local operator StokesDG.
Definition: stokesdgparameter.hh:27
Traits::RangeField getFaceIP(const I &ig, const typename Traits::IntersectionDomain &) const
Interior penalty parameter.
Definition: stokesdgparameter.hh:122
Wrap intersection.
Definition: geometrywrapper.hh:56
const std::string s
Definition: function.hh:1103
Parameter class for local operator NavierStokesDG.
Definition: stokesdgparameter.hh:155
StokesDGParameters(const Dune::ParameterTree &configuration, F &f, B &b, V &v, J &j, IP &ip)
Definition: stokesdgparameter.hh:89
Base::Traits Traits
Traits class.
Definition: stokesdgparameter.hh:79