Merge branch 'feature-wall-functions' into 'develop'

ENH: New wall-function blending approaches

See merge request Development/openfoam!350
This commit is contained in:
Andrew Heather 2020-06-10 10:24:25 +01:00
commit ff568aa67f
20 changed files with 1076 additions and 355 deletions

View File

@ -6,7 +6,7 @@
\\/ M anipulation |
-------------------------------------------------------------------------------
Copyright (C) 2011-2019 OpenFOAM Foundation
Copyright (C) 2017-2019 OpenCFD Ltd.
Copyright (C) 2017-2020 OpenCFD Ltd.
-------------------------------------------------------------------------------
License
This file is part of OpenFOAM.
@ -32,9 +32,20 @@ License
#include "fvMatrix.H"
#include "addToRunTimeSelectionTable.H"
// * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
const Foam::Enum
<
Foam::epsilonWallFunctionFvPatchScalarField::blendingType
>
Foam::epsilonWallFunctionFvPatchScalarField::blendingTypeNames
({
{ blendingType::STEPWISE , "stepwise" },
{ blendingType::MAX , "max" },
{ blendingType::BINOMIAL , "binomial" },
{ blendingType::EXPONENTIAL, "exponential" }
});
Foam::scalar Foam::epsilonWallFunctionFvPatchScalarField::tolerance_ = 1e-5;
// * * * * * * * * * * * * * Protected Member Functions * * * * * * * * * * //
@ -106,17 +117,17 @@ void Foam::epsilonWallFunctionFvPatchScalarField::createAveragingWeights()
epsilonPatches.append(patchi);
const labelUList& faceCells = bf[patchi].patch().faceCells();
forAll(faceCells, i)
for (const auto& faceCell : faceCells)
{
weights[faceCells[i]]++;
++weights[faceCell];
}
}
}
cornerWeights_.setSize(bf.size());
forAll(epsilonPatches, i)
for (const auto& patchi : epsilonPatches)
{
label patchi = epsilonPatches[i];
const fvPatchScalarField& wf = weights.boundaryField()[patchi];
cornerWeights_[patchi] = 1.0/wf.patchInternalField();
}
@ -217,24 +228,71 @@ void Foam::epsilonWallFunctionFvPatchScalarField::calculate
const scalar w = cornerWeights[facei];
// Default high-Re form
scalar epsilonc = w*Cmu75*pow(k[celli], 1.5)/(nutw.kappa()*y[facei]);
scalar Gc =
w
*(nutw[facei] + nuw[facei])
*magGradUw[facei]
*Cmu25*sqrt(k[celli])
/(nutw.kappa()*y[facei]);
scalar epsilonBlended = 0.0;
if (lowReCorrection_ && yPlus < nutw.yPlusLam())
// Contribution from the viscous sublayer
const scalar epsilonVis = w*2.0*k[celli]*nuw[facei]/sqr(y[facei]);
// Contribution from the inertial sublayer
const scalar epsilonLog =
w*Cmu75*pow(k[celli], 1.5)/(nutw.kappa()*y[facei]);
switch (blending_)
{
epsilonc = w*2.0*k[celli]*nuw[facei]/sqr(y[facei]);
Gc = 0;
case blendingType::STEPWISE:
{
if (lowReCorrection_ && yPlus < nutw.yPlusLam())
{
epsilonBlended = epsilonVis;
}
else
{
epsilonBlended = epsilonLog;
}
break;
}
case blendingType::MAX:
{
// (PH:Eq. 27)
epsilonBlended = max(epsilonVis, epsilonLog);
break;
}
case blendingType::BINOMIAL:
{
// (ME:Eqs. 15-16)
epsilonBlended =
pow
(
pow(epsilonVis, n_) + pow(epsilonLog, n_),
1.0/n_
);
break;
}
case blendingType::EXPONENTIAL:
{
// (PH:p. 193)
const scalar Gamma = 0.001*pow4(yPlus)/(1.0 + yPlus);
const scalar invGamma = 1.0/(Gamma + ROOTVSMALL);
epsilonBlended =
epsilonVis*exp(-Gamma) + epsilonLog*exp(-invGamma);
break;
}
}
epsilon0[celli] += epsilonc;
epsilon0[celli] += epsilonBlended;
G0[celli] += Gc;
if (!(lowReCorrection_ && yPlus < nutw.yPlusLam()))
{
G0[celli] +=
w
*(nutw[facei] + nuw[facei])
*magGradUw[facei]
*Cmu25*sqrt(k[celli])
/(nutw.kappa()*y[facei]);
}
}
}
@ -249,6 +307,8 @@ epsilonWallFunctionFvPatchScalarField
)
:
fixedValueFvPatchField<scalar>(p, iF),
blending_(blendingType::STEPWISE),
n_(2.0),
lowReCorrection_(false),
initialised_(false),
master_(-1),
@ -268,6 +328,8 @@ epsilonWallFunctionFvPatchScalarField
)
:
fixedValueFvPatchField<scalar>(ptf, p, iF, mapper),
blending_(ptf.blending_),
n_(ptf.n_),
lowReCorrection_(ptf.lowReCorrection_),
initialised_(false),
master_(-1),
@ -286,6 +348,24 @@ epsilonWallFunctionFvPatchScalarField
)
:
fixedValueFvPatchField<scalar>(p, iF, dict),
blending_
(
blendingTypeNames.getOrDefault
(
"blending",
dict,
blendingType::STEPWISE
)
),
n_
(
dict.getCheckOrDefault<scalar>
(
"n",
2.0,
scalarMinMax::ge(0)
)
),
lowReCorrection_(dict.getOrDefault("lowReCorrection", false)),
initialised_(false),
master_(-1),
@ -305,6 +385,8 @@ epsilonWallFunctionFvPatchScalarField
)
:
fixedValueFvPatchField<scalar>(ewfpsf),
blending_(ewfpsf.blending_),
n_(ewfpsf.n_),
lowReCorrection_(ewfpsf.lowReCorrection_),
initialised_(false),
master_(-1),
@ -322,6 +404,8 @@ epsilonWallFunctionFvPatchScalarField
)
:
fixedValueFvPatchField<scalar>(ewfpsf, iF),
blending_(ewfpsf.blending_),
n_(ewfpsf.n_),
lowReCorrection_(ewfpsf.lowReCorrection_),
initialised_(false),
master_(-1),
@ -459,7 +543,7 @@ void Foam::epsilonWallFunctionFvPatchScalarField::updateWeightedCoeffs
{
const scalar w = weights[facei];
if (tolerance_ < w)
if (w > tolerance_)
{
const label celli = patch().faceCells()[facei];
@ -509,7 +593,7 @@ void Foam::epsilonWallFunctionFvPatchScalarField::manipulateMatrix
forAll(weights, facei)
{
// Only set the values if the weights are > tolerance
if (tolerance_ < weights[facei])
if (weights[facei] > tolerance_)
{
const label celli = faceCells[facei];
@ -538,6 +622,9 @@ void Foam::epsilonWallFunctionFvPatchScalarField::write
) const
{
os.writeEntry("lowReCorrection", lowReCorrection_);
os.writeKeyword("blending") << blendingTypeNames[blending_]
<< token::END_STATEMENT << nl;
os.writeEntry("n", n_);
fixedValueFvPatchField<scalar>::write(os);
}

View File

@ -6,7 +6,7 @@
\\/ M anipulation |
-------------------------------------------------------------------------------
Copyright (C) 2011-2016, 2019 OpenFOAM Foundation
Copyright (C) 2019 OpenCFD Ltd.
Copyright (C) 2019-2020 OpenCFD Ltd.
-------------------------------------------------------------------------------
License
This file is part of OpenFOAM.
@ -31,52 +31,138 @@ Group
grpWallFunctions
Description
This boundary condition provides a wall constraint on turbulent kinetic
energy dissipation rate, i.e. \c epsilon, for low- and high-Reynolds number
turbulence models.
This boundary condition provides a wall constraint on the turbulent kinetic
energy dissipation rate, i.e. \c epsilon, and the turbulent kinetic
energy production contribution, i.e. \c G, for low- and high-Reynolds
number turbulence models.
The condition can be applied to wall boundaries for which it
- calculates \c epsilon and \c G
- specifies the near-wall \c epsilon value
Reference:
\verbatim
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.
where
\vartable
epsilon | turbulent kinetic energy dissipation rate field
G | turbulent kinetic energy production field (divergence-free)
\endvartable
The low-Re correction is activated by setting the entry
\c lowReCorrection to 'on'; in this mode the model switches between
viscous and turbulent functions based on the viscous-to-turbulent
\c y+ value derived from the \c kappa and \c E.
When the \c lowReCorrection is inactive, the wall function operates
in high-Re mode.
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
\endverbatim
Usage
\table
Property | Description | Required | Default value
lowReCorrection | Low-Re correction active | no | off
\endtable
Example of the boundary condition specification:
\verbatim
<patchName>
{
// Mandatory entries
type epsilonWallFunction;
// Mandatory entries (unmodifiable)
type epsilonWallFunction;
// Optional entries
// Optional entries (unmodifiable)
lowReCorrection false;
blending stepwise;
n 2.0;
// Optional (inherited) entries
...
}
\endverbatim
where the entries mean:
\table
Property | Description | Type | Req'd | Dflt
type | Type name: epsilonWallFunction | word | yes | -
lowReCorrection | Flag for low-Re correction | bool | no | false
blending | Viscous/inertial sublayer blending method <!--
--> | word | no | stepwise
n | Binomial blending exponent | scalar | no | 2.0
\endtable
The inherited entries are elaborated in:
- \link fixedValueFvPatchField.H \endlink
- \link nutWallFunctionFvPatchScalarField.H \endlink
Options for the \c blending entry:
\verbatim
stepwise | Stepwise switch (discontinuous)
max | Maximum value switch (discontinuous)
binomial | Binomial blending (smooth)
exponential | Exponential blending (smooth)
\endverbatim
wherein \c epsilon predictions for the viscous and inertial sublayers are
blended according to the following expressions:
- \c stepwise (default):
\f[
\epsilon = \epsilon_{vis} \qquad if \quad y^+ < y^+_{lam}
\f]
\f[
\epsilon = \epsilon_{log} \qquad if \quad y^+ >= y^+_{lam}
\f]
where
\vartable
\epsilon | \f$\epsilon\f$ at \f$y^+\f$
\epsilon_{vis} | \f$\epsilon\f$ computed by viscous subl. assumptions
\epsilon_{log} | \f$\epsilon\f$ 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
\endvartable
- \c max (PH:Eq. 27):
\f[
\epsilon = max(\epsilon_{vis}, \epsilon_{log})
\f]
- \c binomial (ME:Eqs. 15-16):
\f[
\epsilon = ((\epsilon_{vis})^n + (\epsilon_{log})^n)^{1/n}
\f]
where
\vartable
n | Binomial blending exponent
\endvartable
- \c exponential (PH:Eq. 32):
\f[
\epsilon = \epsilon_{vis} \exp[-\Gamma] +\epsilon_{log} \exp[-1/\Gamma]
\f]
where (PH:p. 193)
\vartable
\Gamma_\epsilon | \f$\Gamma = 0.001 (y^+)^4 / (1.0 + y^+)\f$
\Gamma_G | \f$\Gamma = 0.01 (y^+)^4 / (1.0 + 5.0 y^+)\f$
\Gamma_\epsilon | Blending expression for \f$\epsilon\f$
\Gamma_G | Blending expression for \f$G\f$
\endvartable
\c G predictions for the viscous and inertial sublayers are blended
in a stepwise manner, and \c G below \f$y^+_{lam}\f$ (i.e. in the viscous
sublayer) is presumed to be zero.
Note
The coefficients \c Cmu, \c kappa, and \c E are obtained from
- The coefficients \c Cmu, \c kappa, and \c E are obtained from
the specified \c nutWallFunction in order to ensure that each patch
possesses the same set of values for these coefficients.
- \c lowReCorrection operates with only \c stepwise blending treatment to
ensure the backward compatibility.
- If \c lowReCorrection is \c on, \c stepwise blending treatment is fully
active.
- If \c lowReCorrection is \c off, only the inertial sublayer prediction
is used in the wall function, hence high-Re mode operation.
See also
Foam::fixedInternalValueFvPatchField
- Foam::omegaWallFunctionFvPatchScalarField
SourceFiles
epsilonWallFunctionFvPatchScalarField.C
@ -103,6 +189,33 @@ class epsilonWallFunctionFvPatchScalarField
:
public fixedValueFvPatchField<scalar>
{
private:
// Private Enumerations
//- Options for the blending treatment of viscous and inertial sublayers
enum blendingType
{
STEPWISE, //!< "Stepwise switch (discontinuous)"
MAX, //!< "Maximum value switch (discontinuous)"
BINOMIAL, //!< "Binomial blending (smooth)"
EXPONENTIAL //!< "Exponential blending (smooth)"
};
//- Names for blendingType
static const Enum<blendingType> blendingTypeNames;
// Private Data
//- Blending treatment (default = blendingType::STEPWISE)
const enum blendingType blending_;
//- Binomial blending exponent being used when
//- blendingType is blendingType::BINOMIAL (default = 2)
const scalar n_;
protected:
// Protected Data
@ -110,8 +223,8 @@ protected:
//- Tolerance used in weighted calculations
static scalar tolerance_;
//- Apply low-Re correction term; default = no
bool lowReCorrection_;
//- Apply low-Re correction term (default = no)
const bool lowReCorrection_;
//- Initialised flag
bool initialised_;

View File

@ -72,7 +72,15 @@ Foam::kLowReWallFunctionFvPatchScalarField::kLowReWallFunctionFvPatchScalarField
)
:
fixedValueFvPatchField<scalar>(p, iF, dict),
Ceps2_(dict.getOrDefault<scalar>("Ceps2", 1.9)),
Ceps2_
(
dict.getCheckOrDefault<scalar>
(
"Ceps2",
1.9,
scalarMinMax::ge(SMALL)
)
),
Ck_(dict.getOrDefault<scalar>("Ck", -0.416)),
Bk_(dict.getOrDefault<scalar>("Bk", 8.366)),
C_(dict.getOrDefault<scalar>("C", 11.0))
@ -145,12 +153,10 @@ void Foam::kLowReWallFunctionFvPatchScalarField::updateCoeffs()
forAll(kw, facei)
{
const label celli = patch().faceCells()[facei];
const scalar uTau = Cmu25*sqrt(k[celli]);
const scalar yPlus = uTau*y[facei]/nuw[facei];
if (nutw.yPlusLam() < yPlus)
if (yPlus > nutw.yPlusLam())
{
kw[facei] = Ck_/nutw.kappa()*log(yPlus) + Bk_;
}

View File

@ -6,7 +6,7 @@
\\/ M anipulation |
-------------------------------------------------------------------------------
Copyright (C) 2012-2016, 2019 OpenFOAM Foundation
Copyright (C) 2019 OpenCFD Ltd.
Copyright (C) 2019-2020 OpenCFD Ltd.
-------------------------------------------------------------------------------
License
This file is part of OpenFOAM.
@ -31,39 +31,66 @@ Group
grpWallFunctions
Description
This boundary condition provides a wall constraint on turbulent kinetic
This boundary condition provides a wall constraint on the turbulent kinetic
energy, i.e. \c k, for low- and high-Reynolds number turbulence models.
The model operates in two modes, based on the computed viscous-to-turbulent
switch-over \c y+ value derived from \c kappa and \c E.
Usage
\table
Property | Description | Required | Default value
Ceps2 | Model coefficient | no | 1.9
Ck | Model coefficient | no | -0.416
Bk | Model coefficient | no | 8.366
C | Model coefficient | no | 11.0
\endtable
Example of the boundary condition specification:
\verbatim
<patchName>
{
// Mandatory entries
// Mandatory entries (unmodifiable)
type kLowReWallFunction;
// Optional entries
// Optional entries (unmodifiable)
Ceps2 1.9;
Ck -0.416;
Bk 8.366;
C 11.0;
// Optional (inherited) entries
...
}
\endverbatim
where the entries mean:
\table
Property | Description | Type | Req'd | Dflt
type | Type name: kLowReWallFunction | word | yes | -
Ceps2 | Model coefficient | scalar | no | 1.9
Ck | Model coefficient | scalar | no | -0.416
Bk | Model coefficient | scalar | no | 8.366
C | Model coefficient | scalar | no | 11.0
\endtable
The inherited entries are elaborated in:
- \link fixedValueFvPatchField.H \endlink
Viscous and inertial sublayer predictions for \c k are blended in
a stepwise manner:
\f[
k = k_{log} \qquad if \quad y^+ > y^+_{lam}
\f]
\f[
k = k_{vis} \qquad if \quad y^+ <= y^+_{lam}
\f]
where
\vartable
k_{vis} | k prediction in the viscous sublayer
k_{log} | k 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
\endvartable
Note
The coefficients \c Cmu, \c kappa, and \c E are obtained from
the specified \c nutWallFunction in order to ensure that each patch
possesses the same set of values for these coefficients.
See also
Foam::fixedValueFvPatchField
- Foam::fixedValueFvPatchField
- Foam::kqRWallFunctionFvPatchScalarField
SourceFiles
kLowReWallFunctionFvPatchScalarField.C
@ -129,7 +156,7 @@ public:
);
//- Construct by mapping given kLowReWallFunctionFvPatchScalarField
// onto a new patch
//- onto a new patch
kLowReWallFunctionFvPatchScalarField
(
const kLowReWallFunctionFvPatchScalarField&,

View File

@ -6,7 +6,7 @@
\\/ M anipulation |
-------------------------------------------------------------------------------
Copyright (C) 2011-2017, 2019 OpenFOAM Foundation
Copyright (C) 2019 OpenCFD Ltd.
Copyright (C) 2019-2020 OpenCFD Ltd.
-------------------------------------------------------------------------------
License
This file is part of OpenFOAM.
@ -40,13 +40,22 @@ Usage
\verbatim
<patchName>
{
// Mandatory entries
// Mandatory entries (unmodifiable)
type kqRWallFunction;
// Optional (inherited) entries
...
}
\endverbatim
See also
Foam::zeroGradientFvPatchField
where the entries mean:
\table
Property | Description | Type | Req'd | Dflt
type | Type name: kqRWallFunction | word | yes | -
\endtable
The inherited entries are elaborated in:
- \link zeroGradientFvPatchField.H \endlink
SourceFiles
kqRWallFunctionFvPatchField.C

View File

@ -6,7 +6,7 @@
\\/ M anipulation |
-------------------------------------------------------------------------------
Copyright (C) 2011-2016, 2019 OpenFOAM Foundation
Copyright (C) 2019 OpenCFD Ltd.
Copyright (C) 2019-2020 OpenCFD Ltd.
-------------------------------------------------------------------------------
License
This file is part of OpenFOAM.
@ -32,30 +32,30 @@ Group
Description
This boundary condition provides a wall constraint on the turbulent
kinematic viscosity, i.e. \c nut for use with low Reynolds number models.
viscosity, i.e. \c nut for low Reynolds number models.
It sets \c nut to zero, and provides an access function to calculate \c y+.
Usage
\table
Property | Description | Required | Default value
Cmu | Model coefficient | no | 0.09
kappa | Von Karman constant | no | 0.41
E | Model coefficient | no | 9.8
\endtable
Example of the boundary condition specification:
\verbatim
<patchName>
{
// Mandatory entries
// Mandatory entries (unmodifiable)
type nutLowReWallFunction;
// Optional entries
// Optional (inherited) entries
...
}
\endverbatim
See also
Foam::nutWallFunctionFvPatchScalarField
where the entries mean:
\table
Property | Description | Type | Req'd | Dflt
type | Type name: nutLowReWallFunction | word | yes | -
\endtable
The inherited entries are elaborated in:
- \link nutWallFunctionFvPatchScalarField.H \endlink
SourceFiles
nutLowReWallFunctionFvPatchScalarField.C
@ -84,7 +84,7 @@ protected:
// Protected Member Functions
//- Calculate the turbulence viscosity
//- Calculate the turbulent viscosity
virtual tmp<scalarField> calcNut() const;

View File

@ -155,7 +155,7 @@ nutUBlendedWallFunctionFvPatchScalarField
)
:
nutWallFunctionFvPatchScalarField(p, iF, dict),
n_(dict.getOrDefault<scalar>("n", 4))
n_(dict.getOrDefault<scalar>("n", 4.0))
{}
@ -214,8 +214,8 @@ void Foam::nutUBlendedWallFunctionFvPatchScalarField::write
{
fvPatchField<scalar>::write(os);
writeLocalEntries(os);
writeEntry("value", os);
os.writeEntry("n", n_);
writeEntry("value", os);
}

View File

@ -5,7 +5,7 @@
\\ / A nd | www.openfoam.com
\\/ M anipulation |
-------------------------------------------------------------------------------
Copyright (C) 2016-2019 OpenCFD Ltd.
Copyright (C) 2016-2020 OpenCFD Ltd.
-------------------------------------------------------------------------------
License
This file is part of OpenFOAM.
@ -31,8 +31,8 @@ Group
Description
This boundary condition provides a wall constraint on the turbulent
kinematic viscosity, i.e. \c nut, when using wall functions based on
a blending of laminar sub-layer and log region contributions.
viscosity, i.e. \c nut, when using wall functions based on
a blending of viscous and inertial sublayer contributions.
\f[
u_\tau = (u_{\tau,v}^n + u_{\tau,l}^n)^{1/n}
@ -40,9 +40,9 @@ Description
where
\vartable
u_\tau | friction velocity
u_{\tau,v} | friction velocity in the viscous region
u_{\tau,l} | friction velocity in the log region
u_\tau | Friction velocity
u_{\tau,v} | Friction velocity in the viscous sublayer
u_{\tau,l} | Friction velocity in the inertial sublayer
\endvartable
Reference:
@ -55,35 +55,43 @@ Description
\endverbatim
Usage
\table
Property | Description | Required | Default value
Cmu | Model coefficient | no | 0.09
kappa | Von Karman constant | no | 0.41
E | Model coefficient | no | 9.8
\endtable
Example of the boundary condition specification:
\verbatim
<patchName>
{
// Mandatory entries
// Mandatory entries (unmodifiable)
type nutUBlendedWallFunction;
// Optional entries
n 4.0;
// Optional (inherited) entries
...
}
\endverbatim
where the entries mean:
\table
Property | Description | Type | Req'd | Dflt
type | Type name: nutUBlendedWallFunction | word | yes | -
n | Blending factor | scalar | no | 4.0
\endtable
The inherited entries are elaborated in:
- \link nutWallFunctionFvPatchScalarField.H \endlink
Note
The full 'automatic wall treatment' description also requires use of the
Foam::omegaWallFunction with the \c blended flag set to 'on'
- The full 'automatic wall treatment' description also requires use of the
\link omegaWallFunctionFvPatchScalarField.H \endlink with the \c blending
option \c binomial or with the deprecated \c blended flag set to \c on.
- Suffers from non-exact restart since \c correctNut() (called through
\c turbulence->validate) returns a slightly different value every time
it is called.
See \link nutUSpaldingWallFunctionFvPatchScalarField.C \endlink.
Suffers from non-exact restart since correctNut() (called through
turbulence->validate) returns a slightly different value every time
it is called. See nutUSpaldingWallFunctionFvPatchScalarField.C
SeeAlso
Foam::nutWallFunctionFvPatchScalarField
Foam::omegaWallFunctionFvPatchScalarField
See also
- Foam::nutWallFunctionFvPatchScalarField
- Foam::omegaWallFunctionFvPatchScalarField
SourceFiles
nutUBlendedWallFunctionFvPatchScalarField.C
@ -118,7 +126,7 @@ protected:
// Protected Member Functions
//- Calculate the turbulence viscosity
//- Calculate the turbulent viscosity
virtual tmp<scalarField> calcNut() const;
//- Calculate the friction velocity

View File

@ -6,7 +6,7 @@
\\/ M anipulation |
-------------------------------------------------------------------------------
Copyright (C) 2011-2016 OpenFOAM Foundation
Copyright (C) 2019 OpenCFD Ltd.
Copyright (C) 2019-2020 OpenCFD Ltd.
-------------------------------------------------------------------------------
License
This file is part of OpenFOAM.
@ -32,46 +32,51 @@ Group
Description
This boundary condition provides a wall constraint on the turbulent
kinematic viscosity, i.e. \c nut, when using wall functions for rough walls,
viscosity, i.e. \c nut, when using wall functions for rough walls,
based on velocity, \c U.
Usage
\table
Property | Description | Required | Default value
roughnessHeight | Roughness height | yes |
roughnessConstant | Roughness constanr | yes |
roughnessFactor | Scaling factor | yes |
Cmu | Model coefficient | no | 0.09
kappa | Von Karman constant | no | 0.41
E | Model coefficient | no | 9.8
maxIter | Number of Newton-Raphson iterations | no | 10
tolerance | Convergence tolerance | no | 0.0001
\endtable
Example of the boundary condition specification:
\verbatim
<patchName>
{
// Mandatory entries
type nutURoughWallFunction;
// Mandatory entries (unmodifiable)
type nutURoughWallFunction;
roughnessHeight 1e-5;
roughnessConstant 0.5;
roughnessFactor 1;
// Optional entries
// Optional entries (unmodifiable)
maxIter 10;
tolerance 0.0001;
// Optional (inherited) entries
...
}
\endverbatim
Note
Suffers from non-exact restart since correctNut() (called through
turbulence->validate) returns a slightly different value every time
it is called. See nutUSpaldingWallFunctionFvPatchScalarField.C.
Can be avoided by seeding the NR with e.g. the laminar viscosity
or tightening the convergence tolerance
to e.g. 1e-7 and the max number of iterations to 100.
where the entries mean:
\table
Property | Description | Type | Req'd | Dflt
type | Type name: nutURoughWallFunction | word | yes | -
roughnessHeight | Roughness height | scalar | yes | -
roughnessConstant | Roughness constant | scalar | yes | -
roughnessFactor | Scaling factor | scalar | yes | -
maxIter | Number of Newton-Raphson iterations | label | no | 10
tolerance | Convergence tolerance | scalar | no | 0.0001
\endtable
See also
Foam::nutWallFunctionFvPatchScalarField
The inherited entries are elaborated in:
- \link nutWallFunctionFvPatchScalarField.H \endlink
Note
- Suffers from non-exact restart since \c correctNut() (called through
\c turbulence->validate) returns a slightly different value every time
it is called.
See \link nutUSpaldingWallFunctionFvPatchScalarField.C \endlink.
- Can be avoided by seeding the NR with e.g. the laminar viscosity
or tightening the convergence tolerance to e.g. 1e-7 and the max
number of iterations to 100.
SourceFiles
nutURoughWallFunctionFvPatchScalarField.C
@ -214,7 +219,6 @@ public:
return roughnessHeight_;
}
//- Return the roughness constant scale
scalar roughnessConstant() const
{

View File

@ -6,7 +6,7 @@
\\/ M anipulation |
-------------------------------------------------------------------------------
Copyright (C) 2011-2016 OpenFOAM Foundation
Copyright (C) 2019 OpenCFD Ltd.
Copyright (C) 2019-2020 OpenCFD Ltd.
-------------------------------------------------------------------------------
License
This file is part of OpenFOAM.
@ -32,7 +32,7 @@ Group
Description
This boundary condition provides a wall constraint on the turbulent
kinematic viscosity, i.e. \c nut, when using wall functions for rough walls,
viscosity, i.e. \c nut, when using wall functions for rough walls,
based on velocity, \c U, using Spalding's law to give a continuous \c nut
profile to the wall (y+ = 0)
@ -43,52 +43,54 @@ Description
where
\vartable
y^+ | non-dimensional position
u^+ | non-dimensional velocity
\kappa | Von Karman constant
y^+ | Non-dimensional position
u^+ | Non-dimensional velocity
\kappa | von Kármán constant
\endvartable
Usage
\table
Property | Description | Required | Default value
Cmu | Model coefficient | no | 0.09
kappa | Von Karman constant | no | 0.41
E | Model coefficient | no | 9.8
maxIter | Number of Newton-Raphson iterations | no | 10
tolerance | Convergence tolerance | no | 0.0001
\endtable
Example of the boundary condition specification:
\verbatim
<patchName>
{
// Mandatory entries
// Mandatory entries (unmodifiable)
type nutUSpaldingWallFunction;
// Optional entries
// Optional entries (unmodifiable)
maxIter 10;
tolerance 0.0001;
// Optional (inherited) entries
...
}
\endverbatim
See also
Foam::nutWallFunctionFvPatchScalarField
where the entries mean:
\table
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
\endtable
The inherited entries are elaborated in:
- \link nutWallFunctionFvPatchScalarField.H \endlink
Note
Suffers from non-exact restart since correctNut() (called through
turbulence->validate) returns a slightly different value every time
- Suffers from non-exact restart since \c correctNut() (called through
\c 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).
uses the current value of \c *this (\c =nut ).
This can be avoided by overriding the tolerance. This also switches on
- 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, e.g.
a tight tolerance, to make sure to kick every iteration, e.g.
maxIter 100;
tolerance 1e-7;
to make sure to kick every iteration.
SourceFiles
nutUSpaldingWallFunctionFvPatchScalarField.C
@ -131,7 +133,7 @@ protected:
// Protected Member Functions
//- Calculate the turbulence viscosity
//- Calculate the turbulent viscosity
virtual tmp<scalarField> calcNut() const;
//- Calculate the friction velocity

View File

@ -6,7 +6,7 @@
\\/ M anipulation |
-------------------------------------------------------------------------------
Copyright (C) 2011-2016 OpenFOAM Foundation
Copyright (C) 2019 OpenCFD Ltd.
Copyright (C) 2019-2020 OpenCFD Ltd.
-------------------------------------------------------------------------------
License
This file is part of OpenFOAM.
@ -32,41 +32,41 @@ Group
Description
This boundary condition provides a wall constraint on the turbulent
kinematic viscosity, i.e. \c nut, when using wall functions, based on
viscosity, i.e. \c nut, when using wall functions, based on
velocity, i.e. \c U.
As input, the user specifies a look-up table of \c U+ as a function of
near-wall Reynolds number.
The table should be located in the $FOAM_CASE/constant directory.
The table should be located in the \c $FOAM_CASE/constant directory.
Usage
\table
Property | Description | Required | Default value
uPlusTable | U+ as a function of Re table name | yes |
Cmu | Model coefficient | no | 0.09
kappa | Von Karman constant | no | 0.41
E | Model coefficient | no | 9.8
\endtable
Example of the boundary condition specification:
\verbatim
<patchName>
{
// Mandatory entries
// Mandatory entries (unmodifiable)
type nutTabulatedWallFunction;
uPlusTable myUPlusTable;
// Optional entries
// Optional (inherited) entries
...
}
\endverbatim
Note
The tables are not registered since the same table object may be used for
more than one patch.
where the entries mean:
\table
Property | Description | Type | Req'd | Dflt
type | Type name: nutUTabulatedWallFunction | word | yes | -
uPlusTable | U+ as a function of Re table name | word | yes | -
\endtable
See also
Foam::nutWallFunctionFvPatchScalarField
The inherited entries are elaborated in:
- \link nutWallFunctionFvPatchScalarField.H \endlink
Note
- The tables are not registered since the same table object may be used for
more than one patch.
SourceFiles
nutUTabulatedWallFunctionFvPatchScalarField.C
@ -105,7 +105,7 @@ protected:
// Protected Member Functions
//- Calculate the turbulence viscosity
//- Calculate the turbulent viscosity
virtual tmp<scalarField> calcNut() const;
//- Calculate wall u+ from table

View File

@ -61,11 +61,14 @@ Foam::nutUWallFunctionFvPatchScalarField::calcNut() const
forAll(yPlus, facei)
{
if (yPlusLam_ < yPlus[facei])
{
nutw[facei] =
nuw[facei]*(yPlus[facei]*kappa_/log(E_*yPlus[facei]) - 1.0);
}
// Viscous sublayer contribution
const scalar nutVis = 0.0;
// Inertial sublayer contribution
const scalar nutLog =
nuw[facei]*(yPlus[facei]*kappa_/log(E_*yPlus[facei]) - 1.0);
nutw[facei] = blend(nutVis, nutLog, yPlus[facei]);
}
return tnutw;

View File

@ -6,7 +6,7 @@
\\/ M anipulation |
-------------------------------------------------------------------------------
Copyright (C) 2011-2016 OpenFOAM Foundation
Copyright (C) 2019 OpenCFD Ltd.
Copyright (C) 2019-2020 OpenCFD Ltd.
-------------------------------------------------------------------------------
License
This file is part of OpenFOAM.
@ -32,35 +32,38 @@ Group
Description
This boundary condition provides a wall constraint on the turbulent
kinematic viscosity, i.e. \c nut, when using wall functions, based on
viscosity, i.e. \c nut, when using wall functions, based on
velocity, \c U.
Usage
\table
Property | Description | Required | Default value
Cmu | Model coefficient | no | 0.09
kappa | Von Karman constant | no | 0.41
E | Model coefficient | no | 9.8
\endtable
Example of the boundary condition specification:
\verbatim
<patchName>
{
// Mandatory entries
// Mandatory entries (unmodifiable)
type nutUWallFunction;
// Optional entries
// Optional (inherited) entries
...
}
\endverbatim
See also
Foam::nutWallFunctionFvPatchScalarField
where the entries mean:
\table
Property | Description | Type | Req'd | Dflt
type | Type name: nutUWallFunction | word | yes | -
\endtable
The inherited entries are elaborated in:
- \link nutWallFunctionFvPatchScalarField.H \endlink
Note
Suffers from non-exact restart since correctNut() (called through
turbulence->validate) returns a slightly different value every time
it is called. See nutUSpaldingWallFunctionFvPatchScalarField.C
- Suffers from non-exact restart since \c correctNut() (called through
\c turbulence->validate) returns a slightly different value every time
it is called.
See \link nutUSpaldingWallFunctionFvPatchScalarField.C \endlink.
- See \link nutWallFunctionFvPatchScalarField.H \endlink for the wall
function blending treatments.
SourceFiles
nutUWallFunctionFvPatchScalarField.C
@ -92,7 +95,7 @@ protected:
//- Calculate yPlus
virtual tmp<scalarField> calcYPlus(const scalarField& magUp) const;
//- Calculate the turbulence viscosity
//- Calculate the turbulent viscosity
virtual tmp<scalarField> calcNut() const;

View File

@ -40,6 +40,21 @@ namespace Foam
defineTypeNameAndDebug(nutWallFunctionFvPatchScalarField, 0);
}
// * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
const Foam::Enum
<
Foam::nutWallFunctionFvPatchScalarField::blendingType
>
Foam::nutWallFunctionFvPatchScalarField::blendingTypeNames
({
{ blendingType::STEPWISE , "stepwise" },
{ blendingType::MAX , "max" },
{ blendingType::BINOMIAL , "binomial" },
{ blendingType::EXPONENTIAL, "exponential" }
});
// * * * * * * * * * * * * Protected Member Functions * * * * * * * * * * * //
void Foam::nutWallFunctionFvPatchScalarField::checkType()
@ -77,6 +92,14 @@ void Foam::nutWallFunctionFvPatchScalarField::writeLocalEntries
Ostream& os
) const
{
os.writeKeyword("blending") << blendingTypeNames[blending_]
<< token::END_STATEMENT << nl;
if (blending_ == blendingType::BINOMIAL)
{
os.writeEntry("n", n_);
}
if (UName_ != word::null)
{
os.writeEntry("U", UName_);
@ -97,6 +120,8 @@ Foam::nutWallFunctionFvPatchScalarField::nutWallFunctionFvPatchScalarField
)
:
fixedValueFvPatchScalarField(p, iF),
blending_(blendingType::STEPWISE),
n_(4.0),
UName_(word::null),
Cmu_(0.09),
kappa_(0.41),
@ -116,6 +141,8 @@ Foam::nutWallFunctionFvPatchScalarField::nutWallFunctionFvPatchScalarField
)
:
fixedValueFvPatchScalarField(ptf, p, iF, mapper),
blending_(ptf.blending_),
n_(ptf.n_),
UName_(ptf.UName_),
Cmu_(ptf.Cmu_),
kappa_(ptf.kappa_),
@ -134,6 +161,24 @@ Foam::nutWallFunctionFvPatchScalarField::nutWallFunctionFvPatchScalarField
)
:
fixedValueFvPatchScalarField(p, iF, dict),
blending_
(
blendingTypeNames.getOrDefault
(
"blending",
dict,
blendingType::STEPWISE
)
),
n_
(
dict.getCheckOrDefault<scalar>
(
"n",
4.0,
scalarMinMax::ge(0)
)
),
UName_(dict.getOrDefault<word>("U", word::null)),
Cmu_(dict.getOrDefault<scalar>("Cmu", 0.09)),
kappa_(dict.getOrDefault<scalar>("kappa", 0.41)),
@ -150,6 +195,8 @@ Foam::nutWallFunctionFvPatchScalarField::nutWallFunctionFvPatchScalarField
)
:
fixedValueFvPatchScalarField(wfpsf),
blending_(wfpsf.blending_),
n_(wfpsf.n_),
UName_(wfpsf.UName_),
Cmu_(wfpsf.Cmu_),
kappa_(wfpsf.kappa_),
@ -167,6 +214,8 @@ Foam::nutWallFunctionFvPatchScalarField::nutWallFunctionFvPatchScalarField
)
:
fixedValueFvPatchScalarField(wfpsf, iF),
blending_(wfpsf.blending_),
n_(wfpsf.n_),
UName_(wfpsf.UName_),
Cmu_(wfpsf.Cmu_),
kappa_(wfpsf.kappa_),
@ -203,15 +252,73 @@ Foam::scalar Foam::nutWallFunctionFvPatchScalarField::yPlusLam
{
scalar ypl = 11.0;
for (int i = 0; i < 10; ++i)
for (label i = 0; i < 10; ++i)
{
ypl = log(max(E*ypl, 1))/kappa;
ypl = log(max(E*ypl, 1.0))/kappa;
}
return ypl;
}
Foam::scalar Foam::nutWallFunctionFvPatchScalarField::blend
(
const scalar nutVis,
const scalar nutLog,
const scalar yPlus
) const
{
scalar nutw = 0.0;
switch (blending_)
{
case blendingType::STEPWISE:
{
if (yPlus > yPlusLam_)
{
nutw = nutLog;
}
else
{
nutw = nutVis;
}
break;
}
case blendingType::MAX:
{
// (PH:Eq. 27)
nutw = max(nutVis, nutLog);
break;
}
case blendingType::BINOMIAL:
{
// (ME:Eqs. 15-16)
nutw =
pow
(
pow(nutVis, n_) + pow(nutLog, n_),
1.0/n_
);
break;
}
case blendingType::EXPONENTIAL:
{
// (PH:Eq. 31)
const scalar Gamma = 0.01*pow4(yPlus)/(1.0 + 5.0*yPlus);
const scalar invGamma = 1.0/(Gamma + ROOTVSMALL);
nutw = nutVis*exp(-Gamma) + nutLog*exp(-invGamma);
break;
}
}
return nutw;
}
Foam::scalar Foam::nutWallFunctionFvPatchScalarField::yPlusLam() const
{
return yPlusLam_;

View File

@ -33,37 +33,127 @@ Group
Description
The class \c nutWallFunction is a base class that parents the derived
boundary conditions which provide a wall constraint on various fields, such
as turbulence kinematic viscosity, i.e. \c nut, for low- and high-Reynolds
number turbulence models.
as turbulent viscosity, i.e. \c nut, or turbulent kinetic energy dissipation
rate, i.e. \c epsilon, for low- and high-Reynolds number turbulence models.
Reference:
\verbatim
Default model coefficients:
Default model coefficients (tag:VM):
Versteeg, H. K., & Malalasekera, W. (2011).
An introduction to computational fluid dynamics: the finite
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
\endverbatim
Usage
\table
Property | Description | Required | Default value
Cmu | Cmu coefficient | no | 0.09
kappa | Von Karman constant | no | 0.41
E | E coefficient | no | 9.8
\endtable
Examples of the boundary condition specification:
Example of the boundary condition specification:
\verbatim
<patchName>
{
// Mandatory entries
type nutWallFunction;
// Mandatory and other optional entries
...
// Optional entries
// Optional (inherited) entries
Cmu 0.09;
kappa 0.41;
E 9.8;
blending stepwise;
n 4.0;
// Optional (inherited) entries
...
}
\endverbatim
where the entries mean:
\table
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
\endtable
The inherited entries are elaborated in:
- \link fixedValueFvPatchFields.H \endlink
Options for the \c blending entry:
\verbatim
stepwise | Stepwise switch (discontinuous)
max | Maximum value switch (discontinuous)
binomial | Binomial blending (smooth)
exponential | Exponential blending (smooth)
\endverbatim
wherein \c nut predictions for the viscous and inertial sublayers are
blended according to the following expressions:
- \c stepwise (default):
\f[
\nu_t = {\nu_t}_{log} \qquad if \quad y^+ > y^+_{lam}
\f]
\f[
\nu_t = {\nu_t}_{vis} \qquad if \quad y^+ <= y^+_{lam}
\f]
where
\vartable
{\nu_t}_{vis} | \f$\nu_t\f$ prediction in the viscous sublayer
{\nu_t}_{log} | \f$\nu_t\f$ 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
\endvartable
- \c max (PH:Eq. 27):
\f[
\nu_t = max({\nu_t}_{vis}, {\nu_t}_{log})
\f]
- \c binomial (ME:Eqs. 15-16):
\f[
\nu_t = (({\nu_t}_{vis})^n + ({\nu_t}_{log})^n)^{1/n}
\f]
where
\vartable
n | Binomial blending exponent
\endvartable
- \c exponential (PH:Eq. 32):
\f[
\nu_t = {\nu_t}_{vis} \exp[-\Gamma] + {\nu_t}_{log} \exp[-1/\Gamma]
\f]
where (PH:Eq. 31)
\vartable
\Gamma | Blending expression
\Gamma | \f$0.01 (y^+)^4 / (1.0 + 5.0 y^+)\f$
\endvartable
Note
- \c nutWallFunction is not directly usable.
See also
Foam::fixedValueFvPatchField
@ -94,23 +184,46 @@ class nutWallFunctionFvPatchScalarField
{
protected:
// Protected Enumerations
//- Options for the blending treatment of viscous and inertial sublayers
enum blendingType
{
STEPWISE, //!< "Stepwise switch (discontinuous)"
MAX, //!< "Maximum value switch (discontinuous)"
BINOMIAL, //!< "Binomial blending (smooth)"
EXPONENTIAL //!< "Exponential blending (smooth)"
};
//- Names for blendingType
static const Enum<blendingType> blendingTypeNames;
// Protected Data
//- Blending treatment (default = blendingType::STEPWISE)
const enum blendingType blending_;
//- Binomial blending exponent being used when
//- blendingType is blendingType::BINOMIAL (default = 4)
const scalar n_;
//- Name of velocity field
// Default is null (not specified) in which case the velocity is
// retrieved from the turbulence model
word UName_;
//- Cmu coefficient
//- Empirical model coefficient
scalar Cmu_;
//- Von Karman constant
//- von Kármán constant
scalar kappa_;
//- E coefficient
//- Wall roughness parameter
scalar E_;
//- Estimated y+ value at the edge of the viscous sublayer
//- Estimated y+ value at the intersection
//- of the viscous and inertial sublayers
scalar yPlusLam_;
@ -123,7 +236,7 @@ protected:
//- Check the type of the patch
virtual void checkType();
//- Calculate the turbulence viscosity
//- Calculate the turbulent viscosity
virtual tmp<scalarField> calcNut() const = 0;
//- Write local wall function variables
@ -207,13 +320,22 @@ public:
const label patchi
);
//- Calculate the y+ at the edge of the viscous sublayer
//- Estimate the y+ at the intersection of the two sublayers
static scalar yPlusLam(const scalar kappa, const scalar E);
//- Return the y+ at the edge of the viscous sublayer
//- Return the estimated y+ at the two-sublayer intersection
scalar yPlusLam() const;
//- Return the blended nut according to the chosen blending treatment
scalar blend
(
const scalar nutVis,
const scalar nutLog,
const scalar yPlus
) const;
//- Calculate and return the yPlus at the boundary
// yPlus is the first-cell-centre height from boundary in wall units
virtual tmp<scalarField> yPlus() const = 0;

View File

@ -6,7 +6,7 @@
\\/ M anipulation |
-------------------------------------------------------------------------------
Copyright (C) 2011-2016, 2019 OpenFOAM Foundation
Copyright (C) 2019 OpenCFD Ltd.
Copyright (C) 2019-2020 OpenCFD Ltd.
-------------------------------------------------------------------------------
License
This file is part of OpenFOAM.
@ -33,38 +33,38 @@ Group
Description
This boundary condition provides a wall constraint on the turbulent
kinematic viscosity, i.e. \c nut, when using wall functions for rough walls,
based on turbulent kinetic energy, \c k. The condition manipulates the \c E
parameter to account for roughness effects.
based on the turbulent kinetic energy, \c k. The condition manipulates the
wall roughness parameter, i.e. \c E, to account for roughness effects.
Parameter ranges
Parameter ranges:
- roughness height = sand-grain roughness (0 for smooth walls)
- roughness constant = 0.5-1.0
Usage
\table
Property | Description | Required | Default value
Ks | Sand-grain roughness height | yes |
Cs | Roughness constant | yes |
Cmu | Model coefficient | no | 0.09
kappa | Von Karman constant | no | 0.41
E | Model coefficient | no | 9.8
\endtable
Example of the boundary condition specification:
\verbatim
<patchName>
{
// Mandatory entries
// Mandatory entries (unmodifiable)
type nutkRoughWallFunction;
Ks uniform 0;
Cs uniform 0.5;
// Optional entries
// Optional (inherited) entries
...
}
\endverbatim
See also
Foam::nutkRoughWallFunctionFvPatchScalarField
where the entries mean:
\table
Property | Description | Type | Req'd | Dflt
type | Type name: nutkRoughWallFunction | word | yes | -
Ks | Sand-grain roughness height | scalarField | yes | -
Cs | Roughness constant | scalarField | yes | -
\endtable
The inherited entries are elaborated in:
- \link nutkWallFunctionFvPatchScalarField.H \endlink
SourceFiles
nutkRoughWallFunctionFvPatchScalarField.C

View File

@ -67,10 +67,13 @@ calcNut() const
const scalar yPlus = Cmu25*y[facei]*sqrt(k[celli])/nuw[facei];
if (yPlusLam_ < yPlus)
{
nutw[facei] = nuw[facei]*(yPlus*kappa_/log(E_*yPlus) - 1.0);
}
// Viscous sublayer contribution
const scalar nutVis = 0.0;
// Inertial sublayer contribution
const scalar nutLog = nuw[facei]*(yPlus*kappa_/log(E_*yPlus) - 1.0);
nutw[facei] = blend(nutVis, nutLog, yPlus);
}
return tnutw;

View File

@ -6,7 +6,7 @@
\\/ M anipulation |
-------------------------------------------------------------------------------
Copyright (C) 2011-2019 OpenFOAM Foundation
Copyright (C) 2019 OpenCFD Ltd.
Copyright (C) 2019-2020 OpenCFD Ltd.
-------------------------------------------------------------------------------
License
This file is part of OpenFOAM.
@ -32,30 +32,34 @@ Group
Description
This boundary condition provides a wall constraint on the turbulent
kinematic viscosity, i.e. \c nut, when using wall functions,
based on turbulent kinetic energy, \c k.
viscosity, i.e. \c nut, when using wall functions,
based on the turbulent kinetic energy, \c k.
Usage
\table
Property | Description | Required | Default value
Cmu | Model coefficient | no | 0.09
kappa | Von Karman constant | no | 0.41
E | Model coefficient | no | 9.8
\endtable
Example of the boundary condition specification:
\verbatim
<patchName>
{
// Mandatory entries
// Mandatory entries (unmodifiable)
type nutkWallFunction;
// Optional entries
// Optional (inherited) entries
...
}
\endverbatim
See also
Foam::nutWallFunctionFvPatchScalarField
where the entries mean:
\table
Property | Description | Type | Req'd | Dflt
type | Type name: nutkWallFunction | word | yes | -
\endtable
The inherited entries are elaborated in:
- \link nutWallFunctionFvPatchScalarField.H \endlink
Note
- See \link nutWallFunctionFvPatchScalarField.H \endlink for the wall
function blending treatments.
SourceFiles
nutkWallFunctionFvPatchScalarField.C
@ -84,7 +88,7 @@ protected:
// Protected Member Functions
//- Calculate the turbulence viscosity
//- Calculate the turbulent viscosity
virtual tmp<scalarField> calcNut() const;

View File

@ -32,9 +32,22 @@ License
#include "fvMatrix.H"
#include "addToRunTimeSelectionTable.H"
// * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
const Foam::Enum
<
Foam::omegaWallFunctionFvPatchScalarField::blendingType
>
Foam::omegaWallFunctionFvPatchScalarField::blendingTypeNames
({
{ blendingType::STEPWISE , "stepwise" },
{ blendingType::MAX , "max" },
{ blendingType::BINOMIAL2 , "binomial2" },
{ blendingType::BINOMIAL , "binomial" },
{ blendingType::EXPONENTIAL, "exponential" },
{ blendingType::TANH, "tanh" }
});
Foam::scalar Foam::omegaWallFunctionFvPatchScalarField::tolerance_ = 1e-5;
// * * * * * * * * * * * * Protected Member Functions * * * * * * * * * * * //
@ -106,18 +119,16 @@ void Foam::omegaWallFunctionFvPatchScalarField::createAveragingWeights()
omegaPatches.append(patchi);
const labelUList& faceCells = bf[patchi].patch().faceCells();
forAll(faceCells, i)
for (const auto& celli : faceCells)
{
const label celli = faceCells[i];
weights[celli]++;
++weights[celli];
}
}
}
cornerWeights_.setSize(bf.size());
forAll(omegaPatches, i)
for (const auto& patchi : omegaPatches)
{
const label patchi = omegaPatches[i];
const fvPatchScalarField& wf = weights.boundaryField()[patchi];
cornerWeights_[patchi] = 1.0/wf.patchInternalField();
}
@ -212,43 +223,80 @@ void Foam::omegaWallFunctionFvPatchScalarField::calculate
forAll(nutw, facei)
{
const label celli = patch.faceCells()[facei];
const scalar yPlus = Cmu25*y[facei]*sqrt(k[celli])/nuw[facei];
const scalar w = cornerWeights[facei];
const scalar omegaVis = 6*nuw[facei]/(beta1_*sqr(y[facei]));
// Contribution from the viscous sublayer
const scalar omegaVis = 6.0*nuw[facei]/(beta1_*sqr(y[facei]));
// Contribution from the inertial sublayer
const scalar omegaLog = sqrt(k[celli])/(Cmu25*nutw.kappa()*y[facei]);
// Switching between the laminar sub-layer and the log-region rather
// than blending has been found to provide more accurate results over a
// range of near-wall y+.
//
// For backward-compatibility the blending method is provided as an
// option
// Generation contribution is included using the blended option, or
// when using the switching option if operating in the laminar
// sub-layer
bool includeG = true;
if (blended_)
switch (blending_)
{
omega0[celli] += w*sqrt(sqr(omegaVis) + sqr(omegaLog));
}
else
{
if (nutw.yPlusLam() < yPlus)
case blendingType::STEPWISE:
{
omega0[celli] += w*omegaLog;
if (yPlus > nutw.yPlusLam())
{
omega0[celli] += w*omegaLog;
}
else
{
omega0[celli] += w*omegaVis;
}
break;
}
else
case blendingType::MAX:
{
omega0[celli] += w*omegaVis;
includeG = false;
// (PH:Eq. 27)
omega0[celli] += max(omegaVis, omegaLog);
break;
}
case blendingType::BINOMIAL2:
{
// (ME:Eq. 15)
omega0[celli] += w*sqrt(sqr(omegaVis) + sqr(omegaLog));
break;
}
case blendingType::BINOMIAL:
{
omega0[celli] +=
w*pow
(
pow(omegaVis, n_) + pow(omegaLog, n_),
1.0/n_
);
break;
}
case blendingType::EXPONENTIAL:
{
// (PH:Eq. 31)
const scalar Gamma = 0.01*pow4(yPlus)/(1.0 + 5.0*yPlus);
const scalar invGamma = 1.0/(Gamma + ROOTVSMALL);
omega0[celli] +=
w*(omegaVis*exp(-Gamma) + omegaLog*exp(-invGamma));
break;
}
case blendingType::TANH:
{
// (KAS:Eqs. 33-34)
const scalar phiTanh = tanh(pow4(yPlus/10.0));
const scalar omegab1 = omegaVis + omegaLog;
const scalar omegab2 =
pow(pow(omegaVis, 1.2) + pow(omegaLog, 1.2), 1.0/1.2);
omega0[celli] += phiTanh*omegab1 + (1.0 - phiTanh)*omegab2;
break;
}
}
if (includeG)
if (!(blending_ == blendingType::STEPWISE) || yPlus > nutw.yPlusLam())
{
G0[celli] +=
w
@ -270,7 +318,8 @@ Foam::omegaWallFunctionFvPatchScalarField::omegaWallFunctionFvPatchScalarField
)
:
fixedValueFvPatchField<scalar>(p, iF),
blended_(true),
blending_(blendingType::BINOMIAL2),
n_(2.0),
initialised_(false),
master_(-1),
beta1_(0.075),
@ -289,7 +338,8 @@ Foam::omegaWallFunctionFvPatchScalarField::omegaWallFunctionFvPatchScalarField
)
:
fixedValueFvPatchField<scalar>(ptf, p, iF, mapper),
blended_(ptf.blended_),
blending_(ptf.blending_),
n_(ptf.n_),
initialised_(false),
master_(-1),
beta1_(ptf.beta1_),
@ -307,7 +357,24 @@ Foam::omegaWallFunctionFvPatchScalarField::omegaWallFunctionFvPatchScalarField
)
:
fixedValueFvPatchField<scalar>(p, iF, dict),
blended_(dict.getOrDefault("blended", true)),
blending_
(
blendingTypeNames.getOrDefault
(
"blending",
dict,
blendingType::BINOMIAL2
)
),
n_
(
dict.getCheckOrDefault<scalar>
(
"n",
2.0,
scalarMinMax::ge(0)
)
),
initialised_(false),
master_(-1),
beta1_(dict.getOrDefault<scalar>("beta1", 0.075)),
@ -315,6 +382,29 @@ Foam::omegaWallFunctionFvPatchScalarField::omegaWallFunctionFvPatchScalarField
omega_(),
cornerWeights_()
{
// The deprecated 'blended' keyword is superseded by the enum 'blending'
if (dict.found("blended"))
{
IOWarningInFunction(dict)
<< "Using deprecated 'blended' keyword"
<< nl << " Please use either of the below for the same behaviour:"
<< nl << " 'blending binomial2;' for 'blended on;'"
<< nl << " 'blending stepwise;' for 'blended off;'"
<< nl << " OVERWRITING: 'blended' keyword -> 'blending' enum"
<< endl;
bool blended = dict.get<bool>("blended");
if (blended)
{
blending_ = blendingType::BINOMIAL2;
}
else
{
blending_ = blendingType::STEPWISE;
}
}
// apply zero-gradient condition on start-up
this->operator==(patchInternalField());
}
@ -326,7 +416,8 @@ Foam::omegaWallFunctionFvPatchScalarField::omegaWallFunctionFvPatchScalarField
)
:
fixedValueFvPatchField<scalar>(owfpsf),
blended_(owfpsf.blended_),
blending_(owfpsf.blending_),
n_(owfpsf.n_),
initialised_(false),
master_(-1),
beta1_(owfpsf.beta1_),
@ -343,7 +434,8 @@ Foam::omegaWallFunctionFvPatchScalarField::omegaWallFunctionFvPatchScalarField
)
:
fixedValueFvPatchField<scalar>(owfpsf, iF),
blended_(owfpsf.blended_),
blending_(owfpsf.blending_),
n_(owfpsf.n_),
initialised_(false),
master_(-1),
beta1_(owfpsf.beta1_),
@ -481,7 +573,7 @@ void Foam::omegaWallFunctionFvPatchScalarField::updateWeightedCoeffs
{
const scalar w = weights[facei];
if (tolerance_ < w)
if (w > tolerance_)
{
const label celli = patch().faceCells()[facei];
@ -559,7 +651,9 @@ void Foam::omegaWallFunctionFvPatchScalarField::write
Ostream& os
) const
{
os.writeEntry("blended", blended_);
os.writeKeyword("blending") << blendingTypeNames[blending_]
<< token::END_STATEMENT << nl;
os.writeEntry("n", n_);
os.writeEntry("beta1", beta1_);
fixedValueFvPatchField<scalar>::write(os);
}

View File

@ -6,7 +6,7 @@
\\/ M anipulation |
-------------------------------------------------------------------------------
Copyright (C) 2011-2016, 2019 OpenFOAM Foundation
Copyright (C) 2019 OpenCFD Ltd.
Copyright (C) 2019-2020 OpenCFD Ltd.
-------------------------------------------------------------------------------
License
This file is part of OpenFOAM.
@ -31,68 +31,167 @@ Group
grpWallFunctions
Description
This boundary condition provides a wall constraint on specific turbulent
kinetic energy dissipation rate, i.e. \c omega, for low- and high-Reynolds
number turbulence models.
This boundary condition provides a wall constraint on the specific
dissipation rate, i.e. \c omega, and the turbulent kinetic energy
production contribution, i.e. \c G, for low- and high-Reynolds number
turbulence models.
The near-wall omega may be either blended between the viscous region and
logarithmic region values using:
\f[
\omega = sqrt(\omega_{vis}^2 + \omega_{log}^2)
\f]
where
\vartable
\omega_{vis} | omega in viscous region
\omega_{log} | omega in logarithmic region
\endvartable
see eq.(15) of:
Reference:
\verbatim
Menter, F., Esch, T.
"Elements of Industrial Heat Transfer Prediction"
16th Brazilian Congress of Mechanical Engineering (COBEM),
Nov. 2001
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
\endverbatim
or switched between these values based on the viscous-to-turbulent \c y+
value derived from \c kappa and \c E.
Usage
\table
Property | Description | Required | Default value
beta1 | Model coefficient | no | 0.075
blended | Blending switch | no | false
\endtable
Example of the boundary condition specification:
\verbatim
<patchName>
{
// Mandatory entries
// Mandatory entries (unmodifiable)
type omegaWallFunction;
// Optional entries
// Optional entries (unmodifiable)
beta1 0.075;
blending binomial2;
n 2.0;
// Optional (inherited) entries
...
}
\endverbatim
\table
Property | Description | Type | Req'd | Dflt
type | Type name: omegaWallFunction | word | yes | -
blended | Blending switch (Deprecated) | bool | no | true
beta1 | Model coefficient | scalar | no | 0.075
blending | Viscous/inertial sublayer blending method <!--
--> | word | no | binomial2
n | Binomial blending exponent | sclar | no | 2.0
\endtable
The inherited entries are elaborated in:
- \link fixedValueFvPatchField.H \endlink
- \link nutWallFunctionFvPatchScalarField.H \endlink
Options for the \c blending entry:
\verbatim
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)
\endverbatim
wherein \c omega predictions for the viscous and inertial sublayers are
blended according to the following expressions:
- \c stepwise:
\f[
\omega = \omega_{log} \qquad if \quad y^+ > y^+_{lam}
\f]
\f[
\omega = \omega_{vis} \qquad if \quad y^+ <= y^+_{lam}
\f]
where
\vartable
\omega | \f$\omega\f$ at \f$y^+\f$
\omega_{vis} | \f$\omega\f$ computed by using viscous sublayer assumptions
\omega_{log} |\f$\omega\f$ 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
\endvartable
- \c max (PH:Eq. 27):
\f[
\omega = max(\omega_{vis}, \omega_{log})
\f]
- \c binomial2 (ME:Eq. 15) (default):
\f[
\omega = \sqrt{(\omega_{vis})^2 + (\omega_{log})^2}
\f]
- \c binomial:
\f[
\omega = ((\omega_{vis})^n + (\omega_{log})^n)^{1/n}
\f]
where
\vartable
n | Binomial blending exponent
\endvartable
- \c exponential (PH:Eq. 32):
\f[
\omega = \omega_{vis} \exp[-\Gamma] + \omega_{log} \exp[-1/\Gamma]
\f]
where (PH:Eq. 31)
\vartable
\Gamma | Blending expression
\Gamma | \f$0.01 (y^+)^4 / (1.0 + 5.0 y^+)\f$
\endvartable
- \c tanh (KAS:Eqs. 33-34):
\f[
\omega = \phi \omega_{b1} + (1 - \phi)\omega_{b2}
\f]
where
\vartable
\phi | \f$tanh((y^+/10)^4)\f$
\omega_{b1} | \f$\omega_{vis} + \omega_{log}\f$
\omega_{b2} | \f$(\omega_{vis}^{1.2} + \omega_{log}^1.2)^{1/1.2}\f$
\endvartable
\c G predictions for the viscous and inertial sublayers are blended
in a stepwise manner, and \c G below \f$y^+_{lam}\f$ (i.e. in the viscous
sublayer) is presumed to be zero.
Note
The coefficients \c Cmu, \c kappa, and \c E are obtained from
- The coefficients \c Cmu, \c kappa, and \c E are obtained from
the specified \c nutWallFunction in order to ensure that each patch
possesses the same set of values for these coefficients.
Some tests have shown that the stepwise switching method provides
more accurate predictions for 10 < y+ < 30 when used with high Reynolds
number wall-functions.
In addition, the stepwise switching method provides accurate predictions
when used with continuous wall-functions.
- The reason why \c binomial2 and \c binomial blending methods exist at
the same time is to ensure the elementwise backward compatibility with
the previous versions since \c binomial2 and \c 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
Foam::fixedInternalValueFvPatchField
Foam::epsilonWallFunctionFvPatchScalarField
- Foam::epsilonWallFunctionFvPatchScalarField
SourceFiles
omegaWallFunctionFvPatchScalarField.C
@ -119,6 +218,35 @@ class omegaWallFunctionFvPatchScalarField
:
public fixedValueFvPatchField<scalar>
{
private:
// Private Enumerations
//- Options for the blending treatment of viscous and inertial sublayers
enum blendingType
{
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)"
};
//- Names for blendingType
static const Enum<blendingType> blendingTypeNames;
// Private Data
//- Blending treatment (default = blendingType::BINOMIAL2)
enum blendingType blending_;
//- Blending exponent being used when
//- blendingType is blendingType::BINOMIAL (default = 2)
const scalar n_;
protected:
// Protected Data
@ -126,7 +254,8 @@ protected:
//- Tolerance used in weighted calculations
static scalar tolerance_;
//- Blending switch
//- Deprecated(2019-11) Blending switch
// \deprecated(2019-11) - use blending:: options
bool blended_;
//- Initialised flag