particleIO.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) 2011-2017 OpenFOAM Foundation
9  Copyright (C) 2016-2019 OpenCFD Ltd.
10 -------------------------------------------------------------------------------
11 License
12  This file is part of OpenFOAM.
13 
14  OpenFOAM is free software: you can redistribute it and/or modify it
15  under the terms of the GNU General Public License as published by
16  the Free Software Foundation, either version 3 of the License, or
17  (at your option) any later version.
18 
19  OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
20  ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
21  FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
22  for more details.
23 
24  You should have received a copy of the GNU General Public License
25  along with OpenFOAM. If not, see <http://www.gnu.org/licenses/>.
26 
27 \*---------------------------------------------------------------------------*/
28 
29 #include "particle.H"
30 #include "IOstreams.H"
31 
32 // * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
33 
35 
36 const std::size_t Foam::particle::sizeofPosition
37 (
38  offsetof(particle, facei_) - offsetof(particle, coordinates_)
39 );
40 
41 const std::size_t Foam::particle::sizeofFields
42 (
43  sizeof(particle) - offsetof(particle, coordinates_)
44 );
45 
46 
47 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
48 
50 (
51  const polyMesh& mesh,
52  Istream& is,
53  bool readFields,
54  bool newFormat
55 )
56 :
57  mesh_(mesh),
58  coordinates_(),
59  celli_(-1),
60  tetFacei_(-1),
61  tetPti_(-1),
62  facei_(-1),
63  stepFraction_(0.0),
64  origProc_(Pstream::myProcNo()),
65  origId_(-1)
66 {
67  if (newFormat)
68  {
69  if (is.format() == IOstream::ASCII)
70  {
71  is >> coordinates_ >> celli_ >> tetFacei_ >> tetPti_;
72  if (readFields)
73  {
74  is >> facei_ >> stepFraction_ >> origProc_ >> origId_;
75  }
76  }
77  else if (!is.checkLabelSize<>() || !is.checkScalarSize<>())
78  {
79  // Non-native label or scalar size
80 
81  is.beginRawRead();
82 
83  readRawScalar(is, coordinates_.data(), barycentric::nComponents);
84  readRawLabel(is, &celli_);
85  readRawLabel(is, &tetFacei_);
86  readRawLabel(is, &tetPti_);
87 
88  if (readFields)
89  {
90  readRawLabel(is, &facei_);
91  readRawScalar(is, &stepFraction_);
92  readRawLabel(is, &origProc_);
93  readRawLabel(is, &origId_);
94  }
95 
96  is.endRawRead();
97  }
98  else
99  {
100  if (readFields)
101  {
102  is.read(reinterpret_cast<char*>(&coordinates_), sizeofFields);
103  }
104  else
105  {
106  is.read(reinterpret_cast<char*>(&coordinates_), sizeofPosition);
107  }
108  }
109  }
110  else
111  {
113 
114  if (is.format() == IOstream::ASCII)
115  {
116  is >> p.position >> p.celli;
117 
118  if (readFields)
119  {
120  is >> p.facei
121  >> p.stepFraction
122  >> p.tetFacei
123  >> p.tetPti
124  >> p.origProc
125  >> p.origId;
126  }
127  }
128  else if (!is.checkLabelSize<>() || !is.checkScalarSize<>())
129  {
130  // Non-native label or scalar size
131 
132  is.beginRawRead();
133 
134  readRawScalar(is, p.position.data(), vector::nComponents);
135  readRawLabel(is, &p.celli);
136 
137  if (readFields)
138  {
139  readRawLabel(is, &p.facei);
140  readRawScalar(is, &p.stepFraction);
141  readRawLabel(is, &p.tetFacei);
142  readRawLabel(is, &p.tetPti);
143  readRawLabel(is, &p.origProc);
144  readRawLabel(is, &p.origId);
145  }
146 
147  is.endRawRead();
148  }
149  else
150  {
151  if (readFields)
152  {
153  // Read whole struct
154  const size_t s =
155  (
156  sizeof(positionsCompat1706)
157  - offsetof(positionsCompat1706, position)
158  );
159  is.read(reinterpret_cast<char*>(&p.position), s);
160  }
161  else
162  {
163  // Read only position and cell
164  const size_t s =
165  (
166  offsetof(positionsCompat1706, facei)
167  - offsetof(positionsCompat1706, position)
168  );
169  is.read(reinterpret_cast<char*>(&p.position), s);
170  }
171  }
172 
173  if (readFields)
174  {
175  // Note: other position-based properties are set using locate(...)
176  stepFraction_ = p.stepFraction;
177  origProc_ = p.origProc;
178  origId_ = p.origId;
179  }
180 
181  locate
182  (
183  p.position,
184  nullptr,
185  p.celli,
186  false,
187  "Particle initialised with a location outside of the mesh."
188  );
189  }
190 
191  // Check state of Istream
192  is.check(FUNCTION_NAME);
193 }
194 
195 
197 (
198  Ostream& os,
199  const wordRes& filters,
200  const word& delim,
201  const bool namesOnly
202 ) const
203 {
204  #undef writeProp
205  #define writeProp(Name, Value) \
206  particle::writeProperty(os, Name, Value, namesOnly, delim, filters)
207 
208  writeProp("coordinates", coordinates_);
209  writeProp("position", position());
210  writeProp("celli", celli_);
211  writeProp("tetFacei", tetFacei_);
212  writeProp("tetPti", tetPti_);
213  writeProp("facei", facei_);
214  writeProp("stepFraction", stepFraction_);
215  writeProp("origProc", origProc_);
216  writeProp("origId", origId_);
217 
218  #undef writeProp
219 }
220 
221 
223 {
224  if (os.format() == IOstream::ASCII)
225  {
226  os << coordinates_
227  << token::SPACE << celli_
228  << token::SPACE << tetFacei_
229  << token::SPACE << tetPti_;
230  }
231  else
232  {
233  os.write(reinterpret_cast<const char*>(&coordinates_), sizeofPosition);
234  }
235 
236  // Check state of Ostream
237  os.check(FUNCTION_NAME);
238 }
239 
240 
242 {
243  if (os.format() == IOstream::ASCII)
244  {
245  os << position() << token::SPACE << celli_;
246  }
247  else
248  {
250 
251  const size_t s =
252  (
253  offsetof(positionsCompat1706, facei)
254  - offsetof(positionsCompat1706, position)
255  );
256 
257  p.position = position();
258  p.celli = celli_;
259 
260  os.write(reinterpret_cast<const char*>(&p.position), s);
261  }
262 
263  // Check state of Ostream
265 }
266 
267 
269 {
270  if (os.format() == IOstream::ASCII)
271  {
272  os << p.coordinates_
273  << token::SPACE << p.celli_
274  << token::SPACE << p.tetFacei_
275  << token::SPACE << p.tetPti_
276  << token::SPACE << p.facei_
277  << token::SPACE << p.stepFraction_
278  << token::SPACE << p.origProc_
279  << token::SPACE << p.origId_;
280  }
281  else
282  {
283  os.write
284  (
285  reinterpret_cast<const char*>(&p.coordinates_),
286  particle::sizeofFields
287  );
288  }
289 
290  return os;
291 }
292 
293 
294 // ************************************************************************* //
p
volScalarField & p
Definition: createFieldRefs.H:8
Foam::particle::writeCoordinates
void writeCoordinates(Ostream &os) const
Write the particle barycentric coordinates and cell info.
Definition: particleIO.C:222
IOstreams.H
Useful combination of include files which define Sin, Sout and Serr and the use of IO streams general...
Foam::word
A class for handling words, derived from Foam::string.
Definition: word.H:65
s
gmvFile<< "tracers "<< particles.size()<< nl;for(const passiveParticle &p :particles){ gmvFile<< p.position().x()<< " ";}gmvFile<< nl;for(const passiveParticle &p :particles){ gmvFile<< p.position().y()<< " ";}gmvFile<< nl;for(const passiveParticle &p :particles){ gmvFile<< p.position().z()<< " ";}gmvFile<< nl;forAll(lagrangianScalarNames, i){ word name=lagrangianScalarNames[i];IOField< scalar > s(IOobject(name, runTime.timeName(), cloud::prefix, mesh, IOobject::MUST_READ, IOobject::NO_WRITE))
Definition: gmvOutputSpray.H:25
Foam::IOstreamOption::format
streamFormat format() const noexcept
Get the current stream format.
Definition: IOstreamOption.H:286
Foam::string
A class for handling character strings derived from std::string.
Definition: string.H:76
Foam::OBJstream::write
virtual Ostream & write(const char c)
Write character.
Definition: OBJstream.C:78
writeProp
#define writeProp(Name, Value)
Foam::polyMesh
Mesh consisting of general polyhedral cells.
Definition: polyMesh.H:77
Foam::particle::writeProperties
void writeProperties(Ostream &os, const wordRes &filters, const word &delim, const bool namesOnly) const
Write individual particle properties to stream.
Definition: particleIO.C:197
Foam::operator<<
Ostream & operator<<(Ostream &, const boundaryPatch &p)
Write boundaryPatch as dictionary entries (without surrounding braces)
Definition: boundaryPatch.C:83
Foam::IOstream::checkLabelSize
std::enable_if< std::is_integral< T >::value, bool >::type checkLabelSize() const noexcept
Definition: IOstream.H:300
Foam::Istream::endRawRead
virtual bool endRawRead()=0
End of low-level raw binary read.
Foam::Istream
An Istream is an abstract base class for all input systems (streams, files, token lists etc)....
Definition: Istream.H:61
Foam::particle::propertyList
static string propertyList()
Definition: particle.H:349
Foam::particle::writePosition
virtual void writePosition(Ostream &os) const
Write the particle position and cell id.
Definition: particleIO.C:241
Foam::IOstream::checkScalarSize
std::enable_if< std::is_floating_point< T >::value, bool >::type checkScalarSize() const noexcept
Definition: IOstream.H:309
Foam::IOstream::check
virtual bool check(const char *operation) const
Check IOstream status for given operation.
Definition: IOstream.C:58
os
OBJstream os(runTime.globalPath()/outputName)
mesh
dynamicFvMesh & mesh
Definition: createDynamicFvMesh.H:6
Foam::particle::propertyList_
static string propertyList_
String representation of properties.
Definition: particle.H:349
Foam::readFields
void readFields(const typename GeoFieldType::Mesh &mesh, const IOobjectList &objects, const wordHashSet &selectedFields, LIFOStack< regIOobject * > &storedObjects)
Read the selected GeometricFields of the templated type.
Definition: ReadFieldsTemplates.C:312
Foam::particle::positionsCompat1706
Old particle positions content for OpenFOAM-1706 and earlier.
Definition: particle.H:116
Foam::IOstreamOption::ASCII
"ascii" (normal default)
Definition: IOstreamOption.H:72
Foam::Istream::beginRawRead
virtual bool beginRawRead()=0
Start of low-level raw binary read.
Foam::token::SPACE
Space [isspace].
Definition: token.H:125
Foam::particle
Base particle class.
Definition: particle.H:76
Foam::wordRes
A List of wordRe with additional matching capabilities.
Definition: wordRes.H:51
FUNCTION_NAME
#define FUNCTION_NAME
Definition: messageStream.H:295
particle.H
Foam::readRawLabel
label readRawLabel(Istream &is)
Read raw label from binary stream.
Definition: label.C:46
Foam::Ostream
An Ostream is an abstract base class for all output systems (streams, files, token lists,...
Definition: Ostream.H:56
Foam::particle::particle
particle(const polyMesh &mesh, const barycentric &coordinates, const label celli, const label tetFacei, const label tetPti)
Construct from components.
Definition: particle.C:516
Foam::Istream::read
virtual Istream & read(token &)=0
Return next token from stream.