Go to the documentation of this file.
40 namespace incompressible
52 if (!isA<wallFvPatch>(
patch()))
55 <<
"Invalid wall function specification" <<
nl
56 <<
" Patch type for patch " <<
patch().name()
57 <<
" must be wall" <<
nl
58 <<
" Current patch type is " <<
patch().type() <<
nl <<
endl
69 const label patchi =
patch().index();
73 if (isA<nutWallFunctionFvPatchScalarField>(
nut.boundaryField()[patchi]))
78 nut.boundaryField()[patchi]
90 /turbModel.
nu(patchi);
100 return 9.24*(
pow(Prat, 0.75) - 1.0)*(1.0 + 0.28*
exp(-0.007*Prat));
112 for (
int i=0; i<maxIters_; i++)
114 scalar
f = ypt - (
log(E_*ypt)/kappa_ + P)/Prat;
115 scalar df = 1.0 - 1.0/(ypt*kappa_*Prat);
116 scalar yptNew = ypt -
f/df;
122 else if (
mag(yptNew - ypt) < tolerance_)
145 fixedValueFvPatchScalarField(
p, iF),
163 fixedValueFvPatchScalarField(ptf,
p, iF, mapper),
180 fixedValueFvPatchScalarField(
p, iF,
dict),
181 Prt_(
dict.get<scalar>(
"Prt")),
182 kappa_(
dict.getOrDefault<scalar>(
"kappa", 0.41)),
183 E_(
dict.getOrDefault<scalar>(
"E", 9.8))
195 fixedValueFvPatchScalarField(wfpsf),
211 fixedValueFvPatchScalarField(wfpsf, iF),
229 const label patchi =
patch().index();
238 internalField().
group()
261 scalar
yPlus = yPlusp[facei];
273 scalar
nu = nuw[facei];
275 alphatw[facei] =
max(0.0, kt);
279 alphatw[facei] = 0.0;
293 writeEntry(
"value", os);
virtual void checkType()
Check the type of the patch.
IOdictionary is derived from dictionary and IOobject to give the dictionary automatic IO functionalit...
virtual void write(Ostream &) const
Write.
virtual tmp< Field< Type > > snGrad() const
Return patch-normal gradient.
const dimensionSet dimless(0, 0, 0, 0, 0, 0, 0)
Dimensionless.
constexpr const char *const group
Group name for atomic constants.
A class for managing temporary objects.
virtual tmp< volScalarField > nut() const =0
Return the turbulence viscosity.
IOdictionary transportProperties(IOobject("transportProperties", runTime.constant(), mesh, IOobject::MUST_READ, IOobject::NO_WRITE))
static const word propertiesName
Default name of the turbulence properties dictionary.
virtual tmp< scalarField > yPlus() const =0
Calculate and return the yPlus at the boundary.
Ostream & endl(Ostream &os)
Add newline and flush stream.
dimensionedScalar exp(const dimensionedScalar &ds)
const nearWallDist & y() const
Return the near wall distances.
#define forAll(list, i)
Loop across all elements in list.
tmp< scalarField > yPlus(const turbulenceModel &turbModel) const
Return the patch y+.
virtual void updateCoeffs()
Update the coefficients associated with the patch field.
This boundary condition provides a kinematic turbulent thermal conductivity for using wall functions,...
scalar Psmooth(const scalar Prat) const
`P' function
virtual void write(Ostream &) const
Write.
alphatJayatillekeWallFunctionFvPatchScalarField(const fvPatch &, const DimensionedField< scalar, volMesh > &)
Construct from patch and internal field.
scalar Prt_
Turbulent Prandtl number.
A finiteVolume patch using a polyPatch and a fvBoundaryMesh.
dimensionedScalar Pr("Pr", dimless, laminarTransport)
dimensioned< scalar > dimensionedScalar
Dimensioned scalar obtained from generic dimensioned type.
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.
A list of keyword definitions, which are a keyword followed by a number of values (eg,...
scalar yPlusTherm(const scalar P, const scalar Prat) const
Calculate y+ at the edge of the thermal laminar sublayer.
dimensionedScalar log(const dimensionedScalar &ds)
Macros for easy insertion into run-time selection tables.
errorManip< error > abort(error &err)
virtual void updateCoeffs()
Update the coefficients associated with the patch field.
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
scalar kappa_
Von Karman constant.
virtual tmp< volScalarField > nu() const
Return the laminar viscosity.
const std::string patch
OpenFOAM patch number as a std::string.
dimensionedScalar sqrt(const dimensionedScalar &ds)
dimensioned< typename typeOfMag< Type >::type > mag(const dimensioned< Type > &dt)
Ostream & writeEntry(const keyType &key, const T &value)
Write a keyword/value entry.
Templated abstract base class for single-phase incompressible turbulence models.
Foam::fvPatchFieldMapper.
static word groupName(StringType base, const word &group)
Create dot-delimited name.group string.
An Ostream is an abstract base class for all output systems (streams, files, token lists,...
makePatchTypeField(fvPatchScalarField, alphatJayatillekeWallFunctionFvPatchScalarField)
virtual tmp< volScalarField > nuEff() const =0
Return the effective viscosity.
The class nutWallFunction is a base class that parents the derived boundary conditions which provide ...
const Boundary & boundaryField() const
Return const-reference to the boundary field.
Field with dimensions and associated with geometry type GeoMesh which is used to size the field and a...
const volVectorField & U() const
Access function to velocity field.