ReactingMultiphaseParcel.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-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 \*---------------------------------------------------------------------------*/
25 
27 #include "mathematicalConstants.H"
28 
29 using namespace Foam::constant::mathematical;
30 
31 // * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
32 
33 template<class ParcelType>
35 
36 template<class ParcelType>
38 
39 template<class ParcelType>
41 
42 
43 // * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * * //
44 
45 template<class ParcelType>
46 template<class TrackCloudType>
48 (
49  TrackCloudType& cloud,
50  trackingData& td,
51  const scalar p,
52  const scalar T,
53  const label idG,
54  const label idL,
55  const label idS
56 ) const
57 {
58  return
59  this->Y_[GAS]*cloud.composition().Cp(idG, YGas_, p, T)
60  + this->Y_[LIQ]*cloud.composition().Cp(idL, YLiquid_, p, T)
61  + this->Y_[SLD]*cloud.composition().Cp(idS, YSolid_, p, T);
62 }
63 
64 
65 template<class ParcelType>
66 template<class TrackCloudType>
68 (
69  TrackCloudType& cloud,
70  trackingData& td,
71  const scalar p,
72  const scalar T,
73  const label idG,
74  const label idL,
75  const label idS
76 ) const
77 {
78  return
79  this->Y_[GAS]*cloud.composition().Hs(idG, YGas_, p, T)
80  + this->Y_[LIQ]*cloud.composition().Hs(idL, YLiquid_, p, T)
81  + this->Y_[SLD]*cloud.composition().Hs(idS, YSolid_, p, T);
82 }
83 
84 
85 template<class ParcelType>
86 template<class TrackCloudType>
88 (
89  TrackCloudType& cloud,
90  trackingData& td,
91  const scalar p,
92  const scalar T,
93  const label idG,
94  const label idL,
95  const label idS
96 ) const
97 {
98  return
99  this->Y_[GAS]*cloud.composition().L(idG, YGas_, p, T)
100  + this->Y_[LIQ]*cloud.composition().L(idL, YLiquid_, p, T)
101  + this->Y_[SLD]*cloud.composition().L(idS, YSolid_, p, T);
102 }
103 
104 
105 template<class ParcelType>
107 (
108  const scalar mass0,
109  const scalarField& dMassGas,
110  const scalarField& dMassLiquid,
111  const scalarField& dMassSolid
112 )
113 {
114  scalarField& YMix = this->Y_;
115 
116  scalar massGas =
117  this->updateMassFraction(mass0*YMix[GAS], dMassGas, YGas_);
118  scalar massLiquid =
119  this->updateMassFraction(mass0*YMix[LIQ], dMassLiquid, YLiquid_);
120  scalar massSolid =
121  this->updateMassFraction(mass0*YMix[SLD], dMassSolid, YSolid_);
122 
123  scalar massNew = max(massGas + massLiquid + massSolid, rootVSmall);
124 
125  YMix[GAS] = massGas/massNew;
126  YMix[LIQ] = massLiquid/massNew;
127  YMix[SLD] = 1.0 - YMix[GAS] - YMix[LIQ];
128 
129  return massNew;
130 }
131 
132 
133 // * * * * * * * * * * * Protected Member Functions * * * * * * * * * * * * //
134 
135 template<class ParcelType>
136 template<class TrackCloudType>
138 (
139  TrackCloudType& cloud,
140  trackingData& td
141 )
142 {
143  ParcelType::setCellValues(cloud, td);
144 }
145 
146 
147 template<class ParcelType>
148 template<class TrackCloudType>
150 (
151  TrackCloudType& cloud,
152  trackingData& td,
153  const scalar dt
154 )
155 {
156  // Re-use correction from reacting parcel
157  ParcelType::cellValueSourceCorrection(cloud, td, dt);
158 }
159 
160 
161 template<class ParcelType>
162 template<class TrackCloudType>
164 (
165  TrackCloudType& cloud,
166  trackingData& td,
167  const scalar dt
168 )
169 {
170  typedef typename TrackCloudType::reactingCloudType reactingCloudType;
172  cloud.composition();
173 
174 
175  // Define local properties at beginning of timestep
176  // ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
177 
178  const scalar np0 = this->nParticle_;
179  const scalar d0 = this->d_;
180  const vector& U0 = this->U_;
181  const scalar T0 = this->T_;
182  const scalar mass0 = this->mass();
183 
184  const scalar pc = td.pc();
185 
186  const scalarField& YMix = this->Y_;
187  const label idG = composition.idGas();
188  const label idL = composition.idLiquid();
189  const label idS = composition.idSolid();
190 
191 
192  // Calc surface values
193  scalar Ts, rhos, mus, Prs, kappas;
194  this->calcSurfaceValues(cloud, td, T0, Ts, rhos, mus, Prs, kappas);
195  scalar Res = this->Re(rhos, U0, td.Uc(), d0, mus);
196 
197 
198  // Sources
199  //~~~~~~~~
200 
201  // Explicit momentum source for particle
202  vector Su = Zero;
203 
204  // Linearised momentum source coefficient
205  scalar Spu = 0.0;
206 
207  // Momentum transfer from the particle to the carrier phase
208  vector dUTrans = Zero;
209 
210  // Explicit enthalpy source for particle
211  scalar Sh = 0.0;
212 
213  // Linearised enthalpy source coefficient
214  scalar Sph = 0.0;
215 
216  // Sensible enthalpy transfer from the particle to the carrier phase
217  scalar dhsTrans = 0.0;
218 
219 
220  // 1. Compute models that contribute to mass transfer - U, T held constant
221  // ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
222 
223  // Phase change in liquid phase
224  // ~~~~~~~~~~~~~~~~~~~~~~~~~~~~
225 
226  // Mass transfer due to phase change
227  scalarField dMassPC(YLiquid_.size(), 0.0);
228 
229  // Molar flux of species emitted from the particle (kmol/m^2/s)
230  scalar Ne = 0.0;
231 
232  // Sum Ni*Cpi*Wi of emission species
233  scalar NCpW = 0.0;
234 
235  // Surface concentrations of emitted species
236  scalarField Cs(composition.carrier().species().size(), 0.0);
237 
238  // Calc mass and enthalpy transfer due to phase change
239  this->calcPhaseChange
240  (
241  cloud,
242  td,
243  dt,
244  Res,
245  Prs,
246  Ts,
247  mus/rhos,
248  d0,
249  T0,
250  mass0,
251  idL,
252  YMix[LIQ],
253  YLiquid_,
254  dMassPC,
255  Sh,
256  Ne,
257  NCpW,
258  Cs
259  );
260 
261 
262  // Devolatilisation
263  // ~~~~~~~~~~~~~~~~
264 
265  // Mass transfer due to devolatilisation
266  scalarField dMassDV(YGas_.size(), 0.0);
267 
268  // Calc mass and enthalpy transfer due to devolatilisation
269  calcDevolatilisation
270  (
271  cloud,
272  td,
273  dt,
274  this->age_,
275  Ts,
276  d0,
277  T0,
278  mass0,
279  this->mass0_,
280  YMix[GAS]*YGas_,
281  YMix[LIQ]*YLiquid_,
282  YMix[SLD]*YSolid_,
283  canCombust_,
284  dMassDV,
285  Sh,
286  Ne,
287  NCpW,
288  Cs
289  );
290 
291 
292  // Surface reactions
293  // ~~~~~~~~~~~~~~~~~
294 
295  // Change in carrier phase composition due to surface reactions
296  scalarField dMassSRGas(YGas_.size(), 0.0);
297  scalarField dMassSRLiquid(YLiquid_.size(), 0.0);
298  scalarField dMassSRSolid(YSolid_.size(), 0.0);
299  scalarField dMassSRCarrier(composition.carrier().species().size(), 0.0);
300 
301  // Calc mass and enthalpy transfer due to surface reactions
302  calcSurfaceReactions
303  (
304  cloud,
305  td,
306  dt,
307  d0,
308  T0,
309  mass0,
310  canCombust_,
311  Ne,
312  YMix,
313  YGas_,
314  YLiquid_,
315  YSolid_,
316  dMassSRGas,
317  dMassSRLiquid,
318  dMassSRSolid,
319  dMassSRCarrier,
320  Sh,
321  dhsTrans
322  );
323 
324 
325  // 2. Update the parcel properties due to change in mass
326  // ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
327 
328  scalarField dMassGas(dMassDV + dMassSRGas);
329  scalarField dMassLiquid(dMassPC + dMassSRLiquid);
330  scalarField dMassSolid(dMassSRSolid);
331  scalar mass1 =
332  updateMassFractions(mass0, dMassGas, dMassLiquid, dMassSolid);
333 
334  this->Cp_ = CpEff(cloud, td, pc, T0, idG, idL, idS);
335 
336  // Update particle density or diameter
337  if (cloud.constProps().constantVolume())
338  {
339  this->rho_ = mass1/this->volume();
340  }
341  else
342  {
343  this->d_ = cbrt(mass1/this->rho_*6.0/pi);
344  }
345 
346  // Remove the particle when mass falls below minimum threshold
347  if (np0*mass1 < cloud.constProps().minParcelMass())
348  {
349  td.keepParticle = false;
350 
351  if (cloud.solution().coupled())
352  {
353  scalar dm = np0*mass0;
354 
355  // Absorb parcel into carrier phase
356  forAll(YGas_, i)
357  {
358  label gid = composition.localToCarrierId(GAS, i);
359  cloud.rhoTrans(gid)[this->cell()] += dm*YMix[GAS]*YGas_[i];
360  }
361  forAll(YLiquid_, i)
362  {
363  label gid = composition.localToCarrierId(LIQ, i);
364  cloud.rhoTrans(gid)[this->cell()] += dm*YMix[LIQ]*YLiquid_[i];
365  }
366 
367  // No mapping between solid components and carrier phase
368  /*
369  forAll(YSolid_, i)
370  {
371  label gid = composition.localToCarrierId(SLD, i);
372  cloud.rhoTrans(gid)[this->cell()] += dm*YMix[SLD]*YSolid_[i];
373  }
374  */
375 
376  cloud.UTrans()[this->cell()] += dm*U0;
377 
378  cloud.hsTrans()[this->cell()] +=
379  dm*HsEff(cloud, td, pc, T0, idG, idL, idS);
380 
381  cloud.phaseChange().addToPhaseChangeMass(np0*mass1);
382  }
383 
384  return;
385  }
386 
387  // Correct surface values due to emitted species
388  this->correctSurfaceValues(cloud, td, Ts, Cs, rhos, mus, Prs, kappas);
389  Res = this->Re(rhos, U0, td.Uc(), this->d_, mus);
390 
391 
392  // 3. Compute heat- and momentum transfers
393  // ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
394 
395  // Heat transfer
396  // ~~~~~~~~~~~~~
397 
398  // Calculate new particle temperature
399  this->T_ =
400  this->calcHeatTransfer
401  (
402  cloud,
403  td,
404  dt,
405  Res,
406  Prs,
407  kappas,
408  NCpW,
409  Sh,
410  dhsTrans,
411  Sph
412  );
413 
414 
415  this->Cp_ = CpEff(cloud, td, pc, this->T_, idG, idL, idS);
416 
417 
418  // Motion
419  // ~~~~~~
420 
421  // Calculate new particle velocity
422  this->U_ =
423  this->calcVelocity(cloud, td, dt, Res, mus, mass1, Su, dUTrans, Spu);
424 
425 
426  // 4. Accumulate carrier phase source terms
427  // ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
428 
429  if (cloud.solution().coupled())
430  {
431  // Transfer mass lost to carrier mass, momentum and enthalpy sources
432  forAll(YGas_, i)
433  {
434  scalar dm = np0*dMassGas[i];
435  label gid = composition.localToCarrierId(GAS, i);
436  scalar hs = composition.carrier().Hs(gid, pc, T0);
437  cloud.rhoTrans(gid)[this->cell()] += dm;
438  cloud.UTrans()[this->cell()] += dm*U0;
439  cloud.hsTrans()[this->cell()] += dm*hs;
440  }
441  forAll(YLiquid_, i)
442  {
443  scalar dm = np0*dMassLiquid[i];
444  label gid = composition.localToCarrierId(LIQ, i);
445  scalar hs = composition.carrier().Hs(gid, pc, T0);
446  cloud.rhoTrans(gid)[this->cell()] += dm;
447  cloud.UTrans()[this->cell()] += dm*U0;
448  cloud.hsTrans()[this->cell()] += dm*hs;
449  }
450 
451  // No mapping between solid components and carrier phase
452  /*
453  forAll(YSolid_, i)
454  {
455  scalar dm = np0*dMassSolid[i];
456  label gid = composition.localToCarrierId(SLD, i);
457  scalar hs = composition.carrier().Hs(gid, pc, T0);
458  cloud.rhoTrans(gid)[this->cell()] += dm;
459  cloud.UTrans()[this->cell()] += dm*U0;
460  cloud.hsTrans()[this->cell()] += dm*hs;
461  }
462  */
463 
464  forAll(dMassSRCarrier, i)
465  {
466  scalar dm = np0*dMassSRCarrier[i];
467  scalar hs = composition.carrier().Hs(i, pc, T0);
468  cloud.rhoTrans(i)[this->cell()] += dm;
469  cloud.UTrans()[this->cell()] += dm*U0;
470  cloud.hsTrans()[this->cell()] += dm*hs;
471  }
472 
473  // Update momentum transfer
474  cloud.UTrans()[this->cell()] += np0*dUTrans;
475  cloud.UCoeff()[this->cell()] += np0*Spu;
476 
477  // Update sensible enthalpy transfer
478  cloud.hsTrans()[this->cell()] += np0*dhsTrans;
479  cloud.hsCoeff()[this->cell()] += np0*Sph;
480 
481  // Update radiation fields
482  if (cloud.radiation())
483  {
484  const scalar ap = this->areaP();
485  const scalar T4 = pow4(T0);
486  cloud.radAreaP()[this->cell()] += dt*np0*ap;
487  cloud.radT4()[this->cell()] += dt*np0*T4;
488  cloud.radAreaPT4()[this->cell()] += dt*np0*ap*T4;
489  }
490  }
491 }
492 
493 
494 // * * * * * * * * * * * * Protected Member Functions * * * * * * * * * * * //
495 
496 template<class ParcelType>
497 template<class TrackCloudType>
499 (
500  TrackCloudType& cloud,
501  trackingData& td,
502  const scalar dt,
503  const scalar age,
504  const scalar Ts,
505  const scalar d,
506  const scalar T,
507  const scalar mass,
508  const scalar mass0,
509  const scalarField& YGasEff,
510  const scalarField& YLiquidEff,
511  const scalarField& YSolidEff,
512  label& canCombust,
513  scalarField& dMassDV,
514  scalar& Sh,
515  scalar& N,
516  scalar& NCpW,
517  scalarField& Cs
518 ) const
519 {
520  // Check that model is active
521  if (!cloud.devolatilisation().active())
522  {
523  if (canCombust != -1)
524  {
525  canCombust = 1;
526  }
527  return;
528  }
529 
530  // Initialise demand-driven constants
531  (void)cloud.constProps().TDevol();
532  (void)cloud.constProps().LDevol();
533 
534  // Check that the parcel temperature is within necessary limits for
535  // devolatilisation to occur
536  if (T < cloud.constProps().TDevol() || canCombust == -1)
537  {
538  return;
539  }
540 
541  typedef typename TrackCloudType::reactingCloudType reactingCloudType;
543  cloud.composition();
544 
545 
546  // Total mass of volatiles evolved
547  cloud.devolatilisation().calculate
548  (
549  dt,
550  age,
551  mass0,
552  mass,
553  T,
554  YGasEff,
555  YLiquidEff,
556  YSolidEff,
557  canCombust,
558  dMassDV
559  );
560 
561  scalar dMassTot = sum(dMassDV);
562 
563  cloud.devolatilisation().addToDevolatilisationMass
564  (
565  this->nParticle_*dMassTot
566  );
567 
568  Sh -= dMassTot*cloud.constProps().LDevol()/dt;
569 
570  // Update molar emissions
571  if (cloud.heatTransfer().BirdCorrection())
572  {
573  // Molar average molecular weight of carrier mix
574  const scalar Wc = max(small, td.rhoc()*RR*td.Tc()/td.pc());
575 
576  // Note: hardcoded gaseous diffusivities for now
577  // TODO: add to carrier thermo
578  const scalar beta = sqr(cbrt(15.0) + cbrt(15.0));
579 
580  forAll(dMassDV, i)
581  {
582  const label id = composition.localToCarrierId(GAS, i);
583  const scalar Cp = composition.carrier().Cp(id, td.pc(), Ts);
584  const scalar W = composition.carrier().W(id);
585  const scalar Ni = dMassDV[i]/(this->areaS(d)*dt*W);
586 
587  // Dab calc'd using API vapour mass diffusivity function
588  const scalar Dab =
589  3.6059e-3*(pow(1.8*Ts, 1.75))
590  *sqrt(1.0/W + 1.0/Wc)
591  /(td.pc()*beta);
592 
593  N += Ni;
594  NCpW += Ni*Cp*W;
595  Cs[id] += Ni*d/(2.0*Dab);
596  }
597  }
598 }
599 
600 
601 template<class ParcelType>
602 template<class TrackCloudType>
604 (
605  TrackCloudType& cloud,
606  trackingData& td,
607  const scalar dt,
608  const scalar d,
609  const scalar T,
610  const scalar mass,
611  const label canCombust,
612  const scalar N,
613  const scalarField& YMix,
614  const scalarField& YGas,
615  const scalarField& YLiquid,
616  const scalarField& YSolid,
617  scalarField& dMassSRGas,
618  scalarField& dMassSRLiquid,
619  scalarField& dMassSRSolid,
620  scalarField& dMassSRCarrier,
621  scalar& Sh,
622  scalar& dhsTrans
623 ) const
624 {
625  // Check that model is active
626  if (!cloud.surfaceReaction().active())
627  {
628  return;
629  }
630 
631  // Initialise demand-driven constants
632  (void)cloud.constProps().hRetentionCoeff();
633  (void)cloud.constProps().TMax();
634 
635  // Check that model is active
636  if (canCombust != 1)
637  {
638  return;
639  }
640 
641 
642  // Update surface reactions
643  const scalar hReaction = cloud.surfaceReaction().calculate
644  (
645  dt,
646  this->cell(),
647  d,
648  T,
649  td.Tc(),
650  td.pc(),
651  td.rhoc(),
652  mass,
653  YGas,
654  YLiquid,
655  YSolid,
656  YMix,
657  N,
658  dMassSRGas,
659  dMassSRLiquid,
660  dMassSRSolid,
661  dMassSRCarrier
662  );
663 
664  cloud.surfaceReaction().addToSurfaceReactionMass
665  (
666  this->nParticle_
667  *(sum(dMassSRGas) + sum(dMassSRLiquid) + sum(dMassSRSolid))
668  );
669 
670  const scalar xsi = min(T/cloud.constProps().TMax(), 1.0);
671  const scalar coeff =
672  (1.0 - xsi*xsi)*cloud.constProps().hRetentionCoeff();
673 
674  Sh += coeff*hReaction/dt;
675 
676  dhsTrans += (1.0 - coeff)*hReaction;
677 }
678 
679 
680 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
681 
682 template<class ParcelType>
684 (
686 )
687 :
688  ParcelType(p),
689  YGas_(p.YGas_),
690  YLiquid_(p.YLiquid_),
691  YSolid_(p.YSolid_),
692  canCombust_(p.canCombust_)
693 {}
694 
695 
696 template<class ParcelType>
698 (
700  const polyMesh& mesh
701 )
702 :
703  ParcelType(p, mesh),
704  YGas_(p.YGas_),
705  YLiquid_(p.YLiquid_),
706  YSolid_(p.YSolid_),
707  canCombust_(p.canCombust_)
708 {}
709 
710 
711 // * * * * * * * * * * * * * * IOStream operators * * * * * * * * * * * * * //
712 
714 
715 // ************************************************************************* //
#define forAll(list, i)
Loop across all elements in list.
Definition: UList.H:428
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
zeroField Su
Definition: alphaSuSp.H:1
dimensioned< Type > max(const dimensioned< Type > &, const dimensioned< Type > &)
basicSpecieMixture & composition
label canCombust_
Flag to identify if the particle can devolatilise and combust.
dimensionedSymmTensor sqr(const dimensionedVector &dv)
virtual label idSolid() const =0
Solid id.
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.
dimensionedScalar sqrt(const dimensionedScalar &ds)
label localToCarrierId(const label phaseI, const label id, const bool allowNotFound=false) const
Return carrier id of component given local id.
void calc(TrackCloudType &cloud, trackingData &td, const scalar dt)
Update parcel properties over the time interval.
dimensioned< Type > sum(const DimensionedField< Type, GeoMesh > &df)
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 [].
mathematical constants.
dimensionedScalar cbrt(const dimensionedScalar &ds)
static const zero Zero
Definition: zero.H:97
virtual label idLiquid() const =0
Liquid id.
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 [].
dimensioned< Type > min(const dimensioned< Type > &, const dimensioned< Type > &)
scalarField YLiquid_
Mass fractions of liquids [].
dimensionedScalar pow(const dimensionedScalar &ds, const dimensionedScalar &expt)
A cell is defined as a list of faces with extra functionality.
Definition: cell.H:56
const scalarList W(::W(thermo))
dimensionedScalar pow4(const dimensionedScalar &ds)
Multiphase variant of the reacting parcel class with one/two-way coupling with the continuous phase...
dimensionedScalar beta("beta", dimless/dimTemperature, laminarTransport)
label N
Definition: createFields.H:22
const scalar RR
Universal gas constant (default in [J/(kmol K)])
Mesh consisting of general polyhedral cells.
Definition: polyMesh.H:74
virtual label idGas() const =0
Gas id.
Templated reacting parcel composition model class Consists of carrier species (via thermo package)...
Definition: ReactingCloud.H:52
void cellValueSourceCorrection(TrackCloudType &cloud, trackingData &td, const scalar dt)
Correct cell values using latest transfer information.
scalarField Re(const UList< complex > &cf)
Definition: complexFields.C:97
const basicSpecieMixture & carrier() const
Return the carrier components (wrapper function)
void setCellValues(TrackCloudType &cloud, trackingData &td)
Set cell values.