34 namespace compressible
39 scalar alphatJayatillekeWallFunctionFvPatchScalarField::maxExp_ = 50.0;
40 scalar alphatJayatillekeWallFunctionFvPatchScalarField::tolerance_ = 0.01;
41 label alphatJayatillekeWallFunctionFvPatchScalarField::maxIters_ = 10;
54 fixedValueFvPatchScalarField(
p, iF,
dict),
55 Prt_(
dict.lookupOrDefault<scalar>(
"Prt", 0.85))
68 fixedValueFvPatchScalarField(ptf,
p, iF, mapper),
80 fixedValueFvPatchScalarField(awfpsf, iF),
92 return 9.24*(
pow(Prat, 0.75) - 1.0)*(1.0 + 0.28*
exp(-0.007*Prat));
106 const scalar E = nutw.
E();
113 for (
int i=0; i<maxIters_; i++)
116 ypt[facei] - (
log(E*ypt[facei])/
kappa +
P[facei])/Prat[facei];
118 const scalar df = 1.0 - 1.0/(ypt[facei]*nutw.
kappa()*Prat[facei]);
120 const scalar dypt = -
f/df;
124 if (ypt[facei] < vSmall)
130 if (
mag(dypt) < tolerance_)
207 const label celli = nutw.patch().faceCells()[facei];
209 const scalar
uTau = Cmu25*
sqrt(
k[celli]);
211 const scalar
yPlus =
uTau*
y[facei]/nuw[facei];
214 scalar alphaEff = 0.0;
217 const scalar
A = gradHew[facei]*rhow[facei]*
uTau*
y[facei];
219 const scalar
B = gradHew[facei]*Pr[facei]*
yPlus;
221 const scalar
C = Pr[facei]*0.5*rhow[facei]*
uTau*
sqr(
magUp[facei]);
223 alphaEff = (
A -
C)/(
B +
sign(
B)*rootVSmall);
227 const scalar
A = gradHew[facei]*rhow[facei]*
uTau*
y[facei];
239 *(Prt*
sqr(
magUp[facei]) + (Pr[facei] - Prt)*
sqr(magUc));
241 alphaEff = (
A -
C)/(
B +
sign(
B)*rootVSmall);
245 static const scalar alphatwMin = 0;
246 const scalar alphatwMax = great*alphaw[facei]*nutw[facei]/nuw[facei];
250 min(
max(alphaEff - alphaw[facei], alphatwMin), alphatwMax);
267 internalField().group()
static const Foam::dimensionedScalar A("A", Foam::dimPressure, 611.21)
static const Foam::dimensionedScalar B("B", Foam::dimless, 18.678)
static const Foam::dimensionedScalar C("C", Foam::dimTemperature, 234.5)
#define forAll(list, i)
Loop across all elements in list.
Macros for easy insertion into run-time selection tables.
Graphite solid properties.
Field with dimensions and associated with geometry type GeoMesh which is used to size the field and a...
Generic GeometricField class.
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 volScalarField & he()=0
Enthalpy/Internal energy [J/kg].
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.
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 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.
Foam::fvPatchFieldMapper.
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 tmp< Field< Type > > patchInternalField() const
Return internal field next to patch as patch field.
virtual tmp< Field< Type > > snGrad() const
Return patch-normal gradient.
virtual void updateCoeffs()
Update the coefficients associated with the patch field.
A finiteVolume patch using a polyPatch and a fvBoundaryMesh.
virtual tmp< volScalarField > k() const =0
Return the turbulence kinetic energy.
const volVectorField & U() const
Access function to velocity field.
const volScalarField::Boundary & y() const
Return the near wall distance.
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.
scalar Cmu() const
Return Cmu.
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)
dimensionedScalar sign(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.
dimensionedSymmTensor sqr(const dimensionedVector &dv)
tmp< fvMatrix< Type > > operator==(const fvMatrix< Type > &, const fvMatrix< Type > &)
Field< scalar > scalarField
Specialisation of Field<T> for scalar.
dimensionedScalar log(const dimensionedScalar &ds)
layerAndWeight min(const layerAndWeight &a, const layerAndWeight &b)
dimensionedScalar pow(const dimensionedScalar &ds, const dimensionedScalar &expt)
void writeEntry(Ostream &os, const HashTable< T, Key, Hash > &ht)
dimensionedScalar sqrt(const dimensionedScalar &ds)
dimensioned< scalar > mag(const dimensioned< Type > &)
layerAndWeight max(const layerAndWeight &a, const layerAndWeight &b)
dimensionedScalar pow025(const dimensionedScalar &ds)