Base class implementation of the k-omega-SST turbulence model for incompressible and compressible flows. More...
#include <kOmegaSSTBase.H>


Public Types | |
| typedef BasicEddyViscosityModel::alphaField | alphaField |
| typedef BasicEddyViscosityModel::rhoField | rhoField |
| typedef BasicEddyViscosityModel::transportModel | transportModel |
Public Member Functions | |
| kOmegaSSTBase (const kOmegaSSTBase &)=delete | |
| No copy construct. | |
| void | operator= (const kOmegaSSTBase &)=delete |
| No copy assignment. | |
| kOmegaSSTBase (const word &type, const alphaField &alpha, const rhoField &rho, const volVectorField &U, const surfaceScalarField &alphaRhoPhi, const surfaceScalarField &phi, const transportModel &transport, const word &propertiesName=turbulenceModel::propertiesName) | |
| Construct from components. | |
| virtual | ~kOmegaSSTBase ()=default |
| Destructor. | |
| virtual bool | read () |
| Re-read model coefficients if they have changed. | |
| tmp< volScalarField > | DkEff (const volScalarField &F1) const |
| Return the effective diffusivity for k. | |
| tmp< volScalarField > | DomegaEff (const volScalarField &F1) const |
| Return the effective diffusivity for omega. | |
| virtual tmp< volScalarField > | k () const |
| Return the turbulence kinetic energy. | |
| virtual tmp< volScalarField > | omega () const |
| Return the turbulence kinetic energy dissipation rate. | |
| virtual void | correct () |
| Solve the turbulence equations and correct the turbulence viscosity. | |
Protected Attributes | |
| dimensionedScalar | alphaK1_ |
| dimensionedScalar | alphaK2_ |
| dimensionedScalar | alphaOmega1_ |
| dimensionedScalar | alphaOmega2_ |
| dimensionedScalar | gamma1_ |
| dimensionedScalar | gamma2_ |
| dimensionedScalar | beta1_ |
| dimensionedScalar | beta2_ |
| dimensionedScalar | betaStar_ |
| dimensionedScalar | a1_ |
| dimensionedScalar | b1_ |
| dimensionedScalar | c1_ |
| Switch | F3_ |
| Flag to include the F3 term. | |
| const volScalarField & | y_ |
| Wall distance. | |
| volScalarField | k_ |
| Turbulent kinetic energy field [m^2/s^2]. | |
| volScalarField | omega_ |
| Specific dissipation rate field [1/s]. | |
| Switch | decayControl_ |
| Flag to include the decay control. | |
| dimensionedScalar | kInf_ |
| dimensionedScalar | omegaInf_ |
Base class implementation of the k-omega-SST turbulence model for incompressible and compressible flows.
Turbulence model described in:
Menter, F. R. & Esch, T. (2001).
Elements of Industrial Heat Transfer Prediction.
16th Brazilian Congress of Mechanical Engineering (COBEM).
with updated coefficients from
Menter, F. R., Kuntz, M., and Langtry, R. (2003).
Ten Years of Industrial Experience with the SST Turbulence Model.
Turbulence, Heat and Mass Transfer 4, ed: K. Hanjalic, Y. Nagano,
& M. Tummers, Begell House, Inc., 625 - 632.
but with the consistent production terms from the 2001 paper as form in the 2003 paper is a typo, see
http://turbmodels.larc.nasa.gov/sst.html
and the addition of the optional F3 term for rough walls from
Hellsten, A. (1998).
Some Improvements in Menter's k-omega-SST turbulence model
29th AIAA Fluid Dynamics Conference, AIAA-98-2554.
and the optional decay control from:
Spalart, P. R. and Rumsey, C. L. (2007).
Effective Inflow Conditions for Turbulence Models in Aerodynamic
Calculations
AIAA Journal, 45(10), 2544 - 2553.
Note that this implementation is written in terms of alpha diffusion coefficients rather than the more traditional sigma (alpha = 1/sigma) so that the blending can be applied to all coefficients in a consistent manner. The paper suggests that sigma is blended but this would not be consistent with the blending of the k-epsilon and k-omega models.
Also note that the error in the last term of equation (2) relating to sigma has been corrected.
Wall-functions are applied in this implementation by using equations (14) to specify the near-wall omega as appropriate.
The blending functions (15) and (16) are not currently used because of the uncertainty in their origin, range of applicability and that if y+ becomes sufficiently small blending u_tau in this manner clearly becomes nonsense.
The default model coefficients are
kOmegaSSTBaseCoeffs
{
alphaK1 0.85;
alphaK2 1.0;
alphaOmega1 0.5;
alphaOmega2 0.856;
beta1 0.075;
beta2 0.0828;
betaStar 0.09;
gamma1 5/9;
gamma2 0.44;
a1 0.31;
b1 1.0;
c1 10.0;
F3 no;
// Optional decay control
decayControl yes;
kInf \<far-field k value\>;
omegaInf \<far-field omega value\>;
}
Definition at line 125 of file kOmegaSSTBase.H.
| typedef BasicEddyViscosityModel::alphaField alphaField |
Definition at line 324 of file kOmegaSSTBase.H.
| typedef BasicEddyViscosityModel::rhoField rhoField |
Definition at line 325 of file kOmegaSSTBase.H.
| typedef BasicEddyViscosityModel::transportModel transportModel |
Definition at line 326 of file kOmegaSSTBase.H.
|
delete |
No copy construct.
References kOmegaSSTBase().
Referenced by kOmegaSSTBase(), and operator=().


| kOmegaSSTBase | ( | const word & | type, |
| const alphaField & | alpha, | ||
| const rhoField & | rho, | ||
| const volVectorField & | U, | ||
| const surfaceScalarField & | alphaRhoPhi, | ||
| const surfaceScalarField & | phi, | ||
| const transportModel & | transport, | ||
| const word & | propertiesName = turbulenceModel::propertiesName ) |
|
virtualdefault |
Destructor.
|
protected |
Set decay control with kInf and omegaInf.
Definition at line 441 of file kOmegaSSTBase.C.
References decayControl_, dict, Foam::endl(), Foam::Info, kInf_, and omegaInf_.
Referenced by read().


|
protectedvirtual |
Reimplemented in kOmegaSSTLM< BasicTurbulenceModel >, kOmegaSSTLM< EddyDiffusivity< fluidThermoCompressibleTurbulenceModel > >, and kOmegaSSTLM< Foam::EddyDiffusivity< Foam::fluidThermoCompressibleTurbulenceModel > >.
Definition at line 36 of file kOmegaSSTBase.C.
References alphaOmega2_, betaStar_, Foam::dimless, Foam::dimTime, k_, Foam::max(), Foam::min(), mu, omega_, Foam::pow4(), Foam::sqr(), Foam::sqrt(), Foam::tanh(), and y_.
Referenced by alphaK(), alphaOmega(), beta(), blend(), blend(), kOmegaSSTDES< Foam::EddyDiffusivity< Foam::fluidThermoCompressibleTurbulenceModel > >::CDES(), correct(), DkEff(), DomegaEff(), kOmegaSSTDES< Foam::EddyDiffusivity< Foam::fluidThermoCompressibleTurbulenceModel > >::epsilonByk(), and gamma().


|
protectedvirtual |
Definition at line 66 of file kOmegaSSTBase.C.
References betaStar_, k_, Foam::max(), Foam::min(), mu, omega_, Foam::sqr(), Foam::sqrt(), Foam::tanh(), and y_.
Referenced by GbyNu(), kOmegaSSTDDES< Foam::EddyDiffusivity< Foam::fluidThermoCompressibleTurbulenceModel > >::GbyNu(), and kOmegaSSTDES< Foam::EddyDiffusivity< Foam::fluidThermoCompressibleTurbulenceModel > >::GbyNu().


|
protectedvirtual |
Definition at line 83 of file kOmegaSSTBase.C.
References Foam::min(), mu, omega_, Foam::pow4(), Foam::sqr(), Foam::tanh(), and y_.

|
protectedvirtual |
Definition at line 96 of file kOmegaSSTBase.C.
References F2, F3, F3_, and tmp< T >::ref().
Referenced by correct(), and correctNut().


|
inlineprotected |
Return the blended field.
Definition at line 205 of file kOmegaSSTBase.H.
References F1(), psi1, and psi2.
Referenced by alphaK(), alphaOmega(), beta(), and gamma().


|
inlineprotected |
Return the internal blended field.
Definition at line 218 of file kOmegaSSTBase.H.
References F1(), psi1, and psi2.

|
inlineprotected |
|
inlineprotected |
Definition at line 233 of file kOmegaSSTBase.H.
References alphaOmega1_, alphaOmega2_, blend(), and F1().
Referenced by DomegaEff().


|
inlineprotected |
Definition at line 238 of file kOmegaSSTBase.H.
References beta1_, beta2_, blend(), F1(), tmp< T >::New(), IOobject::scopedName(), and Foam::type().
Referenced by correct(), Qsas(), and kOmegaSSTSAS< Foam::EddyDiffusivity< Foam::fluidThermoCompressibleTurbulenceModel > >::Qsas().


|
inlineprotected |
Definition at line 250 of file kOmegaSSTBase.H.
References blend(), F1(), gamma1_, gamma2_, tmp< T >::New(), IOobject::scopedName(), and Foam::type().
Referenced by correct(), Qsas(), and kOmegaSSTSAS< Foam::EddyDiffusivity< Foam::fluidThermoCompressibleTurbulenceModel > >::Qsas().


|
protectedvirtual |
Reimplemented in kOmegaSSTDES< BasicTurbulenceModel >, kOmegaSSTDES< EddyDiffusivity< fluidThermoCompressibleTurbulenceModel > >, kOmegaSSTDES< Foam::EddyDiffusivity< Foam::fluidThermoCompressibleTurbulenceModel > >, kOmegaSST< BasicTurbulenceModel >, kOmegaSST< EddyDiffusivity< fluidThermoCompressibleTurbulenceModel > >, kOmegaSST< Foam::EddyDiffusivity< Foam::fluidThermoCompressibleTurbulenceModel > >, and kOmegaSSTSato< BasicTurbulenceModel >.
Definition at line 110 of file kOmegaSSTBase.C.
References a1_, b1_, optionList::correct(), F23(), k_, Foam::max(), options::New(), omega_, S2(), and Foam::sqrt().
Referenced by correct(), and correctNut().


|
protectedvirtual |
Reimplemented in kOmegaSSTDES< EddyDiffusivity< fluidThermoCompressibleTurbulenceModel > >, kOmegaSSTDES< Foam::EddyDiffusivity< Foam::fluidThermoCompressibleTurbulenceModel > >, kOmegaSST< EddyDiffusivity< fluidThermoCompressibleTurbulenceModel > >, and kOmegaSST< Foam::EddyDiffusivity< Foam::fluidThermoCompressibleTurbulenceModel > >.
Definition at line 123 of file kOmegaSSTBase.C.
References correctNut(), Foam::fvc::grad(), Foam::magSqr(), and Foam::symm().

|
protectedvirtual |
Return square of strain rate.
Reimplemented in kOmegaSSTDDES< BasicTurbulenceModel >, kOmegaSSTDDES< EddyDiffusivity< fluidThermoCompressibleTurbulenceModel > >, kOmegaSSTDDES< Foam::EddyDiffusivity< Foam::fluidThermoCompressibleTurbulenceModel > >, kOmegaSSTDES< BasicTurbulenceModel >, kOmegaSSTDES< EddyDiffusivity< fluidThermoCompressibleTurbulenceModel > >, and kOmegaSSTDES< Foam::EddyDiffusivity< Foam::fluidThermoCompressibleTurbulenceModel > >.
Definition at line 130 of file kOmegaSSTBase.C.
References Foam::magSqr(), and Foam::symm().
Referenced by correct(), correctNut(), kOmegaSST< Foam::EddyDiffusivity< Foam::fluidThermoCompressibleTurbulenceModel > >::correctNut(), GbyNu(), Qsas(), and kOmegaSSTSAS< Foam::EddyDiffusivity< Foam::fluidThermoCompressibleTurbulenceModel > >::Qsas().


|
protectedvirtual |
Return k production rate.
Reimplemented in kOmegaSSTLM< BasicTurbulenceModel >, kOmegaSSTLM< EddyDiffusivity< fluidThermoCompressibleTurbulenceModel > >, and kOmegaSSTLM< Foam::EddyDiffusivity< Foam::fluidThermoCompressibleTurbulenceModel > >.
Definition at line 140 of file kOmegaSSTBase.C.
References betaStar_, c1_, k_, Foam::min(), and omega_.
Referenced by correct().


|
protectedvirtual |
Return epsilon/k which for standard RAS is betaStar*omega.
Reimplemented in kOmegaSSTDES< BasicTurbulenceModel >, kOmegaSSTDES< EddyDiffusivity< fluidThermoCompressibleTurbulenceModel > >, kOmegaSSTDES< Foam::EddyDiffusivity< Foam::fluidThermoCompressibleTurbulenceModel > >, kOmegaSSTLM< BasicTurbulenceModel >, kOmegaSSTLM< EddyDiffusivity< fluidThermoCompressibleTurbulenceModel > >, and kOmegaSSTLM< Foam::EddyDiffusivity< Foam::fluidThermoCompressibleTurbulenceModel > >.
Definition at line 150 of file kOmegaSSTBase.C.
References betaStar_, and omega_.
Referenced by correct().

|
protectedvirtual |
Return (G/nu)_0.
Reimplemented in kOmegaSSTDDES< BasicTurbulenceModel >, kOmegaSSTDDES< EddyDiffusivity< fluidThermoCompressibleTurbulenceModel > >, kOmegaSSTDDES< Foam::EddyDiffusivity< Foam::fluidThermoCompressibleTurbulenceModel > >, kOmegaSSTDES< BasicTurbulenceModel >, kOmegaSSTDES< EddyDiffusivity< fluidThermoCompressibleTurbulenceModel > >, and kOmegaSSTDES< Foam::EddyDiffusivity< Foam::fluidThermoCompressibleTurbulenceModel > >.
Definition at line 161 of file kOmegaSSTBase.C.
References Foam::devTwoSymm(), tmp< T >::New(), IOobject::scopedName(), and Foam::type().
Referenced by correct(), and GbyNu().


|
protectedvirtual |
Return G/nu.
Reimplemented in kOmegaSSTDDES< BasicTurbulenceModel >, kOmegaSSTDDES< EddyDiffusivity< fluidThermoCompressibleTurbulenceModel > >, kOmegaSSTDDES< Foam::EddyDiffusivity< Foam::fluidThermoCompressibleTurbulenceModel > >, kOmegaSSTDES< BasicTurbulenceModel >, kOmegaSSTDES< EddyDiffusivity< fluidThermoCompressibleTurbulenceModel > >, and kOmegaSSTDES< Foam::EddyDiffusivity< Foam::fluidThermoCompressibleTurbulenceModel > >.
Definition at line 176 of file kOmegaSSTBase.C.
References a1_, b1_, betaStar_, c1_, F2(), GbyNu0(), Foam::max(), Foam::min(), omega_, S2(), and Foam::sqrt().
Referenced by correct().


|
protectedvirtual |
Definition at line 192 of file kOmegaSSTBase.C.
References tmp< T >::New().
Referenced by correct().


|
protectedvirtual |
Definition at line 203 of file kOmegaSSTBase.C.
References Foam::dimTime, Foam::dimVolume, tmp< T >::New(), and omega_.
Referenced by correct().


|
protectedvirtual |
Reimplemented in kOmegaSSTSAS< BasicTurbulenceModel >, kOmegaSSTSAS< EddyDiffusivity< fluidThermoCompressibleTurbulenceModel > >, and kOmegaSSTSAS< Foam::EddyDiffusivity< Foam::fluidThermoCompressibleTurbulenceModel > >.
Definition at line 214 of file kOmegaSSTBase.C.
References beta(), Foam::dimTime, Foam::dimVolume, gamma(), tmp< T >::New(), omega_, and S2().
Referenced by correct().


|
delete |
No copy assignment.
References alpha, kOmegaSSTBase(), phi, turbulenceModel::propertiesName, rho, and U.

|
virtual |
Re-read model coefficients if they have changed.
Reimplemented in kOmegaSSTDDES< BasicTurbulenceModel >, kOmegaSSTDDES< EddyDiffusivity< fluidThermoCompressibleTurbulenceModel > >, kOmegaSSTDDES< Foam::EddyDiffusivity< Foam::fluidThermoCompressibleTurbulenceModel > >, kOmegaSSTDES< BasicTurbulenceModel >, kOmegaSSTDES< EddyDiffusivity< fluidThermoCompressibleTurbulenceModel > >, kOmegaSSTDES< Foam::EddyDiffusivity< Foam::fluidThermoCompressibleTurbulenceModel > >, kOmegaSSTIDDES< BasicTurbulenceModel >, kOmegaSSTIDDES< EddyDiffusivity< fluidThermoCompressibleTurbulenceModel > >, kOmegaSSTIDDES< Foam::EddyDiffusivity< Foam::fluidThermoCompressibleTurbulenceModel > >, kOmegaSSTLM< BasicTurbulenceModel >, kOmegaSSTLM< EddyDiffusivity< fluidThermoCompressibleTurbulenceModel > >, kOmegaSSTLM< Foam::EddyDiffusivity< Foam::fluidThermoCompressibleTurbulenceModel > >, kOmegaSSTSAS< BasicTurbulenceModel >, kOmegaSSTSAS< EddyDiffusivity< fluidThermoCompressibleTurbulenceModel > >, kOmegaSSTSAS< Foam::EddyDiffusivity< Foam::fluidThermoCompressibleTurbulenceModel > >, and kOmegaSSTSato< BasicTurbulenceModel >.
Definition at line 465 of file kOmegaSSTBase.C.
References a1_, alphaK1_, alphaK2_, alphaOmega1_, alphaOmega2_, b1_, beta1_, beta2_, betaStar_, c1_, F3_, gamma1_, gamma2_, and setDecayControl().

|
inline |
Return the effective diffusivity for k.
Definition at line 376 of file kOmegaSSTBase.H.
References alphaK(), F1(), tmp< T >::New(), and nu.
Referenced by correct().


|
inline |
Return the effective diffusivity for omega.
Definition at line 388 of file kOmegaSSTBase.H.
References alphaOmega(), F1(), tmp< T >::New(), and nu.
Referenced by correct().


|
inlinevirtual |
Return the turbulence kinetic energy.
Definition at line 400 of file kOmegaSSTBase.H.
References k_.
|
inlinevirtual |
Return the turbulence kinetic energy dissipation rate.
Definition at line 408 of file kOmegaSSTBase.H.
References omega_.
|
virtual |
Solve the turbulence equations and correct the turbulence viscosity.
1. Explicitly swap values on coupled boundary conditions
Update omega and G at the wall
Reimplemented in kOmegaSSTLM< BasicTurbulenceModel >, kOmegaSSTLM< EddyDiffusivity< fluidThermoCompressibleTurbulenceModel > >, kOmegaSSTLM< Foam::EddyDiffusivity< Foam::fluidThermoCompressibleTurbulenceModel > >, and kOmegaSSTSato< BasicTurbulenceModel >.
Definition at line 493 of file kOmegaSSTBase.C.
References Foam::fvc::absolute(), alpha, alphaOmega2_, beta(), beta(), betaStar_, Foam::bound(), tmp< T >::clear(), correctNut(), Foam::fvm::ddt(), Foam::fvc::div(), Foam::fvm::div(), divU, DkEff(), DomegaEff(), epsilonByk(), F1, F1(), F23(), fvOptions, gamma(), gamma, GbyNu(), GbyNu0(), Foam::fvc::grad(), k_, kInf_, kSource(), Foam::fvm::laplacian(), options::New(), nut, omega_, omegaInf_, omegaSource(), phi, Pk(), Qsas(), tmp< T >::ref(), rho, S2(), solve(), Foam::fvm::Sp(), Foam::sqr(), Foam::fvm::SuSp(), and U.

|
protected |
Definition at line 135 of file kOmegaSSTBase.H.
Referenced by alphaK(), kOmegaSSTBase(), and read().
|
protected |
Definition at line 136 of file kOmegaSSTBase.H.
|
protected |
Definition at line 138 of file kOmegaSSTBase.H.
Referenced by alphaOmega(), and read().
|
protected |
Definition at line 139 of file kOmegaSSTBase.H.
Referenced by alphaOmega(), correct(), F1(), and read().
|
protected |
Definition at line 141 of file kOmegaSSTBase.H.
|
protected |
Definition at line 142 of file kOmegaSSTBase.H.
|
protected |
Definition at line 144 of file kOmegaSSTBase.H.
|
protected |
Definition at line 145 of file kOmegaSSTBase.H.
|
protected |
Definition at line 147 of file kOmegaSSTBase.H.
Referenced by correct(), epsilonByk(), F1(), F2(), GbyNu(), Pk(), and read().
|
protected |
Definition at line 149 of file kOmegaSSTBase.H.
Referenced by correctNut(), GbyNu(), and read().
|
protected |
Definition at line 150 of file kOmegaSSTBase.H.
Referenced by correctNut(), GbyNu(), and read().
|
protected |
Definition at line 151 of file kOmegaSSTBase.H.
|
protected |
Flag to include the F3 term.
Definition at line 156 of file kOmegaSSTBase.H.
|
protected |
|
protected |
Turbulent kinetic energy field [m^2/s^2].
Definition at line 172 of file kOmegaSSTBase.H.
Referenced by correct(), correctNut(), F1(), F2(), k(), and Pk().
|
protected |
Specific dissipation rate field [1/s].
Definition at line 177 of file kOmegaSSTBase.H.
Referenced by correct(), correctNut(), epsilonByk(), F1(), F2(), F3(), GbyNu(), omega(), omegaSource(), Pk(), and Qsas().
|
protected |
Flag to include the decay control.
Definition at line 185 of file kOmegaSSTBase.H.
Referenced by setDecayControl().
|
protected |
Definition at line 186 of file kOmegaSSTBase.H.
Referenced by correct(), and setDecayControl().
|
protected |
Definition at line 187 of file kOmegaSSTBase.H.
Referenced by correct(), and setDecayControl().