nutUSpaldingWallFunctionFvPatchScalarField Class Reference

This boundary condition provides a wall constraint on the turbulent viscosity, i.e. nut, based on velocity, i.e. U. Using Spalding's law gives a continuous nut profile to the wall. More...

Inheritance diagram for nutUSpaldingWallFunctionFvPatchScalarField:
[legend]
Collaboration diagram for nutUSpaldingWallFunctionFvPatchScalarField:
[legend]

Public Member Functions

 TypeName ("nutUSpaldingWallFunction")
 Runtime type information. More...
 
 nutUSpaldingWallFunctionFvPatchScalarField (const fvPatch &, const DimensionedField< scalar, volMesh > &)
 Construct from patch and internal field. More...
 
 nutUSpaldingWallFunctionFvPatchScalarField (const fvPatch &, const DimensionedField< scalar, volMesh > &, const dictionary &)
 Construct from patch, internal field and dictionary. More...
 
 nutUSpaldingWallFunctionFvPatchScalarField (const nutUSpaldingWallFunctionFvPatchScalarField &, const fvPatch &, const DimensionedField< scalar, volMesh > &, const fvPatchFieldMapper &)
 
 nutUSpaldingWallFunctionFvPatchScalarField (const nutUSpaldingWallFunctionFvPatchScalarField &)
 Construct as copy. More...
 
virtual tmp< fvPatchScalarFieldclone () const
 Construct and return a clone. More...
 
 nutUSpaldingWallFunctionFvPatchScalarField (const nutUSpaldingWallFunctionFvPatchScalarField &, const DimensionedField< scalar, volMesh > &)
 Construct as copy setting internal field reference. More...
 
virtual tmp< fvPatchScalarFieldclone (const DimensionedField< scalar, volMesh > &iF) const
 Construct and return a clone setting internal field reference. More...
 
virtual ~nutUSpaldingWallFunctionFvPatchScalarField ()
 Destructor. More...
 
virtual tmp< scalarFieldyPlus () const
 Calculate and return the yPlus at the boundary. More...
 
virtual void write (Ostream &os) const
 Write. More...
 
- Public Member Functions inherited from nutWallFunctionFvPatchScalarField
 TypeName ("nutWallFunction")
 Runtime type information. More...
 
 nutWallFunctionFvPatchScalarField (const fvPatch &, const DimensionedField< scalar, volMesh > &)
 Construct from patch and internal field. More...
 
 nutWallFunctionFvPatchScalarField (const fvPatch &, const DimensionedField< scalar, volMesh > &, const dictionary &)
 Construct from patch, internal field and dictionary. More...
 
 nutWallFunctionFvPatchScalarField (const nutWallFunctionFvPatchScalarField &, const fvPatch &, const DimensionedField< scalar, volMesh > &, const fvPatchFieldMapper &)
 
 nutWallFunctionFvPatchScalarField (const nutWallFunctionFvPatchScalarField &)
 Construct as copy. More...
 
 nutWallFunctionFvPatchScalarField (const nutWallFunctionFvPatchScalarField &, const DimensionedField< scalar, volMesh > &)
 Construct as copy setting internal field reference. More...
 
scalar Cmu () const
 Return Cmu. More...
 
scalar kappa () const
 Return kappa. More...
 
scalar E () const
 Return E. More...
 
scalar yPlusLam () const
 Return the estimated y+ at the two-sublayer intersection. More...
 
scalar blend (const scalar nutVis, const scalar nutLog, const scalar yPlus) const
 Return the blended nut according to the chosen blending treatment. More...
 
virtual void updateCoeffs ()
 Update the coefficients associated with the patch field. More...
 

Protected Member Functions

virtual tmp< scalarFieldcalcNut () const
 Uncomment in case of intrumentation. More...
 
virtual tmp< scalarFieldcalcUTau (const scalarField &magGradU) const
 Calculate the friction velocity. More...
 
virtual tmp< scalarFieldcalcUTau (const scalarField &magGradU, const label maxIter, scalarField &err) const
 
virtual void writeLocalEntries (Ostream &) const
 Write local wall function variables. More...
 
- Protected Member Functions inherited from nutWallFunctionFvPatchScalarField
virtual const volVectorFieldU (const turbulenceModel &turb) const
 
virtual void checkType ()
 Check the type of the patch. More...
 

Protected Attributes

const label maxIter_
 Max iterations in calcNut. More...
 
const scalar tolerance_
 Convergence tolerance. More...
 
- Protected Attributes inherited from nutWallFunctionFvPatchScalarField
enum blendingType blending_
 Blending treatment (default = blendingType::STEPWISE) More...
 
const scalar n_
 
word UName_
 Name of velocity field. More...
 
scalar Cmu_
 Empirical model coefficient. More...
 
scalar kappa_
 von Kármán constant More...
 
scalar E_
 Wall roughness parameter. More...
 
scalar yPlusLam_
 

Additional Inherited Members

- Static Public Member Functions inherited from nutWallFunctionFvPatchScalarField
static const nutWallFunctionFvPatchScalarFieldnutw (const turbulenceModel &turbModel, const label patchi)
 Return the nut patchField for the given wall patch. More...
 
static scalar yPlusLam (const scalar kappa, const scalar E)
 Estimate the y+ at the intersection of the two sublayers. More...
 
- Protected Types inherited from nutWallFunctionFvPatchScalarField
enum  blendingType { STEPWISE, MAX, BINOMIAL, EXPONENTIAL }
 Options for the blending treatment of viscous and inertial sublayers. More...
 
- Static Protected Attributes inherited from nutWallFunctionFvPatchScalarField
static const Enum< blendingTypeblendingTypeNames
 Names for blendingType. More...
 

Detailed Description

This boundary condition provides a wall constraint on the turbulent viscosity, i.e. nut, based on velocity, i.e. U. Using Spalding's law gives a continuous nut profile to the wall.

\[ y^+ = u^+ + \frac{1}{E} \left[exp(\kappa u^+) - 1 - \kappa u^+\, - 0.5 (\kappa u^+)^2 - \frac{1}{6} (\kappa u^+)^3\right] \]

where

\( y^+ \) = Non-dimensional position
\( u^+ \) = Non-dimensional velocity
\( \kappa \) = von Kármán constant
Usage
Example of the boundary condition specification:
<patchName>
{
    // Mandatory entries (unmodifiable)
    type            nutUSpaldingWallFunction;

    // Optional entries (unmodifiable)
    maxIter         10;
    tolerance       0.0001;

    // Optional (inherited) entries
    ...
}

where the entries mean:

Property Description Type Req'd Dflt
type Type name: nutUBlendedWallFunction word yes -
maxIter Number of Newton-Raphson iterations label no 10
tolerance Convergence tolerance scalar no 0.0001

The inherited entries are elaborated in:

Note
  • Suffers from non-exact restart since correctNut() (called through turbulence->validate) returns a slightly different value every time it is called. This is since the seed for the Newton-Raphson iteration uses the current value of *this (=nut ).
  • This can be avoided by overriding the tolerance. This also switches on a pre-detection whether the current nut already satisfies the turbulence conditions and if so does not change it at all. This means that the nut only changes if it 'has really changed'. This probably should be used with a tight tolerance, to make sure to kick every iteration, e.g. maxIter 100; tolerance 1e-7;
Source files

Definition at line 145 of file nutUSpaldingWallFunctionFvPatchScalarField.H.

Constructor & Destructor Documentation

◆ nutUSpaldingWallFunctionFvPatchScalarField() [1/5]

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

Construct from patch and internal field.

Definition at line 209 of file nutUSpaldingWallFunctionFvPatchScalarField.C.

Referenced by nutUSpaldingWallFunctionFvPatchScalarField::clone().

Here is the caller graph for this function:

◆ nutUSpaldingWallFunctionFvPatchScalarField() [2/5]

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

Construct from patch, internal field and dictionary.

Definition at line 245 of file nutUSpaldingWallFunctionFvPatchScalarField.C.

◆ nutUSpaldingWallFunctionFvPatchScalarField() [3/5]

Construct by mapping given nutUSpaldingWallFunctionFvPatchScalarField onto a new patch

Definition at line 226 of file nutUSpaldingWallFunctionFvPatchScalarField.C.

◆ nutUSpaldingWallFunctionFvPatchScalarField() [4/5]

◆ nutUSpaldingWallFunctionFvPatchScalarField() [5/5]

Construct as copy setting internal field reference.

Definition at line 279 of file nutUSpaldingWallFunctionFvPatchScalarField.C.

◆ ~nutUSpaldingWallFunctionFvPatchScalarField()

Member Function Documentation

◆ calcNut()

Foam::tmp< Foam::scalarField > calcNut ( ) const
protectedvirtual

Uncomment in case of intrumentation.

mutable uint64_t invocations_; mutable uint64_t nontrivial_; mutable uint64_t nonconvergence_; mutable uint64_t iterations_; Calculate the turbulent viscosity

Implements nutWallFunctionFvPatchScalarField.

Definition at line 39 of file nutUSpaldingWallFunctionFvPatchScalarField.C.

References GeometricField< Type, PatchField, GeoMesh >::boundaryField(), nutUSpaldingWallFunctionFvPatchScalarField::calcUTau(), forAll, Foam::constant::atomic::group, IOobject::groupName(), Foam::mag(), Foam::max(), turbulenceModel::nu(), nutWallFunctionFvPatchScalarField::nutw(), Foam::foamVersion::patch, turbulenceModel::propertiesName, tmp< T >::ref(), Foam::sqr(), nutUSpaldingWallFunctionFvPatchScalarField::tolerance_, and nutWallFunctionFvPatchScalarField::U().

Here is the call graph for this function:

◆ calcUTau() [1/2]

Foam::tmp< Foam::scalarField > calcUTau ( const scalarField magGradU) const
protectedvirtual

Calculate the friction velocity.

Definition at line 94 of file nutUSpaldingWallFunctionFvPatchScalarField.C.

Referenced by nutUSpaldingWallFunctionFvPatchScalarField::calcNut().

Here is the caller graph for this function:

◆ calcUTau() [2/2]

Foam::tmp< Foam::scalarField > calcUTau ( const scalarField magGradU,
const label  maxIter,
scalarField err 
) const
protectedvirtual

Calculate the friction velocity and number of iterations for convergence

Definition at line 105 of file nutUSpaldingWallFunctionFvPatchScalarField.C.

References Foam::exp(), f(), forAll, Foam::constant::atomic::group, IOobject::groupName(), Foam::mag(), magUp, Foam::max(), Foam::min(), turbulenceModel::nu(), Foam::foamVersion::patch, fvPatchField< Type >::patchInternalField(), turbulenceModel::propertiesName, tmp< T >::ref(), Foam::sqr(), Foam::sqrt(), U, uTau, y, turbulenceModel::y(), and Foam::Zero.

Here is the call graph for this function:

◆ writeLocalEntries()

void writeLocalEntries ( Ostream os) const
protectedvirtual

Write local wall function variables.

Reimplemented from nutWallFunctionFvPatchScalarField.

Definition at line 194 of file nutUSpaldingWallFunctionFvPatchScalarField.C.

References os(), Ostream::writeEntryIfDifferent(), and nutWallFunctionFvPatchScalarField::writeLocalEntries().

Here is the call graph for this function:

◆ TypeName()

TypeName ( "nutUSpaldingWallFunction"  )

Runtime type information.

◆ clone() [1/2]

virtual tmp<fvPatchScalarField> clone ( ) const
inlinevirtual

Construct and return a clone.

Definition at line 228 of file nutUSpaldingWallFunctionFvPatchScalarField.H.

References nutUSpaldingWallFunctionFvPatchScalarField::nutUSpaldingWallFunctionFvPatchScalarField().

Here is the call graph for this function:

◆ clone() [2/2]

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

Construct and return a clone setting internal field reference.

Definition at line 245 of file nutUSpaldingWallFunctionFvPatchScalarField.H.

References nutUSpaldingWallFunctionFvPatchScalarField::nutUSpaldingWallFunctionFvPatchScalarField().

Here is the call graph for this function:

◆ yPlus()

Foam::tmp< Foam::scalarField > yPlus ( ) const
virtual

Calculate and return the yPlus at the boundary.

Implements nutWallFunctionFvPatchScalarField.

Definition at line 318 of file nutUSpaldingWallFunctionFvPatchScalarField.C.

References Foam::constant::atomic::group, IOobject::groupName(), Foam::mag(), turbulenceModel::nu(), Foam::foamVersion::patch, turbulenceModel::propertiesName, fvPatchField< Type >::snGrad(), U, y, and turbulenceModel::y().

Here is the call graph for this function:

◆ write()

void write ( Ostream os) const
virtual

Write.

Reimplemented from nutWallFunctionFvPatchScalarField.

Definition at line 340 of file nutUSpaldingWallFunctionFvPatchScalarField.C.

References os(), and fvPatchField< Type >::write().

Here is the call graph for this function:

Member Data Documentation

◆ maxIter_

const label maxIter_
protected

Max iterations in calcNut.

Definition at line 154 of file nutUSpaldingWallFunctionFvPatchScalarField.H.

◆ tolerance_

const scalar tolerance_
protected

Convergence tolerance.

Definition at line 157 of file nutUSpaldingWallFunctionFvPatchScalarField.H.

Referenced by nutUSpaldingWallFunctionFvPatchScalarField::calcNut().


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