39 template<
class ParcelType>
40 template<
class TrackCloudType>
43 TrackCloudType&
cloud,
58 return sum(dVolLiquid);
62 template<
class ParcelType>
63 template<
class TrackCloudType>
66 TrackCloudType&
cloud,
89 auto& phaseChange =
cloud.phaseChange();
92 if (!phaseChange.active())
99 scalar Tvap = phaseChange.Tvap(X);
105 const scalar TMax = phaseChange.TMax(td.pc(), X);
106 const scalar Tdash =
min(
T, TMax);
107 const scalar Tsdash =
min(Ts, TMax);
110 phaseChange.calculate
135 dMassPC[i] =
min(mass*YPhase*
Y[i], dMassPC[i]);
140 const scalar dMassTot =
sum(dMassPC);
143 phaseChange.addToPhaseChangeMass(this->nParticle_*dMassTot);
147 const label cid =
composition.localToCarrierId(idPhase, i);
149 const scalar dh = phaseChange.dh(cid, i, td.pc(), Tdash);
150 Sh -= dMassPC[i]*dh/dt;
155 if (cloud.heatTransfer().BirdCorrection())
158 const scalar Wc = td.rhoc()*
RR*td.Tc()/td.pc();
162 const label cid =
composition.localToCarrierId(idPhase, i);
164 const scalar
Cp =
composition.carrier().Cp(cid, td.pc(), Tsdash);
166 const scalar Ni = dMassPC[i]/(this->areaS(d)*dt*W);
169 composition.liquids().properties()[i].D(td.pc(), Tsdash, Wc);
178 Cs[cid] += Ni*d/(2.0*Dab);
184 template<
class ParcelType>
192 scalar mass1 = mass0 -
sum(dMass);
195 if (mass1 > ROOTVSMALL)
199 Y[i] = (
Y[i]*mass0 - dMass[i])/mass1;
209 template<
class ParcelType>
212 const ReactingParcel<ParcelType>&
p 221 template<
class ParcelType>
224 const ReactingParcel<ParcelType>&
p,
236 template<
class ParcelType>
237 template<
class TrackCloudType>
240 TrackCloudType&
cloud,
244 ParcelType::setCellValues(
cloud, td);
249 this->currentTetIndices()
252 if (td.
pc() <
cloud.constProps().pMin())
257 <<
"Limiting observed pressure in cell " << this->
cell()
258 <<
" to " <<
cloud.constProps().pMin() <<
nl <<
endl;
261 td.
pc() =
cloud.constProps().pMin();
266 template<
class ParcelType>
267 template<
class TrackCloudType>
270 TrackCloudType& cloud,
275 scalar addedMass = 0.0;
276 scalar maxMassI = 0.0;
277 forAll(cloud.rhoTrans(), i)
279 scalar dm = cloud.rhoTrans(i)[this->cell()];
280 maxMassI =
max(maxMassI,
mag(dm));
284 if (maxMassI < ROOTVSMALL)
289 const scalar massCell = this->massCell(td);
291 td.rhoc() += addedMass/cloud.pMesh().cellVolumes()[this->cell()];
293 const scalar massCellNew = massCell + addedMass;
294 td.Uc() = (td.Uc()*massCell + cloud.UTrans()[this->cell()])/massCellNew;
297 forAll(cloud.rhoTrans(), i)
299 scalar
Y = cloud.rhoTrans(i)[this->cell()]/addedMass;
300 CpEff +=
Y*cloud.composition().carrier().Cp(i, td.pc(), td.Tc());
303 const scalar Cpc = td.CpInterp().psi()[this->cell()];
304 td.Cpc() = (massCell*Cpc + addedMass*CpEff)/massCellNew;
306 td.Tc() += cloud.hsTrans()[this->cell()]/(td.Cpc()*massCellNew);
308 if (td.Tc() < cloud.constProps().TMin())
313 <<
"Limiting observed temperature in cell " << this->cell()
314 <<
" to " << cloud.constProps().TMin() <<
nl <<
endl;
317 td.Tc() = cloud.constProps().TMin();
322 template<
class ParcelType>
323 template<
class TrackCloudType>
326 TrackCloudType& cloud,
337 if (!cloud.heatTransfer().BirdCorrection() || (
sum(
Cs) < SMALL))
342 const SLGThermo&
thermo = cloud.thermo();
349 Xinf[i] =
thermo.carrier().Y(i)[this->cell()]/
thermo.carrier().
W(i);
354 const scalar Xsff = 1.0 -
min(
sum(
Cs)*
RR*this->T_/td.pc(), 1.0);
357 const scalar CsTot = td.pc()/(
RR*this->T_);
368 const scalar Csi =
Cs[i] + Xsff*Xinf[i]*CsTot;
370 Xs[i] = (2.0*Csi + Xinf[i]*CsTot)/3.0;
371 Ys[i] = Xs[i]*
thermo.carrier().
W(i);
381 scalar sumYiSqrtW = 0;
382 scalar sumYiCbrtW = 0;
386 const scalar W =
thermo.carrier().
W(i);
387 const scalar sqrtW =
sqrt(W);
388 const scalar cbrtW =
cbrt(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);
395 sumYiSqrtW += Ys[i]*sqrtW;
396 sumYiCbrtW += Ys[i]*cbrtW;
399 Cps =
max(Cps, ROOTVSMALL);
401 rhos *= td.pc()/(
RR*
T);
402 rhos =
max(rhos, ROOTVSMALL);
405 mus =
max(mus, ROOTVSMALL);
407 kappas /= sumYiCbrtW;
408 kappas =
max(kappas, ROOTVSMALL);
410 Prs = Cps*mus/kappas;
414 template<
class ParcelType>
415 template<
class TrackCloudType>
418 TrackCloudType&
cloud,
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_;
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);
462 scalar dhsTrans = 0.0;
513 scalar mass1 = updateMassFraction(mass0, dMass, Y_);
519 cloud.constProps().volUpdateType()
520 == constantProperties::volumeUpdateType::mUndefined
524 if (
cloud.constProps().constantVolume())
526 this->rho_ = mass1/this->
volume();
530 this->d_ =
cbrt(mass1/this->rho_*6/
pi);
535 switch (cloud.constProps().volUpdateType())
537 case constantProperties::volumeUpdateType::mConstRho :
539 this->d_ =
cbrt(mass1/this->rho_*6/
pi);
542 case constantProperties::volumeUpdateType::mConstVol :
544 this->rho_ = mass1/this->
volume();
547 case constantProperties::volumeUpdateType::mUpdateRhoAndVol :
557 this->rho_ = mass1/(this->
volume() - deltaVol);
558 this->d_ =
cbrt(mass1/this->rho_*6/
pi);
566 if (np0*mass1 < cloud.constProps().minParcelMass())
568 td.keepParticle =
false;
570 if (cloud.solution().coupled())
572 scalar dm = np0*mass0;
577 scalar dmi = dm*Y_[i];
581 cloud.rhoTrans(gid)[this->cell()] += dmi;
582 cloud.hsTrans()[this->cell()] += dmi*hs;
584 cloud.UTrans()[this->cell()] += dm*U0;
586 cloud.phaseChange().addToPhaseChangeMass(np0*mass1);
593 correctSurfaceValues(cloud, td, Ts,
Cs, rhos, mus, Prs, kappas);
594 Res = this->
Re(rhos, U0, td.Uc(), this->d(), mus);
605 this->calcHeatTransfer
627 this->calcVelocity(cloud, td, dt, Res, mus, mass1,
Su, dUTrans, Spu);
633 if (cloud.solution().coupled())
638 scalar dm = np0*dMass[i];
642 cloud.rhoTrans(gid)[this->cell()] += dm;
643 cloud.UTrans()[this->cell()] += dm*U0;
644 cloud.hsTrans()[this->cell()] += dm*hs;
648 cloud.UTrans()[this->cell()] += np0*dUTrans;
649 cloud.UCoeff()[this->cell()] += np0*Spu;
652 cloud.hsTrans()[this->cell()] += np0*dhsTrans;
653 cloud.hsCoeff()[this->cell()] += np0*Sph;
656 if (cloud.radiation())
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;
dimensionedScalar Pr("Pr", dimless, laminarTransport)
void size(const label n)
Older name for setAddressableSize.
virtual tmp< volScalarField > kappa() const =0
Thermal diffusivity for temperature of mixture [J/m/s/K].
scalar updatedDeltaVolume(TrackCloudType &cloud, const scalarField &dMass, const scalar p, const scalar T)
Return change of volume due to mass exchange.
dimensioned< typename typeOfMag< Type >::type > mag(const dimensioned< Type > &dt)
basicSpecieMixture & composition
scalar pc() const
Return the continuous phase pressure.
label max(const labelHashSet &set, label maxValue=labelMin)
Find the max value in labelHashSet, optionally limited by second argument.
constexpr char nl
The newline '\n' character (0x0a)
scalar updateMassFraction(const scalar mass0, const scalarField &dMass, scalarField &Y) const
Update mass fraction.
dimensionedScalar sqrt(const dimensionedScalar &ds)
Ostream & endl(Ostream &os)
Add newline and flush stream.
const speciesTable & species() const
Return the table of species.
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.
void cellValueSourceCorrection(TrackCloudType &cloud, trackingData &td, const scalar dt)
Correct cell values using latest transfer information.
dimensioned< Type > sum(const DimensionedField< Type, GeoMesh > &df)
#define forAll(list, i)
Loop across all elements in list.
psiReactionThermo & thermo
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)].
Field< scalar > scalarField
Specialisation of Field<T> for scalar.
A cloud is a registry collection of lagrangian particles.
dimensionedScalar cbrt(const dimensionedScalar &ds)
virtual tmp< volScalarField > mu() const
Dynamic viscosity of mixture [kg/m/s].
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.
constexpr scalar pi(M_PI)
scalarField Re(const UList< complex > &cf)
Extract real component.
label min(const labelHashSet &set, label minValue=labelMax)
Find the min value in labelHashSet, optionally limited by second argument.
const wordList volume
Standard volume field types (scalar, vector, tensor, etc)
A Vector of values with scalar precision, where scalar is float/double depending on the compilation f...
const volScalarField & Cp
int debug
Static debugging option.
void T(FieldField< Field, Type > &f1, const FieldField< Field, Type > &f2)
virtual tmp< volScalarField > Cp() const =0
Heat capacity at constant pressure [J/kg/K].
const Vector< label > N(dict.get< Vector< label >>("N"))
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.
PtrList< volScalarField > & Y
virtual scalar Hs(const label speciei, const scalar p, const scalar T) const =0
Sensible enthalpy [J/kg].
dimensionedScalar pow4(const dimensionedScalar &ds)
void setCellValues(TrackCloudType &cloud, trackingData &td)
Set cell values.
virtual tmp< volScalarField > W() const =0
Molecular weight [kg/kmol].
const scalar RR
Universal gas constant: default in [J/(kmol K)].
void calc(TrackCloudType &cloud, trackingData &td, const scalar dt)
Update parcel properties over the time interval.
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.
static constexpr const zero Zero
Global zero (0)