41 template<
class BasicMomentumTransportModel>
52 const scalar Cpm = pm.
dict().
lookup<scalar>(
C.name());
57 this->mesh_.cellZones()[cellZoneIDs[zonei]];
69 template<
class BasicMomentumTransportModel>
81 const scalar Cpm = pm.
dict().
lookup<scalar>(
C.name());
86 this->mesh_.cellZones()[cellZoneIDs[zonei]];
91 C[celli] = Cpm*Av[celli];
98 template<
class BasicMomentumTransportModel>
109 if (isA<fv::explicitPorositySource>(
fvModels[i]))
112 refCast<const fv::explicitPorositySource>(
fvModels[i]);
114 if (isA<porosityModels::powerLawLopesdaCosta>(eps.
model()))
117 refCast<const porosityModels::powerLawLopesdaCosta>
122 setPorosityCoefficient(Cmu_, pm);
123 setPorosityCoefficient(C1_, pm);
124 setPorosityCoefficient(C2_, pm);
125 setPorosityCoefficient(sigmak_, pm);
126 setPorosityCoefficient(sigmaEps_, pm);
129 setPorosityCoefficient(betap_, pm);
130 setPorosityCoefficient(betad_, pm);
131 setPorosityCoefficient(C4_, pm);
132 setPorosityCoefficient(C5_, pm);
139 template<
class BasicMomentumTransportModel>
142 this->nut_ = Cmu_*
sqr(k_)/epsilon_;
143 this->nut_.correctBoundaryConditions();
148 template<
class BasicMomentumTransportModel>
155 return fvm::Su(CdAv_*(betap_*magU3 - betad_*magU*k_()), k_);
159 template<
class BasicMomentumTransportModel>
170 *(C4_*betap_*epsilon_()/k_()*magU3 - C5_*betad_*magU*epsilon_()),
178 template<
class BasicMomentumTransportModel>
206 this->runTime_.
name(),
222 this->runTime_.
name(),
238 this->runTime_.
name(),
254 this->runTime_.
name(),
270 this->runTime_.
name(),
287 this->runTime_.
name(),
298 this->runTime_.
name(),
309 this->runTime_.
name(),
320 this->runTime_.
name(),
331 this->runTime_.
name(),
343 this->runTime_.
name(),
355 this->runTime_.
name(),
366 if (
type == typeName)
368 this->printCoeffs(
type);
377 template<
class BasicMomentumTransportModel>
391 template<
class BasicMomentumTransportModel>
394 if (!this->turbulence_)
427 epsilon_.boundaryFieldRef().updateCoeffs();
442 + epsilonSource(magU, magU3)
446 epsEqn.
ref().relax();
448 epsEqn.
ref().boundaryManipulate(epsilon_.boundaryFieldRef());
451 bound(epsilon_, this->epsilonMin_);
463 + kSource(magU, magU3)
471 bound(k_, this->kMin_);
#define forAll(list, i)
Loop across all elements in list.
Bound the given scalar field where it is below the specified minimum.
Graphite solid properties.
Field with dimensions and associated with geometry type GeoMesh which is used to size the field and a...
Generic GeometricField class.
IOobject defines the attributes of an object for which implicit objectRegistry management is supporte...
Templated abstract base class for RAS turbulence models.
BasicMomentumTransportModel::alphaField alphaField
virtual tmp< fvScalarMatrix > kSource(const volScalarField::Internal &magU, const volScalarField::Internal &magU3) const
virtual void correct()
Solve the turbulence equations and correct the turbulence viscosity.
void setPorosityCoefficients()
kEpsilonLopesdaCosta(const alphaField &alpha, const rhoField &rho, const volVectorField &U, const surfaceScalarField &alphaRhoPhi, const surfaceScalarField &phi, const viscosity &viscosity, const word &type=typeName)
Construct from components.
void setCdAv(volScalarField::Internal &C, const porosityModels::powerLawLopesdaCosta &pm)
virtual void correctNut()
virtual tmp< fvScalarMatrix > epsilonSource(const volScalarField::Internal &magU, const volScalarField::Internal &magU3) const
void setPorosityCoefficient(volScalarField::Internal &C, const porosityModels::powerLawLopesdaCosta &pm)
virtual bool read()
Re-read model coefficients if they have changed.
BasicMomentumTransportModel::rhoField rhoField
ITstream & lookup(const word &, bool recursive=false, bool patternMatch=true) const
Find and return an entry data stream.
bool found(const word &, bool recursive=false, bool patternMatch=true) const
Search dictionary for given keyword.
static autoPtr< dictionary > New(Istream &)
Construct top-level dictionary on freestore from Istream.
Generic dimensioned Type class.
Eddy viscosity turbulence model base class.
Finite volume constraints.
bool constrain(fvMatrix< Type > &eqn) const
Apply constraints to an equation.
Explicit porosity source.
const porosityModel & model() const
Return the porosity model.
const labelList & cellZoneIDs() const
Return const access to the cell zone IDs.
const dictionary & dict() const
Return dictionary used for model construction.
const scalarField & Av() 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
If object pointer points to valid object:
T & ref() const
Return non-const reference or generate a fatal error.
Abstract base class for all fluid physical properties.
A class for handling words, derived from string.
Foam::fvConstraints & fvConstraints(Foam::fvConstraints::New(mesh))
Foam::fvModels & fvModels(Foam::fvModels::New(mesh))
volScalarField alpha(IOobject("alpha", runTime.name(), mesh, IOobject::READ_IF_PRESENT, IOobject::AUTO_WRITE), lambda *max(Ua &U, zeroSensitivity))
void correct(const RdeltaTType &rDeltaT, const RhoType &rho, volScalarField &psi, const surfaceScalarField &phiCorr, const SpType &Sp, const SuType &Su)
const char *const group
Group name for atomic constants.
const dimensionedScalar G
Newtonian constant of gravitation.
tmp< VolField< typename outerProduct< vector, Type >::type > > grad(const SurfaceField< Type > &ssf)
tmp< VolField< Type > > div(const SurfaceField< Type > &ssf)
tmp< surfaceScalarField > absolute(const tmp< surfaceScalarField > &tphi, const volVectorField &U)
Return the given relative flux in absolute form.
tmp< fvMatrix< Type > > Su(const DimensionedField< Type, volMesh > &, const VolField< Type > &)
tmp< fvMatrix< Type > > laplacian(const VolField< Type > &vf, const word &name)
tmp< fvMatrix< Type > > Sp(const volScalarField::Internal &, const VolField< Type > &)
tmp< fvMatrix< Type > > SuSp(const volScalarField::Internal &, const VolField< Type > &)
tmp< fvMatrix< Type > > div(const surfaceScalarField &flux, const VolField< Type > &vf, const word &name)
tmp< fvMatrix< Type > > ddt(const VolField< Type > &vf)
dimensionedSymmTensor dev(const dimensionedSymmTensor &dt)
bool read(const char *, int32_t &)
intWM_LABEL_SIZE_t label
A label is an int32_t or int64_t as specified by the pre-processor macro WM_LABEL_SIZE.
dimensionedSymmTensor sqr(const dimensionedVector &dv)
dimensionedScalar pow3(const dimensionedScalar &ds)
dimensionedSymmTensor twoSymm(const dimensionedSymmTensor &dt)
const dimensionSet dimless
const dimensionSet dimLength
dimensioned< scalar > mag(const dimensioned< Type > &)
bool bound(volScalarField &, const dimensionedScalar &min)
Bound the given scalar field where it is below the specified min value.
word name(const complex &)
Return a string representation of a complex.
fileType type(const fileName &, const bool checkVariants=true, const bool followLink=true)
Return the file type: directory or file.
SolverPerformance< Type > solve(fvMatrix< Type > &, const word &)
Solve returning the solution statistics given convergence tolerance.