kOmegaSSTSAS.C
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) 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 \*---------------------------------------------------------------------------*/
25 
26 #include "kOmegaSSTSAS.H"
27 
28 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
29 
30 namespace Foam
31 {
32 namespace RASModels
33 {
34 
35 // * * * * * * * * * * * * Protected Member Functions * * * * * * * * * * * //
36 
37 template<class BasicTurbulenceModel>
39 (
40  const volScalarField& S2,
41  const volScalarField& gamma,
42  const volScalarField& beta
43 ) const
44 {
46  (
47  sqrt(this->k_)/(pow025(this->betaStar_)*this->omega_)
48  );
49 
50  volScalarField Lvk
51  (
52  max
53  (
54  kappa_*sqrt(S2)
55  /(
56  mag(fvc::laplacian(this->U_))
58  (
59  "ROOTVSMALL",
60  dimensionSet(0, -1, -1, 0, 0),
61  ROOTVSMALL
62  )
63  ),
64  Cs_*sqrt(kappa_*zeta2_/(beta/this->betaStar_ - gamma))*delta()
65  )
66  );
67 
68  return fvm::Su
69  (
70  this->alpha_*this->rho_
71  *min
72  (
73  max
74  (
75  zeta2_*kappa_*S2*sqr(L/Lvk)
76  - (2*C_/sigmaPhi_)*this->k_
77  *max
78  (
79  magSqr(fvc::grad(this->omega_))/sqr(this->omega_),
80  magSqr(fvc::grad(this->k_))/sqr(this->k_)
81  ),
82  dimensionedScalar("0", dimensionSet(0, 0, -2, 0, 0), 0)
83  ),
84  // Limit SAS production of omega for numerical stability,
85  // particularly during start-up
86  this->omega_/(0.1*this->omega_.time().deltaT())
87  ),
88  this->omega_
89  );
90 }
91 
92 
93 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
94 
95 template<class BasicTurbulenceModel>
97 (
98  const alphaField& alpha,
99  const rhoField& rho,
100  const volVectorField& U,
101  const surfaceScalarField& alphaRhoPhi,
102  const surfaceScalarField& phi,
103  const transportModel& transport,
104  const word& propertiesName,
105  const word& type
106 )
107 :
109  (
110  alpha,
111  rho,
112  U,
113  alphaRhoPhi,
114  phi,
115  transport,
116  propertiesName
117  ),
118 
119  Cs_
120  (
122  (
123  "Cs",
124  this->coeffDict_,
125  0.11
126  )
127  ),
128  kappa_
129  (
131  (
132  "kappa",
133  this->coeffDict_,
134  0.41
135  )
136  ),
137  zeta2_
138  (
140  (
141  "zeta2",
142  this->coeffDict_,
143  3.51
144  )
145  ),
146  sigmaPhi_
147  (
149  (
150  "sigmaPhi",
151  this->coeffDict_,
152  2.0/3.0
153  )
154  ),
155  C_
156  (
158  (
159  "C",
160  this->coeffDict_,
161  2
162  )
163  ),
164 
165  delta_
166  (
168  (
169  IOobject::groupName("delta", U.group()),
170  *this,
171  this->coeffDict_
172  )
173  )
174 {}
175 
176 
177 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
178 
179 template<class BasicTurbulenceModel>
181 {
183  {
184  Cs_.readIfPresent(this->coeffDict());
185  kappa_.readIfPresent(this->coeffDict());
186  sigmaPhi_.readIfPresent(this->coeffDict());
187  zeta2_.readIfPresent(this->coeffDict());
188  C_.readIfPresent(this->coeffDict());
189 
190  return true;
191  }
192  else
193  {
194  return false;
195  }
196 }
197 
198 
199 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
200 
201 } // End namespace RASModels
202 } // End namespace Foam
203 
204 // ************************************************************************* //
scalar delta
surfaceScalarField & phi
tmp< GeometricField< typename outerProduct< vector, Type >::type, fvPatchField, volMesh >> grad(const GeometricField< Type, fvsPatchField, surfaceMesh > &ssf)
Definition: fvcGrad.C:52
BasicTurbulenceModel::transportModel transportModel
Definition: kOmegaSST.H:69
U
Definition: pEqn.H:83
BasicTurbulenceModel::rhoField rhoField
Definition: kOmegaSST.H:68
dimensioned< Type > max(const dimensioned< Type > &, const dimensioned< Type > &)
tmp< fvMatrix< Type > > Su(const DimensionedField< Type, volMesh > &, const GeometricField< Type, fvPatchField, volMesh > &)
dimensionedSymmTensor sqr(const dimensionedVector &dv)
dimensionedScalar sqrt(const dimensionedScalar &ds)
dimensionedScalar pow025(const dimensionedScalar &ds)
tmp< GeometricField< Type, fvPatchField, volMesh > > laplacian(const GeometricField< Type, fvPatchField, volMesh > &vf, const word &name)
Definition: fvcLaplacian.C:45
Dimension set for the base types.
Definition: dimensionSet.H:118
A class for handling words, derived from string.
Definition: word.H:59
static word groupName(Name name, const word &group)
static autoPtr< LESdelta > New(const word &name, const turbulenceModel &turbulence, const dictionary &dict)
Return a reference to the selected LES delta.
Definition: LESdelta.C:66
dimensioned< scalar > magSqr(const dimensioned< Type > &)
dimensioned< Type > min(const dimensioned< Type > &, const dimensioned< Type > &)
virtual bool read()
Re-read model coefficients if they have changed.
Definition: kOmegaSSTSAS.C:180
dimensioned< scalar > dimensionedScalar
Dimensioned scalar obtained from generic dimensioned type.
dimensioned< scalar > mag(const dimensioned< Type > &)
virtual tmp< fvScalarMatrix > Qsas(const volScalarField &S2, const volScalarField &gamma, const volScalarField &beta) const
SAS omega source.
Definition: kOmegaSSTSAS.C:39
const dimensionedScalar alpha
Fine-structure constant: default SI units: [].
BasicTurbulenceModel::alphaField alphaField
Definition: kOmegaSST.H:67
Namespace for OpenFOAM.
Scale-adaptive URAS model based on the k-omega-SST RAS model.
Definition: kOmegaSSTSAS.H:98