Go to the documentation of this file.
56 const auto* filmModelPtr = db().time().findObject
69 const auto& filmModel = *filmModelPtr;
71 const label filmPatchi = filmModel.regionPatchID(patchi);
74 scalarField mDotFilmp = mDotFilm().boundaryField()[filmPatchi];
75 filmModel.toPrimary(filmPatchi, mDotFilmp);
84 internalField().
group()
94 const scalar Cmu25 =
pow(Cmu_, 0.25);
100 scalar ut = Cmu25*
sqrt(
k[faceCelli]);
102 scalar
yPlus =
y[facei]*ut/nuw[facei];
104 scalar mStar = mDotFilmp[facei]/(
y[facei]*ut);
107 if (
yPlus > yPlusCrit_)
109 scalar expTerm =
exp(
min(50.0, B_*mStar));
110 scalar powTerm =
pow(
yPlus, mStar/kappa_);
111 factor = mStar/(expTerm*powTerm - 1.0 + ROOTVSMALL);
115 scalar expTerm =
exp(
min(50.0, mStar));
116 factor = mStar/(expTerm*
yPlus - 1.0 + ROOTVSMALL);
119 uTau[facei] =
sqrt(
max(0, magGradU[facei]*ut*factor));
135 internalField().
group()
147 sqr(
calcUTau(magGradU))/(magGradU + ROOTVSMALL) - nuw
161 filmRegionName_(
"surfaceFilmProperties"),
192 dict.lookupOrDefault<
word>(
"filmRegion",
"surfaceFilmProperties")
194 B_(
dict.lookupOrDefault(
"B", 5.5)),
195 yPlusCrit_(
dict.lookupOrDefault(
"yPlusCrit", 11.05))
235 internalField().
group()
254 "surfaceFilmProperties",
260 writeEntry(
"value", os);
virtual void write(Ostream &) const
Write.
virtual tmp< Field< Type > > snGrad() const
Return patch-normal gradient.
Ostream & writeEntryIfDifferent(const word &key, const T &value1, const T &value2)
Write a keyword/value entry only when the two values differ.
Field< scalar > scalarField
Specialisation of Field<T> for scalar.
virtual void write(Ostream &os) const
Write.
A class for handling words, derived from Foam::string.
scalar yPlusCrit_
y+ value for laminar -> turbulent transition (default = 11.05)
A class for managing temporary objects.
static constexpr const zero Zero
Global zero.
virtual tmp< scalarField > yPlus() const
Calculate and return the yPlus at the boundary.
virtual tmp< scalarField > calcNut() const
Calculate the turbulence viscosity.
scalar B_
B Coefficient (default = 5.5)
makePatchTypeField(fvPatchScalarField, alphatFilmWallFunctionFvPatchScalarField)
static const word propertiesName
Default name of the turbulence properties dictionary.
const char *const group
Group name for atomic constants.
nutkFilmWallFunctionFvPatchScalarField(const fvPatch &, const DimensionedField< scalar, volMesh > &)
Construct from patch and internal field.
dimensionedScalar exp(const dimensionedScalar &ds)
label min(const labelHashSet &set, label minValue=labelMax)
Find the min value in labelHashSet, optionally limited by second argument.
#define forAll(list, i)
Loop across all elements in list.
intWM_LABEL_SIZE_t label
A label is an int32_t or int64_t as specified by the pre-processor macro WM_LABEL_SIZE.
Base class for surface film models.
A finiteVolume patch using a polyPatch and a fvBoundaryMesh.
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.
virtual tmp< scalarField > calcUTau(const scalarField &magGradU) const
Calculate the friction velocity.
word filmRegionName_
Name of film region.
A list of keyword definitions, which are a keyword followed by a number of values (eg,...
Macros for easy insertion into run-time selection tables.
dimensionedSymmTensor sqr(const dimensionedVector &dv)
const std::string patch
OpenFOAM patch number as a std::string.
Templated wrapper class to provide compressible turbulence models thermal diffusivity based thermal t...
dimensionedScalar sqrt(const dimensionedScalar &ds)
static word groupName(StringType name, const word &group)
Create dot-delimited name.group.
dimensioned< typename typeOfMag< Type >::type > mag(const dimensioned< Type > &dt)
label k
Boltzmann constant.
This boundary condition provides a turbulent viscosity condition when using wall functions,...
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,...
virtual void writeLocalEntries(Ostream &) const
Write local wall function variables.
Field with dimensions and associated with geometry type GeoMesh which is used to size the field and a...
This boundary condition provides a wall constraint on the turbulent kinematic viscosity,...