34 namespace compressible
39 const scalar alphatJayatillekeWallFunctionFvPatchScalarField::tolerance_ = 0.01;
40 const label alphatJayatillekeWallFunctionFvPatchScalarField::maxIters_ = 10;
53 fixedValueFvPatchScalarField(
p, iF,
dict),
54 Prt_(
dict.lookupOrDefault<scalar>(
"Prt",
dimless, 0.85))
67 fixedValueFvPatchScalarField(ptf,
p, iF, mapper, false),
70 mapper(*
this, ptf, [&](){
return this->patchInternalField(); });
81 fixedValueFvPatchScalarField(awfpsf, iF),
94 mapper(*
this, ptf, [&](){
return this->patchInternalField(); });
103 return 9.24*(
pow(Prat, 0.75) - 1)*(1 + 0.28*
exp(-0.007*Prat));
117 const scalar E = nutw.
E();
124 for (
int i=0; i<maxIters_; i++)
127 ypt[facei] - (
log(E*ypt[facei])/
kappa +
P[facei])/Prat[facei];
129 const scalar df = 1 - 1/(ypt[facei]*nutw.
kappa()*Prat[facei]);
131 const scalar dypt = -
f/df;
135 if (ypt[facei] < vSmall)
141 if (
mag(dypt) < tolerance_)
165 const scalar E = nutw.
E();
212 const scalar alphaByAlphaEff = Tplus/Pr[facei]/
yPlus[facei];
215 alphaw[facei]*
max(1/alphaByAlphaEff - 1, scalar(0));
233 internalField().group()
#define forAll(list, i)
Loop across all elements in list.
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.
An Ostream is an abstract base class for all output systems (streams, files, token lists,...
virtual const volScalarField & Cp() const =0
Heat capacity at constant pressure [J/kg/K].
virtual const volScalarField & kappa() const =0
Thermal conductivity of mixture [W/m/K].
Base class for single-phase compressible turbulence models.
const rhoField & rho() const
Return the density field.
This boundary condition provides a thermal wall function for turbulent thermal diffusivity (usuallyal...
void write(Ostream &) const
Write.
static tmp< scalarField > alphat(const fluidThermophysicalTransportModel &ttm, const scalar Prt, const label patchi)
Calculate the turbulent thermal diffusivity.
virtual void updateCoeffs()
Update the coefficients associated with the patch field.
virtual void map(const fvPatchScalarField &, const fieldMapper &)
Map the given fvPatchField onto this fvPatchField.
static tmp< scalarField > yPlusTherm(const nutWallFunctionFvPatchScalarField &nutw, const scalarField &P, const scalarField &Prat)
Calculate y+ at the edge of the thermal laminar sublayer.
alphatJayatillekeWallFunctionFvPatchScalarField(const fvPatch &, const DimensionedField< scalar, volMesh > &, const dictionary &)
Construct from patch, internal field and dictionary.
static tmp< scalarField > P(const scalarField &Prat)
Calculate the smoothing function.
A list of keyword definitions, which are a keyword followed by any number of values (e....
Abstract base class for field mapping.
Abstract base class for fluid thermophysical transport models RAS, LES and laminar.
const compressibleMomentumTransportModel & momentumTransport() const
Access function to momentum transport model.
virtual const fluidThermo & thermo() const =0
Access function to fluid thermophysical properties.
Abstract base class with a fat-interface to all derived classes covering all possible ways in which t...
virtual void write(Ostream &) const
Write.
virtual void updateCoeffs()
Update the coefficients associated with the patch field.
A finiteVolume patch using a polyPatch and a fvBoundaryMesh.
virtual tmp< volScalarField > nu() const
Return the laminar viscosity.
This boundary condition provides a turbulent kinematic viscosity condition when using wall functions,...
scalar kappa() const
Return kappa.
scalar E() const
Return E.
static const nutWallFunctionFvPatchScalarField & nutw(const momentumTransportModel &turbModel, const label patchi)
Return the nut patchField for the given wall patch.
virtual tmp< scalarField > yPlus() const =0
Calculate and return the yPlus at the boundary.
A class for managing temporary objects.
T & ref() const
Return non-const reference or generate a fatal error.
makePatchTypeField(fvPatchScalarField, alphatJayatillekeWallFunctionFvPatchScalarField)
const dimensionedScalar kappa
Coulomb constant: default SI units: [N.m2/C2].
dimensionedScalar exp(const dimensionedScalar &ds)
intWM_LABEL_SIZE_t label
A label is an int32_t or int64_t as specified by the pre-processor macro WM_LABEL_SIZE.
tmp< fvMatrix< Type > > operator==(const fvMatrix< Type > &, const fvMatrix< Type > &)
const dimensionSet dimless
Field< scalar > scalarField
Specialisation of Field<T> for scalar.
dimensionedScalar log(const dimensionedScalar &ds)
dimensionedScalar pow(const dimensionedScalar &ds, const dimensionedScalar &expt)
void writeEntry(Ostream &os, const HashTable< T, Key, Hash > &ht)
dimensioned< scalar > mag(const dimensioned< Type > &)
layerAndWeight max(const layerAndWeight &a, const layerAndWeight &b)