Go to the documentation of this file.
39 Foam::omegaWallFunctionFvPatchScalarField::blendingType
41 Foam::omegaWallFunctionFvPatchScalarField::blendingTypeNames
43 { blendingType::STEPWISE ,
"stepwise" },
44 { blendingType::MAX ,
"max" },
45 { blendingType::BINOMIAL2 ,
"binomial2" },
46 { blendingType::BINOMIAL ,
"binomial" },
47 { blendingType::EXPONENTIAL,
"exponential" },
48 { blendingType::TANH,
"tanh" }
65 const volScalarField::Boundary& bf =
omega.boundaryField();
70 if (isA<omegaWallFunctionFvPatchScalarField>(bf[patchi]))
117 if (isA<omegaWallFunctionFvPatchScalarField>(bf[patchi]))
119 omegaPatches.append(patchi);
129 cornerWeights_.setSize(bf.size());
130 for (
const auto& patchi : omegaPatches)
136 G_.setSize(internalField().size(), 0.0);
137 omega_.setSize(internalField().size(), 0.0);
155 refCast<const omegaWallFunctionFvPatchScalarField>(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();
225 const label celli =
patch.faceCells()[facei];
226 const scalar
yPlus = Cmu25*
y[facei]*
sqrt(
k[celli])/nuw[facei];
227 const scalar w = cornerWeights[facei];
230 const scalar omegaVis = 6.0*nuw[facei]/(beta1_*
sqr(
y[facei]));
233 const scalar omegaLog =
sqrt(
k[celli])/(Cmu25*nutw.
kappa()*
y[facei]);
237 case blendingType::STEPWISE:
241 omega0[celli] += w*omegaLog;
245 omega0[celli] += w*omegaVis;
250 case blendingType::MAX:
253 omega0[celli] +=
max(omegaVis, omegaLog);
257 case blendingType::BINOMIAL2:
260 omega0[celli] += w*
sqrt(
sqr(omegaVis) +
sqr(omegaLog));
264 case blendingType::BINOMIAL:
269 pow(omegaVis, n_) +
pow(omegaLog, n_),
275 case blendingType::EXPONENTIAL:
279 const scalar invGamma = 1.0/(Gamma + ROOTVSMALL);
282 w*(omegaVis*
exp(-Gamma) + omegaLog*
exp(-invGamma));
286 case blendingType::TANH:
290 const scalar omegab1 = omegaVis + omegaLog;
291 const scalar omegab2 =
292 pow(
pow(omegaVis, 1.2) +
pow(omegaLog, 1.2), 1.0/1.2);
294 omega0[celli] += phiTanh*omegab1 + (1.0 - phiTanh)*omegab2;
299 if (!(blending_ == blendingType::STEPWISE) ||
yPlus > nutw.
yPlusLam())
303 *(nutw[facei] + nuw[facei])
305 *Cmu25*
sqrt(
k[celli])
321 blending_(blendingType::BINOMIAL2),
341 blending_(ptf.blending_),
362 blendingTypeNames.getOrDefault
366 blendingType::BINOMIAL2
389 <<
"Using deprecated 'blended' keyword"
390 <<
nl <<
" Please use either of the below for the same behaviour:"
391 <<
nl <<
" 'blending binomial2;' for 'blended on;'"
392 <<
nl <<
" 'blending stepwise;' for 'blended off;'"
393 <<
nl <<
" OVERWRITING: 'blended' keyword -> 'blending' enum"
400 blending_ = blendingType::BINOMIAL2;
404 blending_ = blendingType::STEPWISE;
419 blending_(owfpsf.blending_),
437 blending_(owfpsf.blending_),
455 if (
patch().index() == master_)
465 return omegaPatch(master_).G();
474 if (
patch().index() == master_)
484 return omegaPatch(master_).omega(init);
500 internalField().
group()
506 if (
patch().index() == master_)
508 createAveragingWeights();
509 calculateTurbulenceFields(turbModel,
G(
true), omega(
true));
517 FieldType&
G = db().lookupObjectRef<FieldType>(turbModel.
GName());
519 FieldType& omega =
const_cast<FieldType&
>(internalField());
523 const label celli =
patch().faceCells()[facei];
525 G[celli] =
G0[celli];
526 omega[celli] = omega0[celli];
548 internalField().
group()
554 if (
patch().index() == master_)
556 createAveragingWeights();
557 calculateTurbulenceFields(turbModel,
G(
true), omega(
true));
565 FieldType&
G = db().lookupObjectRef<FieldType>(turbModel.
GName());
567 FieldType& omega =
const_cast<FieldType&
>(internalField());
574 const scalar w = weights[facei];
578 const label celli =
patch().faceCells()[facei];
580 G[celli] = (1.0 - w)*
G[celli] + w*
G0[celli];
581 omega[celli] = (1.0 - w)*omega[celli] + w*omega0[celli];
582 omegaf[facei] = omega[celli];
595 if (manipulatedMatrix())
612 if (manipulatedMatrix())
626 if (tolerance_ < weights[facei])
630 constraintCells.append(celli);
631 constraintValues.append(
fld[celli]);
638 <<
": number of constrained cells = " << constraintCells.size()
639 <<
" out of " <<
patch().size()
643 matrix.
setValues(constraintCells, constraintValues);
654 os.
writeEntry(
"blending", blendingTypeNames[blending_]);
int debug
Static debugging option.
virtual void createAveragingWeights()
omegaWallFunctionFvPatchScalarField(const fvPatch &, const DimensionedField< scalar, volMesh > &)
Construct from patch and internal field.
bool changing() const
Is mesh changing (topology changing and/or moving)
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.
virtual void manipulateMatrix(fvMatrix< scalar > &matrix)
Manipulate matrix.
virtual void calculate(const turbulenceModel &turbulence, const List< scalar > &cornerWeights, const fvPatch &patch, scalarField &G, scalarField &omega)
Calculate the omega and G.
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...
const dimensionSet dimless(0, 0, 0, 0, 0, 0, 0)
Dimensionless.
This boundary condition provides a wall constraint on the specific dissipation rate,...
scalar Cmu() const
Return Cmu.
virtual void updateWeightedCoeffs(const scalarField &weights)
Update the coefficients associated with the patch field.
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)
bool found(const word &keyword, enum keyType::option matchOpt=keyType::REGEX) const
Search for an entry (const access) with the given keyword.
static word timeName(const scalar t, const int precision=precision_)
scalar beta1_
beta1 coefficient
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)
T get(const word &keyword, enum keyType::option matchOpt=keyType::REGEX) const
prefixOSstream Pout
An Ostream wrapper for parallel output to std::cout.
label master_
Master patch ID.
const nearWallDist & y() const
Return the near wall distances.
#define forAll(list, i)
Loop across all elements in list.
virtual void write(Ostream &) const
Write.
dimensionedScalar pow025(const dimensionedScalar &ds)
scalarField & omega(bool init=false)
Return non-const access to the master's omega field.
void setValues(const labelUList &cells, const UList< Type > &values)
Set solution in given cells to the specified values.
dimensionedScalar pow4(const dimensionedScalar &ds)
linear/upwind blended differencing scheme.
tmp< faMatrix< Type > > operator==(const faMatrix< Type > &, const faMatrix< Type > &)
dimensionedScalar tanh(const dimensionedScalar &ds)
const dimensionedScalar G0
Conductance quantum: default SI units: [S].
virtual void write(Ostream &) const
Write.
A finiteVolume patch using a polyPatch and a fvBoundaryMesh.
virtual tmp< Field< Type > > patchInternalField() const
Return internal field next to patch as patch field.
dimensioned< scalar > dimensionedScalar
Dimensioned scalar obtained from generic dimensioned type.
virtual omegaWallFunctionFvPatchScalarField & omegaPatch(const label patchi)
Helper function to return non-const access to an omega patch.
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.
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,...
Macros for easy insertion into run-time selection tables.
virtual void calculateTurbulenceFields(const turbulenceModel &turbulence, scalarField &G0, scalarField &omega0)
Main driver to calculate the turbulence fields.
virtual tmp< volScalarField > k() const =0
Return the turbulence kinetic energy.
Mesh data needed to do the Finite Volume discretisation.
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 updateCoeffs()
Update the coefficients associated with the patch field.
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.
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 void updateCoeffs()
Update the coefficients associated with the patch field.
virtual const labelUList & faceCells() const
Return faceCells.
virtual label & master()
Return non-const access to the master patch ID.
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)
#define IOWarningInFunction(ios)
Report an IO warning using Foam::Warning.
static scalar tolerance_
Tolerance used in weighted calculations.
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 volVectorField & U() const
Access function to velocity field.
scalarField & G(bool init=false)
Return non-const access to the master's G field.