SprayParcelIO.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 | Copyright (C) 2011-2016 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 "SprayParcel.H"
27 #include "IOstreams.H"
28 
29 // * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
30 
31 template<class ParcelType>
33 (
34  sizeof(SprayParcel<ParcelType>) - sizeof(ParcelType)
35 );
36 
37 
38 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
39 
40 template<class ParcelType>
42 (
43  const polyMesh& mesh,
44  Istream& is,
45  bool readFields
46 )
47 :
48  ParcelType(mesh, is, readFields),
49  d0_(0.0),
50  position0_(Zero),
51  sigma_(0.0),
52  mu_(0.0),
53  liquidCore_(0.0),
54  KHindex_(0.0),
55  y_(0.0),
56  yDot_(0.0),
57  tc_(0.0),
58  ms_(0.0),
59  injector_(1.0),
60  tMom_(GREAT),
61  user_(0.0)
62 {
63  if (readFields)
64  {
65  if (is.format() == IOstream::ASCII)
66  {
67  d0_ = readScalar(is);
68  is >> position0_;
69  sigma_ = readScalar(is);
70  mu_ = readScalar(is);
71  liquidCore_ = readScalar(is);
72  KHindex_ = readScalar(is);
73  y_ = readScalar(is);
74  yDot_ = readScalar(is);
75  tc_ = readScalar(is);
76  ms_ = readScalar(is);
77  injector_ = readScalar(is);
78  tMom_ = readScalar(is);
79  user_ = readScalar(is);
80  }
81  else
82  {
83  is.read(reinterpret_cast<char*>(&d0_), sizeofFields_);
84  }
85  }
86 
87  // Check state of Istream
88  is.check
89  (
90  "SprayParcel<ParcelType>::SprayParcel"
91  "("
92  "const polyMesh, "
93  "Istream&, "
94  "bool"
95  ")"
96  );
97 }
98 
99 
100 template<class ParcelType>
101 template<class CloudType>
103 {
104  if (!c.size())
105  {
106  return;
107  }
108 
110 }
111 
112 
113 template<class ParcelType>
114 template<class CloudType, class CompositionType>
116 (
117  CloudType& c,
118  const CompositionType& compModel
119 )
120 {
121  if (!c.size())
122  {
123  return;
124  }
125 
126  ParcelType::readFields(c, compModel);
127 
128  IOField<scalar> d0(c.fieldIOobject("d0", IOobject::MUST_READ));
129  c.checkFieldIOobject(c, d0);
130 
131  IOField<vector> position0
132  (
133  c.fieldIOobject("position0", IOobject::MUST_READ)
134  );
135  c.checkFieldIOobject(c, position0);
136 
137  IOField<scalar> sigma(c.fieldIOobject("sigma", IOobject::MUST_READ));
139 
140  IOField<scalar> mu(c.fieldIOobject("mu", IOobject::MUST_READ));
141  c.checkFieldIOobject(c, mu);
142 
143  IOField<scalar> liquidCore(c.fieldIOobject
144  (
145  "liquidCore", IOobject::MUST_READ)
146  );
147  c.checkFieldIOobject(c, liquidCore);
148 
149  IOField<scalar> KHindex(c.fieldIOobject("KHindex", IOobject::MUST_READ));
150  c.checkFieldIOobject(c, KHindex);
151 
152  IOField<scalar> y(c.fieldIOobject("y", IOobject::MUST_READ));
153  c.checkFieldIOobject(c, y);
154 
155  IOField<scalar> yDot(c.fieldIOobject("yDot", IOobject::MUST_READ));
156  c.checkFieldIOobject(c, yDot);
157 
158  IOField<scalar> tc(c.fieldIOobject("tc", IOobject::MUST_READ));
159  c.checkFieldIOobject(c, tc);
160 
161  IOField<scalar> ms(c.fieldIOobject("ms", IOobject::MUST_READ));
162  c.checkFieldIOobject(c, ms);
163 
164  IOField<scalar> injector(c.fieldIOobject("injector", IOobject::MUST_READ));
165  c.checkFieldIOobject(c, injector);
166 
167  IOField<scalar> tMom(c.fieldIOobject("tMom", IOobject::MUST_READ));
168  c.checkFieldIOobject(c, tMom);
169 
170  IOField<scalar> user(c.fieldIOobject("user", IOobject::MUST_READ));
171  c.checkFieldIOobject(c, user);
172 
173  label i = 0;
174  forAllIter(typename Cloud<SprayParcel<ParcelType>>, c, iter)
175  {
176  SprayParcel<ParcelType>& p = iter();
177  p.d0_ = d0[i];
178  p.position0_ = position0[i];
179  p.sigma_ = sigma[i];
180  p.mu_ = mu[i];
181  p.liquidCore_ = liquidCore[i];
182  p.KHindex_ = KHindex[i];
183  p.y_ = y[i];
184  p.yDot_ = yDot[i];
185  p.tc_ = tc[i];
186  p.ms_ = ms[i];
187  p.injector_ = injector[i];
188  p.tMom_ = tMom[i];
189  p.user_ = user[i];
190  i++;
191  }
192 }
193 
194 
195 template<class ParcelType>
196 template<class CloudType>
198 {
199  ParcelType::writeFields(c);
200 }
201 
202 
203 template<class ParcelType>
204 template<class CloudType, class CompositionType>
206 (
207  const CloudType& c,
208  const CompositionType& compModel
209 )
210 {
211  ParcelType::writeFields(c, compModel);
212 
213  label np = c.size();
214 
215  IOField<scalar> d0(c.fieldIOobject("d0", IOobject::NO_READ), np);
216  IOField<vector> position0
217  (
218  c.fieldIOobject("position0", IOobject::NO_READ),
219  np
220  );
221  IOField<scalar> sigma(c.fieldIOobject("sigma", IOobject::NO_READ), np);
222  IOField<scalar> mu(c.fieldIOobject("mu", IOobject::NO_READ), np);
223  IOField<scalar> liquidCore
224  (
225  c.fieldIOobject("liquidCore", IOobject::NO_READ),
226  np
227  );
228  IOField<scalar> KHindex(c.fieldIOobject("KHindex", IOobject::NO_READ), np);
229  IOField<scalar> y(c.fieldIOobject("y", IOobject::NO_READ), np);
230  IOField<scalar> yDot(c.fieldIOobject("yDot", IOobject::NO_READ), np);
231  IOField<scalar> tc(c.fieldIOobject("tc", IOobject::NO_READ), np);
232  IOField<scalar> ms(c.fieldIOobject("ms", IOobject::NO_READ), np);
233  IOField<scalar> injector
234  (
235  c.fieldIOobject("injector", IOobject::NO_READ),
236  np
237  );
238  IOField<scalar> tMom(c.fieldIOobject("tMom", IOobject::NO_READ), np);
239  IOField<scalar> user(c.fieldIOobject("user", IOobject::NO_READ), np);
240 
241  label i = 0;
242  forAllConstIter(typename Cloud<SprayParcel<ParcelType>>, c, iter)
243  {
244  const SprayParcel<ParcelType>& p = iter();
245  d0[i] = p.d0_;
246  position0[i] = p.position0_;
247  sigma[i] = p.sigma_;
248  mu[i] = p.mu_;
249  liquidCore[i] = p.liquidCore_;
250  KHindex[i] = p.KHindex_;
251  y[i] = p.y_;
252  yDot[i] = p.yDot_;
253  tc[i] = p.tc_;
254  ms[i] = p.ms_;
255  injector[i] = p.injector_;
256  tMom[i] = p.tMom_;
257  user[i] = p.user_;
258  i++;
259  }
260 
261  d0.write();
262  position0.write();
263  sigma.write();
264  mu.write();
265  liquidCore.write();
266  KHindex.write();
267  y.write();
268  yDot.write();
269  tc.write();
270  ms.write();
271  injector.write();
272  tMom.write();
273  user.write();
274 }
275 
276 
277 // * * * * * * * * * * * * * * * IOstream Operators * * * * * * * * * * * * //
278 
279 template<class ParcelType>
280 Foam::Ostream& Foam::operator<<
281 (
282  Ostream& os,
284 )
285 {
286  if (os.format() == IOstream::ASCII)
287  {
288  os << static_cast<const ParcelType&>(p)
289  << token::SPACE << p.d0()
290  << token::SPACE << p.position0()
291  << token::SPACE << p.sigma()
292  << token::SPACE << p.mu()
293  << token::SPACE << p.liquidCore()
294  << token::SPACE << p.KHindex()
295  << token::SPACE << p.y()
296  << token::SPACE << p.yDot()
297  << token::SPACE << p.tc()
298  << token::SPACE << p.ms()
299  << token::SPACE << p.injector()
300  << token::SPACE << p.tMom()
301  << token::SPACE << p.user();
302  }
303  else
304  {
305  os << static_cast<const ParcelType&>(p);
306  os.write
307  (
308  reinterpret_cast<const char*>(&p.d0_),
310  );
311  }
312 
313  // Check state of Ostream
314  os.check
315  (
316  "Ostream& operator<<(Ostream&, const SprayParcel<ParcelType>&)"
317  );
318 
319  return os;
320 }
321 
322 
323 // ************************************************************************* //
#define readScalar
Definition: doubleScalar.C:38
streamFormat format() const
Return current stream format.
Definition: IOstream.H:377
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 checkFieldIOobject(const Cloud< ParticleType > &c, const IOField< DataType > &data) const
Check lagrangian data field.
Definition: CloudIO.C:215
scalar tMom_
Momentum relaxation time (needed for calculating parcel acc.)
Definition: SprayParcel.H:170
#define forAllIter(Container, container, iter)
Iterate across all elements in the container object of type.
Definition: UList.H:453
An Istream is an abstract base class for all input systems (streams, files, token lists etc)...
Definition: Istream.H:57
virtual bool check(const char *operation) const
Check IOstream status for given operation.
Definition: IOstream.C:92
const dimensionedScalar sigma
Stefan-Boltzmann constant: default SI units: [W/m2/K4].
scalar liquidCore_
Part of liquid core ( >0.5=liquid, <0.5=droplet )
Definition: SprayParcel.H:148
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
Useful combination of include files which define Sin, Sout and Serr and the use of IO streams general...
scalar y
virtual Istream & read(token &)=0
Return next token from stream.
SprayParcel(const polyMesh &mesh, const vector &position, const label celli, const label tetFacei, const label tetPtI)
Construct from owner, position, and cloud owner.
Definition: SprayParcelI.H:108
static const zero Zero
Definition: zero.H:91
scalar injector_
Injected from injector (needed e.g. for calculating distance.
Definition: SprayParcel.H:167
scalar user_
Passive scalar (extra variable to be defined by user)
Definition: SprayParcel.H:173
static void readFields(CloudType &c, const CompositionType &compModel)
Read.
forAllConstIter(PtrDictionary< phaseModel >, mixture.phases(), phase)
Definition: pEqn.H:29
Base cloud calls templated on particle type.
Definition: Cloud.H:52
An Ostream is an abstract base class for all output systems (streams, files, token lists...
Definition: Ostream.H:53
label size() const
Definition: Cloud.H:175
scalar sigma_
Liquid surface tension [N/m].
Definition: SprayParcel.H:142
static void writeFields(const CloudType &c, const CompositionType &compModel)
Write.
const dimensionedScalar mu
Atomic mass unit.
vector position0_
Injection position.
Definition: SprayParcel.H:139
scalar yDot_
Rate of change of spherical deviation.
Definition: SprayParcel.H:157
scalar ms_
Stripped parcel mass due to breakup.
Definition: SprayParcel.H:163
virtual bool write() const
Write using setting from DB.
scalar d0_
Initial droplet diameter.
Definition: SprayParcel.H:136
Mesh consisting of general polyhedral cells.
Definition: polyMesh.H:74
IOobject fieldIOobject(const word &fieldName, const IOobject::readOption r) const
Helper to construct IOobject for field and current time.
Definition: CloudIO.C:195
scalar tc_
Characteristic time (used in atomization and/or breakup model)
Definition: SprayParcel.H:160
volScalarField & p
scalar mu_
Liquid dynamic viscosity [Pa.s].
Definition: SprayParcel.H:145
scalar y_
Spherical deviation.
Definition: SprayParcel.H:154
scalar KHindex_
Index for KH Breakup.
Definition: SprayParcel.H:151
Templated base class for dsmc cloud.
Definition: DSMCCloud.H:68
A primitive field of type <T> with automated input and output.
Definition: IOField.H:50
Reacing spray parcel, with added functionality for atomization and breakup.
Definition: SprayParcel.H:43