Go to the documentation of this file.
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
73 fixedValueFvPatchScalarField(
p, iF),
93 fixedValueFvPatchScalarField(ptf,
p, iF, mapper),
97 Prt_(ptf.
Prt_.clone(
p.patch())),
98 z0_(ptf.
z0_.clone(
p.patch()))
112 fixedValueFvPatchScalarField(
p, iF,
dict),
145 fixedValueFvPatchScalarField(wfpsf),
149 Prt_(wfpsf.
Prt_.clone(this->patch().patch())),
150 z0_(wfpsf.
z0_.clone(this->patch().patch()))
163 fixedValueFvPatchScalarField(wfpsf, iF),
167 Prt_(wfpsf.
Prt_.clone(this->patch().patch())),
168 z0_(wfpsf.
z0_.clone(this->patch().patch()))
183 const label patchi =
patch().index();
186 const auto& turbModel =
192 internalField().
group()
206 const scalar t = db().time().timeOutputValue();
207 const scalar
Pr =
Pr_->value(t);
213 <<
"Pr cannot be negative or zero. "
214 <<
"Please check input Pr = " <<
Pr
225 if (
Prt[i] < VSMALL || z0[i] < VSMALL)
228 <<
"Elements of input surface fields can only be positive. "
229 <<
"Please check input fields z0 and Prt."
244 const scalar Edash = (
y[facei] + z0[facei])/(z0[facei] + 1
e-4);
254 alphatw =
max(alphatw, scalar(0.01));
265 fixedValueFvPatchScalarField::autoMap(m);
277 fixedValueFvPatchScalarField::rmap(ptf, addr);
280 refCast<const atmAlphatkWallFunctionFvPatchScalarField>(ptf);
282 z0_->rmap(nrwfpsf.
z0_(), addr);
283 Prt_->rmap(nrwfpsf.
Prt_(), addr);
290 os.writeEntry(
"Cmu",
Cmu_);
295 writeEntry(
"value",
os);
atmAlphatkWallFunctionFvPatchScalarField(const fvPatch &, const DimensionedField< scalar, volMesh > &)
Construct from patch and internal field.
virtual void write(Ostream &) const
Write.
virtual void updateCoeffs()
Update the coefficients associated with the patch field.
constexpr const char *const group
Group name for atomic constants.
A class for managing temporary objects.
static const word propertiesName
Default name of the turbulence properties dictionary.
Ostream & endl(Ostream &os)
Add newline and flush stream.
const scalar Cmu_
Empirical model constant.
Top level data entry class for use in dictionaries. Provides a mechanism to specify a variable as a c...
virtual void write(Ostream &) const
Write.
autoPtr< PatchFunction1< scalar > > Prt_
Turbulent Prandtl number field.
#define forAll(list, i)
Loop across all elements in list.
virtual void checkType()
Check the type of the patch.
dimensionedScalar pow025(const dimensionedScalar &ds)
static scalar tolerance_
Solution parameters.
A finiteVolume patch using a polyPatch and a fvBoundaryMesh.
dimensionedScalar Pr("Pr", dimless, laminarTransport)
dimensionedScalar Prt("Prt", dimless, laminarTransport)
T getCheckOrDefault(const word &keyword, const T &deflt, const Predicate &pred, enum keyType::option matchOpt=keyType::REGEX) const
label max(const labelHashSet &set, label maxValue=labelMin)
Find the max value in labelHashSet, optionally limited by second argument.
autoPtr< PatchFunction1< scalar > > z0_
Surface roughness length [m].
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,...
dimensionedScalar log(const dimensionedScalar &ds)
OBJstream os(runTime.globalPath()/outputName)
Macros for easy insertion into run-time selection tables.
Top level data entry class for use in dictionaries. Provides a mechanism to specify a variable as a c...
errorManip< error > abort(error &err)
virtual void rmap(const fvPatchScalarField &, const labelList &)
Reverse map the given fvPatchField onto this fvPatchField.
virtual void updateCoeffs()
Update the coefficients associated with the patch field.
errorManipArg< error, int > exit(error &err, const int errNo=1)
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
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.
dimensionedScalar sqrt(const dimensionedScalar &ds)
autoPtr< Function1< scalar > > Pr_
Molecular Prandtl number.
label k
Boltzmann constant.
const dimensionedScalar e
Elementary charge.
virtual void autoMap(const fvPatchFieldMapper &)
Map (and resize as needed) from self given a mapping object.
This boundary condition provides a wall constraint on the kinematic turbulent thermal conductivity (i...
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, atmBoundaryLayerInletEpsilonFvPatchScalarField)
Smooth ATC in cells next to a set of patches supplied by type.
const scalar kappa_
von Kármán constant
Field with dimensions and associated with geometry type GeoMesh which is used to size the field and a...