34 namespace compressible
39 scalar alphatJayatillekeWallFunctionFvPatchScalarField::maxExp_ = 50.0;
40 scalar alphatJayatillekeWallFunctionFvPatchScalarField::tolerance_ = 0.01;
41 label alphatJayatillekeWallFunctionFvPatchScalarField::maxIters_ = 10;
45 scalar alphatJayatillekeWallFunctionFvPatchScalarField::Psmooth
50 return 9.24*(
pow(Prat, 0.75) - 1.0)*(1.0 + 0.28*
exp(-0.007*Prat));
54 scalar alphatJayatillekeWallFunctionFvPatchScalarField::yPlusTherm
56 const nutWallFunctionFvPatchScalarField& nutw,
63 for (
int i=0; i<maxIters_; i++)
65 scalar f = ypt - (
log(nutw.E()*ypt)/nutw.kappa() + P)/Prat;
66 scalar df = 1.0 - 1.0/(ypt*nutw.kappa()*Prat);
67 scalar yptNew = ypt - f/df;
73 else if (
mag(yptNew - ypt) < tolerance_)
96 fixedValueFvPatchScalarField(p, iF),
110 fixedValueFvPatchScalarField(ptf, p, iF, mapper),
123 fixedValueFvPatchScalarField(p, iF, dict),
134 fixedValueFvPatchScalarField(awfpsf),
146 fixedValueFvPatchScalarField(awfpsf, iF),
168 compressible::turbulenceModel::propertiesName,
169 internalField().
group()
176 const scalar Cmu25 =
pow025(nutw.Cmu());
197 turbModel.transport().he().boundaryField()[
patchi];
208 label celli = patch().faceCells()[facei];
212 scalar
yPlus = uTau*y[facei]/(muw[facei]/rhow[facei]);
215 scalar Pr = muw[facei]/alphaw[facei];
218 scalar Prat = Pr/Prt_;
221 scalar P = Psmooth(Prat);
222 scalar yPlusTherm = this->yPlusTherm(nutw, P, Prat);
225 scalar alphaEff = 0.0;
226 if (yPlus < yPlusTherm)
228 const scalar A = qDot[facei]*rhow[facei]*uTau*y[facei];
230 const scalar B = qDot[facei]*Pr*
yPlus;
232 const scalar
C = Pr*0.5*rhow[facei]*uTau*
sqr(magUp[facei]);
234 alphaEff = A/(B + C + vSmall);
238 const scalar A = qDot[facei]*rhow[facei]*uTau*y[facei];
241 qDot[facei]*Prt_*(1.0/nutw.kappa()*
log(nutw.E()*
yPlus) + P);
244 uTau/nutw.kappa()*
log(nutw.E()*yPlusTherm) -
mag(Uw[facei]);
248 *(Prt_*
sqr(magUp[facei]) + (Pr - Prt_)*
sqr(magUc));
250 alphaEff = A/(B + C + vSmall);
254 alphatw[facei] =
max(0.0, alphaEff - alphaw[facei]);
258 Info<<
" uTau = " << uTau <<
nl 259 <<
" Pr = " << Pr <<
nl 260 <<
" Prt = " << Prt_ <<
nl 261 <<
" qDot = " << qDot[facei] <<
nl 262 <<
" yPlus = " << yPlus <<
nl 263 <<
" yPlusTherm = " << yPlusTherm <<
nl 264 <<
" alphaEff = " << alphaEff <<
nl 265 <<
" alphaw = " << alphaw[facei] <<
nl 266 <<
" alphatw = " << alphatw[facei] <<
nl const char *const group
Group name for atomic constants.
Graphite solid properties.
#define forAll(list, i)
Loop across all elements in list.
intWM_LABEL_SIZE_t label
A label is an int32_t or int64_t as specified by the pre-processor macro WM_LABEL_SIZE.
virtual tmp< Field< Type > > snGrad() const
Return patch-normal gradient.
dimensionedScalar log(const dimensionedScalar &ds)
A list of keyword definitions, which are a keyword followed by any number of values (e...
dimensioned< Type > max(const dimensioned< Type > &, const dimensioned< Type > &)
dimensionedSymmTensor sqr(const dimensionedVector &dv)
This boundary condition provides a thermal wall function for turbulent thermal diffusivity (usuallyal...
dimensionedScalar sqrt(const dimensionedScalar &ds)
Ostream & endl(Ostream &os)
Add newline and flush stream.
A finiteVolume patch using a polyPatch and a fvBoundaryMesh.
dimensionedScalar pow025(const dimensionedScalar &ds)
Abstract base class with a fat-interface to all derived classes covering all possible ways in which t...
void write(Ostream &) const
Write.
label k
Boltzmann constant.
Templated wrapper class to provide compressible turbulence models thermal diffusivity based thermal t...
virtual void write(Ostream &) const
Write.
Macros for easy insertion into run-time selection tables.
virtual tmp< volScalarField > alphaEff() const
Effective thermal turbulent diffusivity of mixture [kg/m/s].
alphatJayatillekeWallFunctionFvPatchScalarField(const fvPatch &, const DimensionedField< scalar, volMesh > &)
Construct from patch and internal field.
This boundary condition provides a turbulent kinematic viscosity condition when using wall functions...
dimensionedScalar exp(const dimensionedScalar &ds)
static word groupName(Name name, const word &group)
Foam::fvPatchFieldMapper.
An Ostream is an abstract base class for all output systems (streams, files, token lists...
virtual tmp< Field< Type > > patchInternalField() const
Return internal field next to patch as patch field.
void writeEntry(Ostream &os, const HashTable< T, Key, Hash > &ht)
dimensionedScalar pow(const dimensionedScalar &ds, const dimensionedScalar &expt)
T lookupOrDefault(const word &, const T &, bool recursive=false, bool patternMatch=true) const
Find and return a T,.
virtual void updateCoeffs()
Update the coefficients associated with the patch field.
Field with dimensions and associated with geometry type GeoMesh which is used to size the field and a...
dimensioned< scalar > mag(const dimensioned< Type > &)
A class for managing temporary objects.
virtual void updateCoeffs()
Update the coefficients associated with the patch field.
static const nutWallFunctionFvPatchScalarField & nutw(const turbulenceModel &turbModel, const label patchi)
Return the nut patchField for the given wall patch.
virtual tmp< volScalarField > alpha() const
Thermal diffusivity for enthalpy of mixture [kg/m/s].
makePatchTypeField(fvPatchScalarField, thermalBaffleFvPatchScalarField)