ReactingParcel.C
Go to the documentation of this file.
1 /*---------------------------------------------------------------------------*\
2  ========= |
3  \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
4  \\ / O peration |
5  \\ / A nd | www.openfoam.com
6  \\/ M anipulation |
7 -------------------------------------------------------------------------------
8  Copyright (C) 2011-2017 OpenFOAM Foundation
9  Copyright (C) 2020 OpenCFD Ltd.
10 -------------------------------------------------------------------------------
11 License
12  This file is part of OpenFOAM.
13 
14  OpenFOAM is free software: you can redistribute it and/or modify it
15  under the terms of the GNU General Public License as published by
16  the Free Software Foundation, either version 3 of the License, or
17  (at your option) any later version.
18 
19  OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
20  ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
21  FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
22  for more details.
23 
24  You should have received a copy of the GNU General Public License
25  along with OpenFOAM. If not, see <http://www.gnu.org/licenses/>.
26 
27 \*---------------------------------------------------------------------------*/
28 
29 #include "ReactingParcel.H"
30 #include "specie.H"
31 #include "CompositionModel.H"
32 #include "PhaseChangeModel.H"
33 #include "mathematicalConstants.H"
34 
35 using namespace Foam::constant::mathematical;
36 
37 // * * * * * * * * * * * Protected Member Functions * * * * * * * * * * * * //
38 
39 template<class ParcelType>
40 template<class TrackCloudType>
42 (
43  TrackCloudType& cloud,
44  const scalarField& dMass,
45  const scalar p,
46  const scalar T
47 )
48 {
49  const auto& composition = cloud.composition();
50 
51  scalarField dVolLiquid(dMass.size(), Zero);
52  forAll(dVolLiquid, i)
53  {
54  dVolLiquid[i] =
55  dMass[i]/composition.liquids().properties()[i].rho(p, T);
56  }
57 
58  return sum(dVolLiquid);
59 }
60 
61 
62 template<class ParcelType>
63 template<class TrackCloudType>
65 (
66  TrackCloudType& cloud,
67  trackingData& td,
68  const scalar dt,
69  const scalar Re,
70  const scalar Pr,
71  const scalar Ts,
72  const scalar nus,
73  const scalar d,
74  const scalar T,
75  const scalar mass,
76  const scalar rho,
77  const label idPhase,
78  const scalar YPhase,
79  const scalarField& Y,
80  const scalarField& Ysol,
81  scalarField& dMassPC,
82  scalar& Sh,
83  scalar& N,
84  scalar& NCpW,
86 )
87 {
88  const auto& composition = cloud.composition();
89  auto& phaseChange = cloud.phaseChange();
90 
91  // Some models allow evaporation and condensation
92  if (!phaseChange.active())
93  {
94  return;
95  }
96 
97  scalarField X(composition.liquids().X(Y));
98 
99  scalar Tvap = phaseChange.Tvap(X);
100  if (T < Tvap)
101  {
102  return;
103  }
104 
105  const scalar TMax = phaseChange.TMax(td.pc(), X);
106  const scalar Tdash = min(T, TMax);
107  const scalar Tsdash = min(Ts, TMax);
108 
109  // Calculate mass transfer due to phase change
110  phaseChange.calculate
111  (
112  dt,
113  this->cell(),
114  Re,
115  Pr,
116  d,
117  nus,
118  rho,
119  Tdash,
120  Tsdash,
121  td.pc(),
122  td.Tc(),
123  X, // components molar fractions purely in the liquid
124  Ysol*mass, // total solid mass
125  YPhase*Y*mass, // total liquid mass
126  dMassPC
127  );
128 
129  // Limit phase change mass by availability of each specie
130  forAll(Y, i)
131  {
132  // evaporation
133  if (dMassPC[i] > 0)
134  {
135  dMassPC[i] = min(mass*YPhase*Y[i], dMassPC[i]);
136  }
137  // condensation Do something?
138  }
139 
140  const scalar dMassTot = sum(dMassPC);
141 
142  // Add to cumulative phase change mass
143  phaseChange.addToPhaseChangeMass(this->nParticle_*dMassTot);
144 
145  forAll(dMassPC, i)
146  {
147  const label cid = composition.localToCarrierId(idPhase, i);
148 
149  const scalar dh = phaseChange.dh(cid, i, td.pc(), Tdash);
150  Sh -= dMassPC[i]*dh/dt;
151  }
152 
153 
154  // Update molar emissions
155  if (cloud.heatTransfer().BirdCorrection())
156  {
157  // Average molecular weight of carrier mix - assumes perfect gas
158  const scalar Wc = td.rhoc()*RR*td.Tc()/td.pc();
159 
160  forAll(dMassPC, i)
161  {
162  const label cid = composition.localToCarrierId(idPhase, i);
163 
164  const scalar Cp = composition.carrier().Cp(cid, td.pc(), Tsdash);
165  const scalar W = composition.carrier().W(cid);
166  const scalar Ni = dMassPC[i]/(this->areaS(d)*dt*W);
167 
168  const scalar Dab =
169  composition.liquids().properties()[i].D(td.pc(), Tsdash, Wc);
170 
171  // Molar flux of species coming from the particle (kmol/m^2/s)
172  N += Ni;
173 
174  // Sum of Ni*Cpi*Wi of emission species
175  NCpW += Ni*Cp*W;
176 
177  // Concentrations of emission species
178  Cs[cid] += Ni*d/(2.0*Dab);
179  }
180  }
181 }
182 
183 
184 template<class ParcelType>
186 (
187  const scalar mass0,
188  const scalarField& dMass,
189  scalarField& Y
190 ) const
191 {
192  scalar mass1 = mass0 - sum(dMass);
193 
194  // only update the mass fractions if the new particle mass is finite
195  if (mass1 > ROOTVSMALL)
196  {
197  forAll(Y, i)
198  {
199  Y[i] = (Y[i]*mass0 - dMass[i])/mass1;
200  }
201  }
202 
203  return mass1;
204 }
205 
206 
207 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
208 
209 template<class ParcelType>
211 (
213 )
214 :
215  ParcelType(p),
216  mass0_(p.mass0_),
217  Y_(p.Y_)
218 {}
219 
220 
221 template<class ParcelType>
223 (
224  const ReactingParcel<ParcelType>& p,
225  const polyMesh& mesh
226 )
227 :
228  ParcelType(p, mesh),
229  mass0_(p.mass0_),
230  Y_(p.Y_)
231 {}
232 
233 
234 // * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * * * //
235 
236 template<class ParcelType>
237 template<class TrackCloudType>
239 (
240  TrackCloudType& cloud,
241  trackingData& td
242 )
243 {
244  ParcelType::setCellValues(cloud, td);
245 
246  td.pc() = td.pInterp().interpolate
247  (
248  this->coordinates(),
249  this->currentTetIndices()
250  );
251 
252  if (td.pc() < cloud.constProps().pMin())
253  {
254  if (debug)
255  {
257  << "Limiting observed pressure in cell " << this->cell()
258  << " to " << cloud.constProps().pMin() << nl << endl;
259  }
260 
261  td.pc() = cloud.constProps().pMin();
262  }
263 }
264 
265 
266 template<class ParcelType>
267 template<class TrackCloudType>
269 (
270  TrackCloudType& cloud,
271  trackingData& td,
272  const scalar dt
273 )
274 {
275  scalar addedMass = 0.0;
276  scalar maxMassI = 0.0;
277  forAll(cloud.rhoTrans(), i)
278  {
279  scalar dm = cloud.rhoTrans(i)[this->cell()];
280  maxMassI = max(maxMassI, mag(dm));
281  addedMass += dm;
282  }
283 
284  if (maxMassI < ROOTVSMALL)
285  {
286  return;
287  }
288 
289  const scalar massCell = this->massCell(td);
290 
291  td.rhoc() += addedMass/cloud.pMesh().cellVolumes()[this->cell()];
292 
293  const scalar massCellNew = massCell + addedMass;
294  td.Uc() = (td.Uc()*massCell + cloud.UTrans()[this->cell()])/massCellNew;
295 
296  scalar CpEff = 0.0;
297  forAll(cloud.rhoTrans(), i)
298  {
299  scalar Y = cloud.rhoTrans(i)[this->cell()]/addedMass;
300  CpEff += Y*cloud.composition().carrier().Cp(i, td.pc(), td.Tc());
301  }
302 
303  const scalar Cpc = td.CpInterp().psi()[this->cell()];
304  td.Cpc() = (massCell*Cpc + addedMass*CpEff)/massCellNew;
305 
306  td.Tc() += cloud.hsTrans()[this->cell()]/(td.Cpc()*massCellNew);
307 
308  if (td.Tc() < cloud.constProps().TMin())
309  {
310  if (debug)
311  {
313  << "Limiting observed temperature in cell " << this->cell()
314  << " to " << cloud.constProps().TMin() << nl << endl;
315  }
316 
317  td.Tc() = cloud.constProps().TMin();
318  }
319 }
320 
321 
322 template<class ParcelType>
323 template<class TrackCloudType>
325 (
326  TrackCloudType& cloud,
327  trackingData& td,
328  const scalar T,
329  const scalarField& Cs,
330  scalar& rhos,
331  scalar& mus,
332  scalar& Prs,
333  scalar& kappas
334 )
335 {
336  // No correction if total concentration of emitted species is small
337  if (!cloud.heatTransfer().BirdCorrection() || (sum(Cs) < SMALL))
338  {
339  return;
340  }
341 
342  const SLGThermo& thermo = cloud.thermo();
343 
344  // Far field carrier molar fractions
345  scalarField Xinf(thermo.carrier().species().size());
346 
347  forAll(Xinf, i)
348  {
349  Xinf[i] = thermo.carrier().Y(i)[this->cell()]/thermo.carrier().W(i);
350  }
351  Xinf /= sum(Xinf);
352 
353  // Molar fraction of far field species at particle surface
354  const scalar Xsff = 1.0 - min(sum(Cs)*RR*this->T_/td.pc(), 1.0);
355 
356  // Surface carrier total molar concentration
357  const scalar CsTot = td.pc()/(RR*this->T_);
358 
359  // Surface carrier composition (molar fraction)
360  scalarField Xs(Xinf.size());
361 
362  // Surface carrier composition (mass fraction)
363  scalarField Ys(Xinf.size());
364 
365  forAll(Xs, i)
366  {
367  // Molar concentration of species at particle surface
368  const scalar Csi = Cs[i] + Xsff*Xinf[i]*CsTot;
369 
370  Xs[i] = (2.0*Csi + Xinf[i]*CsTot)/3.0;
371  Ys[i] = Xs[i]*thermo.carrier().W(i);
372  }
373  Xs /= sum(Xs);
374  Ys /= sum(Ys);
375 
376 
377  rhos = 0;
378  mus = 0;
379  kappas = 0;
380  scalar Cps = 0;
381  scalar sumYiSqrtW = 0;
382  scalar sumYiCbrtW = 0;
383 
384  forAll(Ys, i)
385  {
386  const scalar W = thermo.carrier().W(i);
387  const scalar sqrtW = sqrt(W);
388  const scalar cbrtW = cbrt(W);
389 
390  rhos += Xs[i]*W;
391  mus += Ys[i]*sqrtW*thermo.carrier().mu(i, td.pc(), T);
392  kappas += Ys[i]*cbrtW*thermo.carrier().kappa(i, td.pc(), T);
393  Cps += Xs[i]*thermo.carrier().Cp(i, td.pc(), T);
394 
395  sumYiSqrtW += Ys[i]*sqrtW;
396  sumYiCbrtW += Ys[i]*cbrtW;
397  }
398 
399  Cps = max(Cps, ROOTVSMALL);
400 
401  rhos *= td.pc()/(RR*T);
402  rhos = max(rhos, ROOTVSMALL);
403 
404  mus /= sumYiSqrtW;
405  mus = max(mus, ROOTVSMALL);
406 
407  kappas /= sumYiCbrtW;
408  kappas = max(kappas, ROOTVSMALL);
409 
410  Prs = Cps*mus/kappas;
411 }
412 
413 
414 template<class ParcelType>
415 template<class TrackCloudType>
417 (
418  TrackCloudType& cloud,
419  trackingData& td,
420  const scalar dt
421 )
422 {
423  const auto& composition = cloud.composition();
424 
425 
426  // Define local properties at beginning of time step
427  // ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
428 
429  const scalar np0 = this->nParticle_;
430  const scalar d0 = this->d_;
431  const vector& U0 = this->U_;
432  const scalar T0 = this->T_;
433  const scalar mass0 = this->mass();
434  const scalar rho0 = this->rho_;
435 
436 
437  // Calc surface values
438  scalar Ts, rhos, mus, Prs, kappas;
439  this->calcSurfaceValues(cloud, td, T0, Ts, rhos, mus, Prs, kappas);
440  scalar Res = this->Re(rhos, U0, td.Uc(), d0, mus);
441 
442 
443  // Sources
444  // ~~~~~~~
445 
446  // Explicit momentum source for particle
447  vector Su = Zero;
448 
449  // Linearised momentum source coefficient
450  scalar Spu = 0.0;
451 
452  // Momentum transfer from the particle to the carrier phase
453  vector dUTrans = Zero;
454 
455  // Explicit enthalpy source for particle
456  scalar Sh = 0.0;
457 
458  // Linearised enthalpy source coefficient
459  scalar Sph = 0.0;
460 
461  // Sensible enthalpy transfer from the particle to the carrier phase
462  scalar dhsTrans = 0.0;
463 
464 
465  // 1. Compute models that contribute to mass transfer - U, T held constant
466  // ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
467 
468  // Phase change
469  // ~~~~~~~~~~~~
470 
471  // Mass transfer due to phase change
472  scalarField dMassPC(Y_.size(), Zero);
473 
474  // Molar flux of species emitted from the particle (kmol/m^2/s)
475  scalar Ne = 0.0;
476 
477  // Sum Ni*Cpi*Wi of emission species
478  scalar NCpW = 0.0;
479 
480  // Surface concentrations of emitted species
481  scalarField Cs(composition.carrier().species().size(), Zero);
482 
483  // Calc mass and enthalpy transfer due to phase change
484  calcPhaseChange
485  (
486  cloud,
487  td,
488  dt,
489  Res,
490  Prs,
491  Ts,
492  mus/rhos,
493  d0,
494  T0,
495  mass0,
496  rho0,
497  0,
498  1.0,
499  Y_,
500  scalarField(),
501  dMassPC,
502  Sh,
503  Ne,
504  NCpW,
505  Cs
506  );
507 
508 
509  // 2. Update the parcel properties due to change in mass
510  // ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
511 
512  scalarField dMass(dMassPC);
513  scalar mass1 = updateMassFraction(mass0, dMass, Y_);
514 
515  this->Cp_ = composition.Cp(0, Y_, td.pc(), T0);
516 
517  if
518  (
519  cloud.constProps().volUpdateType()
520  == constantProperties::volumeUpdateType::mUndefined
521  )
522  {
523  // Update particle density or diameter
524  if (cloud.constProps().constantVolume())
525  {
526  this->rho_ = mass1/this->volume();
527  }
528  else
529  {
530  this->d_ = cbrt(mass1/this->rho_*6/pi);
531  }
532  }
533  else
534  {
535  switch (cloud.constProps().volUpdateType())
536  {
537  case constantProperties::volumeUpdateType::mConstRho :
538  {
539  this->d_ = cbrt(mass1/this->rho_*6/pi);
540  break;
541  }
542  case constantProperties::volumeUpdateType::mConstVol :
543  {
544  this->rho_ = mass1/this->volume();
545  break;
546  }
547  case constantProperties::volumeUpdateType::mUpdateRhoAndVol :
548  {
549  scalar deltaVol =
550  updatedDeltaVolume
551  (
552  cloud,
553  dMass,
554  td.pc(),
555  T0
556  );
557  this->rho_ = mass1/(this->volume() - deltaVol);
558  this->d_ = cbrt(mass1/this->rho_*6/pi);
559  break;
560  }
561 
562  }
563  }
564 
565  // Remove the particle when mass falls below minimum threshold
566  if (np0*mass1 < cloud.constProps().minParcelMass())
567  {
568  td.keepParticle = false;
569 
570  if (cloud.solution().coupled())
571  {
572  scalar dm = np0*mass0;
573 
574  // Absorb parcel into carrier phase
575  forAll(Y_, i)
576  {
577  scalar dmi = dm*Y_[i];
578  label gid = composition.localToCarrierId(0, i);
579  scalar hs = composition.carrier().Hs(gid, td.pc(), T0);
580 
581  cloud.rhoTrans(gid)[this->cell()] += dmi;
582  cloud.hsTrans()[this->cell()] += dmi*hs;
583  }
584  cloud.UTrans()[this->cell()] += dm*U0;
585 
586  cloud.phaseChange().addToPhaseChangeMass(np0*mass1);
587  }
588 
589  return;
590  }
591 
592  // Correct surface values due to emitted species
593  correctSurfaceValues(cloud, td, Ts, Cs, rhos, mus, Prs, kappas);
594  Res = this->Re(rhos, U0, td.Uc(), this->d(), mus);
595 
596 
597  // 3. Compute heat- and momentum transfers
598  // ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
599 
600  // Heat transfer
601  // ~~~~~~~~~~~~~
602 
603  // Calculate new particle temperature
604  this->T_ =
605  this->calcHeatTransfer
606  (
607  cloud,
608  td,
609  dt,
610  Res,
611  Prs,
612  kappas,
613  NCpW,
614  Sh,
615  dhsTrans,
616  Sph
617  );
618 
619  this->Cp_ = composition.Cp(0, Y_, td.pc(), T0);
620 
621 
622  // Motion
623  // ~~~~~~
624 
625  // Calculate new particle velocity
626  this->U_ =
627  this->calcVelocity(cloud, td, dt, Res, mus, mass1, Su, dUTrans, Spu);
628 
629 
630  // 4. Accumulate carrier phase source terms
631  // ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
632 
633  if (cloud.solution().coupled())
634  {
635  // Transfer mass lost to carrier mass, momentum and enthalpy sources
636  forAll(dMass, i)
637  {
638  scalar dm = np0*dMass[i];
639  label gid = composition.localToCarrierId(0, i);
640  scalar hs = composition.carrier().Hs(gid, td.pc(), T0);
641 
642  cloud.rhoTrans(gid)[this->cell()] += dm;
643  cloud.UTrans()[this->cell()] += dm*U0;
644  cloud.hsTrans()[this->cell()] += dm*hs;
645  }
646 
647  // Update momentum transfer
648  cloud.UTrans()[this->cell()] += np0*dUTrans;
649  cloud.UCoeff()[this->cell()] += np0*Spu;
650 
651  // Update sensible enthalpy transfer
652  cloud.hsTrans()[this->cell()] += np0*dhsTrans;
653  cloud.hsCoeff()[this->cell()] += np0*Sph;
654 
655  // Update radiation fields
656  if (cloud.radiation())
657  {
658  const scalar ap = this->areaP();
659  const scalar T4 = pow4(T0);
660  cloud.radAreaP()[this->cell()] += dt*np0*ap;
661  cloud.radT4()[this->cell()] += dt*np0*T4;
662  cloud.radAreaPT4()[this->cell()] += dt*np0*ap*T4;
663  }
664  }
665 }
666 
667 
668 // * * * * * * * * * * * * * * IOStream operators * * * * * * * * * * * * * //
669 
670 #include "ReactingParcelIO.C"
671 
672 // ************************************************************************* //
Foam::expressions::patchExpr::debug
int debug
Static debugging option.
Foam::constant::thermodynamic::RR
const scalar RR
Universal gas constant: default in [J/(kmol K)].
Definition: thermodynamicConstants.C:46
Foam::scalarField
Field< scalar > scalarField
Specialisation of Field<T> for scalar.
Definition: primitiveFieldsFwd.H:52
ReactingParcel.H
mathematicalConstants.H
p
volScalarField & p
Definition: createFieldRefs.H:8
Foam::ReactingParcel::trackingData::pInterp
const interpolation< scalar > & pInterp() const
Return const access to the interpolator for continuous phase.
Definition: ReactingParcelTrackingDataI.H:51
Foam::SLGThermo
Thermo package for (S)olids (L)iquids and (G)ases Takes reference to thermo package,...
Definition: SLGThermo.H:64
Foam::Zero
static constexpr const zero Zero
Global zero (0)
Definition: zero.H:131
Foam::basicMultiComponentMixture::species
const speciesTable & species() const
Return the table of species.
Definition: basicMultiComponentMixtureI.H:29
thermo
Basic thermodynamics type based on the use of fitting functions for cp, h, s obtained from the templa...
Foam::ReactingParcel::calc
void calc(TrackCloudType &cloud, trackingData &td, const scalar dt)
Update parcel properties over the time interval.
Definition: ReactingParcel.C:417
Foam::endl
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:369
Foam::fieldTypes::volume
const wordList volume
Standard volume field types (scalar, vector, tensor, etc)
PhaseChangeModel.H
ReactingParcelIO.C
specie.H
rho
rho
Definition: readInitialConditions.H:88
composition
basicSpecieMixture & composition
Definition: createFieldRefs.H:6
Foam::min
label min(const labelHashSet &set, label minValue=labelMax)
Find the min value in labelHashSet, optionally limited by second argument.
Definition: hashSets.C:33
Foam::ReactingParcel::trackingData::pc
scalar pc() const
Return the continuous phase pressure.
Definition: ReactingParcelTrackingDataI.H:58
Foam::basicSpecieMixture::Hs
virtual scalar Hs(const label speciei, const scalar p, const scalar T) const =0
Sensible enthalpy [J/kg].
forAll
#define forAll(list, i)
Loop across all elements in list.
Definition: stdFoam.H:296
Foam::ReactingParcel::calcPhaseChange
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 scalar rho, const label idPhase, const scalar YPhase, const scalarField &YLiq, const scalarField &YSol, scalarField &dMassPC, scalar &Sh, scalar &N, scalar &NCpW, scalarField &Cs)
Calculate Phase change.
Definition: ReactingParcel.C:65
Su
zeroField Su
Definition: alphaSuSp.H:1
Foam::pow4
dimensionedScalar pow4(const dimensionedScalar &ds)
Definition: dimensionedScalar.C:100
Foam::ReactingParcel::updatedDeltaVolume
scalar updatedDeltaVolume(TrackCloudType &cloud, const scalarField &dMass, const scalar p, const scalar T)
Return change of volume due to mass exchange.
Foam::Field< scalar >
Foam::ReactingParcel::setCellValues
void setCellValues(TrackCloudType &cloud, trackingData &td)
Set cell values.
Definition: ReactingParcel.C:239
coordinates
PtrList< coordinateSystem > coordinates(solidRegions.size())
Pr
dimensionedScalar Pr("Pr", dimless, laminarTransport)
Foam::T
void T(FieldField< Field, Type > &f1, const FieldField< Field, Type > &f2)
Definition: FieldFieldFunctions.C:58
Foam::ReactingParcel::trackingData
Definition: ReactingParcel.H:142
Foam::max
label max(const labelHashSet &set, label maxValue=labelMin)
Find the max value in labelHashSet, optionally limited by second argument.
Definition: hashSets.C:47
Foam::ReactingParcel::ReactingParcel
ReactingParcel(const polyMesh &mesh, const barycentric &coordinates, const label celli, const label tetFacei, const label tetPti)
Construct from mesh, coordinates and topology.
Definition: ReactingParcelI.H:108
mesh
dynamicFvMesh & mesh
Definition: createDynamicFvMesh.H:6
Foam::ReactingParcel::cellValueSourceCorrection
void cellValueSourceCorrection(TrackCloudType &cloud, trackingData &td, const scalar dt)
Correct cell values using latest transfer information.
Definition: ReactingParcel.C:269
Y
PtrList< volScalarField > & Y
Definition: createFieldRefs.H:7
Foam::cloud
A cloud is a registry collection of lagrangian particles.
Definition: cloud.H:57
rho0
scalar rho0
Definition: readInitialConditions.H:89
Foam::constant::mathematical
Mathematical constants.
Foam::constant::mathematical::pi
constexpr scalar pi(M_PI)
Foam::nl
constexpr char nl
Definition: Ostream.H:404
Foam::Vector< scalar >
Foam::sqrt
dimensionedScalar sqrt(const dimensionedScalar &ds)
Definition: dimensionedScalar.C:144
CompositionModel.H
Foam::mag
dimensioned< typename typeOfMag< Type >::type > mag(const dimensioned< Type > &dt)
Foam::Re
scalarField Re(const UList< complex > &cf)
Extract real component.
Definition: complexField.C:159
Cp
const volScalarField & Cp
Definition: EEqn.H:7
Foam::sum
dimensioned< Type > sum(const DimensionedField< Type, GeoMesh > &df)
Definition: DimensionedFieldFunctions.C:327
Cs
const scalarField & Cs
Definition: solveBulkSurfactant.H:10
N
const Vector< label > N(dict.get< Vector< label >>("N"))
Foam::ReactingParcel::correctSurfaceValues
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.
Definition: ReactingParcel.C:325
Foam::cbrt
dimensionedScalar cbrt(const dimensionedScalar &ds)
Definition: dimensionedScalar.C:155
Foam::interpolation::interpolate
virtual Type interpolate(const vector &position, const label celli, const label facei=-1) const =0
Interpolate field to the given point in the given cell.
Foam::ReactingParcel
Reacting parcel class with one/two-way coupling with the continuous phase.
Definition: ReactingParcel.H:56
T0
scalar T0
Definition: createFields.H:22
WarningInFunction
#define WarningInFunction
Report a warning using Foam::Warning.
Definition: messageStream.H:328
Foam::cell
A cell is defined as a list of faces with extra functionality.
Definition: cell.H:54
Foam::ReactingParcel::updateMassFraction
scalar updateMassFraction(const scalar mass0, const scalarField &dMass, scalarField &Y) const
Update mass fraction.
Definition: ReactingParcel.C:186
Foam::basicSpecieMixture::Cp
virtual scalar Cp(const label speciei, const scalar p, const scalar T) const =0
Heat capacity at constant pressure [J/(kg K)].