MomentumParcelIO.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-2022 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 "MomentumParcel.H"
27 #include "IOstreams.H"
28 #include "IOField.H"
29 #include "Cloud.H"
30 
31 // * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
32 
33 template<class ParcelType>
36 
37 template<class ParcelType>
39 (
40  sizeof(MomentumParcel<ParcelType>)
41  - offsetof(MomentumParcel<ParcelType>, moving_)
42 );
43 
44 
45 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
46 
47 template<class ParcelType>
49 :
50  ParcelType(is, readFields),
51  moving_(false),
52  typeId_(0),
53  nParticle_(0.0),
54  d_(0.0),
55  dTarget_(0.0),
56  U_(Zero),
57  rho_(0.0),
58  age_(0.0),
59  tTurb_(0.0),
60  UTurb_(Zero)
61 {
62  if (readFields)
63  {
64  if (is.format() == IOstream::ASCII)
65  {
66  moving_ = readBool(is);
67  typeId_ = readLabel(is);
68  nParticle_ = readScalar(is);
69  d_ = readScalar(is);
70  dTarget_ = readScalar(is);
71  is >> U_;
72  rho_ = readScalar(is);
73  age_ = readScalar(is);
74  tTurb_ = readScalar(is);
75  is >> UTurb_;
76  }
77  else
78  {
79  is.read(reinterpret_cast<char*>(&moving_), sizeofFields_);
80  }
81  }
82 
83  // Check state of Istream
84  is.check
85  (
86  "MomentumParcel<ParcelType>::MomentumParcel"
87  "(const polyMesh&, Istream&, bool)"
88  );
89 }
90 
91 
92 template<class ParcelType>
93 template<class CloudType>
95 {
96  bool write = c.size();
97 
99 
100  IOField<label> moving
101  (
102  c.fieldIOobject("active", IOobject::MUST_READ),
103  write
104  );
105  c.checkFieldIOobject(c, moving);
106 
107  IOField<label> typeId
108  (
109  c.fieldIOobject("typeId", IOobject::MUST_READ),
110  write
111  );
112  c.checkFieldIOobject(c, typeId);
113 
114  IOField<scalar> nParticle
115  (
116  c.fieldIOobject("nParticle", IOobject::MUST_READ),
117  write
118  );
119  c.checkFieldIOobject(c, nParticle);
120 
122  (
123  c.fieldIOobject("d", IOobject::MUST_READ),
124  write
125  );
126  c.checkFieldIOobject(c, d);
127 
128  IOField<scalar> dTarget
129  (
130  c.fieldIOobject("dTarget", IOobject::MUST_READ),
131  write
132  );
133  c.checkFieldIOobject(c, dTarget);
134 
136  (
137  c.fieldIOobject("U", IOobject::MUST_READ),
138  write
139  );
140  c.checkFieldIOobject(c, U);
141 
143  (
144  c.fieldIOobject("rho", IOobject::MUST_READ),
145  write
146  );
147  c.checkFieldIOobject(c, rho);
148 
149  IOField<scalar> age
150  (
151  c.fieldIOobject("age", IOobject::MUST_READ),
152  write
153  );
154  c.checkFieldIOobject(c, age);
155 
156  IOField<scalar> tTurb
157  (
158  c.fieldIOobject("tTurb", IOobject::MUST_READ),
159  write
160  );
161  c.checkFieldIOobject(c, tTurb);
162 
163  IOField<vector> UTurb
164  (
165  c.fieldIOobject("UTurb", IOobject::MUST_READ),
166  write
167  );
168  c.checkFieldIOobject(c, UTurb);
169 
170  label i = 0;
171 
172  forAllIter(typename CloudType, c, iter)
173  {
174  MomentumParcel<ParcelType>& p = iter();
175 
176  p.moving_ = moving[i];
177  p.typeId_ = typeId[i];
178  p.nParticle_ = nParticle[i];
179  p.d_ = d[i];
180  p.dTarget_ = dTarget[i];
181  p.U_ = U[i];
182  p.rho_ = rho[i];
183  p.age_ = age[i];
184  p.tTurb_ = tTurb[i];
185  p.UTurb_ = UTurb[i];
186 
187  i++;
188  }
189 }
190 
191 
192 template<class ParcelType>
193 template<class CloudType>
195 {
196  ParcelType::writeFields(c);
197 
198  label np = c.size();
199 
200  IOField<label> moving(c.fieldIOobject("active", IOobject::NO_READ), np);
201  IOField<label> typeId(c.fieldIOobject("typeId", IOobject::NO_READ), np);
202  IOField<scalar> nParticle
203  (
204  c.fieldIOobject("nParticle", IOobject::NO_READ),
205  np
206  );
207  IOField<scalar> d(c.fieldIOobject("d", IOobject::NO_READ), np);
208  IOField<scalar> dTarget(c.fieldIOobject("dTarget", IOobject::NO_READ), np);
209  IOField<vector> U(c.fieldIOobject("U", IOobject::NO_READ), np);
210  IOField<scalar> rho(c.fieldIOobject("rho", IOobject::NO_READ), np);
211  IOField<scalar> age(c.fieldIOobject("age", IOobject::NO_READ), np);
212  IOField<scalar> tTurb(c.fieldIOobject("tTurb", IOobject::NO_READ), np);
213  IOField<vector> UTurb(c.fieldIOobject("UTurb", IOobject::NO_READ), np);
214 
215  label i = 0;
216 
217  forAllConstIter(typename CloudType, c, iter)
218  {
219  const MomentumParcel<ParcelType>& p = iter();
220 
221  moving[i] = p.moving();
222  typeId[i] = p.typeId();
223  nParticle[i] = p.nParticle();
224  d[i] = p.d();
225  dTarget[i] = p.dTarget();
226  U[i] = p.U();
227  rho[i] = p.rho();
228  age[i] = p.age();
229  tTurb[i] = p.tTurb();
230  UTurb[i] = p.UTurb();
231 
232  i++;
233  }
234 
235  const bool write = np > 0;
236 
237  moving.write(write);
238  typeId.write(write);
239  nParticle.write(write);
240  d.write(write);
241  dTarget.write(write);
242  U.write(write);
243  rho.write(write);
244  age.write(write);
245  tTurb.write(write);
246  UTurb.write(write);
247 }
248 
249 
250 // * * * * * * * * * * * * * * * IOstream Operators * * * * * * * * * * * * //
251 
252 template<class ParcelType>
253 Foam::Ostream& Foam::operator<<
254 (
255  Ostream& os,
257 )
258 {
259  if (os.format() == IOstream::ASCII)
260  {
261  os << static_cast<const ParcelType&>(p)
262  << token::SPACE << p.moving()
263  << token::SPACE << p.typeId()
264  << token::SPACE << p.nParticle()
265  << token::SPACE << p.d()
266  << token::SPACE << p.dTarget()
267  << token::SPACE << p.U()
268  << token::SPACE << p.rho()
269  << token::SPACE << p.age()
270  << token::SPACE << p.tTurb()
271  << token::SPACE << p.UTurb();
272  }
273  else
274  {
275  os << static_cast<const ParcelType&>(p);
276  os.write
277  (
278  reinterpret_cast<const char*>(&p.moving_),
279  MomentumParcel<ParcelType>::sizeofFields_
280  );
281  }
282 
283  // Check state of Ostream
284  os.check
285  (
286  "Ostream& operator<<(Ostream&, const MomentumParcel<ParcelType>&)"
287  );
288 
289  return os;
290 }
291 
292 
293 // ************************************************************************* //
Useful combination of include files which define Sin, Sout and Serr and the use of IO streams general...
#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:80
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.
Momentum parcel class with rotational motion (as spherical particles only) and one/two-way coupling w...
label typeId_
Parcel type id.
vector UTurb_
Turbulent velocity fluctuation [m/s].
scalar rho_
Density [kg/m^3].
bool moving_
Moving flag - tracking stopped when moving = false.
static void readFields(TrackCloudType &c)
Read.
scalar tTurb_
Time spent in turbulent eddy [s].
scalar dTarget_
Target diameter [m].
scalar d_
Diameter [m].
MomentumParcel(const polyMesh &mesh, const barycentric &coordinates, const label celli, const label tetFacei, const label tetPti, const label facei)
Construct from mesh, coordinates and topology.
static void writeFields(const TrackCloudType &c)
Write.
scalar age_
Age [s].
vector U_
Velocity of Parcel [m/s].
scalar nParticle_
Number of particles in Parcel.
An Ostream is an abstract base class for all output systems (streams, files, token lists,...
Definition: Ostream.H:57
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
U
Definition: pEqn.H:72
const dimensionedScalar c
Speed of light in a vacuum.
void write(std::ostream &os, const bool binary, List< floatScalar > &fField)
Write floats ascii or binary.
static const zero Zero
Definition: zero.H:97
bool readBool(Istream &)
Definition: boolIO.C:66
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
bool readScalar(const char *buf, doubleScalar &s)
Read whole of buf as a scalar. Return true if successful.
Definition: doubleScalar.H:75
label readLabel(Istream &is)
Definition: label.H:64
volScalarField & p