v2f.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) 2012-2022 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::RASModels::v2f
26 
27 Description
28  Lien and Kalitzin's v2-f turbulence model for incompressible and
29  compressible flows, with a limit imposed on the turbulent viscosity given
30  by Davidson et al.
31 
32  The model solves for turbulence kinetic energy k and turbulence dissipation
33  rate epsilon, with additional equations for the turbulence stress normal to
34  streamlines, v2, and elliptic damping function, f.
35 
36  The variant implemented employs N=6, such that f=0 on walls.
37 
38  Wall boundary conditions are:
39 
40  k = kLowReWallFunction
41  epsilon = epsilonWallFunction
42  v2 = v2WallFunction
43  f = fWallFunction
44 
45  These are applicable to both low- and high-Reynolds number flows.
46 
47  Inlet values can be approximated by:
48 
49  v2 = 2/3 k
50  f = zero-gradient
51 
52  References:
53  \verbatim
54  Lien, F. S., & Kalitzin, G. (2001).
55  Computations of transonic flow with the v2f turbulence model.
56  International Journal of Heat and Fluid Flow, 22(1), 53-61.
57 
58  Davidson, L., Nielsen, P., & Sveningsson, A. (2003).
59  Modifications of the v2-f model for computing the flow in a
60  3D wall jet.
61  Turbulence, Heat and Mass Transfer, 4, 577-584
62  \endverbatim
63 
64  The default model coefficients are
65  \verbatim
66  v2fCoeffs
67  {
68  Cmu 0.22;
69  CmuKEps 0.09;
70  C1 1.4;
71  C2 0.3;
72  CL 0.23;
73  Ceta 70;
74  Ceps2 1.9;
75  Ceps3 -0.33;
76  sigmaEps 1.3;
77  sigmaK 1;
78  }
79  \endverbatim
80 
81  Note:
82  If the kLowReWallFunction is employed, a velocity variant of the
83  turbulent viscosity wall function should be used, e.g. nutUWallFunction.
84  Turbulence k variants (nutk...) for this case will not behave correctly.
85 
86 See also
87  Foam::RASModels::v2fBase
88  Foam::RASModels::kEpsilon
89  Foam::kLowReWallFunctionFvPatchScalarField
90  Foam::epsilonWallFunctionFvPatchScalarField
91  Foam::v2WallFunctionFvPatchScalarField
92  Foam::fWallFunctionFvPatchScalarField
93 
94 SourceFiles
95  v2f.C
96 
97 \*---------------------------------------------------------------------------*/
98 
99 #ifndef v2f_H
100 #define v2f_H
101 
102 #include "v2fBase.H"
103 #include "RASModel.H"
104 #include "eddyViscosity.H"
105 
106 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
107 
108 namespace Foam
109 {
110 namespace RASModels
111 {
112 
113 /*---------------------------------------------------------------------------*\
114  Class v2f Declaration
115 \*---------------------------------------------------------------------------*/
116 
117 template<class BasicMomentumTransportModel>
118 class v2f
119 :
120  public eddyViscosity<RASModel<BasicMomentumTransportModel>>,
121  public v2fBase
122 {
123 
124 protected:
125 
126  // Protected data
127 
128  // Model coefficients
140 
141 
142  // Fields
143 
144  //- Turbulence kinetic energy
146 
147  //- Turbulence dissipation
149 
150  //- Turbulence stress normal to streamlines
152 
153  //- Damping function
155 
156 
157  // Bounding values
161 
162 
163  // Protected Member Functions
164 
165  virtual void correctNut();
166 
167  //- Return time scale, Ts
168  tmp<volScalarField> Ts() const;
169 
170  //- Return length scale, Ls
171  tmp<volScalarField> Ls() const;
172 
173 
174 public:
176  typedef typename BasicMomentumTransportModel::alphaField alphaField;
177  typedef typename BasicMomentumTransportModel::rhoField rhoField;
178 
179 
180  //- Runtime type information
181  TypeName("v2f");
182 
183 
184  // Constructors
185 
186  //- Construct from components
187  v2f
188  (
189  const alphaField& alpha,
190  const rhoField& rho,
191  const volVectorField& U,
192  const surfaceScalarField& alphaRhoPhi,
193  const surfaceScalarField& phi,
194  const viscosity& viscosity,
195  const word& type = typeName
196  );
197 
198 
199  //- Destructor
200  virtual ~v2f()
201  {}
202 
203 
204  // Member Functions
205 
206  //- Read RASProperties dictionary
207  virtual bool read();
208 
209  //- Return the effective diffusivity for k
210  tmp<volScalarField> DkEff() const
211  {
212  return volScalarField::New
213  (
214  "DkEff",
215  this->nut_/sigmaK_ + this->nu()
216  );
217  }
218 
219  //- Return the effective diffusivity for epsilon
221  {
222  return volScalarField::New
223  (
224  "DepsilonEff",
225  this->nut_/sigmaEps_ + this->nu()
226  );
227  }
228 
229  //- Return the turbulence kinetic energy
230  virtual tmp<volScalarField> k() const
231  {
232  return k_;
233  }
234 
235  //- Return the turbulence kinetic energy dissipation rate
236  virtual tmp<volScalarField> epsilon() const
237  {
238  return epsilon_;
239  }
240 
241  //- Return the turbulence specific dissipation rate
242  virtual tmp<volScalarField> omega() const
243  {
244  return volScalarField::New
245  (
246  "omega",
247  epsilon_/(Cmu_*k_),
248  epsilon_.boundaryField().types()
249  );
250  }
251 
252  //- Return turbulence stress normal to streamlines
253  virtual tmp<volScalarField> v2() const
254  {
255  return v2_;
256  }
257 
258  //- Return the damping function
259  virtual tmp<volScalarField> f() const
260  {
261  return f_;
262  }
263 
264  //- Solve the turbulence equations and correct the turbulence viscosity
265  virtual void correct();
266 };
267 
268 
269 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
270 
271 } // End namespace RASModels
272 } // End namespace Foam
273 
274 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
275 
276 #ifdef NoRepository
277  #include "v2f.C"
278 #endif
279 
280 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
281 
282 #endif
283 
284 // ************************************************************************* //
dimensionedScalar v2Min_
Definition: v2f.H:158
BasicMomentumTransportModel::rhoField rhoField
Definition: v2f.H:176
virtual tmp< volScalarField > k() const
Return the turbulence kinetic energy.
Definition: v2f.H:229
dimensionedScalar sigmaEps_
Definition: v2f.H:138
U
Definition: pEqn.H:72
dimensionedScalar Ceta_
Definition: v2f.H:134
const Boundary & boundaryField() const
Return const-reference to the boundary field.
volScalarField alpha(IOobject("alpha", runTime.timeName(), mesh, IOobject::READ_IF_PRESENT, IOobject::AUTO_WRITE), lambda *max(Ua &U, zeroSensitivity))
wordList types() const
Return a list of the patch field types.
dimensionedScalar C1_
Definition: v2f.H:131
static tmp< GeometricField< scalar, fvPatchField, volMesh > > New(const word &name, const Internal &, const PtrList< fvPatchField< scalar >> &)
Return a temporary field constructed from name,.
volScalarField k_
Turbulence kinetic energy.
Definition: v2f.H:144
dimensionedScalar Ceps3_
Definition: v2f.H:136
Eddy viscosity turbulence model base class.
Definition: eddyViscosity.H:49
virtual tmp< volScalarField > epsilon() const
Return the turbulence kinetic energy dissipation rate.
Definition: v2f.H:235
dimensionedScalar CL_
Definition: v2f.H:133
tmp< volScalarField > Ls() const
Return length scale, Ls.
Definition: v2f.C:48
v2f(const alphaField &alpha, const rhoField &rho, const volVectorField &U, const surfaceScalarField &alphaRhoPhi, const surfaceScalarField &phi, const viscosity &viscosity, const word &type=typeName)
Construct from components.
Definition: v2f.C:69
Abstract base-class for v2-f models to provide BCs access to the v2 and f fields. ...
Definition: v2fBase.H:54
volScalarField v2_
Turbulence stress normal to streamlines.
Definition: v2f.H:150
dimensionedScalar fMin_
Definition: v2f.H:159
BasicMomentumTransportModel::alphaField alphaField
Definition: v2f.H:175
A class for handling words, derived from string.
Definition: word.H:59
virtual ~v2f()
Destructor.
Definition: v2f.H:199
virtual tmp< volScalarField > v2() const
Return turbulence stress normal to streamlines.
Definition: v2f.H:252
virtual void correct()
Solve the turbulence equations and correct the turbulence viscosity.
Definition: v2f.C:273
dimensionedScalar CmuKEps_
Definition: v2f.H:130
phi
Definition: correctPhi.H:3
tmp< volScalarField > DkEff() const
Return the effective diffusivity for k.
Definition: v2f.H:209
virtual void correctNut()
Definition: v2f.C:57
Abstract base class for all fluid physical properties.
Definition: viscosity.H:49
volScalarField f_
Damping function.
Definition: v2f.H:153
TypeName("v2f")
Runtime type information.
dimensionedScalar Cmu_
Definition: v2f.H:129
virtual tmp< volScalarField > f() const
Return the damping function.
Definition: v2f.H:258
Lien and Kalitzin&#39;s v2-f turbulence model for incompressible and compressible flows, with a limit imposed on the turbulent viscosity given by Davidson et al.
Definition: v2f.H:117
virtual bool read()
Read RASProperties dictionary.
Definition: v2f.C:248
fileType type(const fileName &, const bool checkVariants=true, const bool followLink=true)
Return the file type: directory or file.
Definition: POSIX.C:488
tmp< volScalarField > DepsilonEff() const
Return the effective diffusivity for epsilon.
Definition: v2f.H:219
A class for managing temporary objects.
Definition: PtrList.H:53
dimensionedScalar C2_
Definition: v2f.H:132
dimensionedScalar Ceps2_
Definition: v2f.H:135
virtual tmp< volScalarField > nu() const
Return the laminar viscosity.
Definition: RASModel.H:208
volScalarField epsilon_
Turbulence dissipation.
Definition: v2f.H:147
virtual tmp< volScalarField > omega() const
Return the turbulence specific dissipation rate.
Definition: v2f.H:241
dimensionedScalar sigmaK_
Definition: v2f.H:137
Namespace for OpenFOAM.
tmp< volScalarField > Ts() const
Return time scale, Ts.
Definition: v2f.C:41