qZeta.H
Go to the documentation of this file.
1 /*---------------------------------------------------------------------------*\
2  ========= |
3  \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
4  \\ / O peration | Website: https://openfoam.org
5  \\ / A nd | Copyright (C) 2011-2018 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 Description
28  Gibson and Dafa'Alla's q-zeta two-equation low-Re turbulence model
29  for incompressible flows
30 
31  This turbulence model is described in:
32  \verbatim
33  Dafa'Alla, A.A., Juntasaro, E. & Gibson, M.M. (1996).
34  Calculation of oscillating boundary layers with the
35  q-zeta turbulence model.
36  Engineering Turbulence Modelling and Experiments 3:
37  Proceedings of the Third International Symposium,
38  Crete, Greece, May 27-29, 141.
39  \endverbatim
40  which is a development of the original q-zeta model described in:
41  \verbatim
42  Gibson, M. M., & Dafa'Alla, A. A. (1995).
43  Two-equation model for turbulent wall flow.
44  AIAA journal, 33(8), 1514-1518.
45  \endverbatim
46 
47 SourceFiles
48  qZeta.C
49 
50 \*---------------------------------------------------------------------------*/
51 
52 #ifndef qZeta_H
53 #define qZeta_H
54 
56 #include "eddyViscosity.H"
57 
58 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
59 
60 namespace Foam
61 {
62 namespace incompressible
63 {
64 namespace RASModels
65 {
66 
67 /*---------------------------------------------------------------------------*\
68  Class qZeta Declaration
69 \*---------------------------------------------------------------------------*/
70 
71 class qZeta
72 :
73  public eddyViscosity<incompressible::RASModel>
74 {
75 
76 protected:
77 
78  // Protected data
79 
80  // Model coefficients
81 
87 
88  //- Lower limit of q
90 
91  //- Lower limit of zeta
93 
94  // Fields
95 
98 
101 
102 
103  // Protected Member Functions
104 
105  tmp<volScalarField> fMu() const;
106  tmp<volScalarField> f2() const;
107  virtual void correctNut();
108 
109 
110 public:
111 
112  //- Runtime type information
113  TypeName("qZeta");
114 
115  // Constructors
116 
117  //- Construct from components
118  qZeta
119  (
120  const geometricOneField& alpha,
121  const geometricOneField& rho,
122  const volVectorField& U,
123  const surfaceScalarField& alphaRhoPhi,
124  const surfaceScalarField& phi,
125  const transportModel& transport,
126  const word& propertiesName = turbulenceModel::propertiesName,
127  const word& type = typeName
128  );
129 
130 
131  //- Destructor
132  virtual ~qZeta()
133  {}
134 
135 
136  // Member Functions
137 
138  //- Read RASProperties dictionary
139  virtual bool read();
140 
141  //- Return the lower allowable limit for q (default: small)
142  const dimensionedScalar& qMin() const
143  {
144  return qMin_;
145  }
146 
147  //- Return the lower allowable limit for zeta (default: small)
148  const dimensionedScalar& zetaMin() const
149  {
150  return zetaMin_;
151  }
152 
153  //- Allow qMin to be changed
155  {
156  return qMin_;
157  }
158 
159  //- Allow zetaMin to be changed
161  {
162  return zetaMin_;
163  }
164 
165  //- Return the effective diffusivity for q
166  tmp<volScalarField> DqEff() const
167  {
168  return tmp<volScalarField>
169  (
170  new volScalarField("DqEff", nut_ + nu())
171  );
172  }
173 
174  //- Return the effective diffusivity for epsilon
176  {
177  return tmp<volScalarField>
178  (
179  new volScalarField("DzetaEff", nut_/sigmaZeta_ + nu())
180  );
181  }
182 
183  //- Return the turbulence kinetic energy
184  virtual tmp<volScalarField> k() const
185  {
186  return k_;
187  }
188 
189  //- Return the turbulence kinetic energy dissipation rate
190  virtual tmp<volScalarField> epsilon() const
191  {
192  return epsilon_;
193  }
195  virtual const volScalarField& q() const
196  {
197  return q_;
198  }
200  virtual const volScalarField& zeta() const
201  {
202  return zeta_;
203  }
204 
205  //- Solve the turbulence equations and correct the turbulence viscosity
206  virtual void correct();
207 };
208 
209 
210 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
211 
212 } // End namespace RASModels
213 } // End namespace incompressible
214 } // End namespace Foam
215 
216 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
217 
218 #endif
219 
220 // ************************************************************************* //
const dimensionedScalar & qMin() const
Return the lower allowable limit for q (default: small)
Definition: qZeta.H:141
surfaceScalarField & phi
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
TypeName("qZeta")
Runtime type information.
virtual ~qZeta()
Destructor.
Definition: qZeta.H:131
Eddy viscosity turbulence model base class.
Definition: eddyViscosity.H:49
tmp< volScalarField > DqEff() const
Return the effective diffusivity for q.
Definition: qZeta.H:165
virtual tmp< volScalarField > epsilon() const
Return the turbulence kinetic energy dissipation rate.
Definition: qZeta.H:189
tmp< volScalarField > DzetaEff() const
Return the effective diffusivity for epsilon.
Definition: qZeta.H:174
virtual const volScalarField & q() const
Definition: qZeta.H:194
GeometricField< scalar, fvPatchField, volMesh > volScalarField
Definition: volFieldsFwd.H:52
tmp< volScalarField > f2() const
Definition: qZeta.C:63
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:70
virtual tmp< volScalarField > k() const
Return the turbulence kinetic energy.
Definition: qZeta.H:183
virtual const volScalarField & zeta() const
Definition: qZeta.H:199
fileName::Type type(const fileName &, const bool followLink=true)
Return the file type: DIRECTORY or FILE.
Definition: POSIX.C:481
U
Definition: pEqn.H:72
tmp< volScalarField > fMu() const
Definition: qZeta.C:46
Base-class for all transport models used by the incompressible turbulence models. ...
A class for managing temporary objects.
Definition: PtrList.H:53
const dimensionedScalar & zetaMin() const
Return the lower allowable limit for zeta (default: small)
Definition: qZeta.H:147
virtual bool read()
Read RASProperties dictionary.
Definition: qZeta.C:217
const dimensionedScalar alpha
Fine-structure constant: default SI units: [].
dimensionedScalar sigmaZeta_
Definition: qZeta.H:84
dimensionedScalar qMin_
Lower limit of q.
Definition: qZeta.H:88
volScalarField & nu
dimensionedScalar zetaMin_
Lower limit of zeta.
Definition: qZeta.H:91
Namespace for OpenFOAM.