kL< BasicTurbulenceModel > Class Template Reference

A one-equation (turbulent kinetic energy k) turbulence closure model for incompressible and compressible geophysical applications. More...

Inheritance diagram for kL< BasicTurbulenceModel >:
[legend]
Collaboration diagram for kL< BasicTurbulenceModel >:
[legend]

Public Types

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

Public Member Functions

 TypeName ("kL")
 Runtime type information. More...
 
 kL (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. More...
 
virtual ~kL ()=default
 Destructor. More...
 
virtual bool read ()
 Re-read model coefficients if they have changed. More...
 
tmp< volScalarFieldDkEff () const
 Return the effective diffusivity for k. More...
 
virtual tmp< volScalarFieldk () const
 Return the turbulent kinetic energy field. More...
 
virtual void correct ()
 Solve the turbulence equations and correct the turbulence viscosity. More...
 
- 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. More...
 
virtual ~eddyViscosity ()=default
 Destructor. More...
 
virtual tmp< volScalarFieldnut () const
 Return the turbulence viscosity. More...
 
virtual tmp< scalarFieldnut (const label patchi) const
 Return the turbulence viscosity on patch. More...
 
virtual tmp< volSymmTensorFieldR () const
 Return the Reynolds stress tensor. More...
 
virtual void validate ()
 Validate the turbulence fields after construction. More...
 
- Public Member Functions inherited from linearViscousStress< RASModel< 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. More...
 
virtual ~linearViscousStress ()=default
 Destructor. More...
 
virtual tmp< volSymmTensorFielddevRhoReff () const
 Return the effective stress tensor. More...
 
virtual tmp< volSymmTensorFielddevRhoReff (const volVectorField &U) const
 Return the effective stress tensor based on a given velocity field. More...
 
virtual tmp< fvVectorMatrixdivDevRhoReff (volVectorField &U) const
 Return the source term for the momentum equation. More...
 
virtual tmp< fvVectorMatrixdivDevRhoReff (const volScalarField &rho, volVectorField &U) const
 Return the source term for the momentum equation. More...
 
- Public Member Functions inherited from RASModel< BasicTurbulenceModel >
 TypeName ("RAS")
 Runtime type information. More...
 
 declareRunTimeSelectionTable (autoPtr, RASModel, dictionary,(const alphaField &alpha, const rhoField &rho, const volVectorField &U, const surfaceScalarField &alphaRhoPhi, const surfaceScalarField &phi, const transportModel &transport, const word &propertiesName),(alpha, rho, U, alphaRhoPhi, phi, transport, propertiesName))
 
 RASModel (const word &type, const alphaField &alpha, const rhoField &rho, const volVectorField &U, const surfaceScalarField &alphaRhoPhi, const surfaceScalarField &phi, const transportModel &transport, const word &propertiesName)
 Construct from components. More...
 
virtual ~RASModel ()=default
 Destructor. More...
 
const dimensionedScalarkMin () const
 Return the lower allowable limit for k (default: SMALL) More...
 
const dimensionedScalarepsilonMin () const
 Return the lower allowable limit for epsilon (default: SMALL) More...
 
const dimensionedScalaromegaMin () const
 Return the lower allowable limit for omega (default: SMALL) More...
 
dimensionedScalarkMin ()
 Allow kMin to be changed. More...
 
dimensionedScalarepsilonMin ()
 Allow epsilonMin to be changed. More...
 
dimensionedScalaromegaMin ()
 Allow omegaMin to be changed. More...
 
virtual const dictionarycoeffDict () const
 Const access to the coefficients dictionary. More...
 
virtual tmp< volScalarFieldnuEff () const
 Return the effective viscosity. More...
 
virtual tmp< scalarFieldnuEff (const label patchi) const
 Return the effective viscosity on patch. More...
 
virtual tmp< volScalarFieldomega () const
 Return the specific dissipation rate. More...
 

Protected Member Functions

virtual void correctNut ()
 Correct the turbulence viscosity. More...
 
virtual tmp< fvScalarMatrixkSource () const
 Add explicit source for turbulent kinetic energy. More...
 
- Protected Member Functions inherited from RASModel< BasicTurbulenceModel >
virtual void printCoeffs (const word &type)
 Print model coefficients. More...
 
 RASModel (const RASModel &)=delete
 No copy construct. More...
 
void operator= (const RASModel &)=delete
 No copy assignment. More...
 

Protected Attributes

dimensionedScalar kappa_
 von Karman constant More...
 
dimensionedScalar sigmak_
 Empirical model coefficient. More...
 
dimensionedScalar beta_
 Thermal expansion coefficient [1/K]. More...
 
dimensionedScalar Cmu0_
 Empirical model coefficient. More...
 
dimensionedScalar Lmax_
 Maximum mixing-length scalar [m]. More...
 
dimensionedScalar CbStable_
 Stable stratification constant. More...
 
dimensionedScalar CbUnstable_
 Unstable stratification constant. More...
 
volScalarField k_
 Turbulent kinetic energy [m2/s2]. More...
 
volScalarField L_
 Characteristic length scale [m]. More...
 
volScalarField Rt_
 Turbulent Richardson number [-]. More...
 
const uniformDimensionedVectorFieldg_
 Gravitational acceleration [m2/s2]. More...
 
const volScalarFieldy_
 Wall distance. More...
 
- Protected Attributes inherited from eddyViscosity< RASModel< BasicTurbulenceModel > >
volScalarField nut_
 
- Protected Attributes inherited from RASModel< BasicTurbulenceModel >
dictionary RASDict_
 RAS coefficients dictionary. More...
 
Switch turbulence_
 Turbulence on/off flag. More...
 
Switch printCoeffs_
 Flag to print the model coeffs at run-time. More...
 
dictionary coeffDict_
 Model coefficients dictionary. More...
 
dimensionedScalar kMin_
 Lower limit of k. More...
 
dimensionedScalar epsilonMin_
 Lower limit of epsilon. More...
 
dimensionedScalar omegaMin_
 Lower limit for omega. More...
 

Additional Inherited Members

- Static Public Member Functions inherited from RASModel< BasicTurbulenceModel >
static autoPtr< RASModelNew (const alphaField &alpha, const rhoField &rho, const volVectorField &U, const surfaceScalarField &alphaRhoPhi, const surfaceScalarField &phi, const transportModel &transport, const word &propertiesName=turbulenceModel::propertiesName)
 Return a reference to the selected RAS model. More...
 

Detailed Description

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

A one-equation (turbulent kinetic energy k) turbulence closure model for incompressible and compressible geophysical applications.

Turbulent kinetic energy (k) is computed with a transport equation and the turbulent length scale (L) is computed with an algebraic expression which depends on the local stratification.

References:

    Standard model (tag:A):
        Axell, L. B., & Liungman, O. (2001).
        A one-equation turbulence model for geophysical applications:
        comparison with data and the k−ε model.
        Environmental Fluid Mechanics, 1(1), 71-106.
        DOI:10.1023/A:1011560202388

    Canopy-related models (tag:W):
        Wilson, J. D., & Flesch, T. K. (1999).
        Wind and remnant tree sway in forest cutblocks.
        III. A windflow model to diagnose spatial variation.
        Agricultural and Forest Meteorology, 93(4), 259-282.
        DOI:10.1016/S0168-1923(98)00121-X
Usage
Example by using constant/turbulenceProperties:
RAS
{
    // Mandatory entries
    RASModel        kL;

    // Optional entries
    kLCoeffs
    {
        kappa       <scalar>;
        sigmak      <scalar>;
        beta        <scalar>;
        Cmu0        <scalar>;
        Lmax        <scalar>;
        CbStable    <scalar>;
        CbUnstable  <scalar>;
    }

    // Inherited entries
    ...
}

where the entries mean:

Property Description Type Reqd Deflt
RASModel Type name: kL word yes -
kappa von Karman constant scalar no 0.41
sigmak Empirical model coefficient scalar no 1.0
beta Thermal expansion coefficient [1/K] scalar no 3.3e-3
Cmu0 Empirical model coefficient scalar no 0.556
Lmax Maximum mixing-length scale [m] scalar no GREAT
CbStable Stable stratification constant scalar no 0.25
CbUnstable Unstable stratification constant scalar no 0.35

The inherited entries are elaborated in:


Input fields (mandatory)

k : Turbulent kinetic energy [m2/s2]
T : Potential temperature [K]


Input fields (optional)

canopyHeight : Canopy height [m]
plantCd : Plant canopy drag coefficient [-]
leafAreaDensity : Leaf area density [1/m]
Rt : Turbulent Richardson number [-]
L : Characteristic length scale [m]
Note
  • Optional input fields can/should be input by using readFields function object.
Source files

Definition at line 222 of file kL.H.

Member Typedef Documentation

◆ alphaField

typedef BasicTurbulenceModel::alphaField alphaField

Definition at line 324 of file kL.H.

◆ rhoField

typedef BasicTurbulenceModel::rhoField rhoField

Definition at line 325 of file kL.H.

◆ transportModel

typedef BasicTurbulenceModel::transportModel transportModel

Definition at line 326 of file kL.H.

Constructor & Destructor Documentation

◆ kL()

kL ( 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.

Definition at line 250 of file kL.C.

References Foam::bound(), and Foam::type().

Here is the call graph for this function:

◆ ~kL()

virtual ~kL ( )
virtualdefault

Destructor.

Member Function Documentation

◆ correctNut()

void correctNut ( )
protectedvirtual

Correct the turbulence viscosity.

Implements eddyViscosity< RASModel< BasicTurbulenceModel > >.

Definition at line 225 of file kL.C.

References optionList::correct(), options::New(), and Foam::sqrt().

Here is the call graph for this function:

◆ kSource()

tmp< fvScalarMatrix > kSource ( ) const
protectedvirtual

Add explicit source for turbulent kinetic energy.

Definition at line 236 of file kL.C.

References Foam::dimTime, Foam::dimVolume, and tmp< T >::New().

Here is the call graph for this function:

◆ TypeName()

TypeName ( "kL< BasicTurbulenceModel >"  )

Runtime type information.

◆ read()

bool read ( )
virtual

Re-read model coefficients if they have changed.

Implements eddyViscosity< RASModel< BasicTurbulenceModel > >.

Definition at line 392 of file kL.C.

References Foam::read().

Here is the call graph for this function:

◆ DkEff()

tmp<volScalarField> DkEff ( ) const
inline

Return the effective diffusivity for k.

Definition at line 359 of file kL.H.

References tmp< T >::New(), nu, and eddyViscosity< RASModel< BasicTurbulenceModel > >::nut_.

Here is the call graph for this function:

◆ k()

virtual tmp<volScalarField> k ( ) const
inlinevirtual

Return the turbulent kinetic energy field.

Implements eddyViscosity< RASModel< BasicTurbulenceModel > >.

Definition at line 369 of file kL.H.

References kL< BasicTurbulenceModel >::k_.

◆ correct()

void correct ( )
virtual

Solve the turbulence equations and correct the turbulence viscosity.

Implements eddyViscosity< RASModel< BasicTurbulenceModel > >.

Definition at line 412 of file kL.C.

References Foam::fvc::absolute(), Foam::constant::atomic::alpha, Foam::bound(), tmp< T >::clear(), correct(), tmp< T >::cref(), Foam::fvm::ddt(), Foam::dev(), Foam::fvm::div(), Foam::fvc::div(), divU, epsilon, fvOptions, Foam::constant::universal::G, Foam::fvc::grad(), Foam::fvm::laplacian(), Foam::magSqr(), options::New(), nut, phi, tmp< T >::ref(), rho, Foam::solve(), Foam::sqr(), Foam::fvm::SuSp(), Foam::symm(), Foam::T(), and U.

Here is the call graph for this function:

Member Data Documentation

◆ kappa_

dimensionedScalar kappa_
protected

von Karman constant

Definition at line 272 of file kL.H.

◆ sigmak_

dimensionedScalar sigmak_
protected

Empirical model coefficient.

Definition at line 275 of file kL.H.

◆ beta_

dimensionedScalar beta_
protected

Thermal expansion coefficient [1/K].

Definition at line 278 of file kL.H.

◆ Cmu0_

dimensionedScalar Cmu0_
protected

Empirical model coefficient.

Definition at line 281 of file kL.H.

◆ Lmax_

dimensionedScalar Lmax_
protected

Maximum mixing-length scalar [m].

Definition at line 284 of file kL.H.

◆ CbStable_

dimensionedScalar CbStable_
protected

Stable stratification constant.

Definition at line 287 of file kL.H.

◆ CbUnstable_

dimensionedScalar CbUnstable_
protected

Unstable stratification constant.

Definition at line 290 of file kL.H.

◆ k_

volScalarField k_
protected

Turbulent kinetic energy [m2/s2].

Definition at line 296 of file kL.H.

Referenced by kL< BasicTurbulenceModel >::k().

◆ L_

volScalarField L_
protected

Characteristic length scale [m].

Definition at line 299 of file kL.H.

◆ Rt_

volScalarField Rt_
protected

Turbulent Richardson number [-].

Definition at line 302 of file kL.H.

◆ g_

const uniformDimensionedVectorField& g_
protected

Gravitational acceleration [m2/s2].

Definition at line 305 of file kL.H.

◆ y_

const volScalarField& y_
protected

Wall distance.

Note: different to wall distance in parent RASModel which is for near-wall cells only

Definition at line 310 of file kL.H.


The documentation for this class was generated from the following files: