dune-pdelab  2.4-dev
dgnavierstokesparameter.hh
Go to the documentation of this file.
1 // -*- tab-width: 2; indent-tabs-mode: nil -*-
2 // vi: set et ts=2 sw=2 sts=2:
3 
4 #ifndef DUNE_PDELAB_LOCALOPERATOR_DGNAVIERSTOKESPARAMETER_HH
5 #define DUNE_PDELAB_LOCALOPERATOR_DGNAVIERSTOKESPARAMETER_HH
6 
7 #include <dune/common/parametertreeparser.hh>
11 
12 namespace Dune {
13  namespace PDELab {
14 
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> >
35  public NavierStokesDefaultParameters<GV,RF,F,B,V,J,navier,full_tensor>
36  {
37 
40 
41  void initFromString(const std::string & method)
42  {
43  std::string s = method;
44  std::transform(s.begin(), s.end(), s.begin(), tolower);
45 
46  // nipg (epsilon=1) 2d p1 -> Klaus sagt sollte auch sigma 1 klappen
47  if (s.find("nipg") != std::string::npos)
48  {
49  _epsilon = 1;
50  return;
51  }
52  // sipg (epsilon=-1) 2d p1 -> Klaus sagt sigma=3.9irgendwas
53  if (s.find("sipg") != std::string::npos)
54  {
55  _epsilon = -1;
56  return;
57  }
58  // obb sigma = 0, epsilon =
59  if (s == "obb")
60  {
61  _epsilon = 1;
62  return;
63  }
64  // extract parameters
65  {
66  double sigma, beta;
67  if (3 == sscanf(s.c_str(), "%d %lg %lg", &_epsilon, &sigma, &beta))
68  return;
69  }
70  DUNE_THROW(Dune::Exception, "Unknown DG type " << method);
71  }
72 
73  public :
74 
76  typedef typename Base::Traits Traits;
77 
79  DGNavierStokesParameters(const std::string& method, const RF mu, const RF rho,
80  F& f, B& b, V& v, J& j, IP& ip) :
81  Base(mu,rho,f,b,v,j) ,
82  _ip(ip)
83  {
84  initFromString(method);
85  }
86 
88  DGNavierStokesParameters(const Dune::ParameterTree& configuration,
89  F& f, B& b, V& v, J& j, IP& ip) :
90  Base(configuration,f,b,v,j) ,
91  _ip(ip) ,
92  _epsilon(configuration.get<int>("epsilon"))
93  {}
94 
95 
97  typename Traits::RangeField
99  {
100  typename Traits::RangeField y(1.0 / dt);
101  return y;
102  }
103 
108  template<typename I>
109  typename Traits::RangeField
110  getFaceIP(const I & ig) const
111  {
112  return _ip.getFaceIP(ig);
113  }
114 
119  template<typename I>
120  typename Traits::RangeField
121  getFaceIP(const I & ig, const typename Traits::IntersectionDomain& ) const
122  {
123  return _ip.getFaceIP(ig);
124  }
125 
128  int
130  {
131  return _epsilon;
132  }
133 
134  private:
135 
136  IP & _ip; // Interior penalty
137  int _epsilon; // IP symmetry factor
138  }; // end class DGNavierStokesParameters
139 
140 
141  namespace NavierStokesDGImp{
158  template< typename PRM, typename Dummy = void>
160 
161  template<typename IntersectionGeometry>
162  static typename PRM::Traits::RangeField
164  (const PRM & ,
165  const IntersectionGeometry& ,
166  const typename PRM::Traits::IntersectionDomain& )
167  {
168  return 1.0;
169  }
170 
171  };
172 
173  template< typename PRM>
175  <PRM,typename Dune::enable_if<PRM::enable_variable_slip>::type>
176  {
177  template<typename IntersectionGeometry>
178  static typename PRM::Traits::RangeField
180  (const PRM & prm,
181  const IntersectionGeometry& ig,
182  const typename PRM::Traits::IntersectionDomain& x)
183  {
184  return prm.boundarySlip(ig,x);
185  }
186 
187  };
189  }
190 
191  } // end namespace PDELab
192 } // end namespace Dune
193 #endif
DGNavierStokesParameters(const Dune::ParameterTree &configuration, F &f, B &b, V &v, J &j, IP &ip)
Constructor.
Definition: dgnavierstokesparameter.hh:88
Definition: stokesparameter.hh:143
Traits::RangeField mu(const EG &e, const typename Traits::Domain &x) const
Dynamic viscosity value from local cell coordinate.
Definition: stokesparameter.hh:205
Traits::VelocityRange j(const IG &ig, const typename Traits::IntersectionDomain &x, const typename Traits::Domain &normal) const
Neumann boundary condition (stress)
const IG & ig
Definition: constraints.hh:147
Traits::RangeField getFaceIP(const I &ig) const
Interior penalty parameter.
Definition: dgnavierstokesparameter.hh:110
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
Definition: stokesparameter.hh:45
Base::Traits Traits
Traits class.
Definition: dgnavierstokesparameter.hh:76
static PRM::Traits::RangeField boundarySlip(const PRM &, const IntersectionGeometry &, const typename PRM::Traits::IntersectionDomain &)
Definition: dgnavierstokesparameter.hh:164
Traits::RangeField getFaceIP(const I &ig, const typename Traits::IntersectionDomain &) const
Interior penalty parameter.
Definition: dgnavierstokesparameter.hh:121
RF RangeField
Export type for range field.
Definition: stokesparameter.hh:66
Parameter class for local operator DGNavierStokes.
Definition: dgnavierstokesparameter.hh:34
Traits::RangeField incompressibilityScaling(typename Traits::RangeField dt) const
Rescaling factor for the incompressibility equation.
Definition: dgnavierstokesparameter.hh:98
Definition: adaptivity.hh:27
Wrap intersection.
Definition: geometrywrapper.hh:56
int epsilonIPSymmetryFactor()
Definition: dgnavierstokesparameter.hh:129
Traits::RangeField rho(const EG &eg, const typename Traits::Domain &x) const
Density value from local cell coordinate.
Definition: stokesparameter.hh:221
const std::string s
Definition: function.hh:1102
Compile-time switch for the boundary slip factor.
Definition: dgnavierstokesparameter.hh:159
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