All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Pages
ReactingParcel.H
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-2020 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::ReactingParcel
26 
27 Description
28  Reacting parcel class with one/two-way coupling with the continuous
29  phase.
30 
31 SourceFiles
32  ReactingParcelI.H
33  ReactingParcel.C
34  ReactingParcelIO.C
35 
36 \*---------------------------------------------------------------------------*/
37 
38 #ifndef ReactingParcel_H
39 #define ReactingParcel_H
40 
41 #include "particle.H"
42 #include "SLGThermo.H"
43 #include "demandDrivenEntry.H"
44 
45 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
46 
47 namespace Foam
48 {
49 
50 template<class ParcelType>
51 class ReactingParcel;
52 
53 template<class ParcelType>
54 Ostream& operator<<
55 (
56  Ostream&,
58 );
59 
60 
61 /*---------------------------------------------------------------------------*\
62  Class ReactingParcel Declaration
63 \*---------------------------------------------------------------------------*/
64 
65 template<class ParcelType>
66 class ReactingParcel
67 :
68  public ParcelType
69 {
70  // Private Data
71 
72  //- Size in bytes of the fields
73  static const std::size_t sizeofFields_;
74 
75 
76 public:
77 
78  //- Class to hold reacting parcel constant properties
79  class constantProperties
80  :
81  public ParcelType::constantProperties
82  {
83  // Private Data
84 
85  //- Minimum pressure [Pa]
87 
88  //- Constant volume flag - e.g. during mass transfer
89  demandDrivenEntry<bool> constantVolume_;
90 
91 
92  public:
93 
94  // Constructors
95 
96  //- Null constructor
98 
99  //- Copy constructor
101 
102  //- Construct from dictionary
103  constantProperties(const dictionary& parentDict);
104 
105 
106  // Access
107 
108  //- Return const access to the minimum pressure
109  inline scalar pMin() const;
110 
111  //- Return const access to the constant volume flag
112  inline bool constantVolume() const;
113  };
114 
116  class trackingData
117  :
118  public ParcelType::trackingData
119  {
120  private:
121 
122  // Private Data
123 
124  // Interpolators for continuous phase fields
125 
126  //- Interpolator for continuous phase pressure field
128 
129 
130  // Cached continuous phase properties
131 
132  //- Pressure [Pa]
133  scalar pc_;
134 
135 
136  public:
137 
138  // Constructors
139 
140  //- Construct from components
141  template<class TrackCloudType>
142  inline trackingData(const TrackCloudType& cloud);
143 
144 
145  // Member Functions
146 
147  //- Return const access to the interpolator for continuous phase
148  // pressure field
149  inline const interpolation<scalar>& pInterp() const;
150 
151  //- Return the continuous phase pressure
152  inline scalar pc() const;
153 
154  //- Access the continuous phase pressure
155  inline scalar& pc();
156  };
157 
158 
159 protected:
160 
161  // Protected data
162 
163  // Parcel properties
164 
165  //- Initial mass [kg]
166  scalar mass0_;
167 
168  //- Mass fractions of mixture []
169  scalarField Y_;
170 
171 
172  // Protected Member Functions
173 
174  //- Calculate Phase change
175  template<class TrackCloudType>
176  void calcPhaseChange
177  (
178  TrackCloudType& cloud,
179  trackingData& td,
180  const scalar dt, // timestep
181  const scalar Re, // Reynolds number
182  const scalar Pr, // Prandtl number
183  const scalar Ts, // Surface temperature
184  const scalar nus, // Surface kinematic viscosity
185  const scalar d, // diameter
186  const scalar T, // temperature
187  const scalar mass, // mass
188  const label idPhase, // id of phase involved in phase change
189  const scalar YPhase, // total mass fraction
190  const scalarField& YComponents, // component mass fractions
191  scalarField& dMassPC, // mass transfer - local to parcel
192  scalar& Sh, // explicit parcel enthalpy source
193  scalar& N, // flux of species emitted from parcel
194  scalar& NCpW, // sum of N*Cp*W of emission species
195  scalarField& Cs // carrier conc. of emission species
196  );
197 
198  //- Update mass fraction
199  scalar updateMassFraction
200  (
201  const scalar mass0,
202  const scalarField& dMass,
203  scalarField& Y
204  ) const;
205 
206 
207 public:
208 
209  // Static Data Members
210 
211  //- Runtime type information
212  TypeName("ReactingParcel");
213 
214  //- String representation of properties
216  (
217  ParcelType,
218  " mass0"
219  + " nPhases(Y1..YN)"
220  );
221 
222 
223  // Constructors
224 
225  //- Construct from mesh, coordinates and topology
226  // Other properties initialised as null
227  inline ReactingParcel
228  (
229  const polyMesh& mesh,
230  const barycentric& coordinates,
231  const label celli,
232  const label tetFacei,
233  const label tetPti
234  );
235 
236  //- Construct from a position and a cell, searching for the rest of the
237  // required topology. Other properties are initialised as null.
238  inline ReactingParcel
239  (
240  const polyMesh& mesh,
241  const vector& position,
242  const label celli
243  );
244 
245  //- Construct from components
246  inline ReactingParcel
247  (
248  const polyMesh& mesh,
249  const barycentric& coordinates,
250  const label celli,
251  const label tetFacei,
252  const label tetPti,
253  const label typeId,
254  const scalar nParticle0,
255  const scalar d0,
256  const scalar dTarget0,
257  const vector& U0,
258  const vector& f0,
259  const vector& angularMomentum0,
260  const vector& torque0,
261  const scalarField& Y0,
262  const constantProperties& constProps
263  );
264 
265  //- Construct from Istream
267  (
268  const polyMesh& mesh,
269  Istream& is,
270  bool readFields = true
271  );
272 
273  //- Construct as a copy
275  (
276  const ReactingParcel& p,
277  const polyMesh& mesh
278  );
279 
280  //- Construct as a copy
282 
283  //- Construct and return a (basic particle) clone
284  virtual autoPtr<particle> clone() const
285  {
287  }
288 
289  //- Construct and return a (basic particle) clone
290  virtual autoPtr<particle> clone(const polyMesh& mesh) const
291  {
292  return autoPtr<particle>
293  (
294  new ReactingParcel<ParcelType>(*this, mesh)
295  );
296  }
297 
298  //- Factory class to read-construct particles used for
299  // parallel transfer
300  class iNew
301  {
302  const polyMesh& mesh_;
303 
304  public:
306  iNew(const polyMesh& mesh)
307  :
308  mesh_(mesh)
309  {}
311  autoPtr<ReactingParcel<ParcelType>> operator()(Istream& is) const
312  {
314  (
315  new ReactingParcel<ParcelType>(mesh_, is, true)
316  );
317  }
318  };
319 
320 
321  // Member Functions
322 
323  // Access
324 
325  //- Return const access to initial mass [kg]
326  inline scalar mass0() const;
327 
328  //- Return const access to mass fractions of mixture []
329  inline const scalarField& Y() const;
330 
331 
332  // Edit
333 
334  //- Return access to initial mass [kg]
335  inline scalar& mass0();
336 
337  //- Return access to mass fractions of mixture []
338  inline scalarField& Y();
339 
340 
341  // Main calculation loop
342 
343  //- Set cell values
344  template<class TrackCloudType>
345  void setCellValues(TrackCloudType& cloud, trackingData& td);
346 
347  //- Correct cell values using latest transfer information
348  template<class TrackCloudType>
350  (
351  TrackCloudType& cloud,
352  trackingData& td,
353  const scalar dt
354  );
355 
356  //- Correct surface values due to emitted species
357  template<class TrackCloudType>
359  (
360  TrackCloudType& cloud,
361  trackingData& td,
362  const scalar T,
363  const scalarField& Cs,
364  scalar& rhos,
365  scalar& mus,
366  scalar& Prs,
367  scalar& kappas
368  );
369 
370  //- Update parcel properties over the time interval
371  template<class TrackCloudType>
372  void calc
373  (
374  TrackCloudType& cloud,
375  trackingData& td,
376  const scalar dt
377  );
378 
379 
380  // I-O
381 
382  //- Read
383  template<class CloudType, class CompositionType>
384  static void readFields
385  (
386  CloudType& c,
387  const CompositionType& compModel
388  );
389 
390  //- Read - no composition
391  template<class CloudType>
392  static void readFields(CloudType& c);
393 
394  //- Write
395  template<class CloudType, class CompositionType>
396  static void writeFields
397  (
398  const CloudType& c,
399  const CompositionType& compModel
400  );
401 
402  //- Write - composition supplied
403  template<class CloudType>
404  static void writeFields(const CloudType& c);
405 
406 
407  // Ostream Operator
408 
409  friend Ostream& operator<< <ParcelType>
410  (
411  Ostream&,
413  );
414 };
415 
416 
417 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
418 
419 } // End namespace Foam
420 
421 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
422 
423 #include "ReactingParcelI.H"
425 
426 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
427 
428 #ifdef NoRepository
429  #include "ReactingParcel.C"
430 #endif
431 
432 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
433 
434 #endif
435 
436 // ************************************************************************* //
Class to hold reacting parcel constant properties.
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:158
void calcPhaseChange(TrackCloudType &cloud, trackingData &td, const scalar dt, const scalar Re, const scalar Pr, const scalar Ts, const scalar nus, const scalar d, const scalar T, const scalar mass, const label idPhase, const scalar YPhase, const scalarField &YComponents, scalarField &dMassPC, scalar &Sh, scalar &N, scalar &NCpW, scalarField &Cs)
Calculate Phase change.
Factory class to read-construct particles used for.
An Istream is an abstract base class for all input systems (streams, files, token lists etc)...
Definition: Istream.H:57
scalar updateMassFraction(const scalar mass0, const scalarField &dMass, scalarField &Y) const
Update mass fraction.
scalarField Y_
Mass fractions of mixture [].
bool cp(const fileName &src, const fileName &dst, const bool followLink=true)
Copy, recursively if necessary, the source to the destination.
Definition: POSIX.C:753
ReactingParcel(const polyMesh &mesh, const barycentric &coordinates, const label celli, const label tetFacei, const label tetPti)
Construct from mesh, coordinates and topology.
bool constantVolume() const
Return const access to the constant volume flag.
static void writeFields(const CloudType &c, const CompositionType &compModel)
Write.
const dimensionedScalar & c
Speed of light in a vacuum.
void cellValueSourceCorrection(TrackCloudType &cloud, trackingData &td, const scalar dt)
Correct cell values using latest transfer information.
scalarList Y0(nSpecie, 0.0)
static void readFields(CloudType &c, const CompositionType &compModel)
Read.
dynamicFvMesh & mesh
scalar pMin() const
Return const access to the minimum pressure.
A cloud is a collection of lagrangian particles.
Definition: cloud.H:51
void correctSurfaceValues(TrackCloudType &cloud, trackingData &td, const scalar T, const scalarField &Cs, scalar &rhos, scalar &mus, scalar &Prs, scalar &kappas)
Correct surface values due to emitted species.
TypeName("ReactingParcel")
Runtime type information.
An Ostream is an abstract base class for all output systems (streams, files, token lists...
Definition: Ostream.H:54
AddToPropertyList(ParcelType, " mass0"+" nPhases(Y1..YN)")
String representation of properties.
void T(FieldField< Field, Type > &f1, const FieldField< Field, Type > &f2)
virtual autoPtr< particle > clone() const
Construct and return a (basic particle) clone.
PtrList< coordinateSystem > coordinates(solidRegions.size())
scalar mass0() const
Return const access to initial mass [kg].
void setCellValues(TrackCloudType &cloud, trackingData &td)
Set cell values.
label N
Definition: createFields.H:22
An auto-pointer similar to the STL auto_ptr but with automatic casting to a reference to the type and...
Definition: PtrList.H:52
Mesh consisting of general polyhedral cells.
Definition: polyMesh.H:74
Reacting parcel class with one/two-way coupling with the continuous phase.
volScalarField & p
scalar mass0_
Initial mass [kg].
Templated base class for dsmc cloud.
Definition: DSMCCloud.H:69
scalarField Re(const UList< complex > &cf)
Definition: complexFields.C:97
const scalarField & Y() const
Return const access to mass fractions of mixture [].
Namespace for OpenFOAM.
void calc(TrackCloudType &cloud, trackingData &td, const scalar dt)
Update parcel properties over the time interval.