ReactingMultiphaseParcel.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-2018 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::ReactingMultiphaseParcel
26 
27 Description
28  Multiphase variant of the reacting parcel class with one/two-way coupling
29  with the continuous phase.
30 
31 SourceFiles
32  ReactingMultiphaseParcelI.H
33  ReactingMultiphaseParcel.C
34  ReactingMultiphaseParcelIO.C
35 
36 \*---------------------------------------------------------------------------*/
37 
38 #ifndef ReactingMultiphaseParcel_H
39 #define ReactingMultiphaseParcel_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>
52 
53 template<class ParcelType>
54 Ostream& operator<<
55 (
56  Ostream&,
58 );
59 
60 /*---------------------------------------------------------------------------*\
61  Class ReactingMultiphaseParcel Declaration
62 \*---------------------------------------------------------------------------*/
63 
64 template<class ParcelType>
66 :
67  public ParcelType
68 {
69  // Private data
70 
71  //- Size in bytes of the fields
72  static const std::size_t sizeofFields_;
73 
74 
75 public:
76 
77  // IDs of phases in ReacingParcel phase list (Y)
78 
79  static const label GAS;
80  static const label LIQ;
81  static const label SLD;
82 
83 
84  //- Class to hold reacting multiphase particle constant properties
85  class constantProperties
86  :
87  public ParcelType::constantProperties
88  {
89  // Private data
90 
91  //- Devolatilisation activation temperature [K]
93 
94  //- Latent heat of devolatilisation [J/kg]
96 
97  //- Fraction of enthalpy retained by parcel due to surface
98  // reactions
99  demandDrivenEntry<scalar> hRetentionCoeff_;
100 
101 
102  public:
103 
104  // Constructors
105 
106  //- Null constructor
108 
109  //- Copy constructor
111 
112  //- Construct from dictionary
113  constantProperties(const dictionary& parentDict);
114 
115 
116  // Access
117 
118  //- Return const access to the devolatilisation temperature
119  inline scalar TDevol() const;
120 
121  //- Return const access to the latent heat of devolatilisation
122  inline scalar LDevol() const;
123 
124  //- Return const access to the fraction of enthalpy retained by
125  // parcel due to surface reactions
126  inline scalar hRetentionCoeff() const;
127  };
128 
129 
130  //- Use base tracking data
131  typedef typename ParcelType::trackingData trackingData;
132 
133 
134 private:
135 
136  // Private Member Functions
137 
138  //- Return the mixture effective specific heat capacity
139  template<class TrackCloudType>
140  scalar CpEff
141  (
142  TrackCloudType& cloud,
143  trackingData& td,
144  const scalar p,
145  const scalar T,
146  const label idG,
147  const label idL,
148  const label idS
149  ) const;
150 
151  //- Return the mixture effective sensible enthalpy
152  template<class TrackCloudType>
153  scalar HsEff
154  (
155  TrackCloudType& cloud,
156  trackingData& td,
157  const scalar p,
158  const scalar T,
159  const label idG,
160  const label idL,
161  const label idS
162  ) const;
163 
164  //- Return the mixture effective latent heat
165  template<class TrackCloudType>
166  scalar LEff
167  (
168  TrackCloudType& cloud,
169  trackingData& td,
170  const scalar p,
171  const scalar T,
172  const label idG,
173  const label idL,
174  const label idS
175  ) const;
176 
177  //- Update the mass fractions (Y, YGas, YLiquid, YSolid)
178  scalar updateMassFractions
179  (
180  const scalar mass0,
181  const scalarField& dMassGas,
182  const scalarField& dMassLiquid,
183  const scalarField& dMassSolid
184  );
185 
186 
187 protected:
188 
189  // Protected data
190 
191  // Parcel properties
192 
193  //- Mass fractions of gases []
195 
196  //- Mass fractions of liquids []
198 
199  //- Mass fractions of solids []
201 
202  //- Flag to identify if the particle can devolatilise and combust
203  // Combustion possible only after volatile content falls below
204  // threshold value. States include:
205  // 0 = can devolatilise, cannot combust but can change
206  // 1 = can devolatilise, can combust
207  // -1 = cannot devolatilise or combust, and cannot change
209 
210 
211  // Protected Member Functions
212 
213  //- Calculate Devolatilisation
214  template<class TrackCloudType>
216  (
217  TrackCloudType& cloud,
218  trackingData& td,
219  const scalar dt, // timestep
220  const scalar age, // age
221  const scalar Ts, // surface temperature
222  const scalar d, // diameter
223  const scalar T, // temperature
224  const scalar mass, // mass
225  const scalar mass0, // mass (initial on injection)
226  const scalarField& YGasEff,// gas component mass fractions
227  const scalarField& YLiquidEff,// liquid component mass fractions
228  const scalarField& YSolidEff,// solid component mass fractions
229  label& canCombust, // 'can combust' flag
230  scalarField& dMassDV, // mass transfer - local to particle
231  scalar& Sh, // explicit particle enthalpy source
232  scalar& N, // flux of species emitted from particle
233  scalar& NCpW, // sum of N*Cp*W of emission species
234  scalarField& Cs // carrier conc. of emission species
235  ) const;
236 
237  //- Calculate surface reactions
238  template<class TrackCloudType>
240  (
241  TrackCloudType& cloud,
242  trackingData& td,
243  const scalar dt, // timestep
244  const scalar d, // diameter
245  const scalar T, // temperature
246  const scalar mass, // mass
247  const label canCombust, // 'can combust' flag
248  const scalar N, // flux of species emitted from particle
249  const scalarField& YMix, // mixture mass fractions
250  const scalarField& YGas, // gas-phase mass fractions
251  const scalarField& YLiquid,// liquid-phase mass fractions
252  const scalarField& YSolid, // solid-phase mass fractions
253  scalarField& dMassSRGas, // gas-phase mass transfer - local
254  scalarField& dMassSRLiquid,// liquid-phase mass transfer - local
255  scalarField& dMassSRSolid, // solid-phase mass transfer - local
256  scalarField& dMassSRCarrier, // carrier phase mass transfer
257  scalar& Sh, // explicit particle enthalpy source
258  scalar& dhsTrans // sensible enthalpy transfer to carrier
259  ) const;
260 
261 
262 public:
263 
264  // Static data members
265 
266  //- Runtime type information
267  TypeName("ReactingMultiphaseParcel");
268 
269  //- String representation of properties
271  (
272  ParcelType,
273  " nGas(Y1..YN)"
274  + " nLiquid(Y1..YN)"
275  + " nSolid(Y1..YN)"
276  );
277 
278 
279  // Constructors
280 
281  //- Construct from mesh, position and topology
282  // Other properties initialised as null
284  (
285  const polyMesh& mesh,
286  const barycentric& coordinates,
287  const label celli,
288  const label tetFacei,
289  const label tetPti
290  );
291 
292  //- Construct from a position and a cell, searching for the rest of the
293  // required topology. Other properties are initialised as null.
295  (
296  const polyMesh& mesh,
297  const vector& position,
298  const label celli
299  );
300 
301  //- Construct from components
303  (
304  const polyMesh& mesh,
305  const barycentric& coordinates,
306  const label celli,
307  const label tetFacei,
308  const label tetPti,
309  const label typeId,
310  const scalar nParticle0,
311  const scalar d0,
312  const scalar dTarget0,
313  const vector& U0,
314  const vector& f0,
315  const vector& angularMomentum0,
316  const vector& torque0,
317  const scalarField& Y0,
318  const scalarField& YGas0,
319  const scalarField& YLiquid0,
320  const scalarField& YSolid0,
321  const constantProperties& constProps
322  );
323 
324  //- Construct from Istream
326  (
327  const polyMesh& mesh,
328  Istream& is,
329  bool readFields = true
330  );
331 
332  //- Construct as a copy
334 
335  //- Construct as a copy
337  (
339  const polyMesh& mesh
340  );
341 
342  //- Construct and return a (basic particle) clone
343  virtual autoPtr<particle> clone() const
344  {
345  return autoPtr<particle>(new ReactingMultiphaseParcel(*this));
346  }
347 
348  //- Construct and return a (basic particle) clone
349  virtual autoPtr<particle> clone(const polyMesh& mesh) const
350  {
351  return autoPtr<particle>(new ReactingMultiphaseParcel(*this, mesh));
352  }
353 
354  //- Factory class to read-construct particles used for
355  // parallel transfer
356  class iNew
357  {
358  const polyMesh& mesh_;
359 
360  public:
362  iNew(const polyMesh& mesh)
363  :
364  mesh_(mesh)
365  {}
366 
368  (
369  Istream& is
370  ) const
371  {
373  (
374  new ReactingMultiphaseParcel<ParcelType>(mesh_, is, true)
375  );
376  }
377  };
378 
379 
380  // Member Functions
381 
382  // Access
383 
384  //- Return const access to mass fractions of gases
385  inline const scalarField& YGas() const;
386 
387  //- Return const access to mass fractions of liquids
388  inline const scalarField& YLiquid() const;
389 
390  //- Return const access to mass fractions of solids
391  inline const scalarField& YSolid() const;
392 
393  //- Return const access to the canCombust flag
394  inline label canCombust() const;
395 
396 
397  // Edit
398 
399  //- Return access to mass fractions of gases
400  inline scalarField& YGas();
401 
402  //- Return access to mass fractions of liquids
403  inline scalarField& YLiquid();
404 
405  //- Return access to mass fractions of solids
406  inline scalarField& YSolid();
407 
408  //- Return access to the canCombust flag
409  inline label& canCombust();
410 
411 
412  // Main calculation loop
413 
414  //- Set cell values
415  template<class TrackCloudType>
416  void setCellValues(TrackCloudType& cloud, trackingData& td);
417 
418  //- Correct cell values using latest transfer information
419  template<class TrackCloudType>
421  (
422  TrackCloudType& cloud,
423  trackingData& td,
424  const scalar dt
425  );
426 
427  //- Update parcel properties over the time interval
428  template<class TrackCloudType>
429  void calc
430  (
431  TrackCloudType& cloud,
432  trackingData& td,
433  const scalar dt
434  );
435 
436 
437  // I-O
438 
439  //- Read
440  template<class CloudType, class CompositionType>
441  static void readFields
442  (
443  CloudType& c,
444  const CompositionType& compModel
445  );
446 
447  //- Read - no composition
448  template<class CloudType>
449  static void readFields(CloudType& c);
450 
451  //- Write
452  template<class CloudType, class CompositionType>
453  static void writeFields
454  (
455  const CloudType& c,
456  const CompositionType& compModel
457  );
458 
459  //- Read - composition supplied
460  template<class CloudType>
461  static void writeFields(const CloudType& c);
462 
463 
464  // Ostream Operator
465 
466  friend Ostream& operator<< <ParcelType>
467  (
468  Ostream&,
470  );
471 };
472 
473 
474 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
475 
476 } // End namespace Foam
477 
478 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
479 
481 
482 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
483 
484 #ifdef NoRepository
485  #include "ReactingMultiphaseParcel.C"
486 #endif
487 
488 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
489 
490 #endif
491 
492 // ************************************************************************* //
label canCombust() const
Return const access to the canCombust flag.
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
static void writeFields(const CloudType &c, const CompositionType &compModel)
Write.
A list of keyword definitions, which are a keyword followed by any number of values (e...
Definition: dictionary.H:137
label canCombust_
Flag to identify if the particle can devolatilise and combust.
An Istream is an abstract base class for all input systems (streams, files, token lists etc)...
Definition: Istream.H:57
void calcDevolatilisation(TrackCloudType &cloud, trackingData &td, const scalar dt, const scalar age, const scalar Ts, const scalar d, const scalar T, const scalar mass, const scalar mass0, const scalarField &YGasEff, const scalarField &YLiquidEff, const scalarField &YSolidEff, label &canCombust, scalarField &dMassDV, scalar &Sh, scalar &N, scalar &NCpW, scalarField &Cs) const
Calculate Devolatilisation.
bool cp(const fileName &src, const fileName &dst, const bool followLink=true)
Copy, recursively if necessary, the source to the destination.
Definition: POSIX.C:737
scalar LDevol() const
Return const access to the latent heat of devolatilisation.
void calc(TrackCloudType &cloud, trackingData &td, const scalar dt)
Update parcel properties over the time interval.
scalarList Y0(nSpecie, 0.0)
virtual autoPtr< particle > clone() const
Construct and return a (basic particle) clone.
static void readFields(CloudType &c, const CompositionType &compModel)
Read.
ReactingMultiphaseParcel(const polyMesh &mesh, const barycentric &coordinates, const label celli, const label tetFacei, const label tetPti)
Construct from mesh, position and topology.
ParcelType::trackingData trackingData
Use base tracking data.
scalarField YGas_
Mass fractions of gases [].
scalar hRetentionCoeff() const
Return const access to the fraction of enthalpy retained by.
dynamicFvMesh & mesh
const scalarField & YGas() const
Return const access to mass fractions of gases.
A cloud is a collection of lagrangian particles.
Definition: cloud.H:51
scalar TDevol() const
Return const access to the devolatilisation temperature.
An Ostream is an abstract base class for all output systems (streams, files, token lists...
Definition: Ostream.H:53
AddToPropertyList(ParcelType, " nGas(Y1..YN)"+" nLiquid(Y1..YN)"+" nSolid(Y1..YN)")
String representation of properties.
void calcSurfaceReactions(TrackCloudType &cloud, trackingData &td, const scalar dt, const scalar d, const scalar T, const scalar mass, const label canCombust, const scalar N, const scalarField &YMix, const scalarField &YGas, const scalarField &YLiquid, const scalarField &YSolid, scalarField &dMassSRGas, scalarField &dMassSRLiquid, scalarField &dMassSRSolid, scalarField &dMassSRCarrier, scalar &Sh, scalar &dhsTrans) const
Calculate surface reactions.
scalarField YSolid_
Mass fractions of solids [].
const scalarField & YLiquid() const
Return const access to mass fractions of liquids.
void T(FieldField< Field, Type > &f1, const FieldField< Field, Type > &f2)
scalarField YLiquid_
Mass fractions of liquids [].
PtrList< coordinateSystem > coordinates(solidRegions.size())
const dimensionedScalar c
Speed of light in a vacuum.
friend Ostream & operator(Ostream &, const ReactingMultiphaseParcel< ParcelType > &)
Multiphase variant of the reacting parcel class with one/two-way coupling with the continuous phase...
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
volScalarField & p
TypeName("ReactingMultiphaseParcel")
Runtime type information.
Factory class to read-construct particles used for.
void cellValueSourceCorrection(TrackCloudType &cloud, trackingData &td, const scalar dt)
Correct cell values using latest transfer information.
Templated base class for dsmc cloud.
Definition: DSMCCloud.H:69
Class to hold reacting multiphase particle constant properties.
Namespace for OpenFOAM.
const scalarField & YSolid() const
Return const access to mass fractions of solids.
void setCellValues(TrackCloudType &cloud, trackingData &td)
Set cell values.