Go to the documentation of this file.
38 template<
class ParcelType>
39 template<
class TrackCloudType>
42 TrackCloudType&
cloud,
63 auto& phaseChange =
cloud.phaseChange();
65 if (!phaseChange.active() || (YPhase < SMALL))
72 scalar Tvap = phaseChange.Tvap(X);
79 const scalar TMax = phaseChange.TMax(td.
pc(), X);
80 const scalar Tdash =
min(
T, TMax);
81 const scalar Tsdash =
min(Ts, TMax);
101 dMassPC =
min(mass*YPhase*
Y, dMassPC);
103 const scalar dMassTot =
sum(dMassPC);
106 phaseChange.addToPhaseChangeMass(this->nParticle_*dMassTot);
110 const label cid =
composition.localToCarrierId(idPhase, i);
112 const scalar dh = phaseChange.dh(cid, i, td.
pc(), Tdash);
113 Sh -= dMassPC[i]*dh/dt;
118 if (
cloud.heatTransfer().BirdCorrection())
121 const scalar Wc = td.rhoc()*
RR*td.Tc()/td.
pc();
125 const label cid =
composition.localToCarrierId(idPhase, i);
129 const scalar Ni = dMassPC[i]/(this->areaS(d)*dt*W);
132 composition.liquids().properties()[i].D(td.
pc(), Tsdash, Wc);
141 Cs[cid] += Ni*d/(2.0*Dab);
147 template<
class ParcelType>
155 scalar mass1 = mass0 -
sum(dMass);
158 if (mass1 > ROOTVSMALL)
162 Y[i] = (
Y[i]*mass0 - dMass[i])/mass1;
172 template<
class ParcelType>
184 template<
class ParcelType>
187 const ReactingParcel<ParcelType>&
p,
199 template<
class ParcelType>
200 template<
class TrackCloudType>
203 TrackCloudType&
cloud,
207 ParcelType::setCellValues(
cloud, td);
212 this->currentTetIndices()
215 if (td.
pc() <
cloud.constProps().pMin())
220 <<
"Limiting observed pressure in cell " << this->
cell()
221 <<
" to " <<
cloud.constProps().pMin() <<
nl <<
endl;
224 td.
pc() =
cloud.constProps().pMin();
229 template<
class ParcelType>
230 template<
class TrackCloudType>
233 TrackCloudType&
cloud,
238 scalar addedMass = 0.0;
239 scalar maxMassI = 0.0;
242 scalar dm =
cloud.rhoTrans(i)[this->
cell()];
243 maxMassI =
max(maxMassI,
mag(dm));
247 if (maxMassI < ROOTVSMALL)
252 const scalar massCell = this->massCell(td);
254 td.rhoc() += addedMass/
cloud.pMesh().cellVolumes()[this->
cell()];
256 const scalar massCellNew = massCell + addedMass;
257 td.Uc() = (td.Uc()*massCell +
cloud.UTrans()[this->
cell()])/massCellNew;
262 scalar
Y =
cloud.rhoTrans(i)[this->
cell()]/addedMass;
263 CpEff +=
Y*
cloud.composition().carrier().Cp(i, td.
pc(), td.Tc());
266 const scalar Cpc = td.CpInterp().psi()[this->
cell()];
267 td.Cpc() = (massCell*Cpc + addedMass*CpEff)/massCellNew;
269 td.Tc() +=
cloud.hsTrans()[this->
cell()]/(td.Cpc()*massCellNew);
271 if (td.Tc() <
cloud.constProps().TMin())
276 <<
"Limiting observed temperature in cell " << this->
cell()
277 <<
" to " <<
cloud.constProps().TMin() <<
nl <<
endl;
280 td.Tc() =
cloud.constProps().TMin();
285 template<
class ParcelType>
286 template<
class TrackCloudType>
289 TrackCloudType&
cloud,
300 if (!
cloud.heatTransfer().BirdCorrection() || (
sum(
Cs) < SMALL))
317 const scalar Xsff = 1.0 -
min(
sum(
Cs)*
RR*this->T_/td.
pc(), 1.0);
320 const scalar CsTot = td.
pc()/(
RR*this->T_);
331 const scalar Csi =
Cs[i] + Xsff*Xinf[i]*CsTot;
333 Xs[i] = (2.0*Csi + Xinf[i]*CsTot)/3.0;
334 Ys[i] = Xs[i]*
thermo.carrier().W(i);
344 scalar sumYiSqrtW = 0;
345 scalar sumYiCbrtW = 0;
349 const scalar W =
thermo.carrier().W(i);
350 const scalar sqrtW =
sqrt(W);
351 const scalar cbrtW =
cbrt(W);
354 mus += Ys[i]*sqrtW*
thermo.carrier().mu(i, td.
pc(),
T);
355 kappas += Ys[i]*cbrtW*
thermo.carrier().kappa(i, td.
pc(),
T);
356 Cps += Xs[i]*
thermo.carrier().Cp(i, td.
pc(),
T);
358 sumYiSqrtW += Ys[i]*sqrtW;
359 sumYiCbrtW += Ys[i]*cbrtW;
362 Cps =
max(Cps, ROOTVSMALL);
364 rhos *= td.
pc()/(
RR*
T);
365 rhos =
max(rhos, ROOTVSMALL);
368 mus =
max(mus, ROOTVSMALL);
370 kappas /= sumYiCbrtW;
371 kappas =
max(kappas, ROOTVSMALL);
373 Prs = Cps*mus/kappas;
377 template<
class ParcelType>
378 template<
class TrackCloudType>
381 TrackCloudType&
cloud,
392 const scalar np0 = this->nParticle_;
393 const scalar d0 = this->d_;
394 const vector& U0 = this->U_;
395 const scalar
T0 = this->T_;
396 const scalar mass0 = this->mass();
400 scalar Ts, rhos, mus, Prs, kappas;
401 this->calcSurfaceValues(
cloud, td,
T0, Ts, rhos, mus, Prs, kappas);
402 scalar Res = this->
Re(rhos, U0, td.Uc(), d0, mus);
424 scalar dhsTrans = 0.0;
473 scalar mass1 = updateMassFraction(mass0, dMass, Y_);
478 if (
cloud.constProps().constantVolume())
480 this->rho_ = mass1/this->
volume();
484 this->d_ =
cbrt(mass1/this->rho_*6.0/
pi);
488 if (np0*mass1 <
cloud.constProps().minParcelMass())
490 td.keepParticle =
false;
492 if (
cloud.solution().coupled())
494 scalar dm = np0*mass0;
499 scalar dmi = dm*Y_[i];
503 cloud.rhoTrans(gid)[this->
cell()] += dmi;
508 cloud.phaseChange().addToPhaseChangeMass(np0*mass1);
515 correctSurfaceValues(
cloud, td, Ts,
Cs, rhos, mus, Prs, kappas);
516 Res = this->
Re(rhos, U0, td.Uc(), this->d(), mus);
527 this->calcHeatTransfer
549 this->calcVelocity(
cloud, td, dt, Res, mus, mass1,
Su, dUTrans, Spu);
555 if (
cloud.solution().coupled())
560 scalar dm = np0*dMass[i];
570 cloud.UTrans()[this->
cell()] += np0*dUTrans;
574 cloud.hsTrans()[this->
cell()] += np0*dhsTrans;
575 cloud.hsCoeff()[this->
cell()] += np0*Sph;
578 if (
cloud.radiation())
580 const scalar ap = this->areaP();
581 const scalar T4 =
pow4(
T0);
582 cloud.radAreaP()[this->
cell()] += dt*np0*ap;
583 cloud.radT4()[this->
cell()] += dt*np0*T4;
584 cloud.radAreaPT4()[this->
cell()] += dt*np0*ap*T4;
int debug
Static debugging option.
const scalar RR
Universal gas constant: default in [J/(kmol K)].
const interpolation< scalar > & pInterp() const
Return const access to the interpolator for continuous phase.
Thermo package for (S)olids (L)iquids and (G)ases Takes reference to thermo package,...
static constexpr const zero Zero
Global zero (0)
const speciesTable & species() const
Return the table of species.
Basic thermodynamics type based on the use of fitting functions for cp, h, s obtained from the templa...
void calc(TrackCloudType &cloud, trackingData &td, const scalar dt)
Update parcel properties over the time interval.
Ostream & endl(Ostream &os)
Add newline and flush stream.
const wordList volume
Standard volume field types (scalar, vector, tensor, etc)
basicSpecieMixture & composition
label min(const labelHashSet &set, label minValue=labelMax)
Find the min value in labelHashSet, optionally limited by second argument.
scalar pc() const
Return the continuous phase pressure.
virtual scalar Hs(const label speciei, const scalar p, const scalar T) const =0
Sensible enthalpy [J/kg].
#define forAll(list, i)
Loop across all elements in list.
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.
dimensionedScalar pow4(const dimensionedScalar &ds)
void setCellValues(TrackCloudType &cloud, trackingData &td)
Set cell values.
PtrList< coordinateSystem > coordinates(solidRegions.size())
dimensionedScalar Pr("Pr", dimless, laminarTransport)
void T(FieldField< Field, Type > &f1, const FieldField< Field, Type > &f2)
label max(const labelHashSet &set, label maxValue=labelMin)
Find the max value in labelHashSet, optionally limited by second argument.
ReactingParcel(const polyMesh &mesh, const barycentric &coordinates, const label celli, const label tetFacei, const label tetPti)
Construct from mesh, coordinates and topology.
void cellValueSourceCorrection(TrackCloudType &cloud, trackingData &td, const scalar dt)
Correct cell values using latest transfer information.
PtrList< volScalarField > & Y
A cloud is a registry collection of lagrangian particles.
constexpr scalar pi(M_PI)
dimensionedScalar sqrt(const dimensionedScalar &ds)
dimensioned< typename typeOfMag< Type >::type > mag(const dimensioned< Type > &dt)
scalarField Re(const UList< complex > &cf)
Extract real component.
const volScalarField & Cp
dimensioned< Type > sum(const DimensionedField< Type, GeoMesh > &df)
const Vector< label > N(dict.get< Vector< label >>("N"))
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.
dimensionedScalar cbrt(const dimensionedScalar &ds)
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.
Reacting parcel class with one/two-way coupling with the continuous phase.
#define WarningInFunction
Report a warning using Foam::Warning.
A cell is defined as a list of faces with extra functionality.
scalar updateMassFraction(const scalar mass0, const scalarField &dMass, scalarField &Y) const
Update mass fraction.
virtual scalar Cp(const label speciei, const scalar p, const scalar T) const =0
Heat capacity at constant pressure [J/(kg K)].