Loading...
Searching...
No Matches
kOmegaSSTLM< BasicTurbulenceModel > Class Template Reference

Langtry-Menter 4-equation transitional SST model based on the k-omega-SST RAS model. More...

#include <kOmegaSSTLM.H>

Inheritance diagram for kOmegaSSTLM< BasicTurbulenceModel >:
Collaboration diagram for kOmegaSSTLM< BasicTurbulenceModel >:

Public Types

typedef BasicTurbulenceModel::alphaField alphaField
typedef BasicTurbulenceModel::rhoField rhoField
typedef BasicTurbulenceModel::transportModel transportModel
Public Types inherited from kOmegaSST< BasicTurbulenceModel >
typedef BasicTurbulenceModel::alphaField alphaField
typedef BasicTurbulenceModel::rhoField rhoField
typedef BasicTurbulenceModel::transportModel transportModel
Public Types inherited from kOmegaSSTBase< eddyViscosity< RASModel< BasicTurbulenceModel > > >
typedef eddyViscosity< RASModel< BasicTurbulenceModel > >::alphaField alphaField
typedef eddyViscosity< RASModel< BasicTurbulenceModel > >::rhoField rhoField
typedef eddyViscosity< RASModel< BasicTurbulenceModel > >::transportModel transportModel
Public Types inherited from eddyViscosity< RASModel< BasicTurbulenceModel > >
typedef BasicTurbulenceModel::alphaField alphaField
typedef BasicTurbulenceModel::rhoField rhoField
typedef BasicTurbulenceModel::transportModel transportModel
Public Types inherited from linearViscousStress< BasicTurbulenceModel >
typedef BasicTurbulenceModel::alphaField alphaField
typedef BasicTurbulenceModel::rhoField rhoField
typedef BasicTurbulenceModel::transportModel transportModel

Public Member Functions

 TypeName ("kOmegaSSTLM")
 Runtime type information.
 kOmegaSSTLM (const alphaField &alpha, const rhoField &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.
virtual ~kOmegaSSTLM ()=default
 Destructor.
virtual bool read ()
 Re-read model coefficients if they have changed.
const volScalarFieldReThetat () const
 Access function transition onset momentum-thickness Reynolds number.
const volScalarFieldgammaInt () const
 Access function to intermittency.
tmp< volScalarFieldDReThetatEff () const
 Return the effective diffusivity for transition onset.
tmp< volScalarFieldDgammaIntEff () const
 Return the effective diffusivity for intermittency.
virtual void correct ()
 Solve the turbulence equations and correct the turbulence viscosity.
Public Member Functions inherited from kOmegaSST< BasicTurbulenceModel >
 TypeName ("kOmegaSST")
 Runtime type information.
 kOmegaSST (const alphaField &alpha, const rhoField &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.
virtual ~kOmegaSST ()=default
 Destructor.
Public Member Functions inherited from kOmegaSSTBase< eddyViscosity< RASModel< BasicTurbulenceModel > > >
 kOmegaSSTBase (const kOmegaSSTBase &)=delete
 No copy construct.
void operator= (const kOmegaSSTBase &)=delete
 No copy assignment.
virtual ~kOmegaSSTBase ()=default
 Destructor.
tmp< volScalarFieldDkEff (const volScalarField &F1) const
 Return the effective diffusivity for k.
tmp< volScalarFieldDomegaEff (const volScalarField &F1) const
 Return the effective diffusivity for omega.
virtual tmp< volScalarFieldk () const
 Return the turbulence kinetic energy.
virtual tmp< volScalarFieldomega () const
 Return the turbulence kinetic energy dissipation rate.
Public Member Functions inherited from eddyViscosity< RASModel< BasicTurbulenceModel > >
 eddyViscosity (const word &modelName, const alphaField &alpha, const rhoField &rho, const volVectorField &U, const surfaceScalarField &alphaRhoPhi, const surfaceScalarField &phi, const transportModel &transport, const word &propertiesName)
 Construct from components.
virtual ~eddyViscosity ()=default
 Destructor.
virtual tmp< volScalarFieldnut () const
 Return the turbulence viscosity.
virtual tmp< volSymmTensorFieldR () const
 Return the Reynolds stress tensor.
virtual void validate ()
 Validate the turbulence fields after construction.
Public Member Functions inherited from linearViscousStress< BasicTurbulenceModel >
 linearViscousStress (const word &modelName, const alphaField &alpha, const rhoField &rho, const volVectorField &U, const surfaceScalarField &alphaRhoPhi, const surfaceScalarField &phi, const transportModel &transport, const word &propertiesName)
 Construct from components.
virtual ~linearViscousStress ()=default
 Destructor.
virtual tmp< volSymmTensorFielddevRhoReff () const
 Return the effective stress tensor.
virtual tmp< volSymmTensorFielddevRhoReff (const volVectorField &U) const
 Return the effective stress tensor based on a given velocity field.
virtual tmp< fvVectorMatrixdivDevRhoReff (volVectorField &U) const
 Return the source term for the momentum equation.
virtual tmp< fvVectorMatrixdivDevRhoReff (const volScalarField &rho, volVectorField &U) const
 Return the source term for the momentum equation.

Protected Member Functions

virtual tmp< volScalarFieldF1 (const volScalarField &CDkOmega) const
 Modified form of the k-omega SST F1 function.
virtual tmp< volScalarField::InternalPk (const volScalarField::Internal &G) const
 Modified form of the k-omega SST k production rate.
virtual tmp< volScalarField::InternalepsilonByk (const volScalarField &F1, const volTensorField &gradU) const
 Modified form of the k-omega SST epsilon/k.
tmp< volScalarField::InternalFthetat (const volScalarField::Internal &Us, const volScalarField::Internal &Omega, const volScalarField::Internal &nu) const
 Freestream blending-function.
tmp< volScalarField::InternalReThetac () const
 Empirical correlation for critical Reynolds number where the.
tmp< volScalarField::InternalFlength (const volScalarField::Internal &nu) const
 Empirical correlation that controls the length of the.
tmp< volScalarField::InternalFonset (const volScalarField::Internal &Rev, const volScalarField::Internal &ReThetac, const volScalarField::Internal &RT) const
 Transition onset location control function.
tmp< volScalarField::InternalReThetat0 (const volScalarField::Internal &Us, const volScalarField::Internal &dUsds, const volScalarField::Internal &nu) const
 Return the transition onset momentum-thickness Reynolds number.
void correctReThetatGammaInt ()
 Solve the turbulence equations and correct the turbulence viscosity.
Protected Member Functions inherited from kOmegaSST< BasicTurbulenceModel >
virtual void correctNut (const volScalarField &S2)
virtual void correctNut ()
Protected Member Functions inherited from kOmegaSSTBase< eddyViscosity< RASModel< BasicTurbulenceModel > > >
void setDecayControl (const dictionary &dict)
 Set decay control with kInf and omegaInf.
virtual tmp< volScalarFieldF2 () const
virtual tmp< volScalarFieldF3 () const
virtual tmp< volScalarFieldF23 () const
tmp< volScalarFieldblend (const volScalarField &F1, const dimensionedScalar &psi1, const dimensionedScalar &psi2) const
 Return the blended field.
tmp< volScalarFieldalphaK (const volScalarField &F1) const
tmp< volScalarFieldalphaOmega (const volScalarField &F1) const
tmp< volScalarField::Internalbeta (const volScalarField::Internal &F1) const
tmp< volScalarField::Internalgamma (const volScalarField::Internal &F1) const
virtual tmp< volScalarFieldS2 (const volTensorField &gradU) const
 Return square of strain rate.
virtual tmp< volScalarField::InternalGbyNu0 (const volTensorField &gradU, const volScalarField &S2) const
 Return (G/nu)_0.
virtual tmp< volScalarField::InternalGbyNu (const volScalarField::Internal &GbyNu0, const volScalarField::Internal &F2, const volScalarField::Internal &S2) const
 Return G/nu.
virtual tmp< fvScalarMatrixkSource () const
virtual tmp< fvScalarMatrixomegaSource () const
virtual tmp< fvScalarMatrixQsas (const volScalarField::Internal &S2, const volScalarField::Internal &gamma, const volScalarField::Internal &beta) const

Protected Attributes

dimensionedScalar ca1_
dimensionedScalar ca2_
dimensionedScalar ce1_
dimensionedScalar ce2_
dimensionedScalar cThetat_
dimensionedScalar sigmaThetat_
scalar lambdaErr_
 Convergence criterion for the lambda/thetat loop.
label maxLambdaIter_
 Maximum number of iterations to converge the lambda/thetat loop.
const dimensionedScalar deltaU_
 Stabilization for division by the magnitude of the velocity.
volScalarField ReThetat_
 Transition onset momentum-thickness Reynolds number.
volScalarField gammaInt_
 Intermittency.
volScalarField::Internal gammaIntEff_
 Effective intermittency.
Protected Attributes inherited from kOmegaSSTBase< eddyViscosity< RASModel< BasicTurbulenceModel > > >
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 volScalarFieldy_
 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_
Protected Attributes inherited from eddyViscosity< RASModel< BasicTurbulenceModel > >
volScalarField nut_

Detailed Description

template<class BasicTurbulenceModel>
class Foam::RASModels::kOmegaSSTLM< BasicTurbulenceModel >

Langtry-Menter 4-equation transitional SST model based on the k-omega-SST RAS model.

References:

    Langtry, R. B., & Menter, F. R. (2009).
    Correlation-based transition modeling for unstructured parallelized
    computational fluid dynamics codes.
    AIAA journal, 47(12), 2894-2906.

    Menter, F. R., Langtry, R., & Volker, S. (2006).
    Transition modelling for general purpose CFD codes.
    Flow, turbulence and combustion, 77(1-4), 277-303.

    Langtry, R. B. (2006).
    A correlation-based transition model using local variables for
    unstructured parallelized CFD codes.
    Phd. Thesis, Universität Stuttgart.

The model coefficients are

    kOmegaSSTCoeffs
    {
        // Default SST coefficients
        alphaK1     0.85;
        alphaK2     1;
        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;
        c1          10;
        F3          no;

        // Default LM coefficients
        ca1         2;
        ca2         0.06;
        ce1         1;
        ce2         50;
        cThetat     0.03;
        sigmaThetat 2;

        lambdaErr   1e-6;
        maxLambdaIter 10;
    }
Source files

Definition at line 103 of file kOmegaSSTLM.H.

Member Typedef Documentation

◆ alphaField

template<class BasicTurbulenceModel>
typedef BasicTurbulenceModel::alphaField alphaField

Definition at line 250 of file kOmegaSSTLM.H.

◆ rhoField

template<class BasicTurbulenceModel>
typedef BasicTurbulenceModel::rhoField rhoField

Definition at line 251 of file kOmegaSSTLM.H.

◆ transportModel

template<class BasicTurbulenceModel>
typedef BasicTurbulenceModel::transportModel transportModel

Definition at line 252 of file kOmegaSSTLM.H.

Constructor & Destructor Documentation

◆ kOmegaSSTLM()

template<class BasicTurbulenceModel>
kOmegaSSTLM ( const alphaField & alpha,
const rhoField & rho,
const volVectorField & U,
const surfaceScalarField & alphaRhoPhi,
const surfaceScalarField & phi,
const transportModel & transport,
const word & propertiesName = turbulenceModel::propertiesName,
const word & type = typeName )

◆ ~kOmegaSSTLM()

template<class BasicTurbulenceModel>
virtual ~kOmegaSSTLM ( )
virtualdefault

Destructor.

Member Function Documentation

◆ F1()

template<class BasicTurbulenceModel>
tmp< volScalarField > F1 ( const volScalarField & CDkOmega) const
protectedvirtual

◆ Pk()

template<class BasicTurbulenceModel>
tmp< volScalarField::Internal > Pk ( const volScalarField::Internal & G) const
protectedvirtual

Modified form of the k-omega SST k production rate.

Reimplemented from kOmegaSSTBase< eddyViscosity< RASModel< BasicTurbulenceModel > > >.

Definition at line 49 of file kOmegaSSTLM.C.

References gammaIntEff_, and kOmegaSSTBase< eddyViscosity< RASModel< BasicTurbulenceModel > > >::Pk().

Here is the call graph for this function:

◆ epsilonByk()

template<class BasicTurbulenceModel>
tmp< volScalarField::Internal > epsilonByk ( const volScalarField & F1,
const volTensorField & gradU ) const
protectedvirtual

Modified form of the k-omega SST epsilon/k.

Reimplemented from kOmegaSSTBase< eddyViscosity< RASModel< BasicTurbulenceModel > > >.

Definition at line 59 of file kOmegaSSTLM.C.

References Foam::clamp(), kOmegaSSTBase< eddyViscosity< RASModel< BasicTurbulenceModel > > >::epsilonByk(), F1(), and gammaIntEff_.

Here is the call graph for this function:

◆ Fthetat()

◆ ReThetac()

template<class BasicTurbulenceModel>
tmp< volScalarField::Internal > ReThetac ( ) const
protected

Empirical correlation for critical Reynolds number where the.

intermittency first starts to increase in the boundary layer

Definition at line 112 of file kOmegaSSTLM.C.

References Foam::dimless, forAll, IOobject::groupName(), GeometricField< scalar, fvPatchField, volMesh >::New(), IOobjectOption::NO_REGISTER, Foam::pow3(), Foam::pow4(), ReThetac(), ReThetat(), ReThetat_, and Foam::sqr().

Referenced by correctReThetatGammaInt(), Fonset(), and ReThetac().

Here is the call graph for this function:
Here is the caller graph for this function:

◆ Flength()

template<class BasicTurbulenceModel>
tmp< volScalarField::Internal > Flength ( const volScalarField::Internal & nu) const
protected

◆ Fonset()

template<class BasicTurbulenceModel>
tmp< volScalarField::Internal > Fonset ( const volScalarField::Internal & Rev,
const volScalarField::Internal & ReThetac,
const volScalarField::Internal & RT ) const
protected

Transition onset location control function.

Definition at line 308 of file kOmegaSSTLM.C.

References IOobject::groupName(), Foam::max(), Foam::min(), GeometricField< scalar, fvPatchField, volMesh >::New(), IOobjectOption::NO_REGISTER, Foam::pow3(), Foam::pow4(), and ReThetac().

Referenced by correctReThetatGammaInt().

Here is the call graph for this function:
Here is the caller graph for this function:

◆ ReThetat0()

template<class BasicTurbulenceModel>
tmp< volScalarField::Internal > ReThetat0 ( const volScalarField::Internal & Us,
const volScalarField::Internal & dUsds,
const volScalarField::Internal & nu ) const
protected

◆ correctReThetatGammaInt()

◆ TypeName()

template<class BasicTurbulenceModel>
TypeName ( "kOmegaSSTLM< BasicTurbulenceModel >" )

Runtime type information.

◆ read()

template<class BasicTurbulenceModel>
bool read ( )
virtual

Re-read model coefficients if they have changed.

Reimplemented from kOmegaSSTBase< eddyViscosity< RASModel< BasicTurbulenceModel > > >.

Definition at line 474 of file kOmegaSSTLM.C.

References ca1_, ca2_, ce1_, ce2_, cThetat_, lambdaErr_, maxLambdaIter_, linearViscousStress< BasicTurbulenceModel >::read(), and sigmaThetat_.

Here is the call graph for this function:

◆ ReThetat()

template<class BasicTurbulenceModel>
const volScalarField & ReThetat ( ) const
inline

Access function transition onset momentum-thickness Reynolds number.

Definition at line 295 of file kOmegaSSTLM.H.

Referenced by Flength(), and ReThetac().

Here is the caller graph for this function:

◆ gammaInt()

template<class BasicTurbulenceModel>
const volScalarField & gammaInt ( ) const
inline

Access function to intermittency.

Definition at line 303 of file kOmegaSSTLM.H.

◆ DReThetatEff()

template<class BasicTurbulenceModel>
tmp< volScalarField > DReThetatEff ( ) const
inline

Return the effective diffusivity for transition onset.

momentum-thickness Reynolds number

Definition at line 313 of file kOmegaSSTLM.H.

Referenced by correctReThetatGammaInt().

Here is the caller graph for this function:

◆ DgammaIntEff()

template<class BasicTurbulenceModel>
tmp< volScalarField > DgammaIntEff ( ) const
inline

Return the effective diffusivity for intermittency.

Definition at line 328 of file kOmegaSSTLM.H.

Referenced by correctReThetatGammaInt().

Here is the caller graph for this function:

◆ correct()

template<class BasicTurbulenceModel>
void correct ( )
virtual

Solve the turbulence equations and correct the turbulence viscosity.

Reimplemented from kOmegaSSTBase< eddyViscosity< RASModel< BasicTurbulenceModel > > >.

Definition at line 593 of file kOmegaSSTLM.C.

References linearViscousStress< BasicTurbulenceModel >::correct(), and correctReThetatGammaInt().

Here is the call graph for this function:

Member Data Documentation

◆ ca1_

template<class BasicTurbulenceModel>
dimensionedScalar ca1_
protected

Definition at line 126 of file kOmegaSSTLM.H.

Referenced by correctReThetatGammaInt(), kOmegaSSTLM(), and read().

◆ ca2_

template<class BasicTurbulenceModel>
dimensionedScalar ca2_
protected

Definition at line 127 of file kOmegaSSTLM.H.

Referenced by correctReThetatGammaInt(), kOmegaSSTLM(), and read().

◆ ce1_

template<class BasicTurbulenceModel>
dimensionedScalar ce1_
protected

Definition at line 129 of file kOmegaSSTLM.H.

Referenced by correctReThetatGammaInt(), kOmegaSSTLM(), and read().

◆ ce2_

template<class BasicTurbulenceModel>
dimensionedScalar ce2_
protected

Definition at line 130 of file kOmegaSSTLM.H.

Referenced by correctReThetatGammaInt(), Fthetat(), kOmegaSSTLM(), and read().

◆ cThetat_

template<class BasicTurbulenceModel>
dimensionedScalar cThetat_
protected

Definition at line 132 of file kOmegaSSTLM.H.

Referenced by correctReThetatGammaInt(), kOmegaSSTLM(), and read().

◆ sigmaThetat_

template<class BasicTurbulenceModel>
dimensionedScalar sigmaThetat_
protected

Definition at line 133 of file kOmegaSSTLM.H.

Referenced by kOmegaSSTLM(), and read().

◆ lambdaErr_

template<class BasicTurbulenceModel>
scalar lambdaErr_
protected

Convergence criterion for the lambda/thetat loop.

Definition at line 138 of file kOmegaSSTLM.H.

Referenced by kOmegaSSTLM(), read(), and ReThetat0().

◆ maxLambdaIter_

template<class BasicTurbulenceModel>
label maxLambdaIter_
protected

Maximum number of iterations to converge the lambda/thetat loop.

Definition at line 143 of file kOmegaSSTLM.H.

Referenced by kOmegaSSTLM(), read(), and ReThetat0().

◆ deltaU_

template<class BasicTurbulenceModel>
const dimensionedScalar deltaU_
protected

Stabilization for division by the magnitude of the velocity.

Definition at line 148 of file kOmegaSSTLM.H.

Referenced by correctReThetatGammaInt(), and kOmegaSSTLM().

◆ ReThetat_

template<class BasicTurbulenceModel>
volScalarField ReThetat_
protected

Transition onset momentum-thickness Reynolds number.

Definition at line 156 of file kOmegaSSTLM.H.

Referenced by correctReThetatGammaInt(), Flength(), Fthetat(), kOmegaSSTLM(), and ReThetac().

◆ gammaInt_

template<class BasicTurbulenceModel>
volScalarField gammaInt_
protected

Intermittency.

Definition at line 161 of file kOmegaSSTLM.H.

Referenced by correctReThetatGammaInt(), Fthetat(), and kOmegaSSTLM().

◆ gammaIntEff_

template<class BasicTurbulenceModel>
volScalarField::Internal gammaIntEff_
protected

Effective intermittency.

Definition at line 166 of file kOmegaSSTLM.H.

Referenced by correctReThetatGammaInt(), epsilonByk(), kOmegaSSTLM(), and Pk().


The documentation for this class was generated from the following files:
  • src/TurbulenceModels/turbulenceModels/RAS/kOmegaSSTLM/kOmegaSSTLM.H
  • src/TurbulenceModels/turbulenceModels/RAS/kOmegaSSTLM/kOmegaSSTLM.C