Go to the documentation of this file.
46 const label patchi =
patch.index();
49 nutWallFunctionFvPatchScalarField::nutw(turbModel, patchi);
64 const scalar Cmu75 =
pow(nutw.
Cmu(), 0.75);
66 const scalar t = db().time().timeOutputValue();
70 for (
const auto& z : z0)
75 <<
"z0 field can only contain positive values. "
76 <<
"Please check input field z0."
89 const scalar
yPlus = Cmu25*
y[facei]*
sqrt(
k[celli])/nuw[facei];
91 const scalar w = cornerWeights[facei];
95 w*Cmu75*
pow(
k[celli], 1.5)/(nutw.
kappa()*(
y[facei] + z0[facei]));
99 *(nutw[facei] + nuw[facei])
101 *Cmu25*
sqrt(
k[celli])
102 /(nutw.
kappa()*(
y[facei] + z0[facei]));
106 epsilonc = w*2.0*
k[celli]*nuw[facei]/
sqr(
y[facei] + z0[facei]);
141 z0_(ptf.
z0_.clone(
p.patch()))
165 z0_(ewfpsf.
z0_.clone(this->patch().patch()))
177 z0_(ewfpsf.
z0_.clone(this->patch().patch()))
188 epsilonWallFunctionFvPatchScalarField::autoMap(m);
199 epsilonWallFunctionFvPatchScalarField::rmap(ptf, addr);
202 refCast<const atmEpsilonWallFunctionFvPatchScalarField>(ptf);
204 z0_->rmap(atmpsf.
z0_(), addr);
virtual tmp< Field< Type > > snGrad() const
Return patch-normal gradient.
virtual tmp< volScalarField > nu() const =0
Return the laminar viscosity.
virtual void autoMap(const fvPatchFieldMapper &)
Map (and resize as needed) from self given a mapping object.
scalar Cmu() const
Return Cmu.
autoPtr< PatchFunction1< scalar > > z0_
Surface roughness length field [m].
A class for managing temporary objects.
const nearWallDist & y() const
Return the near wall distances.
const dimensionedScalar epsilon0
Electric constant: default SI units: [F/m].
#define forAll(list, i)
Loop across all elements in list.
dimensionedScalar pow025(const dimensionedScalar &ds)
const dimensionedScalar G0
Conductance quantum: default SI units: [S].
This boundary condition provides a wall constraint on the turbulent kinetic energy dissipation rate,...
A finiteVolume patch using a polyPatch and a fvBoundaryMesh.
This boundary condition provides a wall constraint on the turbulent kinetic energy dissipation rate,...
dimensionedScalar pow(const dimensionedScalar &ds, const dimensionedScalar &expt)
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,...
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...
virtual tmp< volScalarField > k() const =0
Return the turbulence kinetic energy.
static scalar yPlusLam(const scalar kappa, const scalar E)
Estimate the y+ at the intersection of the two sublayers.
errorManipArg< error, int > exit(error &err, const int errNo=1)
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
dimensionedSymmTensor sqr(const dimensionedVector &dv)
const std::string patch
OpenFOAM patch number as a std::string.
dimensionedScalar sqrt(const dimensionedScalar &ds)
dimensioned< typename typeOfMag< Type >::type > mag(const dimensioned< Type > &dt)
scalar kappa() const
Return kappa.
label k
Boltzmann constant.
void write(vtk::formatter &fmt, const Type &val, const label n=1)
Component-wise write of a value (N times)
Foam::fvPatchFieldMapper.
virtual void rmap(const fvPatchScalarField &, const labelList &)
Reverse map the given fvPatchField onto this fvPatchField.
An Ostream is an abstract base class for all output systems (streams, files, token lists,...
makePatchTypeField(fvPatchScalarField, atmBoundaryLayerInletEpsilonFvPatchScalarField)
The class nutWallFunction is a base class that parents the derived boundary conditions which provide ...
virtual void write(Ostream &) const
Write.
Smooth ATC in cells next to a set of patches supplied by type.
virtual void calculate(const turbulenceModel &turbulence, const List< scalar > &cornerWeights, const fvPatch &patch, scalarField &G, scalarField &epsilon)
Calculate the epsilon and G.
atmEpsilonWallFunctionFvPatchScalarField(const fvPatch &, const DimensionedField< scalar, volMesh > &)
Construct from patch and internal field.
const Boundary & boundaryField() const
Return const-reference to the boundary field.
Field with dimensions and associated with geometry type GeoMesh which is used to size the field and a...
const volVectorField & U() const
Access function to velocity field.