45 scalar alphatJayatillekeWallFunctionFvPatchScalarField::maxExp_ = 50.0;
46 scalar alphatJayatillekeWallFunctionFvPatchScalarField::tolerance_ = 0.01;
47 label alphatJayatillekeWallFunctionFvPatchScalarField::maxIters_ = 10;
51 void alphatJayatillekeWallFunctionFvPatchScalarField::checkType()
53 if (!isA<wallFvPatch>(
patch()))
56 <<
"Patch type for patch " <<
patch().name() <<
" must be wall\n" 57 <<
"Current patch type is " <<
patch().type() <<
nl 63 tmp<scalarField> alphatJayatillekeWallFunctionFvPatchScalarField::yPlus
68 const label patchi =
patch().index();
70 const tmp<volScalarField> tnut = turbModel.nut();
73 if (isA<nutWallFunctionFvPatchScalarField>(
nut.boundaryField()[patchi]))
76 dynamic_cast<const nutWallFunctionFvPatchScalarField&
> 78 nut.boundaryField()[patchi]
89 y*
sqrt(turbModel.nuEff(patchi)*
mag(Uw.snGrad()))
90 /turbModel.nu(patchi);
95 scalar alphatJayatillekeWallFunctionFvPatchScalarField::Psmooth
100 return 9.24*(
pow(Prat, 0.75) - 1.0)*(1.0 + 0.28*
exp(-0.007*Prat));
104 scalar alphatJayatillekeWallFunctionFvPatchScalarField::yPlusTherm
112 for (
int iter = 0; iter < maxIters_; ++iter)
114 const scalar
f = ypt - (
log(E_*ypt)/kappa_ + P)/Prat;
115 const scalar df = 1.0 - 1.0/(ypt*kappa_*Prat);
116 const scalar yptNew = ypt -
f/df;
122 else if (
mag(yptNew - ypt) < tolerance_)
145 fixedValueFvPatchScalarField(
p, iF),
163 fixedValueFvPatchScalarField(ptf,
p, iF, mapper),
178 fixedValueFvPatchScalarField(
p, iF,
dict),
179 Prt_(
dict.getOrDefault<scalar>(
"Prt", 0.85)),
180 kappa_(
dict.getOrDefault<scalar>(
"kappa", 0.41)),
181 E_(
dict.getOrDefault<scalar>(
"E", 9.8))
193 fixedValueFvPatchScalarField(awfpsf),
195 kappa_(awfpsf.kappa_),
209 fixedValueFvPatchScalarField(awfpsf, iF),
211 kappa_(awfpsf.kappa_),
227 const label patchi =
patch().index();
234 compressible::turbulenceModel::propertiesName,
235 internalField().
group()
253 const scalarField& rhow = turbModel.rho().boundaryField()[patchi];
255 turbModel.transport().he().boundaryField()[patchi];
262 turbModel.transport().alphaEff(alphatw, patchi)*hew.
snGrad()
268 const scalar yPlus = yPlusp[facei];
270 const scalar
uTau =
yPlus/
y[facei]*(muw[facei]/rhow[facei]);
273 const scalar
Pr = muw[facei]/alphaw[facei];
276 const scalar Prat =
Pr/Prt_;
279 const scalar P = Psmooth(Prat);
280 const scalar yPlusTherm = this->yPlusTherm(P, Prat);
284 if (yPlus < yPlusTherm)
286 const scalar
A = qDot[facei]*rhow[facei]*
uTau*
y[facei];
287 const scalar
B = qDot[facei]*
Pr*
yPlus;
294 const scalar
A = qDot[facei]*rhow[facei]*
uTau*
y[facei];
295 const scalar
B = qDot[facei]*Prt_*(1.0/kappa_*
log(E_*yPlus) + P);
297 uTau/kappa_*
log(E_*yPlusTherm) -
mag(Uw[facei]);
306 alphatw[facei] =
max(scalar(0),
alphaEff - alphaw[facei]);
311 <<
" Pr = " <<
Pr <<
nl 312 <<
" Prt = " << Prt_ <<
nl 313 <<
" qDot = " << qDot[facei] <<
nl 315 <<
" yPlusTherm = " << yPlusTherm <<
nl 317 <<
" alphaw = " << alphaw[facei] <<
nl 318 <<
" alphatw = " << alphatw[facei] <<
nl 330 os.writeEntryIfDifferent<scalar>(
"Prt", 0.85, Prt_);
331 os.writeEntryIfDifferent<scalar>(
"kappa", 0.41, kappa_);
332 os.writeEntryIfDifferent<scalar>(
"E", 9.8, E_);
342 alphatJayatillekeWallFunctionFvPatchScalarField
dimensionedScalar Pr("Pr", dimless, laminarTransport)
Graphite solid properties.
fvPatchField< vector > fvPatchVectorField
virtual tmp< Field< Type > > snGrad() const
Return patch-normal gradient.
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.
label max(const labelHashSet &set, label maxValue=labelMin)
Find the max value in labelHashSet, optionally limited by second argument.
dimensionedSymmTensor sqr(const dimensionedVector &dv)
constexpr char nl
The newline '\n' character (0x0a)
This boundary condition provides a thermal wall function for turbulent thermal diffusivity (usuallyal...
ThermalDiffusivity< CompressibleTurbulenceModel< fluidThermo > > turbulenceModel
dimensionedScalar sqrt(const dimensionedScalar &ds)
Ostream & endl(Ostream &os)
Add newline and flush stream.
A finiteVolume patch using a polyPatch and a fvBoundaryMesh.
void write(Ostream &) const
Write.
void writeValueEntry(Ostream &os) const
Write *this field as a "value" entry.
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.
#define forAll(list, i)
Loop across all elements in list.
constexpr const char *const group
Group name for atomic constants.
makePatchTypeField(fvPatchScalarField, alphatFixedDmdtWallBoilingWallFunctionFvPatchScalarField)
GeometricField< scalar, fvPatchField, volMesh > volScalarField
alphatJayatillekeWallFunctionFvPatchScalarField(const fvPatch &, const DimensionedField< scalar, volMesh > &)
Construct from patch and internal field.
static word groupName(StringType base, const word &group)
Create dot-delimited name.group string.
dimensionedScalar exp(const dimensionedScalar &ds)
fvPatchField< scalar > fvPatchScalarField
Field< scalar > scalarField
Specialisation of Field<T> for scalar.
A FieldMapper for finite-volume patch fields.
virtual tmp< Field< Type > > patchInternalField() const
Return internal field next to patch.
int debug
Static debugging option.
OBJstream os(runTime.globalPath()/outputName)
dimensionedScalar pow(const dimensionedScalar &ds, const dimensionedScalar &expt)
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...
const std::string patch
OpenFOAM patch number as a std::string.
messageStream Info
Information stream (stdout output on master, null elsewhere)
volScalarField alphaEff("alphaEff", turbulence->nu()/Pr+alphat)
A class for managing temporary objects.
virtual void updateCoeffs()
Update the coefficients associated with the patch field.
static const Foam::dimensionedScalar A("", Foam::dimPressure, 611.21)
static const Foam::dimensionedScalar B("", Foam::dimless, 18.678)