Public Member Functions | Protected Member Functions | Protected Attributes | Static Protected Attributes | List of all members
epsilonWallFunctionFvPatchScalarField Class Reference

This boundary condition provides a wall constraint on the turbulent kinetic energy dissipation rate, i.e. epsilon, and the turbulent kinetic energy production contribution, i.e. G, for low- and high-Reynolds number turbulence models. More...

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

Public Member Functions

 TypeName ("epsilonWallFunction")
 
 epsilonWallFunctionFvPatchScalarField (const fvPatch &, const DimensionedField< scalar, volMesh > &)
 
 epsilonWallFunctionFvPatchScalarField (const fvPatch &, const DimensionedField< scalar, volMesh > &, const dictionary &)
 
 epsilonWallFunctionFvPatchScalarField (const epsilonWallFunctionFvPatchScalarField &, const fvPatch &, const DimensionedField< scalar, volMesh > &, const fvPatchFieldMapper &)
 
 epsilonWallFunctionFvPatchScalarField (const epsilonWallFunctionFvPatchScalarField &)
 
virtual tmp< fvPatchScalarFieldclone () const
 
 epsilonWallFunctionFvPatchScalarField (const epsilonWallFunctionFvPatchScalarField &, const DimensionedField< scalar, volMesh > &)
 
virtual tmp< fvPatchScalarFieldclone (const DimensionedField< scalar, volMesh > &iF) const
 
virtual ~epsilonWallFunctionFvPatchScalarField ()=default
 
scalarFieldG (bool init=false)
 
scalarFieldepsilon (bool init=false)
 
virtual void updateCoeffs ()
 
virtual void updateWeightedCoeffs (const scalarField &weights)
 
virtual void manipulateMatrix (fvMatrix< scalar > &matrix)
 
virtual void manipulateMatrix (fvMatrix< scalar > &matrix, const scalarField &weights)
 
virtual void write (Ostream &) const
 
- Public Member Functions inherited from fixedValueFvPatchField< scalar >
 TypeName ("fixedValue")
 
 fixedValueFvPatchField (const fvPatch &, const DimensionedField< scalar, volMesh > &)
 
 fixedValueFvPatchField (const fvPatch &, const DimensionedField< scalar, volMesh > &, const scalar &value)
 
 fixedValueFvPatchField (const fvPatch &, const DimensionedField< scalar, volMesh > &, const dictionary &, const bool valueRequired=true)
 
 fixedValueFvPatchField (const fixedValueFvPatchField< scalar > &, const fvPatch &, const DimensionedField< scalar, volMesh > &, const fvPatchFieldMapper &)
 
 fixedValueFvPatchField (const fixedValueFvPatchField< scalar > &)
 
 fixedValueFvPatchField (const fixedValueFvPatchField< scalar > &, const DimensionedField< scalar, volMesh > &)
 
virtual tmp< fvPatchField< scalar > > clone () const
 
virtual bool fixesValue () const
 
virtual bool assignable () const
 
virtual tmp< Field< scalar > > valueInternalCoeffs (const tmp< scalarField > &) const
 
virtual tmp< Field< scalar > > valueBoundaryCoeffs (const tmp< scalarField > &) const
 
virtual tmp< Field< scalar > > gradientInternalCoeffs () const
 
virtual tmp< Field< scalar > > gradientBoundaryCoeffs () const
 
virtual void operator= (const UList< scalar > &)
 
virtual void operator= (const fvPatchField< scalar > &)
 
virtual void operator= (const scalar &)
 
virtual void operator+= (const fvPatchField< scalar > &)
 
virtual void operator+= (const Field< scalar > &)
 
virtual void operator+= (const scalar &)
 
virtual void operator-= (const fvPatchField< scalar > &)
 
virtual void operator-= (const Field< scalar > &)
 
virtual void operator-= (const scalar &)
 
virtual void operator*= (const fvPatchField< scalar > &)
 
virtual void operator*= (const Field< scalar > &)
 
virtual void operator*= (const scalar)
 
virtual void operator/= (const fvPatchField< scalar > &)
 
virtual void operator/= (const Field< scalar > &)
 
virtual void operator/= (const scalar)
 
- Public Member Functions inherited from fvPatchField< scalar >
 TypeName ("fvPatchField")
 
 declareRunTimeSelectionTable (tmp, fvPatchField, patch,(const fvPatch &p, const DimensionedField< Type, volMesh > &iF),(p, iF))
 
 declareRunTimeSelectionTable (tmp, fvPatchField, patchMapper,(const fvPatchField< Type > &ptf, const fvPatch &p, const DimensionedField< Type, volMesh > &iF, const fvPatchFieldMapper &m),(dynamic_cast< const fvPatchFieldType & >(ptf), p, iF, m))
 
 declareRunTimeSelectionTable (tmp, fvPatchField, dictionary,(const fvPatch &p, const DimensionedField< Type, volMesh > &iF, const dictionary &dict),(p, iF, dict))
 
 fvPatchField (const fvPatch &, const DimensionedField< Type, volMesh > &)
 
 fvPatchField (const fvPatch &, const DimensionedField< Type, volMesh > &, const Type &value)
 
 fvPatchField (const fvPatch &, const DimensionedField< Type, volMesh > &, const word &patchType)
 
 fvPatchField (const fvPatch &, const DimensionedField< Type, volMesh > &, const Field< Type > &)
 
 fvPatchField (const fvPatch &, const DimensionedField< Type, volMesh > &, const dictionary &, const bool valueRequired=true)
 
 fvPatchField (const fvPatchField< Type > &, const fvPatch &, const DimensionedField< Type, volMesh > &, const fvPatchFieldMapper &)
 
 fvPatchField (const fvPatchField< Type > &)
 
 fvPatchField (const fvPatchField< Type > &, const DimensionedField< Type, volMesh > &)
 
virtual tmp< fvPatchField< Type > > clone (const DimensionedField< Type, volMesh > &iF) const
 
Foam::tmp< Foam::fvPatchField< Type > > NewCalculatedType (const fvPatch &p)
 
Foam::tmp< Foam::fvPatchField< Type > > NewCalculatedType (const fvPatchField< Type2 > &pf)
 
virtual ~fvPatchField ()=default
 
bool useImplicit () const noexcept
 
bool useImplicit (bool on) noexcept
 
virtual bool coupled () const
 
const objectRegistrydb () const
 
const fvPatchpatch () const
 
const DimensionedField< Type, volMesh > & internalField () const
 
const Field< Type > & primitiveField () const
 
const wordpatchType () const
 
wordpatchType ()
 
bool updated () const
 
bool manipulatedMatrix () const
 
virtual void autoMap (const fvPatchFieldMapper &)
 
virtual void rmap (const fvPatchField< Type > &, const labelList &)
 
virtual tmp< Field< Type > > snGrad () const
 
virtual tmp< Field< Type > > snGrad (const scalarField &deltaCoeffs) const
 
virtual tmp< Field< Type > > patchInternalField () const
 
virtual void patchInternalField (Field< Type > &) const
 
virtual tmp< Field< Type > > patchNeighbourField () const
 
virtual void initEvaluate (const Pstream::commsTypes commsType=Pstream::commsTypes::blocking)
 
virtual void evaluate (const Pstream::commsTypes commsType=Pstream::commsTypes::blocking)
 
virtual tmp< Field< Type > > valueInternalCoeffs (const tmp< Field< scalar >> &) const
 
virtual tmp< Field< Type > > valueBoundaryCoeffs (const tmp< Field< scalar >> &) const
 
virtual tmp< Field< Type > > gradientInternalCoeffs (const scalarField &deltaCoeffs) const
 
virtual tmp< Field< Type > > gradientBoundaryCoeffs (const scalarField &deltaCoeffs) const
 
virtual void manipulateMatrix (fvMatrix< Type > &matrix)
 
virtual void manipulateMatrix (fvMatrix< Type > &matrix, const scalarField &weights)
 
virtual void manipulateMatrix (fvMatrix< Type > &matrix, const label iMatrix, const direction cmp)
 
void check (const fvPatchField< Type > &) const
 
virtual void operator= (const UList< Type > &)
 
virtual void operator= (const fvPatchField< Type > &)
 
virtual void operator= (const Type &)
 
virtual void operator+= (const fvPatchField< Type > &)
 
virtual void operator+= (const Field< Type > &)
 
virtual void operator+= (const Type &)
 
virtual void operator-= (const fvPatchField< Type > &)
 
virtual void operator-= (const Field< Type > &)
 
virtual void operator-= (const Type &)
 
virtual void operator== (const fvPatchField< Type > &)
 
virtual void operator== (const Field< Type > &)
 
virtual void operator== (const Type &)
 

Protected Member Functions

virtual void setMaster ()
 
virtual void createAveragingWeights ()
 
virtual epsilonWallFunctionFvPatchScalarFieldepsilonPatch (const label patchi)
 
virtual void calculateTurbulenceFields (const turbulenceModel &turbulence, scalarField &G0, scalarField &epsilon0)
 
virtual void calculate (const turbulenceModel &turbulence, const List< scalar > &cornerWeights, const fvPatch &patch, scalarField &G, scalarField &epsilon)
 
virtual label & master ()
 

Protected Attributes

const bool lowReCorrection_
 
bool initialised_
 
label master_
 
scalarField G_
 
scalarField epsilon_
 
List< List< scalar > > cornerWeights_
 

Static Protected Attributes

static scalar tolerance_ = 1e-5
 

Additional Inherited Members

- Public Types inherited from fvPatchField< scalar >
typedef fvPatch Patch
 
typedef calculatedFvPatchField< Type > Calculated
 
- Static Public Member Functions inherited from fvPatchField< scalar >
static tmp< fvPatchField< Type > > New (const word &, const fvPatch &, const DimensionedField< Type, volMesh > &)
 
static tmp< fvPatchField< Type > > New (const word &, const word &actualPatchType, const fvPatch &, const DimensionedField< Type, volMesh > &)
 
static tmp< fvPatchField< Type > > New (const fvPatchField< Type > &, const fvPatch &, const DimensionedField< Type, volMesh > &, const fvPatchFieldMapper &)
 
static tmp< fvPatchField< Type > > New (const fvPatch &, const DimensionedField< Type, volMesh > &, const dictionary &)
 
static tmp< fvPatchField< Type > > NewCalculatedType (const fvPatch &)
 
static tmp< fvPatchField< Type > > NewCalculatedType (const fvPatchField< Type2 > &)
 
static const wordcalculatedType ()
 
- Static Public Attributes inherited from fvPatchField< scalar >
static int disallowGenericFvPatchField
 

Detailed Description

This boundary condition provides a wall constraint on the turbulent kinetic energy dissipation rate, i.e. epsilon, and the turbulent kinetic energy production contribution, i.e. G, for low- and high-Reynolds number turbulence models.

Reference:

    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 entries (unmodifiable)
    type             epsilonWallFunction;

    // Optional entries (unmodifiable)
    lowReCorrection  false;
    blending         stepwise;
    n                2.0;

    // Optional (inherited) entries
    ...
}

where the entries mean:

Property Description Type Req'd Dflt
type Type name: epsilonWallFunction word yes -
lowReCorrection Flag: apply low-Re correction bool no false
blending Viscous/inertial sublayer blending method word no stepwise
n Binomial blending exponent scalar no 2.0

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 epsilon predictions for the viscous and inertial sublayers are blended according to the following expressions:

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

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

where

$ \epsilon $ = $\epsilon$ at $y^+$
$ \epsilon_{vis} $ = $\epsilon$ computed by viscous subl. assumptions
$ \epsilon_{log} $ = $\epsilon$ computed by inertial subl. assumptions
$ y^+ $ = estimated wall-normal height of the cell centre in wall units
$ y^+_{lam} $ = estimated intersection of the viscous and inertial sublayers

\[ \epsilon = max(\epsilon_{vis}, \epsilon_{log}) \]

\[ \epsilon = ((\epsilon_{vis})^n + (\epsilon_{log})^n)^{1/n} \]

where

$ n $ = Binomial blending exponent

\[ \epsilon = \epsilon_{vis} \exp[-\Gamma] +\epsilon_{log} \exp[-1/\Gamma] \]

where (PH:p. 193)

$ \Gamma_\epsilon $ = $\Gamma = 0.001 (y^+)^4 / (1.0 + y^+)$
$ \Gamma_G $ = $\Gamma = 0.01 (y^+)^4 / (1.0 + 5.0 y^+)$
$ \Gamma_\epsilon $ = Blending expression for $\epsilon$
$ \Gamma_G $ = Blending expression for $G$

G predictions for the viscous and inertial sublayers are blended in a stepwise manner, and G below $y^+_{lam}$ (i.e. in the viscous sublayer) is presumed to be zero.

Note
  • The coefficients Cmu, kappa, and E are obtained from the specified nutWallFunction in order to ensure that each patch possesses the same set of values for these coefficients.
  • lowReCorrection operates with only stepwise blending treatment to ensure the backward compatibility.
  • If lowReCorrection is on, stepwise blending treatment is fully active.
  • If lowReCorrection is off, only the inertial sublayer prediction is used in the wall function, hence high-Re mode operation.
See also
Source files

Definition at line 255 of file epsilonWallFunctionFvPatchScalarField.H.

Constructor & Destructor Documentation

◆ epsilonWallFunctionFvPatchScalarField() [1/5]

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

Definition at line 298 of file epsilonWallFunctionFvPatchScalarField.C.

Referenced by epsilonWallFunctionFvPatchScalarField::clone().

Here is the caller graph for this function:

◆ epsilonWallFunctionFvPatchScalarField() [2/5]

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

Definition at line 338 of file epsilonWallFunctionFvPatchScalarField.C.

References Foam::operator==().

Here is the call graph for this function:

◆ epsilonWallFunctionFvPatchScalarField() [3/5]

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

Definition at line 317 of file epsilonWallFunctionFvPatchScalarField.C.

◆ epsilonWallFunctionFvPatchScalarField() [4/5]

Definition at line 377 of file epsilonWallFunctionFvPatchScalarField.C.

◆ epsilonWallFunctionFvPatchScalarField() [5/5]

Definition at line 395 of file epsilonWallFunctionFvPatchScalarField.C.

◆ ~epsilonWallFunctionFvPatchScalarField()

virtual ~epsilonWallFunctionFvPatchScalarField ( )
virtualdefault

Member Function Documentation

◆ setMaster()

void setMaster ( )
protectedvirtual

◆ createAveragingWeights()

void createAveragingWeights ( )
protectedvirtual

◆ epsilonPatch()

Foam::epsilonWallFunctionFvPatchScalarField & epsilonPatch ( const label  patchi)
protectedvirtual

Definition at line 138 of file epsilonWallFunctionFvPatchScalarField.C.

References epsilon.

Referenced by epsilonWallFunctionFvPatchScalarField::setMaster().

Here is the caller graph for this function:

◆ calculateTurbulenceFields()

void calculateTurbulenceFields ( const turbulenceModel turbulence,
scalarField G0,
scalarField epsilon0 
)
protectedvirtual

◆ calculate()

void calculate ( const turbulenceModel turbulence,
const List< scalar > &  cornerWeights,
const fvPatch patch,
scalarField G,
scalarField epsilon 
)
protectedvirtual

◆ master()

virtual label& master ( )
inlineprotectedvirtual

Definition at line 345 of file epsilonWallFunctionFvPatchScalarField.H.

References epsilonWallFunctionFvPatchScalarField::master_.

Referenced by epsilonWallFunctionFvPatchScalarField::setMaster().

Here is the caller graph for this function:

◆ TypeName()

TypeName ( "epsilonWallFunction"  )

◆ clone() [1/2]

virtual tmp<fvPatchScalarField> clone ( ) const
inlinevirtual

◆ clone() [2/2]

virtual tmp<fvPatchScalarField> clone ( const DimensionedField< scalar, volMesh > &  iF) const
inlinevirtual

◆ G()

Foam::scalarField & G ( bool  init = false)

Definition at line 415 of file epsilonWallFunctionFvPatchScalarField.C.

References init(), and Foam::foamVersion::patch.

Here is the call graph for this function:

◆ epsilon()

Foam::scalarField & epsilon ( bool  init = false)

Definition at line 434 of file epsilonWallFunctionFvPatchScalarField.C.

References init(), and Foam::foamVersion::patch.

Referenced by epsilonWallFunctionFvPatchScalarField::setMaster().

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

◆ updateCoeffs()

void updateCoeffs ( )
virtual

◆ updateWeightedCoeffs()

void updateWeightedCoeffs ( const scalarField weights)
virtual

◆ manipulateMatrix() [1/2]

void manipulateMatrix ( fvMatrix< scalar > &  matrix)
virtual

Definition at line 555 of file epsilonWallFunctionFvPatchScalarField.C.

References fvPatchField::manipulateMatrix(), Foam::foamVersion::patch, and fvMatrix::setValues().

Here is the call graph for this function:

◆ manipulateMatrix() [2/2]

void manipulateMatrix ( fvMatrix< scalar > &  matrix,
const scalarField weights 
)
virtual

◆ write()

void write ( Ostream os) const
virtual

Reimplemented from fixedValueFvPatchField< scalar >.

Reimplemented in atmEpsilonWallFunctionFvPatchScalarField.

Definition at line 614 of file epsilonWallFunctionFvPatchScalarField.C.

References os(), fixedValueFvPatchField< Type >::write(), and Ostream::writeEntry().

Here is the call graph for this function:

Member Data Documentation

◆ tolerance_

Foam::scalar tolerance_ = 1e-5
staticprotected

Definition at line 289 of file epsilonWallFunctionFvPatchScalarField.H.

◆ lowReCorrection_

const bool lowReCorrection_
protected

Definition at line 292 of file epsilonWallFunctionFvPatchScalarField.H.

◆ initialised_

bool initialised_
protected

Definition at line 295 of file epsilonWallFunctionFvPatchScalarField.H.

◆ master_

label master_
protected

◆ G_

scalarField G_
protected

Definition at line 301 of file epsilonWallFunctionFvPatchScalarField.H.

◆ epsilon_

scalarField epsilon_
protected

Definition at line 304 of file epsilonWallFunctionFvPatchScalarField.H.

◆ cornerWeights_

List<List<scalar> > cornerWeights_
protected

Definition at line 307 of file epsilonWallFunctionFvPatchScalarField.H.


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