13#ifndef DUMUX_TURBULENCE_PROPERTIES_HH 
   14#define DUMUX_TURBULENCE_PROPERTIES_HH 
   18#include<dune/common/fvector.hh> 
   25template<
class Scalar, 
unsigned dim, 
bool verbose = false>
 
   34                           const Dune::FieldVector<Scalar, dim> position,
 
   35                           const Scalar kinematicViscosity,
 
   38                           bool print=verbose)
 const 
   43        const Scalar re_x = 
reynoldsNumber(velocity, position[0], kinematicViscosity, 
false);
 
   44        const Scalar c_f = pow((2.0 * log10(re_x) - 0.65), -2.3); 
 
   45        const Scalar wallShearStress = 0.5 * c_f * density * velocity * velocity;
 
   46        const Scalar frictionVelocity = sqrt(wallShearStress / density);
 
   47        const Scalar yPlus = position[yCoordDim] * frictionVelocity / kinematicViscosity;
 
   50            std::cout << 
"turbulence properties at (";
 
   51            for (
unsigned int dimIdx = 0; dimIdx < dim; ++dimIdx)
 
   52                std::cout << position[dimIdx] << 
",";
 
   53            std::cout << 
")" << std::endl;
 
   54            std::cout << 
"estimated Re_x  : " << re_x << 
" [-]" << std::endl;
 
   55            std::cout << 
"estimated c_f   : " << c_f << 
" [-]" << std::endl;
 
   56            std::cout << 
"estimated tau_w : " << wallShearStress << 
" [kg/(m*s^2)]" << std::endl;
 
   57            std::cout << 
"estimated UStar : " << frictionVelocity << 
" [m/s]" << std::endl;
 
   58            std::cout << 
"estimated yPlus : " << yPlus << 
" [-]" << std::endl;
 
   59            std::cout << std::endl;
 
 
   69                          const Scalar kinematicViscosity,
 
   70                          bool print=verbose)
 const 
   77            std::cout << 
"estimated Re_d  : " << re_d << 
" [-]" << std::endl;
 
   78            std::cout << 
"estimated l_ent : " << 
entranceLength << 
" [m]"<< std::endl;
 
   79            std::cout << std::endl;
 
 
   88                          const Scalar charLengthScale,
 
   89                          const Scalar kinematicViscosity,
 
   90                          bool print=verbose)
 const 
   93        return abs(velocity * charLengthScale / kinematicViscosity);
 
 
  101                               bool print=verbose)
 const 
 
  117                                 bool print=verbose)
 const 
 
  133                                  const Scalar kinematicViscosity,
 
  134                                  bool print=verbose)
 const 
  137        const Scalar k = 1.5 * velocity * velocity
 
  141            std::cout << 
"estimated k     : " << k << 
" [m^2/s^2]" << std::endl;
 
 
  152                       const Scalar kinematicViscosity,
 
  153                       bool print=verbose)
 const 
  157        const Scalar factor = 0.1643; 
 
  161            std::cout << 
"estimated eps.  : " << epsilon << 
" [m^2/s^3]" << std::endl;
 
 
  173                           const Scalar kinematicViscosity,
 
  174                           bool print=verbose)
 const 
  178         const Scalar factor = 0.54772; 
 
  180         const Scalar omega = pow(k, 0.5) / (L*factor);
 
  183            std::cout << 
"estimated omega : " << omega << 
" [1/s]" << std::endl;
 
 
  195                          const Scalar kinematicViscosity,
 
  196                          bool print=verbose)
 const 
  206            std::cout << 
"estimated nu~   : " << 
viscosityTilde << 
" [m^2/s]" << std::endl;
 
 
 
This class contains different functions for estimating turbulence properties.
Definition turbulenceproperties.hh:27
Scalar viscosityTilde(const Scalar velocity, const Scalar diameter, const Scalar kinematicViscosity, bool print=verbose) const
Estimates the viscosity tilde based on a formula given in in the ANSYS Fluent user guide ANSYSUserGui...
Definition turbulenceproperties.hh:193
Scalar yPlusEstimation(const Scalar velocity, const Dune::FieldVector< Scalar, dim > position, const Scalar kinematicViscosity, const Scalar density, int yCoordDim=dim-1, bool print=verbose) const
Estimates dimensionless wall distance  based on a formula given in http://www.cfd-online....
Definition turbulenceproperties.hh:33
Scalar entranceLength(const Scalar velocity, const Scalar diameter, const Scalar kinematicViscosity, bool print=verbose) const
Estimates the entrance length for this pipe.
Definition turbulenceproperties.hh:67
Scalar turbulenceLengthScale(const Scalar charLengthScale, bool print=verbose) const
Estimates the turbulence length scale based on a formula given in the ANSYS Fluent user guide ANSYSUs...
Definition turbulenceproperties.hh:116
Scalar turbulentKineticEnergy(const Scalar velocity, const Scalar diameter, const Scalar kinematicViscosity, bool print=verbose) const
Estimates the turbulent kinetic energy based on a formula given in the ANSYS Fluent user guide ANSYSU...
Definition turbulenceproperties.hh:131
Scalar dissipationRate(const Scalar velocity, const Scalar diameter, const Scalar kinematicViscosity, bool print=verbose) const
Estimates the dissipation rate based on a formula given in the ANSYS Fluent user guide ANSYSUserGuide...
Definition turbulenceproperties.hh:171
Scalar turbulenceIntensity(const Scalar reynoldsNumber, bool print=verbose) const
Estimates the turbulence intensity based on a formula given in the ANSYS Fluent user guide ANSYSUserG...
Definition turbulenceproperties.hh:100
Scalar dissipation(const Scalar velocity, const Scalar diameter, const Scalar kinematicViscosity, bool print=verbose) const
Estimates the dissipation based on a formula given in the ANSYS Fluent user guide ANSYSUserGuide12.
Definition turbulenceproperties.hh:150
Scalar reynoldsNumber(const Scalar velocity, const Scalar charLengthScale, const Scalar kinematicViscosity, bool print=verbose) const
Calculates the Reynolds number.
Definition turbulenceproperties.hh:87
Geometry::ctype diameter(const Geometry &geo)
Computes the longest distance between points of a geometry.
Definition diameter.hh:26