Public Member Functions | Static Public Member Functions | Protected Types | Protected Member Functions | Protected Attributes | Static Protected Attributes | List of all members
nutWallFunctionFvPatchScalarField Class Referenceabstract

The class nutWallFunction is a base class that parents the derived boundary conditions which provide a wall constraint on various fields, such as turbulent viscosity, i.e. nut, or turbulent kinetic energy dissipation rate, i.e. epsilon, for low- and high-Reynolds number turbulence models. The class is not an executable itself, yet a provider for common entries to its derived boundary conditions. More...

Inheritance diagram for nutWallFunctionFvPatchScalarField:
Inheritance graph
[legend]
Collaboration diagram for nutWallFunctionFvPatchScalarField:
Collaboration graph
[legend]

Public Member Functions

 TypeName ("nutWallFunction")
 
 nutWallFunctionFvPatchScalarField (const fvPatch &, const DimensionedField< scalar, volMesh > &)
 
 nutWallFunctionFvPatchScalarField (const fvPatch &, const DimensionedField< scalar, volMesh > &, const dictionary &)
 
 nutWallFunctionFvPatchScalarField (const nutWallFunctionFvPatchScalarField &, const fvPatch &, const DimensionedField< scalar, volMesh > &, const fvPatchFieldMapper &)
 
 nutWallFunctionFvPatchScalarField (const nutWallFunctionFvPatchScalarField &)
 
 nutWallFunctionFvPatchScalarField (const nutWallFunctionFvPatchScalarField &, const DimensionedField< scalar, volMesh > &)
 
scalar Cmu () const
 
scalar kappa () const
 
scalar E () const
 
scalar yPlusLam () const
 
scalar blend (const scalar nutVis, const scalar nutLog, const scalar yPlus) const
 
virtual tmp< scalarFieldyPlus () const =0
 
virtual void updateCoeffs ()
 
virtual void write (Ostream &) const
 

Static Public Member Functions

static const nutWallFunctionFvPatchScalarFieldnutw (const turbulenceModel &turbModel, const label patchi)
 
static scalar yPlusLam (const scalar kappa, const scalar E)
 

Protected Types

enum  blendingType { STEPWISE, MAX, BINOMIAL, EXPONENTIAL }
 

Protected Member Functions

virtual const volVectorFieldU (const turbulenceModel &turb) const
 
virtual void checkType ()
 
virtual tmp< scalarFieldcalcNut () const =0
 
virtual void writeLocalEntries (Ostream &) const
 

Protected Attributes

enum blendingType blending_
 
const scalar n_
 
word UName_
 
scalar Cmu_
 
scalar kappa_
 
scalar E_
 
scalar yPlusLam_
 

Static Protected Attributes

static const Enum< blendingTypeblendingTypeNames
 

Detailed Description

The class nutWallFunction is a base class that parents the derived boundary conditions which provide a wall constraint on various fields, such as turbulent viscosity, i.e. nut, or turbulent kinetic energy dissipation rate, i.e. epsilon, for low- and high-Reynolds number turbulence models. The class is not an executable itself, yet a provider for common entries to its derived boundary conditions.

Reference:

    Default model coefficients (tag:VM):
        Versteeg, H. K., & Malalasekera, W. (2011).
        An introduction to computational fluid dynamics: The finite
        volume method. Harlow: Pearson Education.
        Subsection "3.5.2 k-epsilon model".

    Binomial blending of the viscous and inertial sublayers (tag:ME):
        Menter, F., & Esch, T. (2001).
        Elements of industrial heat transfer prediction.
        In Proceedings of the 16th Brazilian Congress of Mechanical
        Engineering (COBEM), November 2001. vol. 20, p. 117-127.

    Exponential/Max blending of the viscous and inertial sublayers (tag:PH):
        Popovac, M., & Hanjalić, K. (2007).
        Compound wall treatment for RANS computation of complex
        turbulent flows and heat transfer.
        Flow, turbulence and combustion, 78(2), 177-202.
        DOI:10.1007/s10494-006-9067-x
Usage
Example of the boundary condition specification:
<patchName>
{
    // Mandatory and other optional entries
    ...

    // Optional (inherited) entries
    Cmu             0.09;
    kappa           0.41;
    E               9.8;
    blending        stepwise;
    n               4.0;
    U               U;

    // Optional (inherited) entries
    ...
}

where the entries mean:

Property Description Type Req'd Dflt
Cmu Empirical model coefficient scalar no 0.09
kappa von Kármán constant scalar no 0.41
E Wall roughness parameter scalar no 9.8
blending Viscous/inertial sublayer blending word no stepwise
n Binomial blending exponent scalar no 2.0
U Name of the velocity field word no U

The inherited entries are elaborated in:

Options for the blending entry:

      stepwise    | Stepwise switch (discontinuous)
      max         | Maximum value switch (discontinuous)
      binomial    | Binomial blending (smooth)
      exponential | Exponential blending (smooth)

wherein nut predictions for the viscous and inertial sublayers are blended according to the following expressions:

\[ \nu_t = {\nu_t}_{log} \qquad if \quad y^+ > y^+_{lam} \]

\[ \nu_t = {\nu_t}_{vis} \qquad if \quad y^+ <= y^+_{lam} \]

where

$ {\nu_t}_{vis} $ = $\nu_t$ prediction in the viscous sublayer
$ {\nu_t}_{log} $ = $\nu_t$ prediction in the inertial sublayer
$ y^+ $ = estimated wall-normal height of the cell centre in wall units
$ y^+_{lam} $ = estimated intersection of the viscous and inertial sublayers

\[ \nu_t = max({\nu_t}_{vis}, {\nu_t}_{log}) \]

\[ \nu_t = (({\nu_t}_{vis})^n + ({\nu_t}_{log})^n)^{1/n} \]

where

$ n $ = Binomial blending exponent

\[ \nu_t = {\nu_t}_{vis} \exp[-\Gamma] + {\nu_t}_{log} \exp[-1/\Gamma] \]

where (PH:Eq. 31)

$ \Gamma $ = Blending expression
$ \Gamma $ = $0.01 (y^+)^4 / (1.0 + 5.0 y^+)$
See also
Source files

Definition at line 247 of file nutWallFunctionFvPatchScalarField.H.

Member Enumeration Documentation

◆ blendingType

enum blendingType
protected
Enumerator
STEPWISE 

"Stepwise switch (discontinuous)"

MAX 

"Maximum value switch (discontinuous)"

BINOMIAL 

"Binomial blending (smooth)"

EXPONENTIAL 

"Exponential blending (smooth)"

Definition at line 256 of file nutWallFunctionFvPatchScalarField.H.

Constructor & Destructor Documentation

◆ nutWallFunctionFvPatchScalarField() [1/5]

nutWallFunctionFvPatchScalarField ( const fvPatch p,
const DimensionedField< scalar, volMesh > &  iF 
)

Definition at line 109 of file nutWallFunctionFvPatchScalarField.C.

◆ nutWallFunctionFvPatchScalarField() [2/5]

nutWallFunctionFvPatchScalarField ( const fvPatch p,
const DimensionedField< scalar, volMesh > &  iF,
const dictionary dict 
)

Definition at line 149 of file nutWallFunctionFvPatchScalarField.C.

◆ nutWallFunctionFvPatchScalarField() [3/5]

nutWallFunctionFvPatchScalarField ( const nutWallFunctionFvPatchScalarField ptf,
const fvPatch p,
const DimensionedField< scalar, volMesh > &  iF,
const fvPatchFieldMapper mapper 
)

Definition at line 128 of file nutWallFunctionFvPatchScalarField.C.

◆ nutWallFunctionFvPatchScalarField() [4/5]

Definition at line 188 of file nutWallFunctionFvPatchScalarField.C.

◆ nutWallFunctionFvPatchScalarField() [5/5]

Definition at line 206 of file nutWallFunctionFvPatchScalarField.C.

Member Function Documentation

◆ U()

const Foam::volVectorField & U ( const turbulenceModel turb) const
protectedvirtual

◆ checkType()

void checkType ( )
protectedvirtual

Definition at line 53 of file nutWallFunctionFvPatchScalarField.C.

References Foam::abort(), Foam::endl(), Foam::FatalError, FatalErrorInFunction, Foam::nl, and Foam::foamVersion::patch.

Here is the call graph for this function:

◆ calcNut()

virtual tmp<scalarField> calcNut ( ) const
protectedpure virtual

◆ writeLocalEntries()

void writeLocalEntries ( Ostream os) const
protectedvirtual

◆ TypeName()

TypeName ( "nutWallFunction"  )

◆ Cmu()

scalar Cmu ( ) const
inline

◆ kappa()

scalar kappa ( ) const
inline

◆ E()

scalar E ( ) const
inline

◆ nutw()

const Foam::nutWallFunctionFvPatchScalarField & nutw ( const turbulenceModel turbModel,
const label  patchi 
)
static

◆ yPlusLam() [1/2]

Foam::scalar yPlusLam ( const scalar  kappa,
const scalar  E 
)
static

◆ yPlusLam() [2/2]

Foam::scalar yPlusLam ( ) const

Definition at line 317 of file nutWallFunctionFvPatchScalarField.C.

◆ blend()

Foam::scalar blend ( const scalar  nutVis,
const scalar  nutLog,
const scalar  yPlus 
) const

Definition at line 260 of file nutWallFunctionFvPatchScalarField.C.

References Foam::exp(), Foam::max(), Foam::pow(), Foam::pow4(), and yPlus.

Referenced by nutkWallFunctionFvPatchScalarField::calcNut(), and nutUWallFunctionFvPatchScalarField::calcNut().

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

◆ yPlus()

virtual tmp<scalarField> yPlus ( ) const
pure virtual

◆ updateCoeffs()

void updateCoeffs ( )
virtual

Definition at line 323 of file nutWallFunctionFvPatchScalarField.C.

References Foam::operator==().

Here is the call graph for this function:

◆ write()

void write ( Ostream os) const
virtual

Member Data Documentation

◆ blendingTypeNames

const Foam::Enum< Foam::nutWallFunctionFvPatchScalarField::blendingType > blendingTypeNames
staticprotected

Definition at line 265 of file nutWallFunctionFvPatchScalarField.H.

◆ blending_

enum blendingType blending_
protected

Definition at line 271 of file nutWallFunctionFvPatchScalarField.H.

◆ n_

const scalar n_
protected

Definition at line 275 of file nutWallFunctionFvPatchScalarField.H.

◆ UName_

word UName_
protected

Definition at line 280 of file nutWallFunctionFvPatchScalarField.H.

◆ Cmu_

scalar Cmu_
protected

◆ kappa_

scalar kappa_
protected

◆ E_

scalar E_
protected

◆ yPlusLam_

scalar yPlusLam_
protected

Definition at line 293 of file nutWallFunctionFvPatchScalarField.H.


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