ejectionModelList.C
Go to the documentation of this file.
1 /*---------------------------------------------------------------------------*\
2  ========= |
3  \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
4  \\ / O peration | Website: https://openfoam.org
5  \\ / A nd | Copyright (C) 2011-2021 OpenFOAM Foundation
6  \\/ M anipulation |
7 -------------------------------------------------------------------------------
8 License
9  This file is part of OpenFOAM.
10 
11  OpenFOAM is free software: you can redistribute it and/or modify it
12  under the terms of the GNU General Public License as published by
13  the Free Software Foundation, either version 3 of the License, or
14  (at your option) any later version.
15 
16  OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
17  ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
18  FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
19  for more details.
20 
21  You should have received a copy of the GNU General Public License
22  along with OpenFOAM. If not, see <http://www.gnu.org/licenses/>.
23 
24 \*---------------------------------------------------------------------------*/
25 
26 #include "ejectionModelList.H"
27 
28 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
29 
30 namespace Foam
31 {
32 namespace regionModels
33 {
34 namespace surfaceFilmModels
35 {
36 
37 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
38 
40 :
42  filmSubModelBase(film)
43 {}
44 
45 
47 (
49  const dictionary& dict
50 )
51 :
54  (
55  "ejectionModelList",
56  film,
57  dict,
58  "ejectionModelList",
59  "ejectionModelList"
60  ),
61  massEjected_(film.intCoupledPatchIDs().size(), 0.0)
62 {
63  Info<< " Selecting film ejection" << endl;
64 
65  if (dict.isDict("ejection"))
66  {
67  const dictionary& ejectionDict(dict.subDict("ejection"));
68  this->setSize(ejectionDict.size());
69 
70  label i = 0;
71  forAllConstIter(dictionary, ejectionDict, iter)
72  {
73  set
74  (
75  i++,
77  (
78  film,
79  ejectionDict.isDict(iter().keyword())
80  ? ejectionDict.subDict(iter().keyword())
82  iter().keyword()
83  )
84  );
85  }
86  }
87  else if (dict.found("ejectionModels"))
88  {
89  const wordList models(dict.lookup("ejectionModels"));
90  this->setSize(models.size());
91 
92  forAll(models, i)
93  {
94  set(i, ejectionModel::New(film, dict, models[i]));
95  }
96  }
97 
98  if (!size())
99  {
100  Info<< " none" << endl;
101  }
102 }
103 
104 
105 // * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
106 
108 {}
109 
110 
111 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
112 
114 (
115  scalarField& availableMass,
116  volScalarField& massToEject,
117  volScalarField& diameterToEject
118 )
119 {
120  // Correct models that accumulate mass and diameter transfers
121  forAll(*this, i)
122  {
123  ejectionModel& im = operator[](i);
124  im.correct(availableMass, massToEject, diameterToEject);
125  }
126 
127  // Push values to boundaries ready for transfer to the primary region
128  massToEject.correctBoundaryConditions();
129  diameterToEject.correctBoundaryConditions();
130 
131  const labelList& patchIDs = film().intCoupledPatchIDs();
132 
133  forAll(patchIDs, i)
134  {
135  label patchi = patchIDs[i];
136  massEjected_[i] =
137  massEjected_[i] + sum(massToEject.boundaryField()[patchi]);
138  }
139 }
140 
141 
143 {
144  const polyBoundaryMesh& pbm = film().regionMesh().boundaryMesh();
145 
146  scalar ejectedMass = 0;
147  scalarField patchEjectedMasses
148  (
150  0
151  );
152 
153  forAll(*this, i)
154  {
155  const ejectionModel& im = operator[](i);
156  ejectedMass += im.ejectedMassTotal();
157  im.patchEjectedMassTotals(patchEjectedMasses);
158  }
159 
160  os << indent << "ejected mass = " << ejectedMass << nl;
161 
162  forAll(patchEjectedMasses, patchi)
163  {
164  if (mag(patchEjectedMasses[patchi]) > vSmall)
165  {
166  os << indent << indent << "from patch " << pbm[patchi].name()
167  << " = " << patchEjectedMasses[patchi] << nl;
168  }
169  }
170 
171  scalarField mass0(massEjected_.size(), 0);
172  this->getBaseProperty("massEjected", mass0);
173 
174  scalarField mass(massEjected_);
176  mass += mass0;
177 
178  const labelList& patchIDs = film().intCoupledPatchIDs();
179 
180  forAll(patchIDs, i)
181  {
182  label patchi = patchIDs[i];
183  Info<< indent << " - patch: " << pbm[patchi].name() << ": "
184  << mass[i] << endl;
185  }
186 
187  if (film().time().writeTime())
188  {
189  setBaseProperty("massEjected", mass);
190  massEjected_ = 0.0;
191  }
192 }
193 
194 
195 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
196 
197 } // End namespace surfaceFilmModels
198 } // End namespace regionModels
199 } // End namespace Foam
200 
201 // ************************************************************************* //
const polyBoundaryMesh & boundaryMesh() const
Return boundary mesh.
Definition: polyMesh.H:434
bool found(const word &, bool recursive=false, bool patternMatch=true) const
Search dictionary for given keyword.
Definition: dictionary.C:643
#define forAll(list, i)
Loop across all elements in list.
Definition: UList.H:434
const labelList & processorPatches() const
Return list of processor patch labels.
static autoPtr< ejectionModel > New(surfaceFilmRegionModel &film, const dictionary &dict, const word &mdoelType)
Return a reference to the selected ejection model.
intWM_LABEL_SIZE_t label
A label is an int32_t or int64_t as specified by the pre-processor macro WM_LABEL_SIZE.
Definition: label.H:59
const word & name() const
Return name.
Definition: IOobject.H:303
Ostream & indent(Ostream &os)
Indent stream.
Definition: Ostream.H:221
Base class for surface film sub-models.
A list of keyword definitions, which are a keyword followed by any number of values (e...
Definition: dictionary.H:156
const Boundary & boundaryField() const
Return const-reference to the boundary field.
Base class for film ejection models, handling mass transfer from the film.
Definition: ejectionModel.H:56
void size(const label)
Override size to be inconsistent with allocated storage.
Definition: ListI.H:164
const T & operator[](const label) const
Return element const reference.
Definition: UPtrListI.H:96
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:251
static void listCombineGather(const List< commsStruct > &comms, List< T > &Value, const CombineOp &cop, const int tag, const label comm)
const dictionary & dict() const
Return const access to the cloud dictionary.
Definition: subModelBase.C:110
bool isDict(const word &) const
Check if entry is a sub-dictionary.
Definition: dictionary.C:936
dimensioned< Type > sum(const DimensionedField< Type, GeoMesh > &df)
const surfaceFilmRegionModel & film() const
Return const access to the film surface film model.
virtual bool writeTime() const
Flag to indicate when to write a property.
const dictionary & subDict(const word &) const
Find and return a sub-dictionary.
Definition: dictionary.C:982
virtual void info(Ostream &os)
Provide some info.
static const dictionary null
Null dictionary.
Definition: dictionary.H:242
virtual scalar ejectedMassTotal() const
Return the total mass ejected.
Definition: ejectionModel.C:91
const fvMesh & regionMesh() const
Return the region mesh database.
Definition: regionModelI.H:55
const globalMeshData & globalData() const
Return parallel info.
Definition: polyMesh.C:1394
forAllConstIter(PtrDictionary< phaseModel >, mixture.phases(), phase)
Definition: pEqn.H:29
void setSize(const label)
Reset size of PtrList. If extending the PtrList, new entries are.
Definition: PtrList.C:131
Foam::polyBoundaryMesh.
An Ostream is an abstract base class for all output systems (streams, files, token lists...
Definition: Ostream.H:54
static const char nl
Definition: Ostream.H:260
label size() const
Return the number of elements in the UPtrList.
Definition: UPtrListI.H:29
label patchi
A templated 1D list of pointers to objects of type <T>, where the size of the array is known and used...
Definition: List.H:70
void setBaseProperty(const word &entryName, const Type &value)
Add generic property to the base model.
virtual void correct(scalarField &availableMass, volScalarField &massToEject, volScalarField &diameterToEject)
Correct.
void correctBoundaryConditions()
Correct boundary field.
messageStream Info
dimensioned< scalar > mag(const dimensioned< Type > &)
Type getBaseProperty(const word &entryName, const Type &defaultValue=pTraits< Type >::zero) const
Retrieve generic property from the base model.
const labelList & intCoupledPatchIDs() const
Return the list of patch IDs internally coupled with the.
Definition: regionModelI.H:173
virtual void patchEjectedMassTotals(scalarField &patchMasses) const
Accumulate the total mass ejected for the patches into the.
ejectionModelList(surfaceFilmRegionModel &film)
Construct null.
Namespace for OpenFOAM.
ITstream & lookup(const word &, bool recursive=false, bool patternMatch=true) const
Find and return an entry data stream.
Definition: dictionary.C:844