zhyang-dev/OpenFOAM-6
1
0
Fork 0
mirror of https://github.com/OpenFOAM/OpenFOAM-6.git synced 2025-12-08 06:57:46 +00:00
Code Issues Actions Packages Projects Releases Wiki Activity

Reduce the stripped parcel mass due to evaporation

Browse Source 
Resolves bug-report http://www.openfoam.org/mantisbt/view.php?id=1295
This commit is contained in:
Henry
2014-12-18 12:50:24 +00:00
parent ca31aaec32
commit 970884a561
1 changed files with 26 additions and 21 deletions
Download Patch File Download Diff File Expand all files Collapse all files

47
src/lagrangian/spray/parcels/Templates/SprayParcel/SprayParcel.C
View File

@ -68,42 +68,47 @@ void Foam::SprayParcel<ParcelType>::calc
const CompositionModel<reactingCloudType>& composition = const CompositionModel<reactingCloudType>& composition =
td.cloud().composition(); td.cloud().composition();
// check if parcel belongs to liquid core // Check if parcel belongs to liquid core
if (liquidCore() > 0.5) if (liquidCore() > 0.5)
{ {
// liquid core parcels should not experience coupled forces // Liquid core parcels should not experience coupled forces
td.cloud().forces().setCalcCoupled(false); td.cloud().forces().setCalcCoupled(false);
} }
// get old mixture composition // Get old mixture composition
const scalarField& Y0(this->Y()); const scalarField& Y0(this->Y());
scalarField X0(composition.liquids().X(Y0)); scalarField X0(composition.liquids().X(Y0));
// check if we have critical or boiling conditions // Check if we have critical or boiling conditions
scalar TMax = composition.liquids().Tc(X0); scalar TMax = composition.liquids().Tc(X0);
const scalar T0 = this->T(); const scalar T0 = this->T();
const scalar pc0 = this->pc_; const scalar pc0 = this->pc_;
if (composition.liquids().pv(pc0, T0, X0) >= pc0*0.999) if (composition.liquids().pv(pc0, T0, X0) >= pc0*0.999)
{ {
// set TMax to boiling temperature // Set TMax to boiling temperature
TMax = composition.liquids().pvInvert(pc0, X0); TMax = composition.liquids().pvInvert(pc0, X0);
} }
// set the maximum temperature limit // Set the maximum temperature limit
td.cloud().constProps().setTMax(TMax); td.cloud().constProps().setTMax(TMax);
// store the parcel properties // Store the parcel properties
this->Cp() = composition.liquids().Cp(pc0, T0, X0); this->Cp() = composition.liquids().Cp(pc0, T0, X0);
sigma_ = composition.liquids().sigma(pc0, T0, X0); sigma_ = composition.liquids().sigma(pc0, T0, X0);
scalar rho0 = composition.liquids().rho(pc0, T0, X0); const scalar rho0 = composition.liquids().rho(pc0, T0, X0);
this->rho() = rho0; this->rho() = rho0;
const scalar mass0 = this->mass();
mu_ = composition.liquids().mu(pc0, T0, X0); mu_ = composition.liquids().mu(pc0, T0, X0);
ParcelType::calc(td, dt, cellI); ParcelType::calc(td, dt, cellI);
if (td.keepParticle) if (td.keepParticle)
{ {
// update Cp, sigma, density and diameter due to change in temperature // Reduce the stripped parcel mass due to evaporation
// assuming the number of particles remains unchanged
this->ms() -= this->ms()*(mass0 - this->mass())/mass0;
// Update Cp, sigma, density and diameter due to change in temperature
// and/or composition // and/or composition
scalar T1 = this->T(); scalar T1 = this->T();
const scalarField& Y1(this->Y()); const scalarField& Y1(this->Y());
@ -125,7 +130,7 @@ void Foam::SprayParcel<ParcelType>::calc
{ {
calcAtomization(td, dt, cellI); calcAtomization(td, dt, cellI);
// preserve the total mass/volume by increasing the number of // Preserve the total mass/volume by increasing the number of
// particles in parcels due to breakup // particles in parcels due to breakup
scalar d2 = this->d(); scalar d2 = this->d();
this->nParticle() *= pow3(d1/d2); this->nParticle() *= pow3(d1/d2);
@ -136,7 +141,7 @@ void Foam::SprayParcel<ParcelType>::calc
} }
} }
// restore coupled forces // Restore coupled forces
td.cloud().forces().setCalcCoupled(true); td.cloud().forces().setCalcCoupled(true);
} }
@ -163,7 +168,7 @@ void Foam::SprayParcel<ParcelType>::calcAtomization
scalar R = specie::RR/Wc; scalar R = specie::RR/Wc;
scalar Tav = atomization.Taverage(this->T(), this->Tc()); scalar Tav = atomization.Taverage(this->T(), this->Tc());
// calculate average gas density based on average temperature // Calculate average gas density based on average temperature
scalar rhoAv = this->pc()/(R*Tav); scalar rhoAv = this->pc()/(R*Tav);
scalar soi = td.cloud().injectors().timeStart(); scalar soi = td.cloud().injectors().timeStart();
@ -171,14 +176,14 @@ void Foam::SprayParcel<ParcelType>::calcAtomization
const vector& pos = this->position(); const vector& pos = this->position();
const vector& injectionPos = this->position0(); const vector& injectionPos = this->position0();
// disregard the continous phase when calculating the relative velocity // Disregard the continous phase when calculating the relative velocity
// (in line with the deactivated coupled assumption) // (in line with the deactivated coupled assumption)
scalar Urel = mag(this->U()); scalar Urel = mag(this->U());
scalar t0 = max(0.0, currentTime - this->age() - soi); scalar t0 = max(0.0, currentTime - this->age() - soi);
scalar t1 = min(t0 + dt, td.cloud().injectors().timeEnd() - soi); scalar t1 = min(t0 + dt, td.cloud().injectors().timeEnd() - soi);
// this should be the vol flow rate from when the parcel was injected // This should be the vol flow rate from when the parcel was injected
scalar volFlowRate = td.cloud().injectors().volumeToInject(t0, t1)/dt; scalar volFlowRate = td.cloud().injectors().volumeToInject(t0, t1)/dt;
scalar chi = 0.0; scalar chi = 0.0;
@ -234,7 +239,7 @@ void Foam::SprayParcel<ParcelType>::calcBreakup
scalar R = specie::RR/Wc; scalar R = specie::RR/Wc;
scalar Tav = td.cloud().atomization().Taverage(this->T(), this->Tc()); scalar Tav = td.cloud().atomization().Taverage(this->T(), this->Tc());
// calculate average gas density based on average temperature // Calculate average gas density based on average temperature
scalar rhoAv = this->pc()/(R*Tav); scalar rhoAv = this->pc()/(R*Tav);
scalar muAv = this->muc(); scalar muAv = this->muc();
vector Urel = this->U() - this->Uc(); vector Urel = this->U() - this->Uc();
@ -316,7 +321,7 @@ Foam::scalar Foam::SprayParcel<ParcelType>::chi
const scalarField& X const scalarField& X
) const ) const
{ {
// modifications to take account of the flash boiling on primary break-up // Modifications to take account of the flash boiling on primary break-up
typedef typename TrackData::cloudType::reactingCloudType reactingCloudType; typedef typename TrackData::cloudType::reactingCloudType reactingCloudType;
const CompositionModel<reactingCloudType>& composition = const CompositionModel<reactingCloudType>& composition =
@ -333,7 +338,7 @@ Foam::scalar Foam::SprayParcel<ParcelType>::chi
{ {
if (pv >= 0.999*pAmb) if (pv >= 0.999*pAmb)
{ {
// liquid is boiling - calc boiling temperature // Liquid is boiling - calc boiling temperature
const liquidProperties& liq = composition.liquids().properties()[i]; const liquidProperties& liq = composition.liquids().properties()[i];
scalar TBoil = liq.pvInvert(p0); scalar TBoil = liq.pvInvert(p0);
@ -368,10 +373,10 @@ void Foam::SprayParcel<ParcelType>::solveTABEq
scalar r2 = r*r; scalar r2 = r*r;
scalar r3 = r*r2; scalar r3 = r*r2;
// inverse of characteristic viscous damping time // Inverse of characteristic viscous damping time
scalar rtd = 0.5*TABCmu*mu_/(this->rho()*r2); scalar rtd = 0.5*TABCmu*mu_/(this->rho()*r2);
// oscillation frequency (squared) // Oscillation frequency (squared)
scalar omega2 = TABComega*sigma_/(this->rho()*r3) - rtd*rtd; scalar omega2 = TABComega*sigma_/(this->rho()*r3) - rtd*rtd;
if (omega2 > 0) if (omega2 > 0)
@ -383,7 +388,7 @@ void Foam::SprayParcel<ParcelType>::solveTABEq
scalar y1 = this->y() - We; scalar y1 = this->y() - We;
scalar y2 = this->yDot()/omega; scalar y2 = this->yDot()/omega;
// update distortion parameters // Update distortion parameters
scalar c = cos(omega*dt); scalar c = cos(omega*dt);
scalar s = sin(omega*dt); scalar s = sin(omega*dt);
scalar e = exp(-rtd*dt); scalar e = exp(-rtd*dt);
@ -402,7 +407,7 @@ void Foam::SprayParcel<ParcelType>::solveTABEq
} }
else else
{ {
// reset distortion parameters // Reset distortion parameters
this->y() = 0; this->y() = 0;
this->yDot() = 0; this->yDot() = 0;
} }