50 if (!isA<wallFvPatch>(
patch()))
53 <<
"Invalid wall function specification" <<
nl 54 <<
" Patch type for patch " <<
patch().name()
55 <<
" must be wall" <<
nl 67 os.writeEntryIfDifferent<scalar>(
"Cmu", 0.09,
Cmu_);
68 os.writeEntryIfDifferent<scalar>(
"kappa", 0.41,
kappa_);
94 fixedValueFvPatchScalarField(
p, iF),
114 fixedValueFvPatchScalarField(ptf,
p, iF, mapper),
117 Pr_(ptf.Pr_.clone()),
118 Prt_(ptf.Prt_.clone(
p.
patch())),
119 z0_(ptf.z0_.clone(
p.
patch()))
133 fixedValueFvPatchScalarField(
p, iF,
dict),
136 dict.getCheckOrDefault<scalar>
145 dict.getCheckOrDefault<scalar>
166 fixedValueFvPatchScalarField(wfpsf),
168 kappa_(wfpsf.kappa_),
169 Pr_(wfpsf.Pr_.clone()),
170 Prt_(wfpsf.Prt_.clone(this->
patch().
patch())),
184 fixedValueFvPatchScalarField(wfpsf, iF),
186 kappa_(wfpsf.kappa_),
187 Pr_(wfpsf.Pr_.clone()),
188 Prt_(wfpsf.Prt_.clone(this->
patch().
patch())),
204 const label patchi =
patch().index();
212 internalField().
group()
218 const tmp<scalarField> tnuw = turbModel.nu(patchi);
221 const tmp<volScalarField> tk = turbModel.k();
226 const scalar t = db().time().timeOutputValue();
227 const scalar
Pr =
Pr_->value(t);
233 <<
"Pr cannot be negative or zero. " 234 <<
"Please check input Pr = " <<
Pr 245 if (
Prt[i] < VSMALL || z0[i] < VSMALL)
248 <<
"Elements of input surface fields can only be positive. " 249 <<
"Please check input fields z0 and Prt." 261 const label celli = faceCells[facei];
264 const scalar Edash = (
y[facei] + z0[facei])/(z0[facei] + 1
e-4);
274 alphatw =
max(alphatw, scalar(0.01));
285 fixedValueFvPatchScalarField::autoMap(m);
304 fixedValueFvPatchScalarField::rmap(ptf, addr);
306 const auto& nrwfpsf =
307 refCast<const atmAlphatkWallFunctionFvPatchScalarField>(ptf);
311 Prt_->rmap(nrwfpsf.Prt_(), addr);
315 z0_->rmap(nrwfpsf.z0_(), addr);
333 atmAlphatkWallFunctionFvPatchScalarField
dimensionedScalar Pr("Pr", dimless, laminarTransport)
Top level data entry class for use in dictionaries. Provides a mechanism to specify a variable as a c...
dimensionedScalar log(const dimensionedScalar &ds)
errorManipArg< error, int > exit(error &err, const int errNo=1)
error FatalError
Error stream (stdout output on all processes), with additional 'FOAM FATAL ERROR' header text and sta...
A list of keyword definitions, which are a keyword followed by a number of values (eg...
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
label max(const labelHashSet &set, label maxValue=labelMin)
Find the max value in labelHashSet, optionally limited by second argument.
virtual void autoMap(const fvPatchFieldMapper &)
Map (and resize as needed) from self given a mapping object.
constexpr char nl
The newline '\n' character (0x0a)
autoPtr< PatchFunction1< scalar > > Prt_
Turbulent Prandtl number field.
ThermalDiffusivity< CompressibleTurbulenceModel< fluidThermo > > turbulenceModel
dimensionedScalar sqrt(const dimensionedScalar &ds)
Ostream & endl(Ostream &os)
Add newline and flush stream.
A finiteVolume patch using a polyPatch and a fvBoundaryMesh.
dimensionedScalar pow025(const dimensionedScalar &ds)
tmp< DimensionedField< TypeR, GeoMesh > > New(const tmp< DimensionedField< TypeR, GeoMesh >> &tf1, const word &name, const dimensionSet &dimensions, const bool initCopy=false)
Global function forwards to reuseTmpDimensionedField::New.
This boundary condition provides a wall constraint on the kinematic turbulent thermal conductivity (i...
void writeValueEntry(Ostream &os) const
Write *this field as a "value" entry.
label k
Boltzmann constant.
virtual void write(Ostream &) const
Write.
virtual void write(Ostream &) const
Write.
Macros for easy insertion into run-time selection tables.
UList< label > labelUList
A UList of labels.
autoPtr< PatchFunction1< scalar > > z0_
Surface roughness length [m].
atmAlphatkWallFunctionFvPatchScalarField(const fvPatch &, const DimensionedField< scalar, volMesh > &)
Construct from patch and internal field.
virtual void checkType()
Check the type of the patch.
const scalar kappa_
von Kármán constant
#define forAll(list, i)
Loop across all elements in list.
constexpr const char *const group
Group name for atomic constants.
GeometricField< scalar, fvPatchField, volMesh > volScalarField
static word groupName(StringType base, const word &group)
Create dot-delimited name.group string.
const dimensionedScalar e
Elementary charge.
fvPatchField< scalar > fvPatchScalarField
static const word propertiesName
Default name of the turbulence properties dictionary.
Field< scalar > scalarField
Specialisation of Field<T> for scalar.
A FieldMapper for finite-volume patch fields.
virtual void rmap(const fvPatchScalarField &, const labelList &)
Reverse map the given fvPatchField onto this fvPatchField.
errorManip< error > abort(error &err)
Top level data entry class for use in dictionaries. Provides a mechanism to specify a variable as a c...
OBJstream os(runTime.globalPath()/outputName)
virtual void updateCoeffs()
Update the coefficients associated with the patch field.
dimensionedScalar Prt("Prt", dimless, laminarTransport)
static scalar tolerance_
Solution parameters.
virtual void updateCoeffs()
Update the coefficients associated with the patch field.
Field with dimensions and associated with geometry type GeoMesh which is used to size the field and a...
const std::string patch
OpenFOAM patch number as a std::string.
List< label > labelList
A List of labels.
makePatchTypeField(fvPatchScalarField, atmBoundaryLayerInletEpsilonFvPatchScalarField)
autoPtr< Function1< scalar > > Pr_
Molecular Prandtl number.
const scalar Cmu_
Empirical model constant.
void writeLocalEntries(Ostream &) const
Write local wall function variables.