72 for (
const auto& z : z0)
77 <<
"z0 field can only contain positive values. "
78 <<
"Please check input field z0."
91 const scalar
yPlus = Cmu25*
y[facei]*
sqrt(
k[celli])/nuw[facei];
93 const scalar w = cornerWeights[facei];
97 w*Cmu75*
pow(
k[celli], 1.5)/(kappa*(
y[facei] + z0[facei]));
101 *(nutw[facei] + nuw[facei])
103 *Cmu25*
sqrt(
k[celli])
104 /(kappa*(
y[facei] + z0[facei]));
108 epsilonc = w*2.0*
k[celli]*nuw[facei]/
sqr(
y[facei] + z0[facei]);
112 epsilon0[celli] += epsilonc;
131 wallCoeffs_.writeEntries(
os);
137Foam::atmEpsilonWallFunctionFvPatchScalarField::
138atmEpsilonWallFunctionFvPatchScalarField
149Foam::atmEpsilonWallFunctionFvPatchScalarField::
150atmEpsilonWallFunctionFvPatchScalarField
159 z0_(ptf.z0_.clone(
p.patch()))
163Foam::atmEpsilonWallFunctionFvPatchScalarField::
164atmEpsilonWallFunctionFvPatchScalarField
176Foam::atmEpsilonWallFunctionFvPatchScalarField::
177atmEpsilonWallFunctionFvPatchScalarField
183 z0_(ewfpsf.z0_.clone(this->patch().patch()))
187Foam::atmEpsilonWallFunctionFvPatchScalarField::
188atmEpsilonWallFunctionFvPatchScalarField
195 z0_(ewfpsf.z0_.clone(this->patch().patch()))
224 refCast<const atmEpsilonWallFunctionFvPatchScalarField>(ptf);
227 z0_->rmap(atmpsf.z0_(), addr);
238 writeLocalEntries(
os);
239 writeEntry(
"value",
os);
Macros for easy insertion into run-time selection tables.
Field with dimensions and associated with geometry type GeoMesh which is used to size the field and a...
const Boundary & boundaryField() const
Return const-reference to the boundary field.
A 1D array of objects of type <T>, where the size of the vector is known and used for subscript bound...
An Ostream is an abstract base class for all output systems (streams, files, token lists,...
Ostream & writeEntryIfDifferent(const word &key, const T &value1, const T &value2)
Write a keyword/value entry only when the two values differ.
Top level data entry class for use in dictionaries. Provides a mechanism to specify a variable as a c...
scalar timeOutputValue() const
Return current time value.
void autoMap(const fvPatchFieldMapper &)
Map (and resize as needed) from self given a mapping object.
void rmap(const atmBoundaryLayer &, const labelList &)
Reverse map the given fvPatchField onto this fvPatchField.
This boundary condition provides a wall constraint on the turbulent kinetic energy dissipation rate (...
virtual void autoMap(const fvPatchFieldMapper &)
Map (and resize as needed) from self given a mapping object.
autoPtr< PatchFunction1< scalar > > z0_
Surface roughness length field [m].
void writeLocalEntries(Ostream &) const
Write local wall function variables.
virtual void rmap(const fvPatchScalarField &, const labelList &)
Reverse map the given fvPatchField onto this fvPatchField.
virtual void calculate(const turbulenceModel &turbulence, const List< scalar > &cornerWeights, const fvPatch &patch, scalarField &G, scalarField &epsilon)
Calculate the epsilon and G.
A list of keyword definitions, which are a keyword followed by a number of values (eg,...
This boundary condition provides wall functions for the turbulent kinetic energy dissipation rate (i....
const bool lowReCorrection_
Apply low-Re correction term (default = no)
wallFunctionCoefficients wallCoeffs_
Wall-function coefficients.
Smooth ATC in cells next to a set of patches supplied by type.
virtual bool write()
Write the output fields.
A FieldMapper for finite-volume patch fields.
virtual tmp< Field< Type > > snGrad() const
Return patch-normal gradient.
const objectRegistry & db() const
Return local objectRegistry.
const fvPatch & patch() const
Return patch.
A finiteVolume patch using a polyPatch and a fvBoundaryMesh.
label index() const
Return the index of this patch in the fvBoundaryMesh.
virtual const labelUList & faceCells() const
Return faceCells.
const Time & time() const noexcept
Return time registry.
A class for managing temporary objects.
Abstract base class for turbulence models (RAS, LES and laminar).
const volVectorField & U() const
Access function to velocity field.
virtual tmp< volScalarField > k() const =0
Return the turbulence kinetic energy.
virtual tmp< volScalarField > nu() const =0
Return the laminar viscosity.
virtual tmp< volScalarField > nut() const =0
Return the turbulence viscosity.
const nearWallDist & y() const
Return the near wall distances.
scalar kappa() const noexcept
Return the object: kappa.
scalar yPlusLam() const noexcept
Return the object: yPlusLam.
scalar Cmu() const noexcept
Return the object: Cmu.
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
OBJstream os(runTime.globalPath()/outputName)
#define makePatchTypeField(PatchTypeField, typePatchTypeField)
dimensionedSymmTensor sqr(const dimensionedVector &dv)
dimensionedScalar pow(const dimensionedScalar &ds, const dimensionedScalar &expt)
dimensionedScalar sqrt(const dimensionedScalar &ds)
dimensioned< typename typeOfMag< Type >::type > mag(const dimensioned< Type > &dt)
tmp< DimensionedField< TypeR, GeoMesh > > New(const tmp< DimensionedField< TypeR, GeoMesh > > &tdf1, const word &name, const dimensionSet &dimensions)
Global function forwards to reuseTmpDimensionedField::New.
errorManipArg< error, int > exit(error &err, const int errNo=1)
dimensionedScalar pow025(const dimensionedScalar &ds)
#define forAll(list, i)
Loop across all elements in list.