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