44 const label patchi =
patch().index();
51 internalField().
group()
66 const scalar t = db().time().timeOutputValue();
70 for (
const scalar z : z0)
75 <<
"z0 field can only contain positive values. " 76 <<
"Please check input field z0." 83 auto&
nutw = tnutw.ref();
87 const scalar Edash = (
y[facei] + z0[facei])/(z0[facei] + 1
e-4);
88 const scalar uStar = magUpn[facei]*
kappa/
log(
max(Edash, 1.0 + 1
e-4));
90 nutw[facei] =
sqr(uStar)/
max(magUpn[facei], 1
e-6)*
y[facei] - nuw[facei];
121 const DimensionedField<scalar, volMesh>& iF
139 boundNut_(ptf.boundNut_),
140 z0_(ptf.z0_.clone(
p.
patch()))
152 boundNut_(
dict.getOrDefault<bool>(
"boundNut", true)),
163 boundNut_(rwfpsf.boundNut_),
175 boundNut_(rwfpsf.boundNut_),
187 nutUWallFunctionFvPatchScalarField::autoMap(m);
202 nutUWallFunctionFvPatchScalarField::rmap(ptf, addr);
205 refCast<const atmNutUWallFunctionFvPatchScalarField>(ptf);
209 z0_->rmap(nrwfpsf.z0_(), addr);
227 atmNutUWallFunctionFvPatchScalarField
virtual tmp< scalarField > calcNut() const
Calculate the turbulent viscosity.
virtual void rmap(const fvPatchScalarField &, const labelList &)
Reverse map the given fvPatchField onto this fvPatchField.
virtual const volVectorField & U(const turbulenceModel &turb) const
Helper to return the velocity field either from the turbulence model (default) or the mesh database...
atmNutUWallFunctionFvPatchScalarField(const fvPatch &, const DimensionedField< scalar, volMesh > &)
Construct from patch and internal field.
dimensionedScalar log(const dimensionedScalar &ds)
errorManipArg< error, int > exit(error &err, const int errNo=1)
dimensioned< typename typeOfMag< Type >::type > mag(const dimensioned< Type > &dt)
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.
This boundary condition provides a wall function for the turbulent viscosity (i.e. nut) based on velocity (i.e. U) for low- and high-Reynolds number applications.
label max(const labelHashSet &set, label maxValue=labelMin)
Find the max value in labelHashSet, optionally limited by second argument.
dimensionedSymmTensor sqr(const dimensionedVector &dv)
A finiteVolume patch using a polyPatch and a fvBoundaryMesh.
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.
void writeValueEntry(Ostream &os) const
Write *this field as a "value" entry.
const dimensionedScalar kappa
Coulomb constant: default SI units: [N.m2/C2].
Abstract base class for turbulence models (RAS, LES and laminar).
Macros for easy insertion into run-time selection tables.
#define forAll(list, i)
Loop across all elements in list.
constexpr const char *const group
Group name for atomic constants.
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.
A FieldMapper for finite-volume patch fields.
virtual void write(Ostream &) const
Write.
static tmp< T > New(Args &&... args)
Construct tmp with forwarding arguments.
void writeLocalEntries(Ostream &) const
Write local wall function variables.
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...
virtual void write(Ostream &) const
Write.
virtual tmp< Field< Type > > patchInternalField() const
Return internal field next to patch.
OBJstream os(runTime.globalPath()/outputName)
virtual void autoMap(const fvPatchFieldMapper &)
Map (and resize as needed) from self given a mapping object.
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.
wallFunctionCoefficients wallCoeffs_
Wall-function coefficients.
List< label > labelList
A List of labels.
scalar kappa() const noexcept
Return the object: kappa.
A class for managing temporary objects.
This boundary condition provides a wall constraint on the turbulent viscosity (i.e. nut) based on velocity (i.e. U) for atmospheric boundary layer modelling. It is designed to be used in conjunction with the atmBoundaryLayerInletVelocity boundary condition.
static const nutWallFunctionFvPatchScalarField & nutw(const turbulenceModel &turbModel, const label patchi)
Return the nut patchField for the given wall patch.
makePatchTypeField(fvPatchScalarField, atmBoundaryLayerInletEpsilonFvPatchScalarField)
const Boundary & boundaryField() const noexcept
Return const-reference to the boundary field.