qZeta.H
Go to the documentation of this file.
1 /*---------------------------------------------------------------------------*\
2  ========= |
3  \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
4  \\ / O peration |
5  \\ / A nd | Copyright (C) 2011-2015 OpenFOAM Foundation
6  \\/ M anipulation |
7 -------------------------------------------------------------------------------
8 License
9  This file is part of OpenFOAM.
10 
11  OpenFOAM is free software: you can redistribute it and/or modify it
12  under the terms of the GNU General Public License as published by
13  the Free Software Foundation, either version 3 of the License, or
14  (at your option) any later version.
15 
16  OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
17  ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
18  FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
19  for more details.
20 
21  You should have received a copy of the GNU General Public License
22  along with OpenFOAM. If not, see <http://www.gnu.org/licenses/>.
23 
24 Class
25  Foam::incompressible::RASModels::qZeta
26 
27 Group
28  grpIcoRASTurbulence
29 
30 Description
31  Gibson and Dafa'Alla's q-zeta two-equation low-Re turbulence model
32  for incompressible flows
33 
34  This turbulence model is described in:
35  \verbatim
36  Dafa'Alla, A.A., Juntasaro, E. & Gibson, M.M. (1996).
37  Calculation of oscillating boundary layers with the
38  q-zeta turbulence model.
39  Engineering Turbulence Modelling and Experiments 3:
40  Proceedings of the Third International Symposium,
41  Crete, Greece, May 27-29, 141.
42  \endverbatim
43  which is a development of the original q-zeta model described in:
44  \verbatim
45  Gibson, M. M., & Dafa'Alla, A. A. (1995).
46  Two-equation model for turbulent wall flow.
47  AIAA journal, 33(8), 1514-1518.
48  \endverbatim
49 
50 SourceFiles
51  qZeta.C
52 
53 \*---------------------------------------------------------------------------*/
54 
55 #ifndef qZeta_H
56 #define qZeta_H
57 
59 #include "eddyViscosity.H"
60 
61 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
62 
63 namespace Foam
64 {
65 namespace incompressible
66 {
67 namespace RASModels
68 {
69 
70 /*---------------------------------------------------------------------------*\
71  Class qZeta Declaration
72 \*---------------------------------------------------------------------------*/
73 
74 class qZeta
75 :
76  public eddyViscosity<incompressible::RASModel>
77 {
78 
79 protected:
80 
81  // Protected data
82 
83  // Model coefficients
84 
90 
91  //- Lower limit of q
93 
94  //- Lower limit of zeta
96 
97  // Fields
98 
104 
105 
106  // Protected Member Functions
107 
108  tmp<volScalarField> fMu() const;
109  tmp<volScalarField> f2() const;
110  virtual void correctNut();
111 
112 
113 public:
114 
115  //- Runtime type information
116  TypeName("qZeta");
117 
118  // Constructors
119 
120  //- Construct from components
121  qZeta
122  (
123  const geometricOneField& alpha,
124  const geometricOneField& rho,
125  const volVectorField& U,
126  const surfaceScalarField& alphaRhoPhi,
127  const surfaceScalarField& phi,
128  const transportModel& transport,
129  const word& propertiesName = turbulenceModel::propertiesName,
130  const word& type = typeName
131  );
132 
133 
134  //- Destructor
135  virtual ~qZeta()
136  {}
137 
138 
139  // Member Functions
140 
141  //- Read RASProperties dictionary
142  virtual bool read();
143 
144  //- Return the lower allowable limit for q (default: SMALL)
145  const dimensionedScalar& qMin() const
146  {
147  return qMin_;
148  }
149 
150  //- Return the lower allowable limit for zeta (default: SMALL)
151  const dimensionedScalar& zetaMin() const
152  {
153  return zetaMin_;
154  }
155 
156  //- Allow qMin to be changed
158  {
159  return qMin_;
160  }
161 
162  //- Allow zetaMin to be changed
164  {
165  return zetaMin_;
166  }
167 
168  //- Return the effective diffusivity for q
169  tmp<volScalarField> DqEff() const
170  {
171  return tmp<volScalarField>
172  (
173  new volScalarField("DqEff", nut_ + nu())
174  );
175  }
176 
177  //- Return the effective diffusivity for epsilon
179  {
180  return tmp<volScalarField>
181  (
182  new volScalarField("DzetaEff", nut_/sigmaZeta_ + nu())
183  );
184  }
185 
186  //- Return the turbulence kinetic energy
187  virtual tmp<volScalarField> k() const
188  {
189  return k_;
190  }
191 
192  //- Return the turbulence kinetic energy dissipation rate
193  virtual tmp<volScalarField> epsilon() const
194  {
195  return epsilon_;
196  }
198  virtual const volScalarField& q() const
199  {
200  return q_;
201  }
203  virtual const volScalarField& zeta() const
204  {
205  return zeta_;
206  }
207 
208  //- Solve the turbulence equations and correct the turbulence viscosity
209  virtual void correct();
210 };
211 
212 
213 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
214 
215 } // End namespace RASModels
216 } // End namespace incompressible
217 } // End namespace Foam
218 
219 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
220 
221 #endif
222 
223 // ************************************************************************* //
surfaceScalarField & phi
U
Definition: pEqn.H:83
virtual void correct()
Solve the turbulence equations and correct the turbulence viscosity.
Definition: qZeta.C:239
A simple wrapper around bool so that it can be read as a word: true/false, on/off, yes/no, y/n, t/f, or none.
Definition: Switch.H:60
tmp< volScalarField > fMu() const
Definition: qZeta.C:46
TypeName("qZeta")
Runtime type information.
const dimensionedScalar & zetaMin() const
Return the lower allowable limit for zeta (default: SMALL)
Definition: qZeta.H:150
virtual ~qZeta()
Destructor.
Definition: qZeta.H:134
Eddy viscosity turbulence model base class.
Definition: eddyViscosity.H:52
const dimensionedScalar & qMin() const
Return the lower allowable limit for q (default: SMALL)
Definition: qZeta.H:144
tmp< volScalarField > f2() const
Definition: qZeta.C:63
incompressible::RASModel::transportModel transportModel
Definition: eddyViscosity.H:75
GeometricField< scalar, fvPatchField, volMesh > volScalarField
Definition: volFieldsFwd.H:52
virtual tmp< volScalarField > epsilon() const
Return the turbulence kinetic energy dissipation rate.
Definition: qZeta.H:192
virtual const volScalarField & zeta() const
Definition: qZeta.H:202
qZeta(const geometricOneField &alpha, const geometricOneField &rho, const volVectorField &U, const surfaceScalarField &alphaRhoPhi, const surfaceScalarField &phi, const transportModel &transport, const word &propertiesName=turbulenceModel::propertiesName, const word &type=typeName)
Construct from components.
Definition: qZeta.C:80
A class representing the concept of a GeometricField of 1 used to avoid unnecessary manipulations for...
static const word propertiesName
Default name of the turbulence properties dictionary.
A class for handling words, derived from string.
Definition: word.H:59
Gibson and Dafa&#39;Alla&#39;s q-zeta two-equation low-Re turbulence model for incompressible flows...
Definition: qZeta.H:73
virtual tmp< volScalarField > k() const
Return the turbulence kinetic energy.
Definition: qZeta.H:186
virtual const volScalarField & q() const
Definition: qZeta.H:197
tmp< volScalarField > DqEff() const
Return the effective diffusivity for q.
Definition: qZeta.H:168
fileName::Type type(const fileName &)
Return the file type: DIRECTORY or FILE.
Definition: POSIX.C:461
tmp< volScalarField > DzetaEff() const
Return the effective diffusivity for epsilon.
Definition: qZeta.H:177
A class for managing temporary objects.
Definition: PtrList.H:54
virtual bool read()
Read RASProperties dictionary.
Definition: qZeta.C:217
const dimensionedScalar alpha
Fine-structure constant: default SI units: [].
dimensionedScalar sigmaZeta_
Definition: qZeta.H:87
dimensionedScalar qMin_
Lower limit of q.
Definition: qZeta.H:91
volScalarField & nu
dimensionedScalar zetaMin_
Lower limit of zeta.
Definition: qZeta.H:94
Namespace for OpenFOAM.