kOmegaSSTLM.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) 2016 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::kOmegaSSTLM
26 
27 Group
28  grpLESTurbulence
29 
30 Description
31  Langtry-Menter 4-equation transitional SST model
32  based on the k-omega-SST RAS model.
33 
34  References:
35  \verbatim
36  Langtry, R. B., & Menter, F. R. (2009).
37  Correlation-based transition modeling for unstructured parallelized
38  computational fluid dynamics codes.
39  AIAA journal, 47(12), 2894-2906.
40 
41  Menter, F. R., Langtry, R., & Volker, S. (2006).
42  Transition modelling for general purpose CFD codes.
43  Flow, turbulence and combustion, 77(1-4), 277-303.
44 
45  Langtry, R. B. (2006).
46  A correlation-based transition model using local variables for
47  unstructured parallelized CFD codes.
48  Phd. Thesis, Universität Stuttgart.
49  \endverbatim
50 
51  The model coefficients are
52  \verbatim
53  kOmegaSSTCoeffs
54  {
55  // Default SST coefficients
56  alphaK1 0.85;
57  alphaK2 1;
58  alphaOmega1 0.5;
59  alphaOmega2 0.856;
60  beta1 0.075;
61  beta2 0.0828;
62  betaStar 0.09;
63  gamma1 5/9;
64  gamma2 0.44;
65  a1 0.31;
66  b1 1;
67  c1 10;
68  F3 no;
69 
70  // Default LM coefficients
71  ca1 2;
72  ca2 0.06;
73  ce1 1;
74  ce2 50;
75  cThetat 0.03;
76  sigmaThetat 2;
77 
78  lambdaErr 1e-6;
79  maxLambdaIter 10;
80  }
81  \endverbatim
82 
83 SourceFiles
84  kOmegaSSTLM.C
85 
86 \*---------------------------------------------------------------------------*/
87 
88 #ifndef kOmegaSSTLM_H
89 #define kOmegaSSTLM_H
90 
91 #include "kOmegaSST.H"
92 
93 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
94 
95 namespace Foam
96 {
97 namespace RASModels
98 {
99 
100 /*---------------------------------------------------------------------------*\
101  Class kOmegaSSTLM Declaration
102 \*---------------------------------------------------------------------------*/
103 
104 template<class BasicTurbulenceModel>
105 class kOmegaSSTLM
106 :
107  public kOmegaSST<BasicTurbulenceModel>
108 {
109  // Private Member Functions
110 
111  // Disallow default bitwise copy construct and assignment
112  kOmegaSSTLM(const kOmegaSSTLM&);
113  void operator=(const kOmegaSSTLM&);
114 
115 
116 protected:
117 
118  // Protected data
119 
120  // Model constants
130 
131  //- Convergence criterion for the lambda/thetat loop
132  scalar lambdaErr_;
133 
134  //- Maximum number of iterations to converge the lambda/thetat loop
136 
137  //- Stabilization for division by the magnitude of the velocity
139 
140 
141  // Fields
142 
143  //- Transition onset momentum-thickness Reynolds number
145 
146  //- Intermittency
148 
149  //- Effective intermittency
151 
152 
153  // Protected Member Functions
154 
155  //- Modified form of the k-omega SST F1 function
156  virtual tmp<volScalarField> F1(const volScalarField& CDkOmega) const;
157 
158  //- Modified form of the k-omega SST k production rate
160  (
162  ) const;
163 
164  //- Modified form of the k-omega SST epsilon/k
166  (
169  ) const;
170 
171  //- Freestream blending-function
173  (
174  const volScalarField::Internal& Us,
175  const volScalarField::Internal& Omega,
177  ) const;
178 
179  //- Empirical correlation for critical Reynolds number where the
180  // intermittency first starts to increase in the boundary layer
182 
183  //- Empirical correlation that controls the length of the
184  // transition region
186  (
188  ) const;
189 
190  //- Transition onset location control function
192  (
193  const volScalarField::Internal& Rev,
195  const volScalarField::Internal& RT
196  ) const;
197 
198  //- Return the transition onset momentum-thickness Reynolds number
199  // (based on freestream conditions)
201  (
202  const volScalarField::Internal& Us,
203  const volScalarField::Internal& dUsds,
205  ) const;
206 
207  //- Solve the turbulence equations and correct the turbulence viscosity
209 
210 
211 public:
213  typedef typename BasicTurbulenceModel::alphaField alphaField;
214  typedef typename BasicTurbulenceModel::rhoField rhoField;
215  typedef typename BasicTurbulenceModel::transportModel transportModel;
216 
217 
218  //- Runtime type information
219  TypeName("kOmegaSSTLM");
220 
221 
222  // Constructors
223 
224  //- Construct from components
226  (
227  const alphaField& alpha,
228  const rhoField& rho,
229  const volVectorField& U,
230  const surfaceScalarField& alphaRhoPhi,
231  const surfaceScalarField& phi,
232  const transportModel& transport,
233  const word& propertiesName = turbulenceModel::propertiesName,
234  const word& type = typeName
235  );
236 
237 
238  //- Destructor
239  virtual ~kOmegaSSTLM()
240  {}
241 
242 
243  // Member Functions
244 
245  //- Re-read model coefficients if they have changed
246  virtual bool read();
247 
248  //- Access function transition onset momentum-thickness Reynolds number
249  const volScalarField& ReThetat() const
250  {
251  return ReThetat_;
252  }
253 
254  //- Access function to intermittency
255  const volScalarField& gammaInt() const
256  {
257  return gammaInt_;
258  }
259 
260  //- Return the effective diffusivity for transition onset
261  // momentum-thickness Reynolds number
263  {
264  return tmp<volScalarField>
265  (
266  new volScalarField
267  (
268  "DReThetatEff",
269  sigmaThetat_*(this->nut_ + this->nu())
270  )
271  );
272  }
273 
274  //- Return the effective diffusivity for intermittency
276  {
277  return tmp<volScalarField>
278  (
279  new volScalarField
280  (
281  "DgammaIntEff",
282  this->nut_ + this->nu()
283  )
284  );
285  }
286 
287  //- Solve the turbulence equations and correct the turbulence viscosity
288  virtual void correct();
289 };
290 
291 
292 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
293 
294 } // End namespace RASModels
295 } // End namespace Foam
296 
297 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
298 
299 #ifdef NoRepository
300  #include "kOmegaSSTLM.C"
301 #endif
302 
303 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
304 
305 #endif
306 
307 // ************************************************************************* //
surfaceScalarField & phi
volScalarField::Internal gammaIntEff_
Effective intermittency.
Definition: kOmegaSSTLM.H:149
U
Definition: pEqn.H:83
intWM_LABEL_SIZE_t label
A label is an int32_t or int64_t as specified by the pre-processor macro WM_LABEL_SIZE.
Definition: label.H:59
tmp< volScalarField::Internal > ReThetat0(const volScalarField::Internal &Us, const volScalarField::Internal &dUsds, const volScalarField::Internal &nu) const
Return the transition onset momentum-thickness Reynolds number.
Definition: kOmegaSSTLM.C:213
dimensionedScalar ca2_
Definition: kOmegaSSTLM.H:122
const dimensionedScalar G
Newtonian constant of gravitation.
tmp< volScalarField > DReThetatEff() const
Return the effective diffusivity for transition onset.
Definition: kOmegaSSTLM.H:261
tmp< volScalarField > DgammaIntEff() const
Return the effective diffusivity for intermittency.
Definition: kOmegaSSTLM.H:274
virtual tmp< volScalarField::Internal > epsilonByk(const volScalarField::Internal &F1, const volScalarField::Internal &F2) const
Modified form of the k-omega SST epsilon/k.
Definition: kOmegaSSTLM.C:62
dimensionedScalar ca1_
Definition: kOmegaSSTLM.H:121
const volScalarField & gammaInt() const
Access function to intermittency.
Definition: kOmegaSSTLM.H:254
virtual ~kOmegaSSTLM()
Destructor.
Definition: kOmegaSSTLM.H:238
volScalarField ReThetat_
Transition onset momentum-thickness Reynolds number.
Definition: kOmegaSSTLM.H:143
GeometricField< scalar, fvPatchField, volMesh > volScalarField
Definition: volFieldsFwd.H:52
tmp< volScalarField::Internal > ReThetac() const
Empirical correlation for critical Reynolds number where the.
Definition: kOmegaSSTLM.C:109
const dimensionedScalar deltaU_
Stabilization for division by the magnitude of the velocity.
Definition: kOmegaSSTLM.H:137
const volScalarField & ReThetat() const
Access function transition onset momentum-thickness Reynolds number.
Definition: kOmegaSSTLM.H:248
dimensionedScalar cThetat_
Definition: kOmegaSSTLM.H:127
BasicTurbulenceModel::alphaField alphaField
Definition: kOmegaSSTLM.H:212
static const word propertiesName
Default name of the turbulence properties dictionary.
A class for handling words, derived from string.
Definition: word.H:59
virtual void correct()
Solve the turbulence equations and correct the turbulence viscosity.
Definition: kOmegaSSTLM.C:607
TypeName("kOmegaSSTLM")
Runtime type information.
virtual bool read()
Re-read model coefficients if they have changed.
Definition: kOmegaSSTLM.C:486
dimensionedScalar sigmaThetat_
Definition: kOmegaSSTLM.H:128
tmp< volScalarField::Internal > Fonset(const volScalarField::Internal &Rev, const volScalarField::Internal &ReThetac, const volScalarField::Internal &RT) const
Transition onset location control function.
Definition: kOmegaSSTLM.C:327
tmp< volScalarField::Internal > Flength(const volScalarField::Internal &nu) const
Empirical correlation that controls the length of the.
Definition: kOmegaSSTLM.C:150
fileName::Type type(const fileName &, const bool followLink=true)
Return the file type: DIRECTORY or FILE.
Definition: POSIX.C:485
dimensionedScalar ce2_
Definition: kOmegaSSTLM.H:125
dimensionedScalar ce1_
Definition: kOmegaSSTLM.H:124
BasicTurbulenceModel::rhoField rhoField
Definition: kOmegaSSTLM.H:213
label maxLambdaIter_
Maximum number of iterations to converge the lambda/thetat loop.
Definition: kOmegaSSTLM.H:134
Langtry-Menter 4-equation transitional SST model based on the k-omega-SST RAS model.
Definition: kOmegaSSTLM.H:104
Field with dimensions and associated with geometry type GeoMesh which is used to size the field and a...
scalar lambdaErr_
Convergence criterion for the lambda/thetat loop.
Definition: kOmegaSSTLM.H:131
void correctReThetatGammaInt()
Solve the turbulence equations and correct the turbulence viscosity.
Definition: kOmegaSSTLM.C:509
volScalarField gammaInt_
Intermittency.
Definition: kOmegaSSTLM.H:146
BasicTurbulenceModel::transportModel transportModel
Definition: kOmegaSSTLM.H:214
tmp< volScalarField::Internal > Fthetat(const volScalarField::Internal &Us, const volScalarField::Internal &Omega, const volScalarField::Internal &nu) const
Freestream blending-function.
Definition: kOmegaSSTLM.C:75
A class for managing temporary objects.
Definition: PtrList.H:53
const dimensionedScalar alpha
Fine-structure constant: default SI units: [].
virtual tmp< volScalarField > F1(const volScalarField &CDkOmega) const
Modified form of the k-omega SST F1 function.
Definition: kOmegaSSTLM.C:39
volScalarField & nu
Namespace for OpenFOAM.
virtual tmp< volScalarField::Internal > Pk(const volScalarField::Internal &G) const
Modified form of the k-omega SST k production rate.
Definition: kOmegaSSTLM.C:52