moleculeIO.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-2024 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 "molecule.H"
27 #include "IOstreams.H"
28 #include "moleculeCloud.H"
29 
30 // * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
31 
32 const std::size_t Foam::molecule::sizeofFields_
33 (
34  offsetof(molecule, siteForces_) - offsetof(molecule, Q_)
35 );
36 
37 
38 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
39 
41 :
42  particle(is, readFields),
43  Q_(Zero),
44  v_(Zero),
45  a_(Zero),
46  pi_(Zero),
47  tau_(Zero),
48  specialPosition_(Zero),
49  potentialEnergy_(0.0),
50  rf_(Zero),
51  special_(0),
52  id_(0),
53  siteForces_(0),
54  sitePositions_(0)
55 {
56  if (readFields)
57  {
58  if (is.format() == IOstream::ASCII)
59  {
60  is >> Q_;
61  is >> v_;
62  is >> a_;
63  is >> pi_;
64  is >> tau_;
65  is >> specialPosition_;
66  potentialEnergy_ = readScalar(is);
67  is >> rf_;
68  special_ = readLabel(is);
69  id_ = readLabel(is);
70  is >> siteForces_;
71  is >> sitePositions_;
72  }
73  else
74  {
75  is.read(reinterpret_cast<char*>(&Q_), sizeofFields_);
76  is >> siteForces_ >> sitePositions_;
77  }
78  }
79 
80  // Check state of Istream
81  is.check
82  (
83  "Foam::molecule::molecule"
84  "(const Cloud<molecule>& cloud, Foam::Istream&), bool"
85  );
86 }
87 
88 
90 {
91  bool write = mC.size();
92 
94 
96  mC.checkFieldIOobject(mC, Q);
97 
99  mC.checkFieldIOobject(mC, v);
100 
102  mC.checkFieldIOobject(mC, a);
103 
105  mC.checkFieldIOobject(mC, pi);
106 
108  mC.checkFieldIOobject(mC, tau);
109 
110  IOField<vector> specialPosition
111  (
112  mC.fieldIOobject("specialPosition", IOobject::MUST_READ),
113  write
114  );
115  mC.checkFieldIOobject(mC, specialPosition);
116 
117  IOField<label> special
118  (
119  mC.fieldIOobject("special", IOobject::MUST_READ),
120  write
121  );
122  mC.checkFieldIOobject(mC, special);
123 
125  mC.checkFieldIOobject(mC, id);
126 
127  label i = 0;
128  forAllIter(moleculeCloud, mC, iter)
129  {
130  molecule& mol = iter();
131 
132  mol.Q_ = Q[i];
133  mol.v_ = v[i];
134  mol.a_ = a[i];
135  mol.pi_ = pi[i];
136  mol.tau_ = tau[i];
137  mol.specialPosition_ = specialPosition[i];
138  mol.special_ = special[i];
139  mol.id_ = id[i];
140  i++;
141  }
142 }
143 
144 
146 {
148 
149  label np = mC.size();
150 
155  IOField<vector> tau(mC.fieldIOobject("tau", IOobject::NO_READ), np);
156  IOField<vector> specialPosition
157  (
158  mC.fieldIOobject("specialPosition", IOobject::NO_READ),
159  np
160  );
161  IOField<label> special(mC.fieldIOobject("special", IOobject::NO_READ), np);
163 
164  // Post processing fields
165 
166  IOField<vector> piGlobal
167  (
168  mC.fieldIOobject("piGlobal", IOobject::NO_READ),
169  np
170  );
171 
172  IOField<vector> tauGlobal
173  (
174  mC.fieldIOobject("tauGlobal", IOobject::NO_READ),
175  np
176  );
177 
178  IOField<vector> orientation1
179  (
180  mC.fieldIOobject("orientation1", IOobject::NO_READ),
181  np
182  );
183 
184  IOField<vector> orientation2
185  (
186  mC.fieldIOobject("orientation2", IOobject::NO_READ),
187  np
188  );
189 
190  IOField<vector> orientation3
191  (
192  mC.fieldIOobject("orientation3", IOobject::NO_READ),
193  np
194  );
195 
196  label i = 0;
197  forAllConstIter(moleculeCloud, mC, iter)
198  {
199  const molecule& mol = iter();
200 
201  Q[i] = mol.Q_;
202  v[i] = mol.v_;
203  a[i] = mol.a_;
204  pi[i] = mol.pi_;
205  tau[i] = mol.tau_;
206  specialPosition[i] = mol.specialPosition_;
207  special[i] = mol.special_;
208  id[i] = mol.id_;
209 
210  piGlobal[i] = mol.Q_ & mol.pi_;
211  tauGlobal[i] = mol.Q_ & mol.tau_;
212 
213  orientation1[i] = mol.Q_ & vector(1,0,0);
214  orientation2[i] = mol.Q_ & vector(0,1,0);
215  orientation3[i] = mol.Q_ & vector(0,0,1);
216 
217  i++;
218  }
219 
220  const bool write = np > 0;
221 
222  Q.write(write);
223  v.write(write);
224  a.write(write);
225  pi.write(write);
226  tau.write(write);
227  specialPosition.write(write);
228  special.write(write);
229  id.write(write);
230 
231  piGlobal.write(write);
232  tauGlobal.write(write);
233 
234  orientation1.write(write);
235  orientation2.write(write);
236  orientation3.write(write);
237 
238  Info<< "writeFields " << mC.name() << endl;
239 
240  if (isA<moleculeCloud>(mC))
241  {
242  const moleculeCloud& m = dynamic_cast<const moleculeCloud&>(mC);
243 
244  m.writeXYZ
245  (
246  m.mesh().time().timePath()/cloud::prefix/"moleculeCloud.xmol"
247  );
248  }
249 }
250 
251 
252 // * * * * * * * * * * * * * * * IOstream Operators * * * * * * * * * * * * //
253 
255 {
256  if (os.format() == IOstream::ASCII)
257  {
258  os << token::SPACE << static_cast<const particle&>(mol)
259  << token::SPACE << mol.Q_
260  << token::SPACE << mol.v_
261  << token::SPACE << mol.a_
262  << token::SPACE << mol.pi_
263  << token::SPACE << mol.tau_
264  << token::SPACE << mol.specialPosition_
265  << token::SPACE << mol.potentialEnergy_
266  << token::SPACE << mol.rf_
267  << token::SPACE << mol.special_
268  << token::SPACE << mol.id_
269  << token::SPACE << mol.siteForces_
270  << token::SPACE << mol.sitePositions_;
271  }
272  else
273  {
274  os << static_cast<const particle&>(mol);
275  os.write
276  (
277  reinterpret_cast<const char*>(&mol.Q_),
278  molecule::sizeofFields_
279  );
280  os << mol.siteForces_ << mol.sitePositions_;
281  }
282 
283  // Check state of Ostream
284  os.check
285  (
286  "Foam::Ostream& Foam::operator<<"
287  "(Foam::Ostream&, const Foam::molecule&)"
288  );
289 
290  return os;
291 }
292 
293 
294 // ************************************************************************* //
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
void checkFieldIOobject(const Cloud< ParticleType > &c, const IOField< DataType > &data) const
Check lagrangian data field.
Definition: CloudIO.C:188
IOobject fieldIOobject(const word &fieldName, const IOobject::readOption r) const
Helper to construct IOobject for field and current time.
Definition: CloudIO.C:168
label size() const
Return the number of particles in the cloud.
Definition: Cloud.H:191
A primitive field of type <Type> with automated input and output.
Definition: IOField.H:53
const word & name() const
Return name.
Definition: IOobject.H:310
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 Ostream & write(const char)=0
Write character.
fileName timePath() const
Return current time path.
Definition: Time.H:281
static const word prefix
The prefix to local: lagrangian.
Definition: cloud.H:63
const polyMesh & mesh() const
void writeXYZ(const fileName &fName) const
Write molecule sites in XYZ format.
Foam::molecule.
Definition: molecule.H:68
molecule(const polyMesh &mesh, const vector &position, const label celli, label &nLocateBoundaryHits, const tensor &Q, const vector &v, const vector &a, const vector &pi, const vector &tau, const vector &specialPosition, const constantProperties &constProps, const label special, const label id)
Construct from a position and a cell, searching for the rest of the.
Definition: moleculeI.H:223
static void writeFields(const Cloud< molecule > &mC)
Definition: moleculeIO.C:145
static void readFields(Cloud< molecule > &mC)
Definition: moleculeIO.C:89
const Time & time() const
Return time.
Base particle class.
Definition: particle.H:83
static void readFields(TrackCloudType &c)
Read the fields associated with the owner cloud.
static void writeFields(const TrackCloudType &c)
Write the fields associated with the owner cloud.
virtual bool write(const bool write=true) const
Write using setting from DB.
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
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:257
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
messageStream Info
Vector< scalar > vector
A scalar version of the templated Vector.
Definition: vector.H:49
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
Ostream & operator<<(Ostream &os, const fvConstraints &constraints)