ReactingMultiphaseParcelIO.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-2023 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 
27 #include "IOstreams.H"
28 
29 // * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
30 
31 template<class ParcelType>
34 
35 template<class ParcelType>
37 (
38  sizeof(scalar)
39 );
40 
41 
42 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
43 
44 template<class ParcelType>
46 (
47  Istream& is,
48  bool readFields
49 )
50 :
51  ParcelType(is, readFields),
52  mass0_(0.0),
53  YGas_(0),
54  YLiquid_(0),
55  YSolid_(0),
56  canCombust_(0)
57 {
58  if (readFields)
59  {
63 
64  if (is.format() == IOstream::ASCII)
65  {
66  is >> mass0_ >> Yg >> Yl >> Ys;
67  }
68  else
69  {
70  is.read(reinterpret_cast<char*>(&mass0_), sizeofFields_);
71  is >> Yg >> Yl >> Ys;
72  }
73 
74  YGas_.transfer(Yg);
75  YLiquid_.transfer(Yl);
76  YSolid_.transfer(Ys);
77 
78  // scale the mass fractions
79  const scalarField& YMix = this->Y_;
80  YGas_ /= YMix[GAS] + rootVSmall;
81  YLiquid_ /= YMix[LIQ] + rootVSmall;
82  YSolid_ /= YMix[SLD] + rootVSmall;
83  }
84 
85  // Check state of Istream
86  is.check
87  (
88  "ReactingMultiphaseParcel<ParcelType>::ReactingMultiphaseParcel"
89  "("
90  "const polyMesh&, "
91  "Istream&, "
92  "bool"
93  ")"
94  );
95 }
96 
97 
98 template<class ParcelType>
99 template<class CloudType>
101 {
103 }
104 
105 
106 template<class ParcelType>
107 template<class CloudType, class CompositionType>
109 (
110  CloudType& c,
111  const CompositionType& compModel
112 )
113 {
114  bool valid = c.size();
115 
116  ParcelType::readFields(c, compModel);
117 
118  IOField<scalar> mass0
119  (
120  c.fieldIOobject("mass0", IOobject::MUST_READ),
121  valid
122  );
123  c.checkFieldIOobject(c, mass0);
124 
125  label i = 0;
126  forAllIter(typename CloudType, c, iter)
127  {
129  p.mass0_ = mass0[i++];
130  }
131 
132  // Get names and sizes for each Y...
133  const label idGas = compModel.idGas();
134  const wordList& gasNames = compModel.componentNames(idGas);
135  const label idLiquid = compModel.idLiquid();
136  const wordList& liquidNames = compModel.componentNames(idLiquid);
137  const label idSolid = compModel.idSolid();
138  const wordList& solidNames = compModel.componentNames(idSolid);
139  const wordList& stateLabels = compModel.stateLabels();
140 
141  // Set storage for each Y... for each parcel
142  forAllIter(typename CloudType, c, iter)
143  {
145  p.YGas_.setSize(gasNames.size(), 0.0);
146  p.YLiquid_.setSize(liquidNames.size(), 0.0);
147  p.YSolid_.setSize(solidNames.size(), 0.0);
148  }
149 
150  // Populate YGas for each parcel
151  forAll(gasNames, j)
152  {
153  IOField<scalar> YGas
154  (
155  c.fieldIOobject
156  (
157  "Y" + gasNames[j] + stateLabels[idGas],
159  ),
160  valid
161  );
162 
163  label i = 0;
164  forAllIter(typename CloudType, c, iter)
165  {
167  p.YGas_[j] = YGas[i++]/(p.Y()[GAS] + rootVSmall);
168  }
169  }
170  // Populate YLiquid for each parcel
171  forAll(liquidNames, j)
172  {
173  IOField<scalar> YLiquid
174  (
175  c.fieldIOobject
176  (
177  "Y" + liquidNames[j] + stateLabels[idLiquid],
179  ),
180  valid
181  );
182 
183  label i = 0;
184  forAllIter(typename CloudType, c, iter)
185  {
187  p.YLiquid_[j] = YLiquid[i++]/(p.Y()[LIQ] + rootVSmall);
188  }
189  }
190  // Populate YSolid for each parcel
191  forAll(solidNames, j)
192  {
193  IOField<scalar> YSolid
194  (
195  c.fieldIOobject
196  (
197  "Y" + solidNames[j] + stateLabels[idSolid],
199  ),
200  valid
201  );
202 
203  label i = 0;
204  forAllIter(typename CloudType, c, iter)
205  {
207  p.YSolid_[j] = YSolid[i++]/(p.Y()[SLD] + rootVSmall);
208  }
209  }
210 }
211 
212 
213 template<class ParcelType>
214 template<class CloudType>
216 {
217  ParcelType::writeFields(c);
218 }
219 
220 
221 template<class ParcelType>
222 template<class CloudType, class CompositionType>
224 (
225  const CloudType& c,
226  const CompositionType& compModel
227 )
228 {
229  ParcelType::writeFields(c, compModel);
230 
231  label np = c.size();
232 
233  {
234  IOField<scalar> mass0(c.fieldIOobject("mass0", IOobject::NO_READ), np);
235 
236  label i = 0;
237  forAllConstIter(typename CloudType, c, iter)
238  {
239  const ReactingMultiphaseParcel<ParcelType>& p = iter();
240  mass0[i++] = p.mass0_;
241  }
242  mass0.write(np > 0);
243 
244  // Write the composition fractions
245  const wordList& stateLabels = compModel.stateLabels();
246 
247  const label idGas = compModel.idGas();
248  const wordList& gasNames = compModel.componentNames(idGas);
249  forAll(gasNames, j)
250  {
251  IOField<scalar> YGas
252  (
253  c.fieldIOobject
254  (
255  "Y" + gasNames[j] + stateLabels[idGas],
257  ),
258  np
259  );
260 
261  label i = 0;
262  forAllConstIter(typename CloudType, c, iter)
263  {
264  const ReactingMultiphaseParcel<ParcelType>& p0 = iter();
265  YGas[i++] = p0.YGas()[j]*p0.Y()[GAS];
266  }
267 
268  YGas.write(np > 0);
269  }
270 
271  const label idLiquid = compModel.idLiquid();
272  const wordList& liquidNames = compModel.componentNames(idLiquid);
273  forAll(liquidNames, j)
274  {
275  IOField<scalar> YLiquid
276  (
277  c.fieldIOobject
278  (
279  "Y" + liquidNames[j] + stateLabels[idLiquid],
281  ),
282  np
283  );
284 
285  label i = 0;
286  forAllConstIter(typename CloudType, c, iter)
287  {
288  const ReactingMultiphaseParcel<ParcelType>& p0 = iter();
289  YLiquid[i++] = p0.YLiquid()[j]*p0.Y()[LIQ];
290  }
291 
292  YLiquid.write(np > 0);
293  }
294 
295  const label idSolid = compModel.idSolid();
296  const wordList& solidNames = compModel.componentNames(idSolid);
297  forAll(solidNames, j)
298  {
299  IOField<scalar> YSolid
300  (
301  c.fieldIOobject
302  (
303  "Y" + solidNames[j] + stateLabels[idSolid],
305  ),
306  np
307  );
308 
309  label i = 0;
310  forAllConstIter(typename CloudType, c, iter)
311  {
312  const ReactingMultiphaseParcel<ParcelType>& p0 = iter();
313  YSolid[i++] = p0.YSolid()[j]*p0.Y()[SLD];
314  }
315 
316  YSolid.write(np > 0);
317  }
318  }
319 }
320 
321 
322 // * * * * * * * * * * * * * * * IOstream Operators * * * * * * * * * * * * //
323 
324 template<class ParcelType>
325 Foam::Ostream& Foam::operator<<
326 (
327  Ostream& os,
329 )
330 {
331  scalarField YGasLoc(p.YGas()*p.Y()[0]);
332  scalarField YLiquidLoc(p.YLiquid()*p.Y()[1]);
333  scalarField YSolidLoc(p.YSolid()*p.Y()[2]);
334  if (os.format() == IOstream::ASCII)
335  {
336  os << static_cast<const ParcelType&>(p)
337  << token::SPACE << p.mass0()
338  << token::SPACE << YGasLoc
339  << token::SPACE << YLiquidLoc
340  << token::SPACE << YSolidLoc;
341  }
342  else
343  {
344  os << static_cast<const ParcelType&>(p);
345  os.write
346  (
347  reinterpret_cast<const char*>(&p.mass0_),
348  ReactingMultiphaseParcel<ParcelType>::sizeofFields_
349  );
350  os << YGasLoc << YLiquidLoc << YSolidLoc;
351  }
352 
353  // Check state of Ostream
354  os.check
355  (
356  "Ostream& operator<<"
357  "("
358  "Ostream&, "
359  "const ReactingMultiphaseParcel<ParcelType>&"
360  ")"
361  );
362 
363  return os;
364 }
365 
366 
367 // ************************************************************************* //
Useful combination of include files which define Sin, Sout and Serr and the use of IO streams general...
#define forAll(list, i)
Loop across all elements in list.
Definition: UList.H:434
#define forAllIter(Container, container, iter)
Iterate across all elements in the container object of type.
Definition: UList.H:459
#define forAllConstIter(Container, container, iter)
Iterate across all elements in the container object of type.
Definition: UList.H:477
Templated base class for dsmc cloud.
Definition: DSMCCloud.H:79
A primitive field of type <Type> with automated input and output.
Definition: IOField.H:53
streamFormat format() const
Return current stream format.
Definition: IOstream.H:374
virtual bool check(const char *operation) const
Check IOstream status for given operation.
Definition: IOstream.C:92
An Istream is an abstract base class for all input systems (streams, files, token lists etc)....
Definition: Istream.H:60
virtual Istream & read(token &)=0
Return next token from stream.
void transfer(List< T > &)
Transfer the contents of the argument List into this list.
Definition: List.C:342
void size(const label)
Override size to be inconsistent with allocated storage.
Definition: ListI.H:164
void setSize(const label)
Reset size of List.
Definition: List.C:281
An Ostream is an abstract base class for all output systems (streams, files, token lists,...
Definition: Ostream.H:57
Multiphase variant of the reacting parcel class with one/two-way coupling with the continuous phase.
scalarField YLiquid_
Mass fractions of liquids [].
scalarField YSolid_
Mass fractions of solids [].
static void writeFields(const CloudType &c, const CompositionType &compModel)
Write.
scalarField YGas_
Mass fractions of gases [].
const scalarField & YGas() const
Return const access to mass fractions of gases.
ReactingMultiphaseParcel(const polyMesh &mesh, const barycentric &coordinates, const label celli, const label tetFacei, const label tetPti, const label facei)
Construct from mesh, coordinates and topology.
const scalarField & YSolid() const
Return const access to mass fractions of solids.
const scalarField & YLiquid() const
Return const access to mass fractions of liquids.
static void readFields(CloudType &c, const CompositionType &compModel)
Read.
virtual bool write(const bool write=true) const
Write using setting from DB.
A class for handling character strings derived from std::string.
Definition: string.H:79
bool valid(const PtrList< ModelType > &l)
const dimensionedScalar c
Speed of light in a vacuum.
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
void readFields(const typename GeoFieldType::Mesh &mesh, const IOobjectList &objects, const HashSet< word > &selectedFields, LIFOStack< regIOobject * > &storedObjects)
Read the selected GeometricFields of the specified type.
Definition: ReadFields.C:244
volScalarField & p