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

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

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

Public Member Functions

 TypeName ("omegaWallFunction")
 
 omegaWallFunctionFvPatchScalarField (const fvPatch &, const DimensionedField< scalar, volMesh > &)
 
 omegaWallFunctionFvPatchScalarField (const fvPatch &, const DimensionedField< scalar, volMesh > &, const dictionary &)
 
 omegaWallFunctionFvPatchScalarField (const omegaWallFunctionFvPatchScalarField &, const fvPatch &, const DimensionedField< scalar, volMesh > &, const fvPatchFieldMapper &)
 
 omegaWallFunctionFvPatchScalarField (const omegaWallFunctionFvPatchScalarField &)
 
virtual tmp< fvPatchScalarFieldclone () const
 
 omegaWallFunctionFvPatchScalarField (const omegaWallFunctionFvPatchScalarField &, const DimensionedField< scalar, volMesh > &)
 
virtual tmp< fvPatchScalarFieldclone (const DimensionedField< scalar, volMesh > &iF) const
 
virtual ~omegaWallFunctionFvPatchScalarField ()=default
 
scalarFieldG (bool init=false)
 
scalarFieldomega (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 omegaWallFunctionFvPatchScalarFieldomegaPatch (const label patchi)
 
virtual void calculateTurbulenceFields (const turbulenceModel &turbulence, scalarField &G0, scalarField &omega0)
 
virtual void calculate (const turbulenceModel &turbulence, const List< scalar > &cornerWeights, const fvPatch &patch, scalarField &G, scalarField &omega)
 
virtual label & master ()
 

Protected Attributes

bool blended_
 
bool initialised_
 
label master_
 
scalar beta1_
 
scalarField G_
 
scalarField omega_
 
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 specific dissipation rate, i.e. omega, 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

    Tanh blending of the viscous and inertial sublayers (tag:KAS):
        Knopp, T., Alrutz, T., & Schwamborn, D. (2006).
        A grid and flow adaptive wall-function method for RANS
        turbulence modelling.
        Journal of Computational Physics, 220(1), 19-40.
        DOI:10.1016/j.jcp.2006.05.003
Usage
Example of the boundary condition specification:
<patchName>
{
    // Mandatory entries (unmodifiable)
    type            omegaWallFunction;

    // Optional entries (unmodifiable)
    beta1           0.075;
    blending        binomial2;
    n               2.0;

    // Optional (inherited) entries
    ...
}
Property Description Type Req'd Dflt
type Type name: omegaWallFunction word yes -
beta1 Model coefficient scalar no 0.075
blending Viscous/inertial sublayer blending method word no binomial2
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)
      binomial2   | Binomial blending (smooth) n = 2
      binomial    | Binomial blending (smooth)
      exponential | Exponential blending (smooth)
      tanh        | Tanh blending (smooth)

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

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

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

where

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

\[ \omega = max(\omega_{vis}, \omega_{log}) \]

\[ \omega = \sqrt{(\omega_{vis})^2 + (\omega_{log})^2} \]

\[ \omega = ((\omega_{vis})^n + (\omega_{log})^n)^{1/n} \]

where

$ n $ = Binomial blending exponent

\[ \omega = \omega_{vis} \exp[-\Gamma] + \omega_{log} \exp[-1/\Gamma] \]

where (PH:Eq. 31)

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

\[ \omega = \phi \omega_{b1} + (1 - \phi)\omega_{b2} \]

where

$ \phi $ = $tanh((y^+/10)^4)$
$ \omega_{b1} $ = $\omega_{vis} + \omega_{log}$
$ \omega_{b2} $ = $(\omega_{vis}^{1.2} + \omega_{log}^1.2)^{1/1.2}$

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.
  • The reason why binomial2 and binomial blending methods exist at the same time is to ensure the bitwise regression with the previous versions since binomial2 and binomial with n=2 will yield slightly different output due to the miniscule differences in the implementation of the basic functions (i.e. pow, sqrt, sqr).
See also
Source files

Definition at line 287 of file omegaWallFunctionFvPatchScalarField.H.

Constructor & Destructor Documentation

◆ omegaWallFunctionFvPatchScalarField() [1/5]

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

Definition at line 308 of file omegaWallFunctionFvPatchScalarField.C.

Referenced by omegaWallFunctionFvPatchScalarField::clone().

Here is the caller graph for this function:

◆ omegaWallFunctionFvPatchScalarField() [2/5]

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

Definition at line 346 of file omegaWallFunctionFvPatchScalarField.C.

References dict, Foam::endl(), dictionary::found(), dictionary::get(), IOWarningInFunction, Foam::nl, and Foam::operator==().

Here is the call graph for this function:

◆ omegaWallFunctionFvPatchScalarField() [3/5]

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

Definition at line 326 of file omegaWallFunctionFvPatchScalarField.C.

◆ omegaWallFunctionFvPatchScalarField() [4/5]

Definition at line 407 of file omegaWallFunctionFvPatchScalarField.C.

◆ omegaWallFunctionFvPatchScalarField() [5/5]

Definition at line 424 of file omegaWallFunctionFvPatchScalarField.C.

◆ ~omegaWallFunctionFvPatchScalarField()

virtual ~omegaWallFunctionFvPatchScalarField ( )
virtualdefault

Member Function Documentation

◆ setMaster()

void setMaster ( )
protectedvirtual

◆ createAveragingWeights()

void createAveragingWeights ( )
protectedvirtual

◆ omegaPatch()

Foam::omegaWallFunctionFvPatchScalarField & omegaPatch ( const label  patchi)
protectedvirtual

Definition at line 138 of file omegaWallFunctionFvPatchScalarField.C.

References GeometricField::boundaryField().

Referenced by omegaWallFunctionFvPatchScalarField::setMaster().

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

◆ calculateTurbulenceFields()

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

◆ calculate()

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

◆ master()

virtual label& master ( )
inlineprotectedvirtual

Definition at line 383 of file omegaWallFunctionFvPatchScalarField.H.

References omegaWallFunctionFvPatchScalarField::master_.

Referenced by omegaWallFunctionFvPatchScalarField::setMaster().

Here is the caller graph for this function:

◆ TypeName()

TypeName ( "omegaWallFunction"  )

◆ clone() [1/2]

virtual tmp<fvPatchScalarField> clone ( ) const
inlinevirtual

Reimplemented in atmOmegaWallFunctionFvPatchScalarField.

Definition at line 430 of file omegaWallFunctionFvPatchScalarField.H.

References omegaWallFunctionFvPatchScalarField::omegaWallFunctionFvPatchScalarField().

Here is the call graph for this function:

◆ clone() [2/2]

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

Reimplemented from fixedValueFvPatchField< scalar >.

Reimplemented in atmOmegaWallFunctionFvPatchScalarField.

Definition at line 447 of file omegaWallFunctionFvPatchScalarField.H.

References omegaWallFunctionFvPatchScalarField::omegaWallFunctionFvPatchScalarField().

Here is the call graph for this function:

◆ G()

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

Definition at line 444 of file omegaWallFunctionFvPatchScalarField.C.

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

Here is the call graph for this function:

◆ omega()

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

Definition at line 463 of file omegaWallFunctionFvPatchScalarField.C.

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

Referenced by omegaWallFunctionFvPatchScalarField::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 584 of file omegaWallFunctionFvPatchScalarField.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 atmOmegaWallFunctionFvPatchScalarField.

Definition at line 643 of file omegaWallFunctionFvPatchScalarField.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 323 of file omegaWallFunctionFvPatchScalarField.H.

◆ blended_

bool blended_
protected

Definition at line 327 of file omegaWallFunctionFvPatchScalarField.H.

◆ initialised_

bool initialised_
protected

Definition at line 330 of file omegaWallFunctionFvPatchScalarField.H.

◆ master_

label master_
protected

◆ beta1_

scalar beta1_
protected

Definition at line 336 of file omegaWallFunctionFvPatchScalarField.H.

◆ G_

scalarField G_
protected

Definition at line 339 of file omegaWallFunctionFvPatchScalarField.H.

◆ omega_

scalarField omega_
protected

Definition at line 342 of file omegaWallFunctionFvPatchScalarField.H.

◆ cornerWeights_

List<List<scalar> > cornerWeights_
protected

Definition at line 345 of file omegaWallFunctionFvPatchScalarField.H.


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