Go to the documentation of this file.
43 template<
class BasicTurbulenceModel>
54 const scalar Cpm = pm.
dict().
get<scalar>(
C.name());
56 for (
const label zonei : cellZoneIDs)
60 for (
const label celli :
cells)
69 template<
class BasicTurbulenceModel>
81 const scalar Cpm = pm.
dict().
get<scalar>(
C.name());
83 for (
const label zonei : cellZoneIDs)
89 const label celli =
cells[i];
90 C[celli] = Cpm*Sigma[i];
97 template<
class BasicTurbulenceModel>
104 if (isA<fv::explicitPorositySource>(
fvOptions[i]))
107 refCast<const fv::explicitPorositySource>(
fvOptions[i]);
109 if (isA<porosityModels::powerLawLopesdaCosta>(eps.
model()))
112 refCast<const porosityModels::powerLawLopesdaCosta>
117 setPorosityCoefficient(Cmu_, pm);
118 Cmu_.correctBoundaryConditions();
119 setPorosityCoefficient(C1_, pm);
120 setPorosityCoefficient(C2_, pm);
121 setPorosityCoefficient(sigmak_, pm);
122 setPorosityCoefficient(sigmaEps_, pm);
123 sigmak_.correctBoundaryConditions();
124 sigmaEps_.correctBoundaryConditions();
126 setCdSigma(CdSigma_, pm);
127 setPorosityCoefficient(betap_, pm);
128 setPorosityCoefficient(betad_, pm);
129 setPorosityCoefficient(C4_, pm);
130 setPorosityCoefficient(C5_, pm);
137 template<
class BasicTurbulenceModel>
140 this->nut_ = Cmu_*
sqr(k_)/epsilon_;
141 this->nut_.correctBoundaryConditions();
144 BasicTurbulenceModel::correctNut();
148 template<
class BasicTurbulenceModel>
155 return fvm::Su(CdSigma_*(betap_*magU3 - betad_*magU*k_()), k_);
159 template<
class BasicTurbulenceModel>
170 *(C4_*betap_*epsilon_()/k_()*magU3 - C5_*betad_*magU*epsilon_()),
178 template<
class BasicTurbulenceModel>
187 const word& propertiesName,
208 this->runTime_.timeName(),
224 this->runTime_.timeName(),
240 this->runTime_.timeName(),
256 this->runTime_.timeName(),
272 this->runTime_.timeName(),
289 this->runTime_.timeName(),
300 this->runTime_.timeName(),
311 this->runTime_.timeName(),
322 this->runTime_.timeName(),
333 this->runTime_.timeName(),
345 this->runTime_.timeName(),
357 this->runTime_.timeName(),
365 bound(k_, this->kMin_);
366 bound(epsilon_, this->epsilonMin_);
368 if (
type == typeName)
370 this->printCoeffs(
type);
373 setPorosityCoefficients();
379 template<
class BasicTurbulenceModel>
391 template<
class BasicTurbulenceModel>
394 if (!this->turbulence_)
423 epsilon_.boundaryFieldRef().updateCoeffs();
438 + epsilonSource(magU, magU3)
442 epsEqn.
ref().relax();
444 epsEqn.
ref().boundaryManipulate(epsilon_.boundaryFieldRef());
447 bound(epsilon_, this->epsilonMin_);
459 + kSource(magU, magU3)
467 bound(k_, this->kMin_);
Defines the attributes of an object for which implicit objectRegistry management is supported,...
virtual bool read()
Re-read model coefficients if they have changed.
tmp< GeometricField< typename outerProduct< vector, Type >::type, fvPatchField, volMesh >> grad(const GeometricField< Type, fvsPatchField, surfaceMesh > &ssf)
A class for handling words, derived from Foam::string.
const labelList & cellZoneIDs() const
Return const access to the cell zone IDs.
void correct(GeometricField< Type, fvPatchField, volMesh > &field)
Apply correction to field.
void clear() const noexcept
const dimensionSet dimLength(0, 1, 0, 0, 0, 0, 0)
const dimensionedScalar G
Newtonian constant of gravitation.
A class for managing temporary objects.
bool found(const word &keyword, enum keyType::option matchOpt=keyType::REGEX) const
Search for an entry (const access) with the given keyword.
optionList(const optionList &)=delete
No copy construct.
Templated abstract base class for RAS turbulence models.
static options & New(const fvMesh &mesh)
Construct fvOptions and register to database if not present.
const dimensionedScalar alpha
Fine-structure constant: default SI units: [].
volScalarField & bound(volScalarField &, const dimensionedScalar &lowerBound)
Bound the given scalar field if it has gone unbounded.
Variant of the standard k-epsilon turbulence model with additional source terms to handle the changes...
bool read(const char *buf, int32_t &val)
Same as readInt32.
tmp< GeometricField< Type, fvPatchField, volMesh > > div(const GeometricField< Type, fvsPatchField, surfaceMesh > &ssf)
void setCdSigma(volScalarField::Internal &C, const porosityModels::powerLawLopesdaCosta &pm)
BasicTurbulenceModel::transportModel transportModel
virtual tmp< fvScalarMatrix > epsilonSource(const volScalarField::Internal &magU, const volScalarField::Internal &magU3) const
T get(const word &keyword, enum keyType::option matchOpt=keyType::REGEX) const
virtual void correct()
Solve the turbulence equations and correct the turbulence viscosity.
tmp< fvMatrix< Type > > Su(const DimensionedField< Type, volMesh > &, const GeometricField< Type, fvPatchField, volMesh > &)
#define forAll(list, i)
Loop across all elements in list.
virtual tmp< fvScalarMatrix > kSource(const volScalarField::Internal &magU, const volScalarField::Internal &magU3) const
tmp< fvMatrix< Type > > SuSp(const volScalarField::Internal &, const GeometricField< Type, fvPatchField, volMesh > &)
Applies an explicit porosity source to the momentum equation within a specified region.
tmp< fvMatrix< Type > > Sp(const volScalarField::Internal &, const GeometricField< Type, fvPatchField, volMesh > &)
dimensionedScalar pow3(const dimensionedScalar &ds)
tmp< fvMatrix< Type > > laplacian(const GeometricField< Type, fvPatchField, volMesh > &vf, const word &name)
DimensionedField< scalar, volMesh > Internal
Type of the internal field from which this GeometricField is derived.
virtual void correctNut()
Bound the given scalar field if it has gone unbounded.
SolverPerformance< Type > solve(faMatrix< Type > &, Istream &)
Solve returning the solution statistics given convergence tolerance.
dimensioned< scalar > dimensionedScalar
Dimensioned scalar obtained from generic dimensioned type.
void setPorosityCoefficients()
const porosityModel & model() const
Access to the porosityModel.
const scalarField & Sigma() const
Return the porosity surface area per unit volume zone field.
Generic dimensioned Type class.
BasicTurbulenceModel::alphaField alphaField
tmp< fvMatrix< Type > > ddt(const GeometricField< Type, fvPatchField, volMesh > &vf)
dimensionedSymmTensor sqr(const dimensionedVector &dv)
void setPorosityCoefficient(volScalarField::Internal &C, const porosityModels::powerLawLopesdaCosta &pm)
dimensioned< typename typeOfMag< Type >::type > mag(const dimensioned< Type > &dt)
fileName::Type type(const fileName &name, const bool followLink=true)
Return the file type: DIRECTORY or FILE, normally following symbolic links.
Eddy viscosity turbulence model base class.
static word groupName(StringType base, const word &group)
Create dot-delimited name.group string.
Graphite solid properties.
tmp< fvMatrix< Type > > div(const surfaceScalarField &flux, const GeometricField< Type, fvPatchField, volMesh > &vf, const word &name)
dimensionedSymmTensor twoSymm(const dimensionedSymmTensor &dt)
const dictionary & dict() const
Return dictionary used for model construction.
Variant of the power law porosity model with spatially varying drag coefficient.
tmp< surfaceScalarField > absolute(const tmp< surfaceScalarField > &tphi, const volVectorField &U)
Return the given relative flux in absolute form.
BasicTurbulenceModel::rhoField rhoField
const dimensionSet dimless
Dimensionless.
dimensionedSymmTensor dev(const dimensionedSymmTensor &dt)