multiphaseMangrovesSource.C
Go to the documentation of this file.
1/*---------------------------------------------------------------------------*\
2 ========= |
3 \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
4 \\ / O peration |
5 \\ / A nd | www.openfoam.com
6 \\/ M anipulation |
7-------------------------------------------------------------------------------
8 Copyright (C) 2017 IH-Cantabria
9 Copyright (C) 2017-2021 OpenCFD Ltd.
10-------------------------------------------------------------------------------
11License
12 This file is part of OpenFOAM.
13
14 OpenFOAM is free software: you can redistribute it and/or modify it
15 under the terms of the GNU General Public License as published by
16 the Free Software Foundation, either version 3 of the License, or
17 (at your option) any later version.
18
19 OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
20 ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
21 FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
22 for more details.
23
24 You should have received a copy of the GNU General Public License
25 along with OpenFOAM. If not, see <http://www.gnu.org/licenses/>.
26
27\*---------------------------------------------------------------------------*/
28
31#include "fvMesh.H"
32#include "fvMatrices.H"
33#include "fvmDdt.H"
34#include "fvmSup.H"
36
37// * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
38
39namespace Foam
40{
41namespace fv
42{
45 (
46 option,
49 );
50}
51}
52
53
54// * * * * * * * * * * * * Protected Member Functions * * * * * * * * * * * //
55
58{
59 auto tdragCoeff = tmp<volScalarField>::New
60 (
62 (
63 typeName + ":dragCoeff",
65 mesh_.time(),
68 ),
69 mesh_,
71 );
72 auto& dragCoeff = tdragCoeff.ref();
73
75 {
76 const scalar a = aZone_[i];
77 const scalar N = NZone_[i];
78 const scalar Cd = CdZone_[i];
79
80 const labelList& zones = zoneIDs_[i];
81
82 for (label zonei : zones)
83 {
84 const cellZone& cz = mesh_.cellZones()[zonei];
85
86 for (label celli : cz)
87 {
88 const vector& Uc = U[celli];
89
90 dragCoeff[celli] = 0.5*Cd*a*N*mag(Uc);
91 }
92 }
93 }
94
95 dragCoeff.correctBoundaryConditions();
96
97 return tdragCoeff;
98}
99
100
103{
104 auto tinertiaCoeff = tmp<volScalarField>::New
105 (
107 (
108 typeName + ":inertiaCoeff",
109 mesh_.time().timeName(),
110 mesh_.time(),
113 ),
114 mesh_,
116 );
117 auto& inertiaCoeff = tinertiaCoeff.ref();
118
119 const scalar pi = constant::mathematical::pi;
120
121 forAll(zoneIDs_, i)
122 {
123 const scalar a = aZone_[i];
124 const scalar N = NZone_[i];
125 const scalar Cm = CmZone_[i];
126
127 const labelList& zones = zoneIDs_[i];
128
129 for (label zonei : zones)
130 {
131 const cellZone& cz = mesh_.cellZones()[zonei];
132
133 for (label celli : cz)
134 {
135 inertiaCoeff[celli] = 0.25*(Cm+1)*pi*a*a*N;
136 }
137 }
138 }
139
140 inertiaCoeff.correctBoundaryConditions();
141
142 return tinertiaCoeff;
143}
144
145
146// * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
147
149(
150 const word& name,
151 const word& modelType,
152 const dictionary& dict,
153 const fvMesh& mesh
154)
155:
156 fv::option(name, modelType, dict, mesh),
157 aZone_(),
158 NZone_(),
159 CmZone_(),
160 CdZone_()
161{
162 read(dict);
163}
164
165
166// * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
167
169(
170 fvMatrix<vector>& eqn,
171 const label fieldi
172)
173{
174 const volVectorField& U = eqn.psi();
175
176 fvMatrix<vector> mangrovesEqn
177 (
178 - fvm::Sp(dragCoeff(U), U)
179 - inertiaCoeff()*fvm::ddt(U)
180 );
181
182 // Contributions are added to RHS of momentum equation
183 eqn += mangrovesEqn;
184}
185
186
188(
189 const volScalarField& rho,
190 fvMatrix<vector>& eqn,
191 const label fieldi
192)
193{
194 const volVectorField& U = eqn.psi();
195
196 fvMatrix<vector> mangrovesEqn
197 (
198 - fvm::Sp(rho*dragCoeff(U), U)
199 - rho*inertiaCoeff()*fvm::ddt(U)
200 );
201
202 // Contributions are added to RHS of momentum equation
203 eqn += mangrovesEqn;
204}
205
206
208{
210 {
211 if (!coeffs_.readIfPresent("UNames", fieldNames_))
212 {
213 fieldNames_.resize(1);
214 fieldNames_.first() = coeffs_.getOrDefault<word>("U", "U");
215 }
217
218 // Create the Mangroves models - 1 per region
219 const dictionary& regionsDict(coeffs_.subDict("regions"));
220 const wordList regionNames(regionsDict.toc());
221 aZone_.setSize(regionNames.size(), 1);
222 NZone_.setSize(regionNames.size(), 1);
223 CmZone_.setSize(regionNames.size(), 1);
224 CdZone_.setSize(regionNames.size(), 1);
225 zoneIDs_.setSize(regionNames.size());
226
227 forAll(zoneIDs_, i)
228 {
229 const word& regionName = regionNames[i];
230 const dictionary& modelDict = regionsDict.subDict(regionName);
231
232 const word zoneName(modelDict.get<word>("cellZone"));
233
234 zoneIDs_[i] = mesh_.cellZones().indices(zoneName);
235 if (zoneIDs_[i].empty())
236 {
238 << "Unable to find cellZone " << zoneName << nl
239 << "Valid cellZones are:" << mesh_.cellZones().names()
240 << exit(FatalError);
241 }
242
243 modelDict.readEntry("a", aZone_[i]);
244 modelDict.readEntry("N", NZone_[i]);
245 modelDict.readEntry("Cm", CmZone_[i]);
246 modelDict.readEntry("Cd", CdZone_[i]);
247 }
248
249 return true;
250 }
251
252 return false;
253}
254
255
256// ************************************************************************* //
Macros for easy insertion into run-time selection tables.
#define addToRunTimeSelectionTable(baseType, thisType, argNames)
Add to construction table with typeName as the key.
Defines the attributes of an object for which implicit objectRegistry management is supported,...
Definition: IOobject.H:170
virtual bool read()
Re-read model coefficients if they have changed.
static autoPtr< Time > New()
Construct (dummy) Time - no functionObjects or libraries.
Definition: Time.C:717
static word timeName(const scalar t, const int precision=precision_)
Definition: Time.C:780
A subset of mesh cells.
Definition: cellZone.H:65
A list of keyword definitions, which are a keyword followed by a number of values (eg,...
Definition: dictionary.H:126
T get(const word &keyword, enum keyType::option matchOpt=keyType::REGEX) const
const dictionary & subDict(const word &keyword, enum keyType::option matchOpt=keyType::REGEX) const
Find and return a sub-dictionary.
Definition: dictionary.C:460
wordList toc() const
Return the table of contents.
Definition: dictionary.C:602
bool readEntry(const word &keyword, T &val, enum keyType::option matchOpt=keyType::REGEX, bool mandatory=true) const
A special matrix type and solver, designed for finite volume solutions of scalar equations....
Definition: fvMatrix.H:121
const GeometricField< Type, fvPatchField, volMesh > & psi(const label i=0) const
Return psi.
Definition: fvMatrix.H:412
Mesh data needed to do the Finite Volume discretisation.
Definition: fvMesh.H:91
const Time & time() const
Return the top-level database.
Definition: fvMesh.H:290
labelListList zoneIDs_
Porosity coefficient.
virtual bool read(const dictionary &dict)
Read dictionary.
tmp< volScalarField > dragCoeff(const volVectorField &U) const
Return the drag force coefficient.
scalarList aZone_
Width of the vegetation element.
tmp< volScalarField > inertiaCoeff() const
Return the inertia force coefficient.
scalarList NZone_
Number of plants per unit of area.
virtual void addSup(fvMatrix< vector > &eqn, const label fieldi)
Add implicit contribution to momentum equation.
Base abstract class for handling finite volume options (i.e. fvOption).
Definition: fvOption.H:127
const fvMesh & mesh_
Reference to the mesh database.
Definition: fvOption.H:139
void resetApplied()
Resize/reset applied flag list for all fieldNames_ entries.
Definition: fvOption.C:48
const cellZoneMesh & cellZones() const noexcept
Return cell zone mesh.
Definition: polyMesh.H:504
A class for managing temporary objects.
Definition: tmp.H:65
A class for handling words, derived from Foam::string.
Definition: word.H:68
#define defineTypeNameAndDebug(Type, DebugSwitch)
Define the typeName and debug information.
Definition: className.H:121
U
Definition: pEqn.H:72
dynamicFvMesh & mesh
Foam::word regionName(Foam::polyMesh::defaultRegion)
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
Definition: error.H:453
A special matrix type and solver, designed for finite volume solutions of scalar equations.
Calculate the matrix for the first temporal derivative.
Calculate the finiteVolume matrix for implicit and explicit sources.
wordList regionNames
constexpr scalar pi(M_PI)
tmp< fvMatrix< Type > > ddt(const GeometricField< Type, fvPatchField, volMesh > &vf)
Definition: fvmDdt.C:48
zeroField Sp(const Foam::zero, const GeometricField< Type, fvPatchField, volMesh > &)
A no-op source.
Namespace for OpenFOAM.
const dimensionSet dimless
Dimensionless.
const dimensionSet dimTime(0, 0, 1, 0, 0, 0, 0)
Definition: dimensionSets.H:53
dimensioned< typename typeOfMag< Type >::type > mag(const dimensioned< Type > &dt)
static constexpr const zero Zero
Global zero (0)
Definition: zero.H:131
error FatalError
word name(const expressions::valueTypeCode typeCode)
A word representation of a valueTypeCode. Empty for INVALID.
Definition: exprTraits.C:59
errorManipArg< error, int > exit(error &err, const int errNo=1)
Definition: errorManip.H:130
constexpr char nl
The newline '\n' character (0x0a)
Definition: Ostream.H:53
labelList fv(nPoints)
dictionary dict
#define forAll(list, i)
Loop across all elements in list.
Definition: stdFoam.H:333
const Vector< label > N(dict.get< Vector< label > >("N"))