Go to the documentation of this file.
47 multiphaseMangrovesSource,
63 typeName +
":dragCoeff",
76 const scalar a =
aZone_[i];
82 for (label zonei : zones)
86 for (label celli : cz)
108 typeName +
":inertiaCoeff",
109 mesh_.time().timeName(),
117 auto& inertiaCoeff = tinertiaCoeff.ref();
123 const scalar a = aZone_[i];
124 const scalar
N = NZone_[i];
125 const scalar Cm = CmZone_[i];
129 for (label zonei : zones)
131 const cellZone& cz = mesh_.cellZones()[zonei];
133 for (label celli : cz)
135 inertiaCoeff[celli] = 0.25*(Cm+1)*
pi*a*a*
N;
140 inertiaCoeff.correctBoundaryConditions();
142 return tinertiaCoeff;
148 Foam::fv::multiphaseMangrovesSource::multiphaseMangrovesSource
151 const word& modelType,
211 if (!coeffs_.readIfPresent(
"UNames", fieldNames_))
213 fieldNames_.resize(1);
214 fieldNames_.first() = coeffs_.getOrDefault<
word>(
"U",
"U");
217 applied_.setSize(fieldNames_.size(),
false);
220 const dictionary& regionsDict(coeffs_.subDict(
"regions"));
221 const wordList regionNames(regionsDict.toc());
222 aZone_.
setSize(regionNames.size(), 1);
223 NZone_.setSize(regionNames.size(), 1);
224 CmZone_.setSize(regionNames.size(), 1);
225 CdZone_.setSize(regionNames.size(), 1);
226 zoneIDs_.setSize(regionNames.size());
233 const word zoneName(modelDict.
get<
word>(
"cellZone"));
235 zoneIDs_[i] = mesh_.cellZones().indices(zoneName);
236 if (zoneIDs_[i].empty())
239 <<
"Unable to find cellZone " << zoneName <<
nl
240 <<
"Valid cellZones are:" << mesh_.cellZones().names()
Defines the attributes of an object for which implicit objectRegistry management is supported,...
const dimensionSet dimless(0, 0, 0, 0, 0, 0, 0)
Dimensionless.
A class for handling words, derived from Foam::string.
A class for managing temporary objects.
static constexpr const zero Zero
Global zero (0)
bool read(const char *buf, int32_t &val)
Same as readInt32.
static word timeName(const scalar t, const int precision=precision_)
const cellZoneMesh & cellZones() const
Return cell zone mesh.
T get(const word &keyword, enum keyType::option matchOpt=keyType::REGEX) const
scalarList CdZone_
Drag coefficient.
virtual bool read(const dictionary &dict)
Read dictionary.
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)
tmp< fvMatrix< Type > > Sp(const volScalarField::Internal &, const GeometricField< Type, fvPatchField, volMesh > &)
Finite volume options abstract base class. Provides a base set of controls, e.g.:
word name(const complex &c)
Return string representation of complex.
scalarList NZone_
Number of plants per unit of area.
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.
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.
Calculate the matrix for implicit and explicit sources.
errorManipArg< error, int > exit(error &err, const int errNo=1)
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
tmp< fvMatrix< Type > > ddt(const GeometricField< Type, fvPatchField, volMesh > &vf)
constexpr scalar pi(M_PI)
dimensioned< typename typeOfMag< Type >::type > mag(const dimensioned< Type > &dt)
tmp< volScalarField > inertiaCoeff() const
Return the inertia force coefficient.
static tmp< T > New(Args &&... args)
Construct tmp of T with forwarding arguments.
A special matrix type and solver, designed for finite volume solutions of scalar equations....
tmp< volScalarField > dragCoeff(const volVectorField &U) const
Return the drag force coefficient.
const Time & time() const
Return the top-level database.
defineTypeNameAndDebug(atmAmbientTurbSource, 0)
virtual void addSup(fvMatrix< vector > &eqn, const label fieldi)
Add implicit contribution to momentum equation.
addToRunTimeSelectionTable(option, atmAmbientTurbSource, dictionary)
const Vector< label > N(dict.get< Vector< label >>("N"))
labelListList zoneIDs_
Porosity coefficient.
void setSize(const label newSize)
Alias for resize(const label)
scalarList aZone_
Width of the vegetation element.
Calculate the matrix for the first temporal derivative.