SprayParcel.H
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-2017 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 Class
25  Foam::SprayParcel
26 
27 Description
28  Reacing spray parcel, with added functionality for atomization and breakup
29 
30 \*---------------------------------------------------------------------------*/
31 
32 #ifndef SprayParcel_H
33 #define SprayParcel_H
34 
35 #include "particle.H"
36 #include "demandDrivenEntry.H"
37 
38 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
39 
40 namespace Foam
41 {
42 
43 template<class ParcelType>
44 class SprayParcel;
45 
46 template<class ParcelType>
47 Ostream& operator<<
48 (
49  Ostream&,
51 );
52 
53 /*---------------------------------------------------------------------------*\
54  Class SprayParcel Declaration
55 \*---------------------------------------------------------------------------*/
56 
57 template<class ParcelType>
58 class SprayParcel
59 :
60  public ParcelType
61 {
62  // Private data
63 
64  //- Size in bytes of the fields
65  static const std::size_t sizeofFields_;
66 
67 
68 public:
69 
70  //- Class to hold reacting particle constant properties
71  class constantProperties
72  :
73  public ParcelType::constantProperties
74  {
75  // Private data
76 
77  //- Particle initial surface tension [N/m]
79 
80  //- Particle initial dynamic viscosity [Pa.s]
82 
83 
84  public:
85 
86  // Constructors
87 
88  //- Null constructor
90 
91  //- Copy constructor
93 
94  //- Construct from dictionary
95  constantProperties(const dictionary& parentDict);
96 
97  //- Construct from components
99  (
100  const label parcelTypeId,
101  const scalar rhoMin,
102  const scalar rho0,
103  const scalar minParcelMass,
104  const scalar youngsModulus,
105  const scalar poissonsRatio,
106  const scalar T0,
107  const scalar TMin,
108  const scalar TMax,
109  const scalar Cp0,
110  const scalar epsilon0,
111  const scalar f0,
112  const scalar Pr,
113  const scalar pMin,
114  const Switch& constantVolume,
115  const scalar sigma0,
116  const scalar mu0
117  );
118 
119 
120  // Access
121 
122  //- Return const access to the initial surface tension
123  inline scalar sigma0() const;
124 
125  //- Return const access to the initial dynamic viscosity
126  inline scalar mu0() const;
127  };
128 
129 
130 protected:
131 
132  // Protected data
133 
134  // Spray parcel properties
135 
136  //- Initial droplet diameter
137  scalar d0_;
138 
139  //- Injection position
141 
142  //- Liquid surface tension [N/m]
143  scalar sigma_;
144 
145  //- Liquid dynamic viscosity [Pa.s]
146  scalar mu_;
147 
148  //- Part of liquid core ( >0.5=liquid, <0.5=droplet )
149  scalar liquidCore_;
150 
151  //- Index for KH Breakup
152  scalar KHindex_;
153 
154  //- Spherical deviation
155  scalar y_;
156 
157  //- Rate of change of spherical deviation
158  scalar yDot_;
159 
160  //- Characteristic time (used in atomization and/or breakup model)
161  scalar tc_;
162 
163  //- Stripped parcel mass due to breakup
164  scalar ms_;
165 
166  //- Injected from injector (needed e.g. for calculating distance
167  // from injector)
168  scalar injector_;
169 
170  //- Momentum relaxation time (needed for calculating parcel acc.)
171  scalar tMom_;
172 
173  //- Passive scalar (extra variable to be defined by user)
174  scalar user_;
175 
176 
177 public:
178 
179  // Static data members
180 
181  //- Runtime type information
182  TypeName("SprayParcel");
183 
184 
185  // Constructors
186 
187  //- Construct from mesh, coordinates and topology
188  // Other properties initialised as null
189  inline SprayParcel
190  (
191  const polyMesh& mesh,
192  const barycentric& coordinates,
193  const label celli,
194  const label tetFacei,
195  const label tetPti
196  );
197 
198  //- Construct from a position and a cell, searching for the rest of the
199  // required topology. Other properties are initialised as null.
200  inline SprayParcel
201  (
202  const polyMesh& mesh,
203  const vector& position,
204  const label celli
205  );
206 
207  //- Construct from components
208  inline SprayParcel
209  (
210  const polyMesh& mesh,
211  const barycentric& coordinates,
212  const label celli,
213  const label tetFacei,
214  const label tetPti,
215  const label typeId,
216  const scalar nParticle0,
217  const scalar d0,
218  const scalar dTarget0,
219  const vector& U0,
220  const vector& f0,
221  const vector& angularMomentum0,
222  const vector& torque0,
223  const scalarField& Y0,
224  const scalar liquidCore,
225  const scalar KHindex,
226  const scalar y,
227  const scalar yDot,
228  const scalar tc,
229  const scalar ms,
230  const scalar injector,
231  const scalar tMom,
232  const scalar user,
233  const typename ParcelType::constantProperties& constProps
234  );
235 
236  //- Construct from Istream
238  (
239  const polyMesh& mesh,
240  Istream& is,
241  bool readFields = true
242  );
243 
244  //- Construct as a copy
246  (
247  const SprayParcel& p,
248  const polyMesh& mesh
249  );
250 
251  //- Construct as a copy
252  SprayParcel(const SprayParcel& p);
253 
254  //- Construct and return a (basic particle) clone
255  virtual autoPtr<particle> clone() const
256  {
257  return autoPtr<particle>(new SprayParcel<ParcelType>(*this));
258  }
259 
260  //- Construct and return a (basic particle) clone
261  virtual autoPtr<particle> clone(const polyMesh& mesh) const
262  {
263  return autoPtr<particle>
264  (
265  new SprayParcel<ParcelType>(*this, mesh)
266  );
267  }
268 
269  //- Factory class to read-construct particles used for
270  // parallel transfer
271  class iNew
272  {
273  const polyMesh& mesh_;
274 
275  public:
277  iNew(const polyMesh& mesh)
278  :
279  mesh_(mesh)
280  {}
282  autoPtr<SprayParcel<ParcelType>> operator()(Istream& is) const
283  {
285  (
286  new SprayParcel<ParcelType>(mesh_, is, true)
287  );
288  }
289  };
290 
291 
292  // Member Functions
293 
294  // Access
295 
296  //- Return const access to initial droplet diameter
297  inline scalar d0() const;
298 
299  //- Return const access to initial droplet position
300  inline const vector& position0() const;
301 
302  //- Return const access to the liquid surface tension
303  inline scalar sigma() const;
304 
305  //- Return const access to the liquid dynamic viscosity
306  inline scalar mu() const;
307 
308  //- Return const access to liquid core
309  inline scalar liquidCore() const;
310 
311  //- Return const access to Kelvin-Helmholtz breakup index
312  inline scalar KHindex() const;
313 
314  //- Return const access to spherical deviation
315  inline scalar y() const;
316 
317  //- Return const access to rate of change of spherical deviation
318  inline scalar yDot() const;
319 
320  //- Return const access to atomization characteristic time
321  inline scalar tc() const;
322 
323  //- Return const access to stripped parcel mass
324  inline scalar ms() const;
325 
326  //- Return const access to injector id
327  inline scalar injector() const;
328 
329  //- Return const access to momentum relaxation time
330  inline scalar tMom() const;
331 
332  //- Return const access to passive user scalar
333  inline scalar user() const;
334 
335 
336  // Edit
337 
338  //- Return access to initial droplet diameter
339  inline scalar& d0();
340 
341  //- Return access to initial droplet position
342  inline vector& position0();
343 
344  //- Return access to the liquid surface tension
345  inline scalar& sigma();
346 
347  //- Return access to the liquid dynamic viscosity
348  inline scalar& mu();
349 
350  //- Return access to liquid core
351  inline scalar& liquidCore();
352 
353  //- Return access to Kelvin-Helmholtz breakup index
354  inline scalar& KHindex();
355 
356  //- Return access to spherical deviation
357  inline scalar& y();
358 
359  //- Return access to rate of change of spherical deviation
360  inline scalar& yDot();
361 
362  //- Return access to atomization characteristic time
363  inline scalar& tc();
364 
365  //- Return access to stripped parcel mass
366  inline scalar& ms();
367 
368  //- Return access to injector id
369  inline scalar& injector();
370 
371  //- Return access to momentum relaxation time
372  inline scalar& tMom();
373 
374  //- Return access to passive user scalar
375  inline scalar& user();
376 
377 
378  // Main calculation loop
379 
380  //- Set cell values
381  template<class TrackData>
382  void setCellValues
383  (
384  TrackData& td,
385  const scalar dt,
386  const label celli
387  );
388 
389  //- Correct parcel properties according to atomization model
390  template<class TrackData>
391  void calcAtomization
392  (
393  TrackData& td,
394  const scalar dt,
395  const label celli
396  );
397 
398  //- Correct parcel properties according to breakup model
399  template<class TrackData>
400  void calcBreakup
401  (
402  TrackData& td,
403  const scalar dt,
404  const label celli
405  );
406 
407  //- Correct cell values using latest transfer information
408  template<class TrackData>
410  (
411  TrackData& td,
412  const scalar dt,
413  const label celli
414  );
415 
416  //- Correct surface values due to emitted species
417  template<class TrackData>
419  (
420  TrackData& td,
421  const label celli,
422  const scalar T,
423  const scalarField& Cs,
424  scalar& rhos,
425  scalar& mus,
426  scalar& Pr,
427  scalar& kappa
428  );
429 
430  //- Update parcel properties over the time interval
431  template<class TrackData>
432  void calc
433  (
434  TrackData& td,
435  const scalar dt,
436  const label celli
437  );
438 
439  //- Calculate the chi-factor for flash-boiling for the
440  // atomization model
441  template<class TrackData>
442  scalar chi
443  (
444  TrackData& td,
445  const scalarField& X
446  ) const;
447 
448  //- Solve the TAB equation
449  template<class TrackData>
450  void solveTABEq
451  (
452  TrackData& td,
453  const scalar dt
454  );
455 
456 
457  // I-O
458 
459  //- Read
460  template<class CloudType, class CompositionType>
461  static void readFields
462  (
463  CloudType& c,
464  const CompositionType& compModel
465  );
466 
467  //- Read - no composition
468  template<class CloudType>
469  static void readFields(CloudType& c);
470 
471  //- Write
472  template<class CloudType, class CompositionType>
473  static void writeFields
474  (
475  const CloudType& c,
476  const CompositionType& compModel
477  );
478 
479  //- Write - composition supplied
480  template<class CloudType>
481  static void writeFields(const CloudType& c);
482 
483 
484  // Ostream Operator
485 
486  friend Ostream& operator<< <ParcelType>
487  (
488  Ostream&,
490  );
491 };
492 
493 
494 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
495 
496 } // End namespace Foam
497 
498 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
499 
500 #include "SprayParcelI.H"
501 
502 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
503 
504 #ifdef NoRepository
505  #include "SprayParcel.C"
506 #endif
507 
508 
509 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
510 
511 #endif
512 
513 // ************************************************************************* //
dimensionedScalar Pr("Pr", dimless, laminarTransport)
scalar mu0() const
Return const access to the initial dynamic viscosity.
Definition: SprayParcelI.H:233
Factory class to read-construct particles used for.
Definition: SprayParcel.H:270
void solveTABEq(TrackData &td, const scalar dt)
Solve the TAB equation.
Definition: SprayParcel.C:364
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
A list of keyword definitions, which are a keyword followed by any number of values (e...
Definition: dictionary.H:137
scalar tMom_
Momentum relaxation time (needed for calculating parcel acc.)
Definition: SprayParcel.H:170
scalar KHindex() const
Return const access to Kelvin-Helmholtz breakup index.
Definition: SprayParcelI.H:277
scalar tMom() const
Return const access to momentum relaxation time.
Definition: SprayParcelI.H:319
An Istream is an abstract base class for all input systems (streams, files, token lists etc)...
Definition: Istream.H:57
scalar sigma0() const
Return const access to the initial surface tension.
Definition: SprayParcelI.H:225
PtrList< dimensionedScalar > rhoMin(fluidRegions.size())
scalar liquidCore() const
Return const access to liquid core.
Definition: SprayParcelI.H:270
void calc(TrackData &td, const scalar dt, const label celli)
Update parcel properties over the time interval.
Definition: SprayParcel.C:61
bool cp(const fileName &src, const fileName &dst, const bool followLink=true)
Copy, recursively if necessary, the source to the destination.
Definition: POSIX.C:739
A simple wrapper around bool so that it can be read as a word: true/false, on/off, yes/no, y/n, t/f, or none.
Definition: Switch.H:60
scalar rho0
scalar chi(TrackData &td, const scalarField &X) const
Calculate the chi-factor for flash-boiling for the.
scalar ms() const
Return const access to stripped parcel mass.
Definition: SprayParcelI.H:305
const dimensionedScalar epsilon0
Electric constant: default SI units: [F/m].
const dimensionedScalar kappa
Coulomb constant: default SI units: [N.m2/C2].
void calcAtomization(TrackData &td, const scalar dt, const label celli)
Correct parcel properties according to atomization model.
Definition: SprayParcel.C:150
scalar liquidCore_
Part of liquid core ( >0.5=liquid, <0.5=droplet )
Definition: SprayParcel.H:148
SprayParcel(const polyMesh &mesh, const barycentric &coordinates, const label celli, const label tetFacei, const label tetPti)
Construct from mesh, coordinates and topology.
Definition: SprayParcelI.H:108
scalar y() const
Return const access to spherical deviation.
Definition: SprayParcelI.H:284
scalarList Y0(nSpecie, 0.0)
void correctSurfaceValues(TrackData &td, const label celli, const scalar T, const scalarField &Cs, scalar &rhos, scalar &mus, scalar &Pr, scalar &kappa)
Correct surface values due to emitted species.
scalar d0() const
Return const access to initial droplet diameter.
Definition: SprayParcelI.H:242
dynamicFvMesh & mesh
scalar mu() const
Return const access to the liquid dynamic viscosity.
Definition: SprayParcelI.H:263
scalar tc() const
Return const access to atomization characteristic time.
Definition: SprayParcelI.H:298
scalar injector_
Injected from injector (needed e.g. for calculating distance.
Definition: SprayParcel.H:167
Class to hold reacting particle constant properties.
Definition: SprayParcel.H:70
scalar user_
Passive scalar (extra variable to be defined by user)
Definition: SprayParcel.H:173
static void readFields(CloudType &c, const CompositionType &compModel)
Read.
void cellValueSourceCorrection(TrackData &td, const scalar dt, const label celli)
Correct cell values using latest transfer information.
Definition: SprayParcel.C:48
scalar yDot() const
Return const access to rate of change of spherical deviation.
Definition: SprayParcelI.H:291
An Ostream is an abstract base class for all output systems (streams, files, token lists...
Definition: Ostream.H:53
scalar user() const
Return const access to passive user scalar.
Definition: SprayParcelI.H:326
constantProperties()
Null constructor.
Definition: SprayParcelI.H:29
scalar sigma_
Liquid surface tension [N/m].
Definition: SprayParcel.H:142
void T(FieldField< Field, Type > &f1, const FieldField< Field, Type > &f2)
dimensionedScalar pMin("pMin", dimPressure, fluid)
static void writeFields(const CloudType &c, const CompositionType &compModel)
Write.
vector position0_
Injection position.
Definition: SprayParcel.H:139
PtrList< coordinateSystem > coordinates(solidRegions.size())
void calcBreakup(TrackData &td, const scalar dt, const label celli)
Correct parcel properties according to breakup model.
Definition: SprayParcel.C:217
scalar injector() const
Return const access to injector id.
Definition: SprayParcelI.H:312
const vector & position0() const
Return const access to initial droplet position.
Definition: SprayParcelI.H:249
scalar yDot_
Rate of change of spherical deviation.
Definition: SprayParcel.H:157
const dimensionedScalar c
Speed of light in a vacuum.
scalar ms_
Stripped parcel mass due to breakup.
Definition: SprayParcel.H:163
scalar d0_
Initial droplet diameter.
Definition: SprayParcel.H:136
An auto-pointer similar to the STL auto_ptr but with automatic casting to a reference to the type and...
Definition: PtrList.H:52
virtual autoPtr< particle > clone() const
Construct and return a (basic particle) clone.
Definition: SprayParcel.H:254
Mesh consisting of general polyhedral cells.
Definition: polyMesh.H:74
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
TypeName("SprayParcel")
Runtime type information.
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:69
scalar sigma() const
Return const access to the liquid surface tension.
Definition: SprayParcelI.H:256
void setCellValues(TrackData &td, const scalar dt, const label celli)
Set cell values.
Definition: SprayParcel.C:35
Namespace for OpenFOAM.
Reacing spray parcel, with added functionality for atomization and breakup.
Definition: SprayParcel.H:43