36 template<
class ParcelType>
37 template<
class TrackCloudType>
40 TrackCloudType&
cloud,
60 typedef typename TrackCloudType::thermoCloudType thermoCloudType;
64 typedef typename TrackCloudType::reactingCloudType reactingCloudType;
74 scalar Tvap = phaseChange.
Tvap(X);
81 const scalar TMax = phaseChange.
TMax(td.pc(), X);
82 const scalar Tdash =
min(
T, TMax);
83 const scalar Tsdash =
min(Ts, TMax);
85 const typename TrackCloudType::parcelType&
p =
86 static_cast<const typename TrackCloudType::parcelType&
>(*this);
87 typename TrackCloudType::parcelType::trackingData& ttd =
88 static_cast<typename TrackCloudType::parcelType::trackingData&
>(td);
109 dMassPC =
min(mass*YPhase*
Y, dMassPC);
111 const scalar dMassTot =
sum(dMassPC);
120 const scalar dh = phaseChange.
dh(cid, i, td.pc(), Tdash);
121 Sh -= dMassPC[i]*dh/dt;
126 if (
cloud.heatTransfer().BirdCorrection())
129 const scalar Wc = td.rhoc()*
RR*td.Tc()/td.pc();
135 const scalar
Cp =
composition.carrier().Cp(cid, td.pc(), Tsdash);
137 const scalar Ni = dMassPC[i]/(this->areaS(d)*dt*
W);
140 composition.liquids().properties()[i].D(td.pc(), Tsdash, Wc);
149 Cs[cid] += Ni*d/(2.0*Dab);
155 template<
class ParcelType>
163 scalar mass1 = mass0 -
sum(dMass);
166 if (mass1 > rootVSmall)
170 Y[i] = (
Y[i]*mass0 - dMass[i])/mass1;
180 template<
class ParcelType>
193 template<
class ParcelType>
194 template<
class TrackCloudType>
197 TrackCloudType&
cloud,
201 ParcelType::setCellValues(
cloud, td);
203 td.pc() = td.pInterp().interpolate
206 this->currentTetIndices(td.mesh)
209 if (td.pc() <
cloud.constProps().pMin())
214 <<
"Limiting observed pressure in cell " << this->
cell()
215 <<
" to " <<
cloud.constProps().pMin() <<
nl <<
endl;
218 td.pc() =
cloud.constProps().pMin();
223 template<
class ParcelType>
224 template<
class TrackCloudType>
227 TrackCloudType&
cloud,
232 scalar addedMass = 0.0;
233 scalar maxMassI = 0.0;
236 scalar dm =
cloud.rhoTrans(i)[this->
cell()];
237 maxMassI =
max(maxMassI,
mag(dm));
241 if (maxMassI < rootVSmall)
246 const scalar massCell = this->massCell(td);
248 td.rhoc() += addedMass/
cloud.pMesh().cellVolumes()[this->
cell()];
250 const scalar massCellNew = massCell + addedMass;
251 td.Uc() = (td.Uc()*massCell +
cloud.UTransRef()[this->
cell()])/massCellNew;
256 scalar
Y =
cloud.rhoTrans(i)[this->
cell()]/addedMass;
257 CpEff +=
Y*
cloud.composition().carrier().Cp(i, td.pc(), td.Tc());
260 const scalar Cpc = td.CpInterp().psi()[this->
cell()];
261 td.Cpc() = (massCell*Cpc + addedMass*CpEff)/massCellNew;
263 td.Tc() +=
cloud.hsTransRef()[this->
cell()]/(td.Cpc()*massCellNew);
265 if (td.Tc() <
cloud.constProps().TMin())
270 <<
"Limiting observed temperature in cell " << this->
cell()
271 <<
" to " <<
cloud.constProps().TMin() <<
nl <<
endl;
274 td.Tc() =
cloud.constProps().TMin();
279 template<
class ParcelType>
280 template<
class TrackCloudType>
283 TrackCloudType&
cloud,
294 if (!
cloud.heatTransfer().BirdCorrection() || (
sum(Cs) < small))
306 Xinf[i] = carrier.
Y(i)[this->
cell()]/carrier.
Wi(i);
311 const scalar Xsff = 1.0 -
min(
sum(Cs)*
RR*this->T_/td.pc(), 1.0);
314 const scalar CsTot = td.pc()/(
RR*this->T_);
325 const scalar Csi = Cs[i] + Xsff*Xinf[i]*CsTot;
327 Xs[i] = (2.0*Csi + Xinf[i]*CsTot)/3.0;
328 Ys[i] = Xs[i]*carrier.
Wi(i);
338 scalar sumYiSqrtW = 0;
339 scalar sumYiCbrtW = 0;
343 const scalar
W = carrier.
Wi(i);
344 const scalar sqrtW =
sqrt(
W);
345 const scalar cbrtW =
cbrt(
W);
348 mus += Ys[i]*sqrtW*carrier.
mu(i, td.pc(),
T);
349 kappas += Ys[i]*cbrtW*carrier.
kappa(i, td.pc(),
T);
350 Cps += Xs[i]*carrier.
Cp(i, td.pc(),
T);
352 sumYiSqrtW += Ys[i]*sqrtW;
353 sumYiCbrtW += Ys[i]*cbrtW;
356 Cps =
max(Cps, rootVSmall);
358 rhos *= td.pc()/(
RR*
T);
359 rhos =
max(rhos, rootVSmall);
362 mus =
max(mus, rootVSmall);
364 kappas /= sumYiCbrtW;
365 kappas =
max(kappas, rootVSmall);
367 Prs = Cps*mus/kappas;
371 template<
class ParcelType>
372 template<
class TrackCloudType>
375 TrackCloudType&
cloud,
380 typedef typename TrackCloudType::thermoCloudType thermoCloudType;
388 const scalar np0 = this->nParticle_;
389 const scalar d0 = this->d_;
390 const vector& U0 = this->U_;
391 const scalar
T0 = this->T_;
392 const scalar mass0 = this->mass();
396 scalar Ts, rhos, mus, Prs, kappas;
397 this->calcSurfaceValues(
cloud, td,
T0, Ts, rhos, mus, Prs, kappas);
398 scalar Res = this->
Re(rhos, U0, td.Uc(), d0, mus);
420 scalar dhsTrans = 0.0;
469 scalar mass1 = updateMassFraction(mass0, dMass, Y_);
474 if (
cloud.constProps().constantVolume())
476 this->rho_ = mass1/this->volume();
480 this->d_ =
cbrt(mass1/this->rho_*6.0/
pi);
484 if (np0*mass1 <
cloud.constProps().minParcelMass())
486 td.keepParticle =
false;
488 if (
cloud.solution().coupled())
490 scalar dm = np0*mass0;
495 scalar dmi = dm*Y_[i];
499 cloud.rhoTrans(gid)[this->
cell()] += dmi;
500 cloud.hsTransRef()[this->
cell()] += dmi*hs;
502 cloud.UTransRef()[this->
cell()] += dm*U0;
504 cloud.phaseChange().addToPhaseChangeMass(np0*mass1);
511 correctSurfaceValues(
cloud, td, Ts, Cs, rhos, mus, Prs, kappas);
512 Res = this->
Re(rhos, U0, td.Uc(), this->d(), mus);
523 this->calcHeatTransfer
545 this->calcVelocity(
cloud, td, dt, Res, mus, mass1,
Su, dUTrans, Spu);
551 if (
cloud.solution().coupled())
556 scalar dm = np0*dMass[i];
561 cloud.UTransRef()[this->
cell()] += dm*U0;
562 cloud.hsTransRef()[this->
cell()] += dm*hs;
566 cloud.UTransRef()[this->
cell()] += np0*dUTrans;
567 cloud.UCoeffRef()[this->
cell()] += np0*Spu;
570 cloud.hsTransRef()[this->
cell()] += np0*dhsTrans;
571 cloud.hsCoeffRef()[this->
cell()] += np0*Sph;
574 if (
cloud.radiation())
576 const scalar ap = this->areaP();
577 const scalar T4 =
pow4(
T0);
578 cloud.radAreaP()[this->
cell()] += dt*np0*ap;
579 cloud.radT4()[this->
cell()] += dt*np0*T4;
580 cloud.radAreaPT4()[this->
cell()] += dt*np0*ap*T4;
scalar Cp(const scalar p, const scalar T) const
#define forAll(list, i)
Loop across all elements in list.
Templated reacting parcel composition model class Consists of carrier species (via thermo package),...
void size(const label)
Override size to be inconsistent with allocated storage.
Templated phase change model class.
virtual scalar dh(const label idc, const label idl, const scalar p, const scalar T) const
Return the enthalpy per unit mass.
virtual void calculate(const typename CloudType::parcelType &p, const typename CloudType::parcelType::trackingData &td, const scalar dt, 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.
void addToPhaseChangeMass(const scalar dMass)
Add to phase change mass.
virtual scalar Tvap(const scalarField &X) const
Return vapourisation temperature.
virtual scalar TMax(const scalar p, const scalarField &X) const
Return maximum/limiting temperature.
Reacting parcel class with one/two-way coupling with the continuous phase.
void cellValueSourceCorrection(TrackCloudType &cloud, trackingData &td, const scalar dt)
Correct cell values using latest transfer information.
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.
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.
ParcelType::trackingData trackingData
Use base tracking data.
scalar updateMassFraction(const scalar mass0, const scalarField &dMass, scalarField &Y) const
Update mass fraction.
void setCellValues(TrackCloudType &cloud, trackingData &td)
Set cell values.
ReactingParcel(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 calc(TrackCloudType &cloud, trackingData &td, const scalar dt)
Update parcel properties over the time interval.
Specialisation of basicMixture for a mixture consisting of a number for molecular species.
virtual scalar Hs(const label speciei, const scalar p, const scalar T) const =0
Sensible enthalpy [J/kg].
virtual scalar Wi(const label speciei) const =0
Molecular weight of the given specie [kg/kmol].
const speciesTable & species() const
Return the table of species.
virtual scalar Cp(const label speciei, const scalar p, const scalar T) const =0
Heat capacity at constant pressure [J/kg/K].
PtrList< volScalarField > & Y()
Return the mass-fraction fields.
virtual scalar mu(const label speciei, const scalar p, const scalar T) const =0
Dynamic viscosity [kg/m/s].
virtual scalar kappa(const label speciei, const scalar p, const scalar T) const =0
Thermal conductivity [W/m/K].
A cell is defined as a list of faces with extra functionality.
A cloud is a collection of lagrangian particles.
#define WarningInFunction
Report a warning using Foam::Warning.
const dimensionedScalar RR
Universal gas constant: default SI units: [J/kmol/K].
tmp< VolField< Type > > Su(const VolField< Type > &su, const VolField< Type > &vf)
intWM_LABEL_SIZE_t label
A label is an int32_t or int64_t as specified by the pre-processor macro WM_LABEL_SIZE.
Ostream & endl(Ostream &os)
Add newline and flush stream.
dimensioned< Type > sum(const DimensionedField< Type, GeoMesh > &df)
layerAndWeight min(const layerAndWeight &a, const layerAndWeight &b)
scalarList W(const fluidMulticomponentThermo &thermo)
dimensionedScalar sqrt(const dimensionedScalar &ds)
dimensioned< scalar > mag(const dimensioned< Type > &)
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)
basicSpecieMixture & composition
PtrList< volScalarField > & Y