CollidingParcelIO.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 "CollidingParcel.H"
27 #include "IOstreams.H"
28 #include "IOField.H"
29 
30 // * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
31 
32 template<class ParcelType>
35 
36 template<class ParcelType>
38 (
39  offsetof(CollidingParcel<ParcelType>, collisionRecords_)
40  - offsetof(CollidingParcel<ParcelType>, f_)
41 );
42 
43 
44 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
45 
46 template<class ParcelType>
48 :
49  ParcelType(is, readFields),
50  f_(Zero),
51  angularMomentum_(Zero),
52  torque_(Zero),
53  collisionRecords_()
54 {
55  if (readFields)
56  {
57  if (is.format() == IOstream::ASCII)
58  {
59  is >> f_;
60  is >> angularMomentum_;
61  is >> torque_;
62  }
63  else
64  {
65  is.read(reinterpret_cast<char*>(&f_), sizeofFields_);
66  }
67 
68  is >> collisionRecords_;
69  }
70 
71  // Check state of Istream
72  is.check
73  (
74  "CollidingParcel<ParcelType>::Collisions"
75  "(const polyMesh&, Istream&, bool)"
76  );
77 }
78 
79 
80 template<class ParcelType>
81 template<class CloudType>
83 {
84  bool write = c.size();
85 
87 
88  IOField<vector> f(c.fieldIOobject("f", IOobject::MUST_READ), write);
89  c.checkFieldIOobject(c, f);
90 
91  IOField<vector> angularMomentum
92  (
93  c.fieldIOobject("angularMomentum", IOobject::MUST_READ),
94  write
95  );
96  c.checkFieldIOobject(c, angularMomentum);
97 
98  IOField<vector> torque
99  (
100  c.fieldIOobject("torque", IOobject::MUST_READ),
101  write
102  );
103  c.checkFieldIOobject(c, torque);
104 
105  labelFieldCompactIOField collisionRecordsPairAccessed
106  (
107  c.fieldIOobject("collisionRecordsPairAccessed", IOobject::MUST_READ),
108  write
109  );
110  c.checkFieldFieldIOobject(c, collisionRecordsPairAccessed);
111 
112  labelFieldCompactIOField collisionRecordsPairOrigProcOfOther
113  (
114  c.fieldIOobject
115  (
116  "collisionRecordsPairOrigProcOfOther",
118  ),
119  write
120  );
121  c.checkFieldFieldIOobject(c, collisionRecordsPairOrigProcOfOther);
122 
123  labelFieldCompactIOField collisionRecordsPairOrigIdOfOther
124  (
125  c.fieldIOobject
126  (
127  "collisionRecordsPairOrigIdOfOther",
129  ),
130  write
131  );
132  c.checkFieldFieldIOobject(c, collisionRecordsPairOrigProcOfOther);
133 
134  pairDataFieldCompactIOField collisionRecordsPairData
135  (
136  c.fieldIOobject("collisionRecordsPairData", IOobject::MUST_READ),
137  write
138  );
139  c.checkFieldFieldIOobject(c, collisionRecordsPairData);
140 
141  labelFieldCompactIOField collisionRecordsWallAccessed
142  (
143  c.fieldIOobject("collisionRecordsWallAccessed", IOobject::MUST_READ),
144  write
145  );
146  c.checkFieldFieldIOobject(c, collisionRecordsWallAccessed);
147 
148  vectorFieldCompactIOField collisionRecordsWallPRel
149  (
150  c.fieldIOobject("collisionRecordsWallPRel", IOobject::MUST_READ),
151  write
152  );
153  c.checkFieldFieldIOobject(c, collisionRecordsWallPRel);
154 
155  wallDataFieldCompactIOField collisionRecordsWallData
156  (
157  c.fieldIOobject("collisionRecordsWallData", IOobject::MUST_READ),
158  write
159  );
160  c.checkFieldFieldIOobject(c, collisionRecordsWallData);
161 
162  label i = 0;
163 
164  forAllIter(typename CloudType, c, iter)
165  {
166  CollidingParcel<ParcelType>& p = iter();
167 
168  p.f_ = f[i];
169  p.angularMomentum_ = angularMomentum[i];
170  p.torque_ = torque[i];
171 
172  p.collisionRecords_ = collisionRecordList
173  (
174  collisionRecordsPairAccessed[i],
175  collisionRecordsPairOrigProcOfOther[i],
176  collisionRecordsPairOrigIdOfOther[i],
177  collisionRecordsPairData[i],
178  collisionRecordsWallAccessed[i],
179  collisionRecordsWallPRel[i],
180  collisionRecordsWallData[i]
181  );
182 
183  i++;
184  }
185 }
186 
187 
188 template<class ParcelType>
189 template<class CloudType>
191 {
192  ParcelType::writeFields(c);
193 
194  label np = c.size();
195 
196  IOField<vector> f(c.fieldIOobject("f", IOobject::NO_READ), np);
197  IOField<vector> angularMomentum
198  (
199  c.fieldIOobject("angularMomentum", IOobject::NO_READ),
200  np
201  );
202  IOField<vector> torque(c.fieldIOobject("torque", IOobject::NO_READ), np);
203 
204  labelFieldCompactIOField collisionRecordsPairAccessed
205  (
206  c.fieldIOobject("collisionRecordsPairAccessed", IOobject::NO_READ),
207  np
208  );
209  labelFieldCompactIOField collisionRecordsPairOrigProcOfOther
210  (
211  c.fieldIOobject
212  (
213  "collisionRecordsPairOrigProcOfOther",
215  ),
216  np
217  );
218  labelFieldCompactIOField collisionRecordsPairOrigIdOfOther
219  (
220  c.fieldIOobject("collisionRecordsPairOrigIdOfOther", IOobject::NO_READ),
221  np
222  );
223  pairDataFieldCompactIOField collisionRecordsPairData
224  (
225  c.fieldIOobject("collisionRecordsPairData", IOobject::NO_READ),
226  np
227  );
228  labelFieldCompactIOField collisionRecordsWallAccessed
229  (
230  c.fieldIOobject("collisionRecordsWallAccessed", IOobject::NO_READ),
231  np
232  );
233  vectorFieldCompactIOField collisionRecordsWallPRel
234  (
235  c.fieldIOobject("collisionRecordsWallPRel", IOobject::NO_READ),
236  np
237  );
238  wallDataFieldCompactIOField collisionRecordsWallData
239  (
240  c.fieldIOobject("collisionRecordsWallData", IOobject::NO_READ),
241  np
242  );
243 
244  label i = 0;
245 
246  forAllConstIter(typename CloudType, c, iter)
247  {
248  const CollidingParcel<ParcelType>& p = iter();
249 
250  f[i] = p.f();
251  angularMomentum[i] = p.angularMomentum();
252  torque[i] = p.torque();
253 
254  collisionRecordsPairAccessed[i] = p.collisionRecords().pairAccessed();
255  collisionRecordsPairOrigProcOfOther[i] =
256  p.collisionRecords().pairOrigProcOfOther();
257  collisionRecordsPairOrigIdOfOther[i] =
258  p.collisionRecords().pairOrigIdOfOther();
259  collisionRecordsPairData[i] = p.collisionRecords().pairData();
260  collisionRecordsWallAccessed[i] = p.collisionRecords().wallAccessed();
261  collisionRecordsWallPRel[i] = p.collisionRecords().wallPRel();
262  collisionRecordsWallData[i] = p.collisionRecords().wallData();
263 
264  i++;
265  }
266 
267  const bool write = (np > 0);
268 
269  f.write(write);
270  angularMomentum.write(write);
271  torque.write(write);
272 
273  collisionRecordsPairAccessed.write(write);
274  collisionRecordsPairOrigProcOfOther.write(write);
275  collisionRecordsPairOrigIdOfOther.write(write);
276  collisionRecordsPairData.write(write);
277  collisionRecordsWallAccessed.write(write);
278  collisionRecordsWallPRel.write(write);
279  collisionRecordsWallData.write(write);
280 }
281 
282 
283 // * * * * * * * * * * * * * * * IOstream Operators * * * * * * * * * * * * //
284 
285 template<class ParcelType>
286 Foam::Ostream& Foam::operator<<
287 (
288  Ostream& os,
290 )
291 {
292  if (os.format() == IOstream::ASCII)
293  {
294  os << static_cast<const ParcelType&>(p)
295  << token::SPACE << p.f_
296  << token::SPACE << p.angularMomentum_
297  << token::SPACE << p.torque_
298  << token::SPACE << p.collisionRecords_;
299  }
300  else
301  {
302  os << static_cast<const ParcelType&>(p);
303  os.write
304  (
305  reinterpret_cast<const char*>(&p.f_),
306  CollidingParcel<ParcelType>::sizeofFields_
307  );
308  os << p.collisionRecords();
309  }
310 
311  // Check state of Ostream
312  os.check
313  (
314  "Ostream& operator<<(Ostream&, const CollidingParcel<ParcelType>&)"
315  );
316 
317  return os;
318 }
319 
320 
321 // ************************************************************************* //
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
Wrapper around parcel types to add collision modelling.
vector f_
Force on particle due to collisions [N].
CollidingParcel(const polyMesh &mesh, const barycentric &coordinates, const label celli, const label tetFacei, const label tetPti, const label facei)
Construct from mesh, coordinates and topology.
collisionRecordList collisionRecords_
Particle collision records.
static void writeFields(const CloudType &c)
Write.
vector angularMomentum_
Angular momentum of Parcel in global reference frame [kg m2/s].
static void readFields(CloudType &c)
Read.
vector torque_
Torque on particle due to collisions in global.
A Field of objects of type <Type> with automated input and output using a compact storage....
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.
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
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
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
CollisionRecordList< vector, vector > collisionRecordList
labelList f(nPoints)
volScalarField & p