Go to the documentation of this file.
39 Foam::nutURoughWallFunctionFvPatchScalarField::calcNut()
const
41 const label patchi =
patch().index();
48 internalField().
group()
53 const tmp<scalarField> tnuw = turbModel.nu(patchi);
59 tmp<scalarField> tyPlus = calcYPlus(
magUp);
69 const scalar
Re =
magUp[facei]*
y[facei]/nuw[facei] + ROOTVSMALL;
79 Foam::nutURoughWallFunctionFvPatchScalarField::calcYPlus
84 const label patchi =
patch().index();
91 internalField().
group()
95 const tmp<scalarField> tnuw = turbModel.nu(patchi);
101 if (0.0 < roughnessHeight_)
104 const scalar c_1 = 1/(90 - 2.25) + roughnessConstant_;
105 static const scalar c_2 = 2.25/(90 - 2.25);
106 static const scalar c_3 = 2.0*
atan(1.0)/
log(90/2.25);
107 static const scalar c_4 = c_3*
log(2.25);
115 const scalar magUpara =
magUp[facei];
116 const scalar
Re = magUpara*
y[facei]/nuw[facei];
117 const scalar kappaRe = kappa_*
Re;
119 scalar yp = yPlusLam_;
120 const scalar ryPlusLam = 1.0/yp;
123 scalar yPlusLast = 0.0;
124 scalar dKsPlusdYPlus = roughnessHeight_/
y[facei];
127 dKsPlusdYPlus *= roughnessFactor_;
134 scalar KsPlus = yp*dKsPlusdYPlus;
139 scalar yPlusGPrime = 0.0;
143 const scalar t_1 = 1 + roughnessConstant_*KsPlus;
145 yPlusGPrime = roughnessConstant_*KsPlus/t_1;
147 else if (KsPlus > 2.25)
149 const scalar t_1 = c_1*KsPlus - c_2;
150 const scalar t_2 = c_3*
log(KsPlus) - c_4;
151 const scalar sint_2 =
sin(t_2);
152 const scalar logt_1 =
log(t_1);
155 (c_1*sint_2*KsPlus/t_1) + (c_3*logt_1*
cos(t_2));
158 scalar denom = 1.0 +
log(E_*yp) -
G - yPlusGPrime;
159 if (
mag(denom) > VSMALL)
161 yp = (kappaRe + yp*(1 - yPlusGPrime))/denom;
165 mag(ryPlusLam*(yp - yPlusLast)) > tolerance_
179 const scalar magUpara =
magUp[facei];
180 const scalar
Re = magUpara*
y[facei]/nuw[facei];
181 const scalar kappaRe = kappa_*
Re;
183 scalar yp = yPlusLam_;
184 const scalar ryPlusLam = 1.0/yp;
187 scalar yPlusLast = 0.0;
192 yp = (kappaRe + yp)/(1.0 +
log(E_*yp));
197 mag(ryPlusLam*(yp - yPlusLast)) > tolerance_
209 void Foam::nutURoughWallFunctionFvPatchScalarField::writeLocalEntries
215 os.writeEntry(
"roughnessHeight", roughnessHeight_);
216 os.writeEntry(
"roughnessConstant", roughnessConstant_);
217 os.writeEntry(
"roughnessFactor", roughnessFactor_);
218 os.writeEntryIfDifferent<label>(
"maxIter", 10, maxIter_);
219 os.writeEntryIfDifferent<scalar>(
"tolerance", 0.0001, tolerance_);
233 roughnessHeight_(
Zero),
234 roughnessConstant_(
Zero),
235 roughnessFactor_(
Zero),
251 roughnessHeight_(ptf.roughnessHeight_),
252 roughnessConstant_(ptf.roughnessConstant_),
253 roughnessFactor_(ptf.roughnessFactor_),
254 maxIter_(ptf.maxIter_),
255 tolerance_(ptf.tolerance_)
268 roughnessHeight_(
dict.
get<scalar>(
"roughnessHeight")),
269 roughnessConstant_(
dict.
get<scalar>(
"roughnessConstant")),
270 roughnessFactor_(
dict.
get<scalar>(
"roughnessFactor")),
283 roughnessHeight_(rwfpsf.roughnessHeight_),
284 roughnessConstant_(rwfpsf.roughnessConstant_),
285 roughnessFactor_(rwfpsf.roughnessFactor_),
286 maxIter_(rwfpsf.maxIter_),
287 tolerance_(rwfpsf.tolerance_)
299 roughnessHeight_(rwfpsf.roughnessHeight_),
300 roughnessConstant_(rwfpsf.roughnessConstant_),
301 roughnessFactor_(rwfpsf.roughnessFactor_),
302 maxIter_(rwfpsf.maxIter_),
303 tolerance_(rwfpsf.tolerance_)
312 const label patchi =
patch().index();
319 internalField().
group()
325 return calcYPlus(
magUp());
335 writeLocalEntries(os);
336 writeEntry(
"value", os);
virtual void write(Ostream &) const
Write.
Field< scalar > scalarField
Specialisation of Field<T> for scalar.
constexpr const char *const group
Group name for atomic constants.
const dimensionedScalar G
Newtonian constant of gravitation.
A class for managing temporary objects.
static constexpr const zero Zero
Global zero (0)
dimensionedScalar sin(const dimensionedScalar &ds)
This boundary condition provides a wall constraint on the turbulent viscosity, i.e....
static const word propertiesName
Default name of the turbulence properties dictionary.
virtual const volVectorField & U(const turbulenceModel &turb) const
T get(const word &keyword, enum keyType::option matchOpt=keyType::REGEX) const
#define forAll(list, i)
Loop across all elements in list.
ThermalDiffusivity< CompressibleTurbulenceModel< fluidThermo > > turbulenceModel
A finiteVolume patch using a polyPatch and a fvBoundaryMesh.
label max(const labelHashSet &set, label maxValue=labelMin)
Find the max value in labelHashSet, optionally limited by second argument.
virtual void write(Ostream &os) const
Write.
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,...
nutURoughWallFunctionFvPatchScalarField(const fvPatch &, const DimensionedField< scalar, volMesh > &)
Construct from patch and internal field.
dimensionedScalar log(const dimensionedScalar &ds)
Macros for easy insertion into run-time selection tables.
virtual tmp< scalarField > yPlus() const
Calculate and return the yPlus at the boundary.
fvPatchField< vector > fvPatchVectorField
dimensionedSymmTensor sqr(const dimensionedVector &dv)
const std::string patch
OpenFOAM patch number as a std::string.
static const nutWallFunctionFvPatchScalarField & nutw(const turbulenceModel &turbModel, const label patchi)
Return the nut patchField for the given wall patch.
dimensioned< typename typeOfMag< Type >::type > mag(const dimensioned< Type > &dt)
scalarField Re(const UList< complex > &cf)
Extract real component.
dimensionedScalar atan(const dimensionedScalar &ds)
Foam::fvPatchFieldMapper.
static word groupName(StringType base, const word &group)
Create dot-delimited name.group string.
An Ostream is an abstract base class for all output systems (streams, files, token lists,...
makePatchTypeField(fvPatchScalarField, atmBoundaryLayerInletEpsilonFvPatchScalarField)
T getOrDefault(const word &keyword, const T &deflt, enum keyType::option matchOpt=keyType::REGEX) const
The class nutWallFunction is a base class that parents the derived boundary conditions which provide ...
virtual void writeLocalEntries(Ostream &) const
Write local wall function variables.
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...
dimensionedScalar cos(const dimensionedScalar &ds)