36 template<
class ParcelType>
37 template<
class TrackCloudType>
40 TrackCloudType& cloud,
60 typedef typename TrackCloudType::reactingCloudType reactingCloudType;
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);
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();
127 const scalar
Cp = composition.
carrier().Cp(cid, td.
pc(), Tsdash);
128 const scalar
W = composition.
carrier().Wi(cid);
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>
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;
240 forAll(cloud.rhoTrans(), i)
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;
260 forAll(cloud.rhoTrans(), i)
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;
344 scalar sumYiSqrtW = 0;
345 scalar sumYiCbrtW = 0;
350 const scalar sqrtW =
sqrt(W);
351 const scalar cbrtW =
cbrt(W);
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,
386 typedef typename TrackCloudType::reactingCloudType reactingCloudType;
394 const scalar np0 = this->nParticle_;
395 const scalar d0 = this->d_;
396 const vector& U0 = this->U_;
397 const scalar
T0 = this->T_;
398 const scalar mass0 = this->mass();
402 scalar Ts, rhos, mus, Prs, kappas;
403 this->calcSurfaceValues(cloud, td, T0, Ts, rhos, mus, Prs, kappas);
404 scalar Res = this->
Re(rhos, U0, td.Uc(), d0, mus);
426 scalar dhsTrans = 0.0;
475 scalar mass1 = updateMassFraction(mass0, dMass, Y_);
477 this->Cp_ = composition.
Cp(0, Y_, td.
pc(),
T0);
480 if (cloud.constProps().constantVolume())
482 this->rho_ = mass1/this->volume();
486 this->d_ =
cbrt(mass1/this->rho_*6.0/
pi);
490 if (np0*mass1 < cloud.constProps().minParcelMass())
492 td.keepParticle =
false;
494 if (cloud.solution().coupled())
496 scalar dm = np0*mass0;
501 scalar dmi = dm*Y_[i];
503 scalar hs = composition.
carrier().Hs(gid, td.
pc(),
T0);
505 cloud.rhoTrans(gid)[this->
cell()] += dmi;
506 cloud.hsTrans()[this->
cell()] += dmi*hs;
508 cloud.UTrans()[this->
cell()] += dm*U0;
510 cloud.phaseChange().addToPhaseChangeMass(np0*mass1);
517 correctSurfaceValues(cloud, td, Ts, Cs, rhos, mus, Prs, kappas);
518 Res = this->
Re(rhos, U0, td.Uc(), this->d(), mus);
529 this->calcHeatTransfer
543 this->Cp_ = composition.
Cp(0, Y_, td.
pc(),
T0);
551 this->calcVelocity(cloud, td, dt, Res, mus, mass1, Su, dUTrans, Spu);
557 if (cloud.solution().coupled())
562 scalar dm = np0*dMass[i];
564 scalar hs = composition.
carrier().Hs(gid, td.
pc(),
T0);
566 cloud.rhoTrans(gid)[this->
cell()] += dm;
567 cloud.UTrans()[this->
cell()] += dm*U0;
568 cloud.hsTrans()[this->
cell()] += dm*hs;
572 cloud.UTrans()[this->
cell()] += np0*dUTrans;
573 cloud.UCoeff()[this->
cell()] += np0*Spu;
576 cloud.hsTrans()[this->
cell()] += np0*dhsTrans;
577 cloud.hsCoeff()[this->
cell()] += np0*Sph;
580 if (cloud.radiation())
582 const scalar ap = this->areaP();
583 const scalar T4 =
pow4(T0);
584 cloud.radAreaP()[this->
cell()] += dt*np0*ap;
585 cloud.radT4()[this->
cell()] += dt*np0*T4;
586 cloud.radAreaPT4()[this->
cell()] += dt*np0*ap*T4;
virtual scalar TMax(const scalar p, const scalarField &X) const
Return maximum/limiting temperature.
#define forAll(list, i)
Loop across all elements in list.
intWM_LABEL_SIZE_t label
A label is an int32_t or int64_t as specified by the pre-processor macro WM_LABEL_SIZE.
Templated phase change model class.
dimensioned< Type > max(const dimensioned< Type > &, const dimensioned< Type > &)
virtual scalar kappa(const label speciei, const scalar p, const scalar T) const =0
Thermal conductivity [W/m/K].
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.
basicSpecieMixture & composition
scalar pc() const
Return the continuous phase pressure.
scalar updateMassFraction(const scalar mass0, const scalarField &dMass, scalarField &Y) const
Update mass fraction.
void size(const label)
Override size to be inconsistent with allocated storage.
scalarField Y_
Mass fractions of mixture [].
dimensionedScalar sqrt(const dimensionedScalar &ds)
Ostream & endl(Ostream &os)
Add newline and flush stream.
const interpolation< scalar > & pInterp() const
Return const access to the interpolator for continuous phase.
ReactingParcel(const polyMesh &mesh, const barycentric &coordinates, const label celli, const label tetFacei, const label tetPti)
Construct from mesh, coordinates and topology.
virtual scalar Wi(const label speciei) const =0
Molecular weight of the given specie [kg/kmol].
label localToCarrierId(const label phaseI, const label id, const bool allowNotFound=false) const
Return carrier id of component given local id.
void addToPhaseChangeMass(const scalar dMass)
Add to phase change mass.
rhoReactionThermo & thermo
virtual scalar Cp(const label phaseI, const scalarField &Y, const scalar p, const scalar T) const
Return specific heat capacity for the phase phaseI.
void cellValueSourceCorrection(TrackCloudType &cloud, trackingData &td, const scalar dt)
Correct cell values using latest transfer information.
dimensioned< Type > sum(const DimensionedField< Type, GeoMesh > &df)
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.
virtual scalar Cp(const label speciei, const scalar p, const scalar T) const =0
Heat capacity at constant pressure [J/kg/K].
virtual scalar mu(const label speciei, const scalar p, const scalar T) const =0
Dynamic viscosity [kg/m/s].
dimensionedScalar cbrt(const dimensionedScalar &ds)
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.
Thermo package for (S)olids (L)iquids and (G)ases Takes reference to thermo package, and provides:
void T(FieldField< Field, Type > &f1, const FieldField< Field, Type > &f2)
dimensioned< Type > min(const dimensioned< Type > &, const dimensioned< Type > &)
const basicSpecieMixture & carrier() const
Return reference to the gaseous components.
virtual scalar dh(const label idc, const label idl, const scalar p, const scalar T) const
Return the enthalpy per unit mass.
PtrList< coordinateSystem > coordinates(solidRegions.size())
#define WarningInFunction
Report a warning using Foam::Warning.
A cell is defined as a list of faces with extra functionality.
const scalarList W(::W(thermo))
dimensionedScalar pow4(const dimensionedScalar &ds)
scalar Cp(const scalar p, const scalar T) const
void setCellValues(TrackCloudType &cloud, trackingData &td)
Set cell values.
dimensioned< scalar > mag(const dimensioned< Type > &)
virtual void calculate(const scalar dt, const label celli, const scalar Re, const scalar Pr, const scalar d, const scalar nu, const scalar T, const scalar Ts, const scalar pc, const scalar Tc, const scalarField &X, scalarField &dMassPC) const =0
Update model.
const scalar RR
Universal gas constant (default in [J/kmol/K])
Mesh consisting of general polyhedral cells.
Reacting parcel class with one/two-way coupling with the continuous phase.
scalar mass0_
Initial mass [kg].
PtrList< volScalarField > & Y()
Return the mass-fraction fields.
Templated reacting parcel composition model class Consists of carrier species (via thermo package)...
virtual scalar Tvap(const scalarField &X) const
Return vapourisation temperature.
virtual bool active() const
Return the model 'active' status - default active = true.
scalarField Re(const UList< complex > &cf)
const basicSpecieMixture & carrier() const
Return the carrier components (wrapper function)
void calc(TrackCloudType &cloud, trackingData &td, const scalar dt)
Update parcel properties over the time interval.
const speciesTable & species() const
Return the table of species.
const liquidMixtureProperties & liquids() const
Return the global (additional) liquids.