Commit 458738f0 authored by Kutalmis Bercin's avatar Kutalmis Bercin Committed by Andrew Heather
Browse files

ENH: Subtle renovations in wall function documentation and input

  - Improves header file documentations
  - Adds const specifier to various objects
  - Allows various model constants to be user-defined
  - Changes from lookupOrDefault() to getOrDefault()
  - Consistent namespace usage:
    - If the WF belongs to only Foam:: namespace, use Foam:: explicitly
    - If the WF belongs to more than one namespaces, put the WF in namespace
      parentheses
  - Adds the missing dashes in comments required by Doxygen
  - Corrects capitalisation in comments
  - The terminology 'viscous sublayer' was preferred over 'laminar sublayer'
    due to the fact that the sublayer is not laminar.
  - The order of initialisation of fields is changed to be the same across all.
parent ce9ea1ca
......@@ -2,10 +2,10 @@
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration |
\\ / A nd | Copyright (C) 2017 OpenCFD Ltd.
\\ / A nd | Copyright (C) 2017-2019 OpenCFD Ltd.
\\/ M anipulation |
-------------------------------------------------------------------------------
| Copyright (C) 2011-2016, 2019 OpenFOAM Foundation
| Copyright (C) 2011-2019 OpenFOAM Foundation
-------------------------------------------------------------------------------
License
This file is part of OpenFOAM.
......@@ -31,6 +31,7 @@ License
#include "fvMatrix.H"
#include "addToRunTimeSelectionTable.H"
// * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
Foam::scalar Foam::epsilonWallFunctionFvPatchScalarField::tolerance_ = 1e-5;
......@@ -127,7 +128,10 @@ void Foam::epsilonWallFunctionFvPatchScalarField::createAveragingWeights()
Foam::epsilonWallFunctionFvPatchScalarField&
Foam::epsilonWallFunctionFvPatchScalarField::epsilonPatch(const label patchi)
Foam::epsilonWallFunctionFvPatchScalarField::epsilonPatch
(
const label patchi
)
{
const volScalarField& epsilon =
static_cast<const volScalarField&>(this->internalField());
......@@ -244,11 +248,11 @@ epsilonWallFunctionFvPatchScalarField
)
:
fixedValueFvPatchField<scalar>(p, iF),
G_(),
epsilon_(),
lowReCorrection_(false),
initialised_(false),
master_(-1),
G_(),
epsilon_(),
cornerWeights_()
{}
......@@ -263,11 +267,11 @@ epsilonWallFunctionFvPatchScalarField
)
:
fixedValueFvPatchField<scalar>(ptf, p, iF, mapper),
G_(),
epsilon_(),
lowReCorrection_(ptf.lowReCorrection_),
initialised_(false),
master_(-1),
G_(),
epsilon_(),
cornerWeights_()
{}
......@@ -281,11 +285,11 @@ epsilonWallFunctionFvPatchScalarField
)
:
fixedValueFvPatchField<scalar>(p, iF, dict),
G_(),
epsilon_(),
lowReCorrection_(dict.lookupOrDefault("lowReCorrection", false)),
lowReCorrection_(dict.getOrDefault("lowReCorrection", false)),
initialised_(false),
master_(-1),
G_(),
epsilon_(),
cornerWeights_()
{
// Apply zero-gradient condition on start-up
......@@ -300,11 +304,11 @@ epsilonWallFunctionFvPatchScalarField
)
:
fixedValueFvPatchField<scalar>(ewfpsf),
G_(),
epsilon_(),
lowReCorrection_(ewfpsf.lowReCorrection_),
initialised_(false),
master_(-1),
G_(),
epsilon_(),
cornerWeights_()
{}
......@@ -317,18 +321,21 @@ epsilonWallFunctionFvPatchScalarField
)
:
fixedValueFvPatchField<scalar>(ewfpsf, iF),
G_(),
epsilon_(),
lowReCorrection_(ewfpsf.lowReCorrection_),
initialised_(false),
master_(-1),
G_(),
epsilon_(),
cornerWeights_()
{}
// * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
Foam::scalarField& Foam::epsilonWallFunctionFvPatchScalarField::G(bool init)
Foam::scalarField& Foam::epsilonWallFunctionFvPatchScalarField::G
(
bool init
)
{
if (patch().index() == master_)
{
......@@ -398,7 +405,7 @@ void Foam::epsilonWallFunctionFvPatchScalarField::updateCoeffs()
forAll(*this, facei)
{
label celli = patch().faceCells()[facei];
const label celli = patch().faceCells()[facei];
G[celli] = G0[celli];
epsilon[celli] = epsilon0[celli];
......@@ -449,11 +456,11 @@ void Foam::epsilonWallFunctionFvPatchScalarField::updateWeightedCoeffs
// Only set the values if the weights are > tolerance
forAll(weights, facei)
{
scalar w = weights[facei];
const scalar w = weights[facei];
if (w > tolerance_)
if (tolerance_ < w)
{
label celli = patch().faceCells()[facei];
const label celli = patch().faceCells()[facei];
G[celli] = (1.0 - w)*G[celli] + w*G0[celli];
epsilon[celli] = (1.0 - w)*epsilon[celli] + w*epsilon0[celli];
......@@ -501,7 +508,7 @@ void Foam::epsilonWallFunctionFvPatchScalarField::manipulateMatrix
forAll(weights, facei)
{
// Only set the values if the weights are > tolerance
if (weights[facei] > tolerance_)
if (tolerance_ < weights[facei])
{
const label celli = faceCells[facei];
......@@ -524,6 +531,16 @@ void Foam::epsilonWallFunctionFvPatchScalarField::manipulateMatrix
}
void Foam::epsilonWallFunctionFvPatchScalarField::write
(
Ostream& os
) const
{
os.writeEntry("lowReCorrection", lowReCorrection_);
fixedValueFvPatchField<scalar>::write(os);
}
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
namespace Foam
......
......@@ -30,32 +30,31 @@ Group
grpWallFunctions
Description
This boundary condition provides a turbulence dissipation wall constraint
for low- and high-Reynolds number turbulence models.
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.
The condition can be applied to wall boundaries for which it
- calculates \c epsilon and \c G
- specifies the near-wall epsilon value
- specifies the near-wall \c epsilon value
where
\vartable
epsilon | turblence dissipation field
G | turblence generation field
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
laminar and turbulent functions based on the laminar-to-turbulent y+ value
derived from the kappa and E specified in the corresponding nutWallFunction. When the \c lowReCorrection is inactive, the
wall function operates in high-Re mode.
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.
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
lowReCorrection | Low-Re correction active | no | off
\endtable
......@@ -63,13 +62,20 @@ Usage
\verbatim
<patchName>
{
// Mandatory entries
type epsilonWallFunction;
// Optional entries
}
\endverbatim
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::fixedInternalValueFvPatchField
Foam::omegaWallFunctionFvPatchScalarField
SourceFiles
epsilonWallFunctionFvPatchScalarField.C
......@@ -98,17 +104,11 @@ class epsilonWallFunctionFvPatchScalarField
{
protected:
// Protected data
// Protected Data
//- Tolerance used in weighted calculations
static scalar tolerance_;
//- Local copy of turbulence G field
scalarField G_;
//- Local copy of turbulence epsilon field
scalarField epsilon_;
//- Apply low-Re correction term; default = no
bool lowReCorrection_;
......@@ -118,6 +118,12 @@ protected:
//- Master patch ID
label master_;
//- Local copy of turbulence G field
scalarField G_;
//- Local copy of turbulence epsilon field
scalarField epsilon_;
//- List of averaging corner weights
List<List<scalar>> cornerWeights_;
......@@ -125,11 +131,11 @@ protected:
// Protected Member Functions
//- Set the master patch - master is responsible for updating all
// wall function patches
//- wall function patches
virtual void setMaster();
//- Create the averaging weights for cells which are bounded by
// multiple wall function faces
//- multiple wall function faces
virtual void createAveragingWeights();
//- Helper function to return non-const access to an epsilon patch
......@@ -187,8 +193,8 @@ public:
);
//- Construct by mapping given
// epsilonWallFunctionFvPatchScalarField
// onto a new patch
//- epsilonWallFunctionFvPatchScalarField
//- onto a new patch
epsilonWallFunctionFvPatchScalarField
(
const epsilonWallFunctionFvPatchScalarField&,
......@@ -231,6 +237,7 @@ public:
);
}
//- Destructor
virtual ~epsilonWallFunctionFvPatchScalarField() = default;
......@@ -263,6 +270,12 @@ public:
fvMatrix<scalar>& matrix,
const scalarField& weights
);
// I-O
//- Write
virtual void write(Ostream&) const;
};
......
......@@ -2,7 +2,7 @@
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration |
\\ / A nd |
\\ / A nd | Copyright (C) 2019 OpenCFD Ltd.
\\/ M anipulation |
-------------------------------------------------------------------------------
| Copyright (C) 2012-2016, 2019 OpenFOAM Foundation
......@@ -117,6 +117,9 @@ void fWallFunctionFvPatchScalarField::updateCoeffs()
const scalarField& y = turbModel.y()[patchi];
const tmp<scalarField> tnuw = turbModel.nu(patchi);
const scalarField& nuw = tnuw();
const tmp<volScalarField> tk = turbModel.k();
const volScalarField& k = tk();
......@@ -126,28 +129,25 @@ void fWallFunctionFvPatchScalarField::updateCoeffs()
const tmp<volScalarField> tv2 = v2fModel.v2();
const volScalarField& v2 = tv2();
const tmp<scalarField> tnuw = turbModel.nu(patchi);
const scalarField& nuw = tnuw();
const scalar Cmu25 = pow025(nutw.Cmu());
const scalar N = 6.0;
scalarField& f = *this;
// Set f wall values
forAll(f, facei)
{
label celli = patch().faceCells()[facei];
const label celli = patch().faceCells()[facei];
scalar uTau = Cmu25*sqrt(k[celli]);
const scalar uTau = Cmu25*sqrt(k[celli]);
scalar yPlus = uTau*y[facei]/nuw[facei];
const scalar yPlus = uTau*y[facei]/nuw[facei];
if (yPlus > nutw.yPlusLam())
if (nutw.yPlusLam() < yPlus)
{
scalar N = 6.0;
scalar v2c = v2[celli];
scalar epsc = epsilon[celli];
scalar kc = k[celli];
const scalar v2c = v2[celli];
const scalar epsc = epsilon[celli];
const scalar kc = k[celli];
f[facei] = N*v2c*epsc/(sqr(kc) + ROOTVSMALL);
f[facei] /= sqr(uTau) + ROOTVSMALL;
......
......@@ -2,7 +2,7 @@
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration |
\\ / A nd |
\\ / A nd | Copyright (C) 2019 OpenCFD Ltd.
\\/ M anipulation |
-------------------------------------------------------------------------------
| Copyright (C) 2012-2016, 2019 OpenFOAM Foundation
......@@ -30,29 +30,44 @@ Group
grpWallFunctions
Description
This boundary condition provides a turbulence damping function, f, wall
function condition for low- and high Reynolds number, turbulent flow cases
The model operates in two modes, based on the computed laminar-to-turbulent
switch-over y+ value derived from kappa and E specified in the corresponding
nutWallFunction.
This boundary condition provides a wall constraint on the elliptic
relaxation factor, \c f, which is executed in the \c v2-f eddy viscosity
turbulence model. The condition is applicable for low- and high-Reynolds
number turbulent flow cases.
For \c f, the viscous sublayer and log-law region blending approaches are
claimed to be inviable (Popovac and Hanjalić (2007), p. 194). Therefore,
the only boundary condition blending mode is the stepwise mode
where the viscous sublayer and log-law region contributions switch over
\c y+ value derived from the \c kappa and \c E.
Reference:
\verbatim
Remark on the blending approach:
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 | 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
type fWallFunction;
// No optional entry
}
\endverbatim
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
......@@ -81,7 +96,6 @@ class fWallFunctionFvPatchScalarField
:
public fixedValueFvPatchField<scalar>
{
public:
//- Runtime type information
......@@ -106,7 +120,7 @@ public:
);
//- Construct by mapping given fWallFunctionFvPatchScalarField
// onto a new patch
//- onto a new patch
fWallFunctionFvPatchScalarField
(
const fWallFunctionFvPatchScalarField&,
......
......@@ -2,7 +2,7 @@
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration |
\\ / A nd |
\\ / A nd | Copyright (C) 2019 OpenCFD Ltd.
\\/ M anipulation |
-------------------------------------------------------------------------------
| Copyright (C) 2012-2016, 2019 OpenFOAM Foundation
......@@ -30,25 +30,24 @@ License
#include "turbulenceModel.H"
#include "addToRunTimeSelectionTable.H"
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
namespace Foam
{
// * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
kLowReWallFunctionFvPatchScalarField::kLowReWallFunctionFvPatchScalarField
Foam::kLowReWallFunctionFvPatchScalarField::kLowReWallFunctionFvPatchScalarField
(
const fvPatch& p,
const DimensionedField<scalar, volMesh>& iF
)
:
fixedValueFvPatchField<scalar>(p, iF),
Ceps2_(1.9)
Ceps2_(1.9),
Ck_(-0.416),
Bk_(8.366),
C_(11.0)
{}
kLowReWallFunctionFvPatchScalarField::kLowReWallFunctionFvPatchScalarField
Foam::kLowReWallFunctionFvPatchScalarField::kLowReWallFunctionFvPatchScalarField
(
const kLowReWallFunctionFvPatchScalarField& ptf,
const fvPatch& p,
......@@ -57,11 +56,14 @@ kLowReWallFunctionFvPatchScalarField::kLowReWallFunctionFvPatchScalarField
)
:
fixedValueFvPatchField<scalar>(ptf, p, iF, mapper),
Ceps2_(ptf.Ceps2_)
Ceps2_(ptf.Ceps2_),
Ck_(ptf.Ck_),
Bk_(ptf.Bk_),
C_(ptf.C_)
{}
kLowReWallFunctionFvPatchScalarField::kLowReWallFunctionFvPatchScalarField
Foam::kLowReWallFunctionFvPatchScalarField::kLowReWallFunctionFvPatchScalarField
(
const fvPatch& p,
const DimensionedField<scalar, volMesh>& iF,
......@@ -69,34 +71,43 @@ kLowReWallFunctionFvPatchScalarField::kLowReWallFunctionFvPatchScalarField
)
:
fixedValueFvPatchField<scalar>(p, iF, dict),
Ceps2_(dict.lookupOrDefault<scalar>("Ceps2", 1.9))
Ceps2_(dict.getOrDefault<scalar>("Ceps2", 1.9)),
Ck_(dict.getOrDefault<scalar>("Ck", -0.416)),
Bk_(dict.getOrDefault<scalar>("Bk", 8.366)),
C_(dict.getOrDefault<scalar>("C", 11.0))
{}
kLowReWallFunctionFvPatchScalarField::kLowReWallFunctionFvPatchScalarField
Foam::kLowReWallFunctionFvPatchScalarField::kLowReWallFunctionFvPatchScalarField
(
const kLowReWallFunctionFvPatchScalarField& kwfpsf
)
:
fixedValueFvPatchField<scalar>(kwfpsf),
Ceps2_(kwfpsf.Ceps2_)
Ceps2_(kwfpsf.Ceps2_),
Ck_(kwfpsf.Ck_),
Bk_(kwfpsf.Bk_),
C_(kwfpsf.C_)
{}
kLowReWallFunctionFvPatchScalarField::kLowReWallFunctionFvPatchScalarField
Foam::kLowReWallFunctionFvPatchScalarField::kLowReWallFunctionFvPatchScalarField
(
const kLowReWallFunctionFvPatchScalarField& kwfpsf,
const DimensionedField<scalar, volMesh>& iF
)
:
fixedValueFvPatchField<scalar>(kwfpsf, iF),
Ceps2_(kwfpsf.Ceps2_)
Ceps2_(kwfpsf.Ceps2_),
Ck_(kwfpsf.Ck_),
Bk_(kwfpsf.Bk_),
C_(kwfpsf.C_)
{}
// * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //