Go to the documentation of this file.
46 multiphaseMangrovesTurbulenceModel,
66 mesh_.time().timeName(),
79 const scalar a = aZone_[i];
80 const scalar
N = NZone_[i];
81 const scalar Ckp = CkpZone_[i];
82 const scalar Cd = CdZone_[i];
84 for (label zonei : zoneIDs_[i])
86 const cellZone& cz = mesh_.cellZones()[zonei];
88 for (label celli : cz)
90 kCoeff[celli] = Ckp*Cd*a*
N*
mag(
U[celli]);
95 kCoeff.correctBoundaryConditions();
111 typeName +
":epsilonCoeff",
112 mesh_.time().timeName(),
125 const scalar a = aZone_[i];
126 const scalar
N = NZone_[i];
127 const scalar Cep = CepZone_[i];
128 const scalar Cd = CdZone_[i];
130 for (label zonei : zoneIDs_[i])
132 const cellZone& cz = mesh_.cellZones()[zonei];
134 for (label celli : cz)
136 epsilonCoeff[celli] = Cep*Cd*a*
N*
mag(
U[celli]);
141 epsilonCoeff.correctBoundaryConditions();
143 return tepsilonCoeff;
149 Foam::fv::multiphaseMangrovesTurbulenceModel::multiphaseMangrovesTurbulenceModel
152 const word& modelType,
165 epsilonName_(
"epsilon")
181 if (eqn.
psi().name() == epsilonName_)
189 else if (eqn.
psi().name() == kName_)
209 if (eqn.
psi().name() == epsilonName_)
217 else if (eqn.
psi().name() == kName_)
251 const wordList regionNames(regionsDict.toc());
264 const word zoneName(modelDict.
get<
word>(
"cellZone"));
270 <<
"Unable to find cellZone " << zoneName <<
nl
Defines the attributes of an object for which implicit objectRegistry management is supported,...
scalarList aZone_
Width of the vegetation element.
const dimensionSet dimless(0, 0, 0, 0, 0, 0, 0)
Dimensionless.
A class for handling words, derived from Foam::string.
tmp< volScalarField > kCoeff(const volVectorField &U) const
Return the k coefficient.
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.
const cellZoneMesh & cellZones() const
Return cell zone mesh.
scalarList CdZone_
Drag coefficient.
virtual void addSup(fvMatrix< scalar > &eqn, const label fieldi)
Add implicit contribution to momentum equation.
T get(const word &keyword, enum keyType::option matchOpt=keyType::REGEX) const
const GeometricField< Type, fvPatchField, volMesh > & psi() const
const fvMesh & mesh_
Reference to the mesh database.
#define forAll(list, i)
Loop across all elements in list.
A special matrix type and solver, designed for finite volume solutions of scalar equations.
const dimensionSet dimTime(0, 0, 1, 0, 0, 0, 0)
Base abstract class for handling finite volume options (i.e. fvOption).
word name(const complex &c)
Return string representation of complex.
dictionary coeffs_
Dictionary containing source coefficients.
wordList fieldNames_
Field names to apply source to - populated by derived models.
virtual bool read(const dictionary &dict)
Read dictionary.
const dictionary & subDict(const word &keyword, enum keyType::option matchOpt=keyType::REGEX) const
Find and return a sub-dictionary.
bool readEntry(const word &keyword, T &val, enum keyType::option matchOpt=keyType::REGEX, bool mandatory=true) const
dimensioned< scalar > dimensionedScalar
Dimensioned scalar obtained from generic dimensioned type.
void resize(const label newSize)
Adjust allocated size of list.
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.
virtual bool read(const dictionary &dict)
Read source dictionary.
tmp< volScalarField > epsilonCoeff(const volVectorField &U) const
Return the epsilon coefficient.
wordList names() const
A list of the zone names.
Calculate the matrix for implicit and explicit sources.
errorManipArg< error, int > exit(error &err, const int errNo=1)
tmp< DimensionedField< TypeR, GeoMesh > > New(const tmp< DimensionedField< TypeR, GeoMesh >> &tdf1, const word &name, const dimensionSet &dimensions)
Global function forwards to reuseTmpDimensionedField::New.
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
labelListList zoneIDs_
Zone indices.
Internal & ref(const bool updateAccessTime=true)
Return a reference to the dimensioned internal field.
dimensioned< typename typeOfMag< Type >::type > mag(const dimensioned< Type > &dt)
A special matrix type and solver, designed for finite volume solutions of scalar equations....
scalarList NZone_
Number of plants per unit of area.
defineTypeNameAndDebug(atmAmbientTurbSource, 0)
addToRunTimeSelectionTable(option, atmAmbientTurbSource, dictionary)
const Vector< label > N(dict.get< Vector< label >>("N"))
List< bool > applied_
Applied flag list - corresponds to each fieldNames_ entry.
labelList indices(const keyType &key) const
Return zone indices for all matches.
void setSize(const label newSize)
Alias for resize(const label)
bool readIfPresent(const word &keyword, T &val, enum keyType::option matchOpt=keyType::REGEX) const