4 #ifndef DUNE_PDELAB_LOCALOPERATOR_DGNAVIERSTOKESPARAMETER_HH 5 #define DUNE_PDELAB_LOCALOPERATOR_DGNAVIERSTOKESPARAMETER_HH 7 #include <dune/common/parametertreeparser.hh> 32 template<
typename GV,
typename RF,
typename F,
typename B,
typename V,
typename J,
33 bool navier =
false,
bool full_tensor =
false,
typename IP = DefaultInteriorPenalty<typename V::Traits::RangeFieldType> >
41 void initFromString(
const std::string & method)
43 std::string
s = method;
44 std::transform(s.begin(), s.end(), s.begin(), tolower);
47 if (s.find(
"nipg") != std::string::npos)
53 if (s.find(
"sipg") != std::string::npos)
67 if (3 == sscanf(s.c_str(),
"%d %lg %lg", &_epsilon, &sigma, &beta))
70 DUNE_THROW(Dune::Exception,
"Unknown DG type " << method);
80 F&
f, B& b, V& v, J&
j, IP& ip) :
81 Base(mu,rho,f,b,v,j) ,
84 initFromString(method);
89 F&
f, B& b, V& v, J&
j, IP& ip) :
90 Base(configuration,f,b,v,j) ,
92 _epsilon(configuration.get<int>(
"epsilon"))
112 return _ip.getFaceIP(ig);
123 return _ip.getFaceIP(ig);
141 namespace NavierStokesDGImp{
158 template<
typename PRM,
typename Dummy =
void>
161 template<
typename IntersectionGeometry>
162 static typename PRM::Traits::RangeField
166 const typename PRM::Traits::IntersectionDomain& )
173 template<
typename PRM>
175 <PRM,typename
Dune::enable_if<PRM::enable_variable_slip>::type>
177 template<
typename IntersectionGeometry>
178 static typename PRM::Traits::RangeField
182 const typename PRM::Traits::IntersectionDomain& x)
184 return prm.boundarySlip(ig,x);
int epsilonIPSymmetryFactor()
Definition: dgnavierstokesparameter.hh:129
Dune::FieldVector< DomainField, dimDomain-1 > IntersectionDomain
domain type
Definition: stokesparameter.hh:63
Base::Traits Traits
Traits class.
Definition: dgnavierstokesparameter.hh:76
Compile-time switch for the boundary slip factor.
Definition: dgnavierstokesparameter.hh:159
RF RangeField
Export type for range field.
Definition: stokesparameter.hh:66
const IG & ig
Definition: constraints.hh:148
Traits::RangeField getFaceIP(const I &ig, const typename Traits::IntersectionDomain &) const
Interior penalty parameter.
Definition: dgnavierstokesparameter.hh:121
Definition: adaptivity.hh:27
Traits::RangeField getFaceIP(const I &ig) const
Interior penalty parameter.
Definition: dgnavierstokesparameter.hh:110
Definition: stokesparameter.hh:45
DGNavierStokesParameters(const Dune::ParameterTree &configuration, F &f, B &b, V &v, J &j, IP &ip)
Constructor.
Definition: dgnavierstokesparameter.hh:88
Traits::VelocityRange f(const EG &e, const typename Traits::Domain &x) const
source term
Definition: stokesparameter.hh:185
Parameter class for local operator DGNavierStokes.
Definition: dgnavierstokesparameter.hh:34
Traits::VelocityRange j(const IG &ig, const typename Traits::IntersectionDomain &x, const typename Traits::Domain &normal) const
Neumann boundary condition (stress)
Definition: stokesparameter.hh:143
Traits::RangeField rho(const EG &eg, const typename Traits::Domain &x) const
Density value from local cell coordinate.
Definition: stokesparameter.hh:221
Wrap intersection.
Definition: geometrywrapper.hh:56
const std::string s
Definition: function.hh:1102
Traits::RangeField incompressibilityScaling(typename Traits::RangeField dt) const
Rescaling factor for the incompressibility equation.
Definition: dgnavierstokesparameter.hh:98
Traits::RangeField mu(const EG &e, const typename Traits::Domain &x) const
Dynamic viscosity value from local cell coordinate.
Definition: stokesparameter.hh:205
DGNavierStokesParameters(const std::string &method, const RF mu, const RF rho, F &f, B &b, V &v, J &j, IP &ip)
Constructor.
Definition: dgnavierstokesparameter.hh:79