Go to the documentation of this file.
39 Foam::epsilonWallFunctionFvPatchScalarField::blendingType
41 Foam::epsilonWallFunctionFvPatchScalarField::blendingTypeNames
43 { blendingType::STEPWISE ,
"stepwise" },
44 { blendingType::MAX ,
"max" },
45 { blendingType::BINOMIAL ,
"binomial" },
46 { blendingType::EXPONENTIAL,
"exponential" }
64 const volScalarField::Boundary& bf =
epsilon.boundaryField();
69 if (isA<epsilonWallFunctionFvPatchScalarField>(bf[patchi]))
89 const volScalarField::Boundary& bf =
epsilon.boundaryField();
116 if (isA<epsilonWallFunctionFvPatchScalarField>(bf[patchi]))
118 epsilonPatches.append(patchi);
128 cornerWeights_.setSize(bf.size());
130 for (
const auto& patchi : epsilonPatches)
136 G_.setSize(internalField().size(),
Zero);
137 epsilon_.setSize(internalField().size(),
Zero);
152 const volScalarField::Boundary& bf =
epsilon.boundaryField();
155 refCast<const epsilonWallFunctionFvPatchScalarField>(bf[patchi]);
169 forAll(cornerWeights_, patchi)
171 if (!cornerWeights_[patchi].empty())
182 forAll(cornerWeights_, patchi)
184 if (!cornerWeights_[patchi].empty())
203 const label patchi =
patch.index();
221 const scalar Cmu75 =
pow(nutw.
Cmu(), 0.75);
226 const label celli =
patch.faceCells()[facei];
228 const scalar
yPlus = Cmu25*
y[facei]*
sqrt(
k[celli])/nuw[facei];
230 const scalar w = cornerWeights[facei];
232 scalar epsilonBlended = 0.0;
235 const scalar epsilonVis = w*2.0*
k[celli]*nuw[facei]/
sqr(
y[facei]);
238 const scalar epsilonLog =
239 w*Cmu75*
pow(
k[celli], 1.5)/(nutw.
kappa()*
y[facei]);
243 case blendingType::STEPWISE:
247 epsilonBlended = epsilonVis;
251 epsilonBlended = epsilonLog;
256 case blendingType::MAX:
259 epsilonBlended =
max(epsilonVis, epsilonLog);
263 case blendingType::BINOMIAL:
269 pow(epsilonVis, n_) +
pow(epsilonLog, n_),
275 case blendingType::EXPONENTIAL:
279 const scalar invGamma = 1.0/(Gamma + ROOTVSMALL);
281 epsilonVis*
exp(-Gamma) + epsilonLog*
exp(-invGamma);
292 *(nutw[facei] + nuw[facei])
294 *Cmu25*
sqrt(
k[celli])
311 blending_(blendingType::STEPWISE),
313 lowReCorrection_(
false),
332 blending_(ptf.blending_),
354 blendingTypeNames.getOrDefault
358 blendingType::STEPWISE
389 blending_(ewfpsf.blending_),
408 blending_(ewfpsf.blending_),
426 if (
patch().index() == master_)
436 return epsilonPatch(master_).G();
445 if (
patch().index() == master_)
455 return epsilonPatch(master_).epsilon(
init);
471 internalField().
group()
477 if (
patch().index() == master_)
479 createAveragingWeights();
480 calculateTurbulenceFields(turbModel,
G(
true),
epsilon(
true));
488 FieldType&
G = db().lookupObjectRef<FieldType>(turbModel.
GName());
490 FieldType&
epsilon =
const_cast<FieldType&
>(internalField());
494 const label celli =
patch().faceCells()[facei];
496 G[celli] =
G0[celli];
519 internalField().
group()
525 if (
patch().index() == master_)
527 createAveragingWeights();
528 calculateTurbulenceFields(turbModel,
G(
true),
epsilon(
true));
536 FieldType&
G = db().lookupObjectRef<FieldType>(turbModel.
GName());
538 FieldType&
epsilon =
const_cast<FieldType&
>(internalField());
545 const scalar w = weights[facei];
549 const label celli =
patch().faceCells()[facei];
551 G[celli] = (1.0 - w)*
G[celli] + w*
G0[celli];
553 epsilonf[facei] =
epsilon[celli];
566 if (manipulatedMatrix())
583 if (manipulatedMatrix())
597 if (weights[facei] > tolerance_)
601 constraintCells.append(celli);
602 constraintValues.append(
fld[celli]);
609 <<
": number of constrained cells = " << constraintCells.size()
610 <<
" out of " <<
patch().size()
614 matrix.
setValues(constraintCells, constraintValues);
626 os.
writeEntry(
"blending", blendingTypeNames[blending_]);
int debug
Static debugging option.
const DimensionedField< scalar, volMesh > & internalField() const
Return dimensioned internal field reference.
virtual tmp< Field< Type > > snGrad() const
Return patch-normal gradient.
Field< scalar > scalarField
Specialisation of Field<T> for scalar.
Defines the attributes of an object for which implicit objectRegistry management is supported,...
virtual tmp< volScalarField > nu() const =0
Return the laminar viscosity.
Enum is a wrapper around a list of names/values that represent particular enumeration (or int) values...
virtual label & master()
Return non-const access to the master patch ID.
scalar Cmu() const
Return Cmu.
Info<< "Reading field U\n"<< endl;volVectorField U(IOobject("U", runTime.timeName(), mesh, IOobject::MUST_READ, IOobject::AUTO_WRITE), mesh);volScalarField rho(IOobject("rho", runTime.timeName(), mesh, IOobject::NO_READ, IOobject::AUTO_WRITE), thermo.rho());volVectorField rhoU(IOobject("rhoU", runTime.timeName(), mesh, IOobject::NO_READ, IOobject::NO_WRITE), rho *U);volScalarField rhoE(IOobject("rhoE", runTime.timeName(), mesh, IOobject::NO_READ, IOobject::NO_WRITE), rho *(e+0.5 *magSqr(U)));surfaceScalarField pos(IOobject("pos", runTime.timeName(), mesh), mesh, dimensionedScalar("pos", dimless, 1.0));surfaceScalarField neg(IOobject("neg", runTime.timeName(), mesh), mesh, dimensionedScalar("neg", dimless, -1.0));surfaceScalarField phi("phi", fvc::flux(rhoU));Info<< "Creating turbulence model\n"<< endl;autoPtr< compressible::turbulenceModel > turbulence(compressible::turbulenceModel::New(rho, U, phi, thermo))
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)
virtual void updateWeightedCoeffs(const scalarField &weights)
Update the coefficients associated with the patch field.
static scalar tolerance_
Tolerance used in weighted calculations.
static word timeName(const scalar t, const int precision=precision_)
static const word propertiesName
Default name of the turbulence properties dictionary.
Ostream & endl(Ostream &os)
Add newline and flush stream.
dimensionedScalar exp(const dimensionedScalar &ds)
prefixOSstream Pout
OSstream wrapped stdout (std::cout) with parallel prefix.
const nearWallDist & y() const
Return the near wall distances.
virtual void createAveragingWeights()
epsilonWallFunctionFvPatchScalarField(const fvPatch &, const DimensionedField< scalar, volMesh > &)
Construct from patch and internal field.
const dimensionedScalar epsilon0
Electric constant: default SI units: [F/m].
#define forAll(list, i)
Loop across all elements in list.
virtual void write(Ostream &) const
Write.
dimensionedScalar pow025(const dimensionedScalar &ds)
dimensionedScalar pow4(const dimensionedScalar &ds)
scalarField & epsilon(bool init=false)
Return non-const access to the master's epsilon field.
tmp< faMatrix< Type > > operator==(const faMatrix< Type > &, const faMatrix< Type > &)
const dimensionedScalar G0
Conductance quantum: default SI units: [S].
A finiteVolume patch using a polyPatch and a fvBoundaryMesh.
This boundary condition provides a wall constraint on the turbulent kinetic energy dissipation rate,...
virtual tmp< Field< Type > > patchInternalField() const
Return internal field next to patch as patch field.
virtual void calculate(const turbulenceModel &turbulence, const List< scalar > &cornerWeights, const fvPatch &patch, scalarField &G, scalarField &epsilon)
Calculate the epsilon and G.
const bool lowReCorrection_
Apply low-Re correction term (default = no)
dimensioned< scalar > dimensionedScalar
Dimensioned scalar obtained from generic dimensioned type.
virtual void write(Ostream &) const
Write.
T getCheckOrDefault(const word &keyword, const T &deflt, const Predicate &pred, enum keyType::option matchOpt=keyType::REGEX) const
dimensionedScalar pow(const dimensionedScalar &ds, const dimensionedScalar &expt)
label max(const labelHashSet &set, label maxValue=labelMin)
Find the max value in labelHashSet, optionally limited by second argument.
gmvFile<< "tracers "<< particles.size()<< nl;for(const passiveParticle &p :particles){ gmvFile<< p.position().x()<< ' ';}gmvFile<< nl;for(const passiveParticle &p :particles){ gmvFile<< p.position().y()<< ' ';}gmvFile<< nl;for(const passiveParticle &p :particles){ gmvFile<< p.position().z()<< ' ';}gmvFile<< nl;for(const word &name :lagrangianScalarNames){ IOField< scalar > fld(IOobject(name, runTime.timeName(), cloud::prefix, mesh, IOobject::MUST_READ, IOobject::NO_WRITE))
word GName() const
Helper function to return the name of the turbulence G field.
scalarField & G(bool init=false)
Return non-const access to the master's G field.
Abstract base class for turbulence models (RAS, LES and laminar).
virtual void updateCoeffs()
Update the coefficients associated with the patch field.
A list of keyword definitions, which are a keyword followed by a number of values (eg,...
OBJstream os(runTime.globalPath()/outputName)
Macros for easy insertion into run-time selection tables.
virtual tmp< volScalarField > k() const =0
Return the turbulence kinetic energy.
Mesh data needed to do the Finite Volume discretisation.
virtual void manipulateMatrix(fvMatrix< scalar > &matrix)
Manipulate matrix.
static scalar yPlusLam(const scalar kappa, const scalar E)
Estimate the y+ at the intersection of the two sublayers.
virtual void manipulateMatrix(fvMatrix< Type > &matrix)
Manipulate matrix.
virtual void calculateTurbulenceFields(const turbulenceModel &turbulence, scalarField &G0, scalarField &epsilon0)
Main driver to calculate the turbulence fields.
virtual void updateCoeffs()
Update the coefficients associated with the patch field.
virtual epsilonWallFunctionFvPatchScalarField & epsilonPatch(const label patchi)
Helper function to return non-const access to an epsilon patch.
void setValues(const labelUList &cellLabels, const Type &value)
bool changing() const noexcept
Is mesh changing (topology changing and/or moving)
dimensionedSymmTensor sqr(const dimensionedVector &dv)
const std::string patch
OpenFOAM patch number as a std::string.
static MinMax< T > ge(const T &minVal)
A semi-infinite range from minVal to the type max.
static const nutWallFunctionFvPatchScalarField & nutw(const turbulenceModel &turbModel, const label patchi)
Return the nut patchField for the given wall patch.
dimensionedScalar sqrt(const dimensionedScalar &ds)
dimensioned< typename typeOfMag< Type >::type > mag(const dimensioned< Type > &dt)
scalar kappa() const
Return kappa.
label k
Boltzmann constant.
const dimensionedScalar e
Elementary charge.
label master_
Master patch ID.
A special matrix type and solver, designed for finite volume solutions of scalar equations....
Ostream & writeEntry(const keyType &key, const T &value)
Write a keyword/value entry.
virtual const labelUList & faceCells() const
Return faceCells.
const Time & time() const
Return the top-level database.
Foam::fvPatchFieldMapper.
const fvPatch & patch() const
Return patch.
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 ...
Smooth ATC in cells next to a set of patches supplied by type.
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...
const dimensionSet dimless
Dimensionless.
const volVectorField & U() const
Access function to velocity field.