ReactingHeterogeneousParcelIO.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) 2018-2019 OpenCFD Ltd.
9-------------------------------------------------------------------------------
10License
11 This file is part of OpenFOAM.
12
13 OpenFOAM is free software: you can redistribute it and/or modify it
14 under the terms of the GNU General Public License as published by
15 the Free Software Foundation, either version 3 of the License, or
16 (at your option) any later version.
17
18 OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
19 ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
20 FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
21 for more details.
22
23 You should have received a copy of the GNU General Public License
24 along with OpenFOAM. If not, see <http://www.gnu.org/licenses/>.
25
26\*---------------------------------------------------------------------------*/
27
29#include "IOstreams.H"
31
32// * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
33
34template<class ParcelType>
37
38
39template<class ParcelType>
41(
42 sizeof(scalar)
43);
44
45
46// * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
47
48template<class ParcelType>
50(
51 const polyMesh& mesh,
52 Istream& is,
53 bool readFields,
54 bool newFormat
55)
56:
57 ParcelType(mesh, is, readFields, newFormat),
58 F_(0),
59 canCombust_(1)
60{
61 if (readFields)
62 {
64
65 is >> F;
66
67 F_.transfer(F);
68 }
69
71}
72
73
74template<class ParcelType>
75template<class CloudType>
77{
79}
80
81
82template<class ParcelType>
83template<class CloudType, class CompositionType>
85(
86 CloudType& c,
87 const CompositionType& compModel
88)
89{
90 const bool valid = c.size();
91
93
94 IOField<scalar> mass0
95 (
96 c.fieldIOobject("mass0", IOobject::MUST_READ),
97 valid
98 );
99 c.checkFieldIOobject(c, mass0);
100
101 label i = 0;
103 {
104 p.mass0() = mass0[i];
105 ++i;
106 }
107
108 const label idSolid = compModel.idSolid();
109 const wordList& solidNames = compModel.componentNames(idSolid);
110
111 // WIP until find a way to get info from Reacting Heterogeneous model
112 label nF(1);
113
114 // Set storage for each Y... for each parcel
116 {
117 p.Y().setSize(solidNames.size(), Zero);
118 p.F_.setSize(nF, Zero);
119 }
120
121 for (label i = 0; i < nF; i++)
122 {
123 // Read F
125 (
126 c.fieldIOobject
127 (
128 "F" + name(i),
130 ),
131 valid
132 );
133
134 label j = 0;
136 {
137 p.F_[i] = F[j];
138 ++j;
139 }
140 }
141
142
144 {
146 (
147 c.fieldIOobject
148 (
149 "Y" + solidNames[j],
151 ),
152 valid
153 );
154
155 label i = 0;
157 {
158 p.Y()[j] = Y[i];
159 ++i;
160 }
161 }
162}
163
164
165template<class ParcelType>
166template<class CloudType>
168(
169 const CloudType& c
170)
171{
173}
174
175
176template<class ParcelType>
177template<class CloudType, class CompositionType>
179(
180 const CloudType& c,
181 const CompositionType& compModel
182)
183{
184 // Skip Reacting layer. This class takes charge of
185 // writing Ysolids and F
187
188 const label np = c.size();
189 const bool valid = np;
190
191 IOField<scalar> mass0(c.fieldIOobject("mass0", IOobject::NO_READ), np);
192
193 label nF = 0;
194 label i = 0;
196 {
197 mass0[i] = p.mass0_;
198 if (!i)
199 {
200 // Assume same size throughout
201 nF = p.F().size();
202 }
203 ++i;
204 }
205 mass0.write(valid);
206
207 for (label i = 0; i < nF; i++)
208 {
210 (
211 c.fieldIOobject
212 (
213 "F" + name(i),
215 ),
216 np
217 );
218
220 {
221 F = p0.F()[i];
222 }
223
224 F.write(valid);
225 }
226
227 const label idSolid = compModel.idSolid();
228 const wordList& solidNames = compModel.componentNames(idSolid);
229
231 {
233 (
234 c.fieldIOobject
235 (
236 "Y" + solidNames[j],
238 ),
239 np
240 );
241
242 label i = 0;
244 {
245 Y[i] = p0.Y()[j];
246 ++i;
247 }
248
249 Y.write(valid);
250 }
251}
252
253
254template<class ParcelType>
256(
257 Ostream& os,
258 const wordRes& filters,
259 const word& delim,
260 const bool namesOnly
261) const
262{
263 ParcelType::writeProperties(os, filters, delim, namesOnly);
264
265 #undef writeProp
266 #define writeProp(Name, Value) \
267 ParcelType::writeProperty(os, Name, Value, namesOnly, delim, filters)
268
269 writeProp("F", F_);
270 writeProp("canCombust", canCombust_);
271
272 #undef writeProp
273}
274
275
276template<class ParcelType>
277template<class CloudType>
279(
280 CloudType& c,
281 const objectRegistry& obr
282)
283{
285}
286
287
288template<class ParcelType>
289template<class CloudType>
291(
292 const CloudType& c,
293 objectRegistry& obr
294)
295{
297}
298
299
300template<class ParcelType>
301template<class CloudType, class CompositionType>
303(
304 CloudType& c,
305 const CompositionType& compModel,
306 const objectRegistry& obr
307)
308{
309 //ParcelType::readObjects(c, obr);
310 // Skip Reacting layer
311 ThermoParcel<KinematicParcel<particle>>::readObjects(c, obr);
312
313 // const label np = c.size();
314
316 << "Reading of objects is still a work-in-progress" << nl;
317
318}
319
320
321template<class ParcelType>
322template<class CloudType, class CompositionType>
324(
325 const CloudType& c,
326 const CompositionType& compModel,
327 objectRegistry& obr
328)
329{
330 //ParcelType::writeObjects(c, obr);
331 // Skip Reacting layer
332 ThermoParcel<KinematicParcel<particle>>::writeObjects(c, obr);
333
334 const label np = c.size();
335
336 // WIP
337 label nF = 0;
339 {
340 nF = p0.F().size();
341 break;
342 }
343
344 if (np > 0)
345 {
346 for (label i = 0; i < nF; i++)
347 {
348 const word fieldName = "F" + name(i);
349 auto& F = cloud::createIOField<scalar>(fieldName, np, obr);
350
351 label j = 0;
353 {
354 F[j] = p0.F()[i];
355 ++j;
356 }
357 }
358
359 const label idSolid = compModel.idSolid();
360 const wordList& solidNames = compModel.componentNames(idSolid);
362 {
363 const word fieldName = "Y" + solidNames[j];
364 auto& Y = cloud::createIOField<scalar>(fieldName, np, obr);
365
366 label i = 0;
368 {
369 Y[i] = p0.Y()[j];
370 ++i;
371 }
372 }
373 }
374}
375
376
377// * * * * * * * * * * * * * * * IOstream Operators * * * * * * * * * * * * //
378
379template<class ParcelType>
380Foam::Ostream& Foam::operator<<
381(
382 Ostream& os,
384)
385{
386 scalarField F(p.F());
387 if (os.format() == IOstream::ASCII)
388 {
389 os << static_cast<const ParcelType&>(p)
390 << token::SPACE << F;
391 }
392 else
393 {
394 os << static_cast<const ParcelType&>(p);
395 os << F ;
396 }
397
399 return os;
400}
401
402
403// ************************************************************************* //
Useful combination of include files which define Sin, Sout and Serr and the use of IO streams general...
Templated base class for dsmc cloud.
Definition: DSMCCloud.H:75
A 1D vector of objects of type <T> that resizes itself as necessary to accept the new objects.
Definition: DynamicList.H:72
A primitive field of type <T> with automated input and output.
Definition: IOField.H:58
streamFormat format() const noexcept
Get the current stream format.
virtual bool check(const char *operation) const
Check IOstream status for given operation.
Definition: IOstream.C:58
An Istream is an abstract base class for all input systems (streams, files, token lists etc)....
Definition: Istream.H:64
void transfer(List< T > &list)
Definition: List.C:447
void setSize(const label n)
Alias for resize()
Definition: List.H:218
An Ostream is an abstract base class for all output systems (streams, files, token lists,...
Definition: Ostream.H:62
static const std::size_t sizeofFields
Size in bytes of the fields.
scalarField F_
Progress variables for reactions.
static void readObjects(CloudType &c, const objectRegistry &obr)
Read particle fields as objects from the obr registry.
const scalarField & F() const
Return const access to F.
void writeProperties(Ostream &os, const wordRes &filters, const word &delim, const bool namesOnly) const
Write individual parcel properties to stream.
static void readFields(CloudType &c, const CompositionType &compModel)
Read - composition supplied.
Thermodynamic parcel class with one/two-way coupling with the continuous phase. Includes Kinematic pa...
Definition: ThermoParcel.H:74
void size(const label n)
Older name for setAddressableSize.
Definition: UList.H:114
Reads fields from the time directories and adds them to the mesh database for further post-processing...
Definition: readFields.H:158
Allows specification of different writing frequency of objects registered to the database.
Definition: writeObjects.H:142
static void readObjects(Cloud< injectedParticle > &c, const objectRegistry &obr)
Read particle fields as objects from the obr registry.
void writeProperties(Ostream &os, const wordRes &filters, const word &delim, const bool namesOnly) const
Write individual parcel properties to stream.
Registry of regIOobjects.
static string propertyList()
Definition: particle.H:372
Mesh consisting of general polyhedral cells.
Definition: polyMesh.H:81
virtual bool write(const bool valid=true) const
Write using setting from DB.
A class for handling character strings derived from std::string.
Definition: string.H:79
@ SPACE
Space [isspace].
Definition: token.H:125
A List of wordRe with additional matching capabilities.
Definition: wordRes.H:54
A class for handling words, derived from Foam::string.
Definition: word.H:68
volScalarField & p
PtrList< volScalarField > & Y
const volScalarField & p0
Definition: EEqn.H:36
dynamicFvMesh & mesh
const wordList solidNames(rp["solid"])
OBJstream os(runTime.globalPath()/outputName)
#define writeProp(Name, Value)
#define WarningInFunction
Report a warning using Foam::Warning.
#define FUNCTION_NAME
volVectorField F(fluid.F())
void readFields(const typename GeoFieldType::Mesh &mesh, const IOobjectList &objects, const wordHashSet &selectedFields, LIFOStack< regIOobject * > &storedObjects)
Read the selected GeometricFields of the templated type.
void writeFields(const fvMesh &mesh, const wordHashSet &selectedFields, const bool writeFaceFields)
static constexpr const zero Zero
Global zero (0)
Definition: zero.H:131
word name(const expressions::valueTypeCode typeCode)
A word representation of a valueTypeCode. Empty for INVALID.
Definition: exprTraits.C:59
constexpr char nl
The newline '\n' character (0x0a)
Definition: Ostream.H:53
#define forAll(list, i)
Loop across all elements in list.
Definition: stdFoam.H:333