Go to the documentation of this file.
57 void Foam::fv::buoyancyTurbSource::buoyancyTurbSourceEpsilon
67 const dictionary& turbDict = turbPtr->coeffDict();
77 const vector gHat(g_.value()/
mag(g_.value()));
91 void Foam::fv::buoyancyTurbSource::buoyancyTurbSourceOmega
106 word(turbPtr->type() +
":gamma")
114 void Foam::fv::buoyancyTurbSource::buoyancyTurbSourceK
116 fvMatrix<scalar>& eqn
131 const word& sourceName,
132 const word& modelType,
139 rhoName_(coeffs_.getOrDefault<
word>(
"rho",
"rho")),
140 alphatName_(coeffs_.getOrDefault<
word>(
"alphat",
"alphat")),
141 Tname_(coeffs_.getOrDefault<
word>(
"T",
"T")),
147 coeffs_.getCheckOrDefault<scalar>
151 [&](
const scalar
x){ return x > SMALL; }
161 if (
mag(g_.value()) < SMALL)
164 <<
"Gravitational field cannot be equal to or less than zero"
168 const auto* turbPtr =
177 <<
"Unable to find a turbulence model."
181 fieldNames_.resize(2);
186 if (!tepsilon.
isTmp())
189 fieldNames_[0] = tepsilon().name();
191 else if (!tomega.
isTmp())
194 fieldNames_[0] = tomega().name();
199 <<
"Unable to find an omega or epsilon field." <<
nl
200 <<
"buoyancyTurbSource needs an omega- or epsilon-based model."
204 fieldNames_[1] = turbPtr->k()().
name();
208 Log <<
" Applying buoyancyTurbSource to: "
209 << fieldNames_[0] <<
" and " << fieldNames_[1]
224 buoyancyTurbSourceK(eqn);
230 buoyancyTurbSourceEpsilon(eqn);
234 buoyancyTurbSourceOmega(eqn);
264 buoyancyTurbSourceK(
alpha,
rho, eqn, fieldi);
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.
Intermediate abstract class for handling cell-set options for the derived fvOptions.
const dimensionSet dimLength(0, 1, 0, 0, 0, 0, 0)
A class for managing temporary objects.
const dimensionSet dimVelocity
void resetApplied()
Resize/reset applied flag list for all fieldNames_ entries.
const dimensionedScalar alpha
Fine-structure constant: default SI units: [].
virtual void addSup(fvMatrix< scalar > &eqn, const label fieldi)
static const word propertiesName
Default name of the turbulence properties dictionary.
static const Foam::dimensionedScalar B("", Foam::dimless, 18.678)
A class representing the concept of a GeometricField of 1 used to avoid unnecessary manipulations for...
bool isTmp() const noexcept
Identical to is_pointer()
Ostream & endl(Ostream &os)
Add newline and flush stream.
const fvMesh & mesh_
Reference to the mesh database.
const dimensionSet dimTime(0, 0, 1, 0, 0, 0, 0)
tmp< fvMatrix< Type > > SuSp(const volScalarField::Internal &, const GeometricField< Type, fvPatchField, volMesh > &)
tmp< fvMatrix< Type > > Sp(const volScalarField::Internal &, const GeometricField< Type, fvPatchField, volMesh > &)
ThermalDiffusivity< CompressibleTurbulenceModel< fluidThermo > > turbulenceModel
dimensionedScalar tanh(const dimensionedScalar &ds)
DimensionedField< scalar, volMesh > Internal
Type of the internal field from which this GeometricField is derived.
void T(FieldField< Field, Type > &f1, const FieldField< Field, Type > &f2)
dimensioned< scalar > dimensionedScalar
Dimensioned scalar obtained from generic dimensioned type.
GeometricField< scalar, fvPatchField, volMesh > volScalarField
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.
Mesh data needed to do the Finite Volume discretisation.
Vector< scalar > vector
A scalar version of the templated Vector.
errorManipArg< error, int > exit(error &err, const int errNo=1)
buoyancyTurbSource(const word &sourceName, const word &modelType, const dictionary &dict, const fvMesh &mesh)
Construct from explicit source name and mesh.
Type & lookupObjectRef(const word &name, const bool recursive=false) const
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
dimensionedSymmTensor sqr(const dimensionedVector &dv)
dimensioned< typename typeOfMag< Type >::type > mag(const dimensioned< Type > &dt)
label k
Boltzmann constant.
A special matrix type and solver, designed for finite volume solutions of scalar equations....
const dimensionSet dimTemperature(0, 0, 0, 1, 0, 0, 0)
word name(const expressions::valueTypeCode typeCode)
A word representation of a valueTypeCode. Empty for INVALID.
defineTypeNameAndDebug(atmAmbientTurbSource, 0)
addToRunTimeSelectionTable(option, atmAmbientTurbSource, dictionary)
static const gravity & New(const Time &runTime)
Construct on Time.
const dimensionSet dimless
Dimensionless.