43template<
class BasicTurbulenceModel>
54 const scalar Cpm = pm.
dict().
get<scalar>(
C.name());
56 for (
const label zonei : cellZoneIDs)
60 for (
const label celli :
cells)
69template<
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];
97template<
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);
137template<
class BasicTurbulenceModel>
140 this->nut_ = Cmu_*
sqr(k_)/epsilon_;
141 this->nut_.correctBoundaryConditions();
144 BasicTurbulenceModel::correctNut();
148template<
class BasicTurbulenceModel>
155 return fvm::Su(CdSigma_*(betap_*magU3 - betad_*magU*k_()), k_);
159template<
class BasicTurbulenceModel>
170 *(C4_*betap_*epsilon_()/k_()*magU3 - C5_*betad_*magU*epsilon_()),
178template<
class BasicTurbulenceModel>
187 const word& propertiesName,
344 IOobject::groupName(
"k", alphaRhoPhi.group()),
356 IOobject::groupName(
"epsilon", alphaRhoPhi.group()),
368 if (
type == typeName)
370 this->printCoeffs(
type);
379template<
class BasicTurbulenceModel>
391template<
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_);
Bound the given scalar field if it has gone unbounded.
Graphite solid properties.
Field with dimensions and associated with geometry type GeoMesh which is used to size the field and a...
Defines the attributes of an object for which implicit objectRegistry management is supported,...
Templated abstract base class for RAS turbulence models.
virtual tmp< fvScalarMatrix > epsilonSource() const
Variant of the standard k-epsilon turbulence model with additional source terms to handle the changes...
BasicTurbulenceModel::alphaField alphaField
BasicTurbulenceModel::rhoField rhoField
virtual void correct()
Solve the turbulence equations and correct the turbulence viscosity.
void setCdSigma(volScalarField::Internal &C, const porosityModels::powerLawLopesdaCosta &pm)
void setPorosityCoefficients()
virtual void correctNut()
BasicTurbulenceModel::transportModel transportModel
void setPorosityCoefficient(volScalarField::Internal &C, const porosityModels::powerLawLopesdaCosta &pm)
virtual bool read()
Re-read model coefficients if they have changed.
virtual tmp< fvScalarMatrix > kSource() const
Add explicit source for turbulent kinetic energy.
static autoPtr< Time > New()
Construct (dummy) Time - no functionObjects or libraries.
T get(const word &keyword, enum keyType::option matchOpt=keyType::REGEX) const
bool found(const word &keyword, enum keyType::option matchOpt=keyType::REGEX) const
Search for an entry (const access) with the given keyword.
Generic dimensioned Type class.
Eddy viscosity turbulence model base class.
List of finite volume options.
Applies an explicit porosity source to the momentum equation within a specified region.
const porosityModel & model() const noexcept
Access to the porosityModel.
const labelList & cellZoneIDs() const
Return const access to the cell zone IDs.
const dictionary & dict() const
Return dictionary used for model construction.
const scalarField & Sigma() const
Return the porosity surface area per unit volume zone field.
Variant of the power law porosity model with spatially varying drag coefficient.
A class for managing temporary objects.
void clear() const noexcept
A class for handling words, derived from Foam::string.
tmp< GeometricField< typename outerProduct< vector, Type >::type, fvPatchField, volMesh > > grad(const GeometricField< Type, fvsPatchField, surfaceMesh > &ssf)
tmp< GeometricField< Type, fvPatchField, volMesh > > div(const GeometricField< Type, fvsPatchField, surfaceMesh > &ssf)
tmp< surfaceScalarField > absolute(const tmp< surfaceScalarField > &tphi, const volVectorField &U)
Return the given relative flux in absolute form.
tmp< fvMatrix< Type > > laplacian(const GeometricField< Type, fvPatchField, volMesh > &vf, const word &name)
tmp< fvMatrix< Type > > div(const surfaceScalarField &flux, const GeometricField< Type, fvPatchField, volMesh > &vf, const word &name)
zeroField Su(const Foam::zero, const GeometricField< Type, fvPatchField, volMesh > &)
A no-op source.
tmp< fvMatrix< Type > > ddt(const GeometricField< Type, fvPatchField, volMesh > &vf)
zeroField SuSp(const Foam::zero, const GeometricField< Type, fvPatchField, volMesh > &)
A no-op source.
zeroField Sp(const Foam::zero, const GeometricField< Type, fvPatchField, volMesh > &)
A no-op source.
dimensionedSymmTensor dev(const dimensionedSymmTensor &dt)
const dimensionSet dimless
Dimensionless.
bool read(const char *buf, int32_t &val)
Same as readInt32.
dimensionedSymmTensor sqr(const dimensionedVector &dv)
const dimensionSet dimLength(0, 1, 0, 0, 0, 0, 0)
dimensionedScalar pow3(const dimensionedScalar &ds)
volScalarField & bound(volScalarField &, const dimensionedScalar &lowerBound)
Bound the given scalar field if it has gone unbounded.
dimensionedSymmTensor twoSymm(const dimensionedSymmTensor &dt)
fileName::Type type(const fileName &name, const bool followLink=true)
Return the file type: DIRECTORY or FILE, normally following symbolic links.
dimensioned< typename typeOfMag< Type >::type > mag(const dimensioned< Type > &dt)
#define forAll(list, i)
Loop across all elements in list.