injectionModelList.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 | Copyright (C) 2011-2017 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 "injectionModelList.H"
27 
28 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
29 
30 namespace Foam
31 {
32 namespace regionModels
33 {
34 namespace surfaceFilmModels
35 {
36 
37 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
38 
39 injectionModelList::injectionModelList(surfaceFilmModel& film)
40 :
42  filmSubModelBase(film)
43 {}
44 
45 
46 injectionModelList::injectionModelList
47 (
49  const dictionary& dict
50 )
51 :
54  (
55  "injectionModelList",
56  film,
57  dict,
58  "injectionModelList",
59  "injectionModelList"
60  ),
61  massInjected_(film.intCoupledPatchIDs().size(), 0.0)
62 {
63  const wordList activeModels(dict.lookup("injectionModels"));
64 
65  wordHashSet models;
66  forAll(activeModels, i)
67  {
68  models.insert(activeModels[i]);
69  }
70 
71  Info<< " Selecting film injection models" << endl;
72  if (models.size() > 0)
73  {
74  this->setSize(models.size());
75 
76  label i = 0;
77  forAllConstIter(wordHashSet, models, iter)
78  {
79  const word& model = iter.key();
80  set(i, injectionModel::New(film, dict, model));
81  i++;
82  }
83  }
84  else
85  {
86  Info<< " none" << endl;
87  }
88 }
89 
90 
91 // * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
92 
94 {}
95 
96 
97 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
98 
100 (
101  scalarField& availableMass,
102  volScalarField& massToInject,
103  volScalarField& diameterToInject
104 )
105 {
106  // Correct models that accumulate mass and diameter transfers
107  forAll(*this, i)
108  {
109  injectionModel& im = operator[](i);
110  im.correct(availableMass, massToInject, diameterToInject);
111  }
112 
113  // Push values to boundaries ready for transfer to the primary region
114  massToInject.correctBoundaryConditions();
115  diameterToInject.correctBoundaryConditions();
116 
117  const labelList& patchIDs = film().intCoupledPatchIDs();
118 
119  forAll(patchIDs, i)
120  {
121  label patchi = patchIDs[i];
122  massInjected_[i] =
123  massInjected_[i] + sum(massToInject.boundaryField()[patchi]);
124  }
125 }
126 
127 
129 {
130  const polyBoundaryMesh& pbm = film().regionMesh().boundaryMesh();
131 
132  scalar injectedMass = 0;
133  scalarField patchInjectedMasses
134  (
136  0
137  );
138 
139  forAll(*this, i)
140  {
141  const injectionModel& im = operator[](i);
142  injectedMass += im.injectedMassTotal();
143  im.patchInjectedMassTotals(patchInjectedMasses);
144  }
145 
146  os << indent << "injected mass = " << injectedMass << nl;
147 
148  forAll(patchInjectedMasses, patchi)
149  {
150  if (mag(patchInjectedMasses[patchi]) > VSMALL)
151  {
152  os << indent << indent << "from patch " << pbm[patchi].name()
153  << " = " << patchInjectedMasses[patchi] << nl;
154  }
155  }
156 
157  scalarField mass0(massInjected_.size(), 0);
158  this->getBaseProperty("massInjected", mass0);
159 
160  scalarField mass(massInjected_);
162  mass += mass0;
163 
164  const labelList& patchIDs = film().intCoupledPatchIDs();
165 
166  forAll(patchIDs, i)
167  {
168  label patchi = patchIDs[i];
169  Info<< indent << " - patch: " << pbm[patchi].name() << ": "
170  << mass[i] << endl;
171  }
172 
173  if (film().time().writeTime())
174  {
175  setBaseProperty("massInjected", mass);
176  massInjected_ = 0.0;
177  }
178 }
179 
180 
181 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
182 
183 } // End namespace surfaceFilmModels
184 } // End namespace regionModels
185 } // End namespace Foam
186 
187 // ************************************************************************* //
const polyBoundaryMesh & boundaryMesh() const
Return boundary mesh.
Definition: polyMesh.H:424
const surfaceFilmModel & film() const
Return const access to the film surface film model.
A HashTable with keys but without contents.
Definition: HashSet.H:59
#define forAll(list, i)
Loop across all elements in list.
Definition: UList.H:428
const labelList & processorPatches() const
Return list of processor patch labels.
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:291
virtual void correct(scalarField &availableMass, volScalarField &massToInject, volScalarField &diameterToInject)
Correct.
Ostream & indent(Ostream &os)
Indent stream.
Definition: Ostream.H:223
A list of keyword definitions, which are a keyword followed by any number of values (e...
Definition: dictionary.H:137
const Boundary & boundaryField() const
Return const-reference to the boundary field.
void size(const label)
Override size to be inconsistent with allocated storage.
Definition: ListI.H:163
static autoPtr< injectionModel > New(surfaceFilmModel &film, const dictionary &dict, const word &mdoelType)
Return a reference to the selected injection model.
const T & operator[](const label) const
Return element const reference.
Definition: UPtrListI.H:103
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:253
virtual void info(Ostream &os)
Provide some info.
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
dimensioned< Type > sum(const DimensionedField< Type, GeoMesh > &df)
virtual bool writeTime() const
Flag to indicate when to write a property.
virtual void patchInjectedMassTotals(scalarField &patchMasses) const
Accumulate the total mass injected for the patches into the.
A class for handling words, derived from string.
Definition: word.H:59
const fvMesh & regionMesh() const
Return the region mesh database.
Definition: regionModelI.H:61
const globalMeshData & globalData() const
Return parallel info.
Definition: polyMesh.C:1182
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:53
static const char nl
Definition: Ostream.H:262
Base class for film injection models, handling mass transfer from the film.
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:63
void setBaseProperty(const word &entryName, const Type &value)
Add generic property to the base model.
virtual scalar injectedMassTotal() const
Return the total mass injected.
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:179
Namespace for OpenFOAM.
ITstream & lookup(const word &, bool recursive=false, bool patternMatch=true) const
Find and return an entry data stream.
Definition: dictionary.C:576