Go to the documentation of this file.
51 { blendingType::STEPWISE ,
"stepwise" },
52 { blendingType::MAX ,
"max" },
53 { blendingType::BINOMIAL ,
"binomial" },
54 { blendingType::EXPONENTIAL,
"exponential" }
62 if (!isA<wallFvPatch>(
patch()))
65 <<
"Invalid wall function specification" <<
nl
66 <<
" Patch type for patch " <<
patch().name()
67 <<
" must be wall" <<
nl
68 <<
" Current patch type is " <<
patch().type() <<
nl <<
endl
97 if (blending_ == blendingType::BINOMIAL)
121 fixedValueFvPatchScalarField(
p, iF),
122 blending_(blendingType::STEPWISE),
128 yPlusLam_(yPlusLam(kappa_, E_))
142 fixedValueFvPatchScalarField(ptf,
p, iF, mapper),
162 fixedValueFvPatchScalarField(
p, iF,
dict),
165 blendingTypeNames.getOrDefault
169 blendingType::STEPWISE
188 yPlusLam_(yPlusLam(kappa_, E_))
199 fixedValueFvPatchScalarField(wfpsf),
218 fixedValueFvPatchScalarField(wfpsf, iF),
241 refCast<const nutWallFunctionFvPatchScalarField>
243 turbModel.
nut()().boundaryField()[patchi],
257 for (label i = 0; i < 10; ++i)
277 case blendingType::STEPWISE:
279 if (
yPlus > yPlusLam_)
290 case blendingType::MAX:
293 nutw =
max(nutVis, nutLog);
297 case blendingType::BINOMIAL:
303 pow(nutVis, n_) +
pow(nutLog, n_),
309 case blendingType::EXPONENTIAL:
313 const scalar invGamma = 1.0/(Gamma + ROOTVSMALL);
315 nutw = nutVis*
exp(-Gamma) + nutLog*
exp(-invGamma);
339 fixedValueFvPatchScalarField::updateCoeffs();
349 writeLocalEntries(
os);
350 writeEntry(
"value",
os);
nutWallFunctionFvPatchScalarField(const fvPatch &, const DimensionedField< scalar, volMesh > &)
Construct from patch and internal field.
virtual void write(Ostream &) const
Write.
blendingType
Options for the blending treatment of viscous and inertial sublayers.
Enum is a wrapper around a list of names/values that represent particular enumeration (or int) values...
A class for handling words, derived from Foam::string.
virtual void updateCoeffs()
Update the coefficients associated with the patch field.
virtual tmp< volScalarField > nut() const =0
Return the turbulence viscosity.
virtual const volVectorField & U(const turbulenceModel &turb) const
Ostream & endl(Ostream &os)
Add newline and flush stream.
dimensionedScalar exp(const dimensionedScalar &ds)
const dimensionedScalar kappa
Coulomb constant: default SI units: [N.m2/C2].
scalar kappa_
von Kármán constant
dimensionedScalar pow4(const dimensionedScalar &ds)
tmp< faMatrix< Type > > operator==(const faMatrix< Type > &, const faMatrix< Type > &)
virtual void write(Ostream &) const
Write.
A finiteVolume patch using a polyPatch and a fvBoundaryMesh.
enum blendingType blending_
Blending treatment (default = blendingType::STEPWISE)
T getCheckOrDefault(const word &keyword, const T &deflt, const Predicate &pred, enum keyType::option matchOpt=keyType::REGEX) const
dimensionedScalar pow(const dimensionedScalar &ds, const dimensionedScalar &expt)
label max(const labelHashSet &set, label maxValue=labelMin)
Find the max value in labelHashSet, optionally limited by second argument.
Abstract base class for turbulence models (RAS, LES and laminar).
A list of keyword definitions, which are a keyword followed by a number of values (eg,...
scalar E_
Wall roughness parameter.
dimensionedScalar log(const dimensionedScalar &ds)
OBJstream os(runTime.globalPath()/outputName)
Macros for easy insertion into run-time selection tables.
scalar Cmu_
Empirical model coefficient.
word UName_
Name of velocity field.
errorManip< error > abort(error &err)
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
scalar blend(const scalar nutVis, const scalar nutLog, const scalar yPlus) const
Return the blended nut according to the chosen blending treatment.
static const Enum< blendingType > blendingTypeNames
Names for blendingType.
const std::string patch
OpenFOAM patch number as a std::string.
static MinMax< T > ge(const T &minVal)
A semi-infinite range from minVal to the type max.
static const nutWallFunctionFvPatchScalarField & nutw(const turbulenceModel &turbModel, const label patchi)
Return the nut patchField for the given wall patch.
static const word null
An empty word.
Ostream & writeEntry(const keyType &key, const T &value)
Write a keyword/value entry.
Foam::fvPatchFieldMapper.
An Ostream is an abstract base class for all output systems (streams, files, token lists,...
T getOrDefault(const word &keyword, const T &deflt, enum keyType::option matchOpt=keyType::REGEX) const
The class nutWallFunction is a base class that parents the derived boundary conditions which provide ...
compressible::turbulenceModel & turb
defineTypeNameAndDebug(combustionModel, 0)
virtual void writeLocalEntries(Ostream &) const
Write local wall function variables.
scalar yPlusLam() const
Return the estimated y+ at the two-sublayer intersection.
Field with dimensions and associated with geometry type GeoMesh which is used to size the field and a...
virtual void checkType()
Check the type of the patch.