From 2ac008348bab9270f13d9685adf162ae7db1180c Mon Sep 17 00:00:00 2001 From: andy <andy> Date: Fri, 6 May 2011 12:29:43 +0100 Subject: [PATCH] ENH: Updated cloud source term computations --- .../KinematicParcel/KinematicParcel.C | 36 +--- .../KinematicParcel/KinematicParcel.H | 3 - .../ReactingMultiphaseParcel.C | 181 ++++++++---------- .../Templates/ReactingParcel/ReactingParcel.C | 156 +++++++-------- .../Templates/ThermoParcel/ThermoParcel.C | 69 +++---- .../Templates/ThermoParcel/ThermoParcel.H | 4 - 6 files changed, 174 insertions(+), 275 deletions(-) diff --git a/src/lagrangian/intermediate/parcels/Templates/KinematicParcel/KinematicParcel.C b/src/lagrangian/intermediate/parcels/Templates/KinematicParcel/KinematicParcel.C index b75b20f8456..9d1a918bbc4 100644 --- a/src/lagrangian/intermediate/parcels/Templates/KinematicParcel/KinematicParcel.C +++ b/src/lagrangian/intermediate/parcels/Templates/KinematicParcel/KinematicParcel.C @@ -104,13 +104,10 @@ void Foam::KinematicParcel<ParcelType>::calc // Define local properties at beginning of time step // ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~ const scalar np0 = nParticle_; - const scalar d0 = d_; - const vector U0 = U_; - const scalar rho0 = rho_; const scalar mass0 = mass(); // Reynolds number - const scalar Re = this->Re(U0, d0, rhoc_, muc_); + const scalar Re = this->Re(U_, d_, rhoc_, muc_); // Sources @@ -119,6 +116,9 @@ void Foam::KinematicParcel<ParcelType>::calc // Explicit momentum source for particle vector Su = vector::zero; + // Linearised momentum source coefficient + scalar Spu = 0.0; + // Momentum transfer from the particle to the carrier phase vector dUTrans = vector::zero; @@ -127,23 +127,7 @@ void Foam::KinematicParcel<ParcelType>::calc // ~~~~~~ // Calculate new particle velocity - scalar Spu = 0.0; - vector U1 = - calcVelocity - ( - td, - dt, - cellI, - Re, - muc_, - d0, - U0, - rho0, - mass0, - Su, - dUTrans, - Spu - ); + this->U_ = calcVelocity(td, dt, cellI, Re, muc_, mass0, Su, dUTrans, Spu); // Accumulate carrier phase source terms @@ -156,11 +140,6 @@ void Foam::KinematicParcel<ParcelType>::calc // Update momentum transfer coefficient td.cloud().UCoeff()[cellI] += np0*Spu; } - - - // Set new particle properties - // ~~~~~~~~~~~~~~~~~~~~~~~~~~~ - U_ = U1; } @@ -173,9 +152,6 @@ const Foam::vector Foam::KinematicParcel<ParcelType>::calcVelocity const label cellI, const scalar Re, const scalar mu, - const scalar d, - const vector& U, - const scalar rho, const scalar mass, const vector& Su, vector& dUTrans, @@ -205,7 +181,7 @@ const Foam::vector Foam::KinematicParcel<ParcelType>::calcVelocity Spu = dt*Feff.Sp(); IntegrationScheme<vector>::integrationResult Ures = - td.cloud().UIntegrator().integrate(U, dt, abp, bp); + td.cloud().UIntegrator().integrate(U_, dt, abp, bp); vector Unew = Ures.value(); diff --git a/src/lagrangian/intermediate/parcels/Templates/KinematicParcel/KinematicParcel.H b/src/lagrangian/intermediate/parcels/Templates/KinematicParcel/KinematicParcel.H index a937ba3f7a3..19b953b7e78 100644 --- a/src/lagrangian/intermediate/parcels/Templates/KinematicParcel/KinematicParcel.H +++ b/src/lagrangian/intermediate/parcels/Templates/KinematicParcel/KinematicParcel.H @@ -297,9 +297,6 @@ protected: const label cellI, // owner cell const scalar Re, // Reynolds number const scalar mu, // local carrier viscosity - const scalar d, // diameter - const vector& U, // velocity - const scalar rho, // density const scalar mass, // mass const vector& Su, // explicit particle momentum source vector& dUTrans, // momentum transfer to carrier diff --git a/src/lagrangian/intermediate/parcels/Templates/ReactingMultiphaseParcel/ReactingMultiphaseParcel.C b/src/lagrangian/intermediate/parcels/Templates/ReactingMultiphaseParcel/ReactingMultiphaseParcel.C index 664fac97810..ae95fa21ef2 100644 --- a/src/lagrangian/intermediate/parcels/Templates/ReactingMultiphaseParcel/ReactingMultiphaseParcel.C +++ b/src/lagrangian/intermediate/parcels/Templates/ReactingMultiphaseParcel/ReactingMultiphaseParcel.C @@ -172,12 +172,11 @@ void Foam::ReactingMultiphaseParcel<ParcelType>::calc // Define local properties at beginning of timestep // ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~ + const scalar np0 = this->nParticle_; const scalar d0 = this->d_; const vector& U0 = this->U_; - const scalar rho0 = this->rho_; const scalar T0 = this->T_; - const scalar Cp0 = this->Cp_; const scalar mass0 = this->mass(); const scalar pc = this->pc_; @@ -189,11 +188,8 @@ void Foam::ReactingMultiphaseParcel<ParcelType>::calc // Calc surface values - // ~~~~~~~~~~~~~~~~~~~ scalar Ts, rhos, mus, Prs, kappas; this->calcSurfaceValues(td, cellI, T0, Ts, rhos, mus, Prs, kappas); - - // Reynolds number scalar Res = this->Re(U0, d0, rhos, mus); @@ -203,16 +199,25 @@ void Foam::ReactingMultiphaseParcel<ParcelType>::calc // Explicit momentum source for particle vector Su = vector::zero; + // Linearised momentum source coefficient + scalar Spu = 0.0; + // Momentum transfer from the particle to the carrier phase vector dUTrans = vector::zero; // Explicit enthalpy source for particle scalar Sh = 0.0; + // Linearised enthalpy source coefficient + scalar Sph = 0.0; + // Sensible enthalpy transfer from the particle to the carrier phase scalar dhsTrans = 0.0; + // 1. Compute models that contribute to mass transfer - U, T held constant + // ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~ + // Phase change in liquid phase // ~~~~~~~~~~~~~~~~~~~~~~~~~~~~ @@ -310,28 +315,78 @@ void Foam::ReactingMultiphaseParcel<ParcelType>::calc ); - // Correct surface values due to emitted species - this->correctSurfaceValues(td, cellI, Ts, Cs, rhos, mus, Prs, kappas); - Res = this->Re(U0, d0, rhos, mus); - - - // Update component mass fractions - // ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~ + // 2. Update the parcel properties due to change in mass + // ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~ scalarField dMassGas(dMassDV + dMassSRGas); scalarField dMassLiquid(dMassPC + dMassSRLiquid); scalarField dMassSolid(dMassSRSolid); - scalar mass1 = updateMassFractions(mass0, dMassGas, dMassLiquid, dMassSolid); + this->Cp_ = CpEff(td, pc, T0, idG, idL, idS); + + // Update particle density or diameter + if (td.cloud().constProps().constantVolume()) + { + this->rho_ = mass1/this->volume(); + } + else + { + this->d_ = cbrt(mass1/this->rho_*6.0/pi); + } + + // Remove the particle when mass falls below minimum threshold + if (np0*mass1 < td.cloud().constProps().minParticleMass()) + { + td.keepParticle = false; + + if (td.cloud().solution().coupled()) + { + scalar dm = np0*mass1; + + // Absorb parcel into carrier phase + forAll(YGas_, i) + { + label gid = composition.localToGlobalCarrierId(GAS, i); + td.cloud().rhoTrans(gid)[cellI] += dm*YMix[GAS]*YGas_[i]; + } + forAll(YLiquid_, i) + { + label gid = composition.localToGlobalCarrierId(LIQ, i); + td.cloud().rhoTrans(gid)[cellI] += dm*YMix[LIQ]*YLiquid_[i]; + } +/* + // No mapping between solid components and carrier phase + forAll(YSolid_, i) + { + label gid = composition.localToGlobalCarrierId(SLD, i); + td.cloud().rhoTrans(gid)[cellI] += dm*YMix[SLD]*YSolid_[i]; + } +*/ + td.cloud().UTrans()[cellI] += dm*U0; + + td.cloud().hsTrans()[cellI] += dm*HsEff(td, pc, T0, idG, idL, idS); + + td.cloud().addToMassPhaseChange(dm); + } + + return; + } + + // Correct surface values due to emitted species + this->correctSurfaceValues(td, cellI, Ts, Cs, rhos, mus, Prs, kappas); + Res = this->Re(U0, this->d_, rhos, mus); + + + // 3. Compute heat- and momentum transfers + // ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~ // Heat transfer // ~~~~~~~~~~~~~ // Calculate new particle temperature - scalar Sph = 0.0; - scalar T1 = + this->T_ = this->calcHeatTransfer ( td, @@ -340,10 +395,6 @@ void Foam::ReactingMultiphaseParcel<ParcelType>::calc Res, Prs, kappas, - d0, - rho0, - T0, - Cp0, NCpW, Sh, dhsTrans, @@ -351,35 +402,23 @@ void Foam::ReactingMultiphaseParcel<ParcelType>::calc ); + this->Cp_ = CpEff(td, pc, this->T_, idG, idL, idS); + + // Motion // ~~~~~~ // Calculate new particle velocity - scalar Spu = 0; - vector U1 = - this->calcVelocity - ( - td, - dt, - cellI, - Res, - mus, - d0, - U0, - rho0, - mass0, - Su, - dUTrans, - Spu - ); + this->U_ = + this->calcVelocity(td, dt, cellI, Res, mus, mass1, Su, dUTrans, Spu); - // Accumulate carrier phase source terms - // ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~ + // 4. Accumulate carrier phase source terms + // ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~ if (td.cloud().solution().coupled()) { - // Transfer mass lost from particle to carrier mass source + // Transfer mass lost to carrier mass, momentum and enthalpy sources forAll(YGas_, i) { scalar dm = np0*dMassGas[i]; @@ -421,76 +460,12 @@ void Foam::ReactingMultiphaseParcel<ParcelType>::calc // Update momentum transfer td.cloud().UTrans()[cellI] += np0*dUTrans; - - // Update momentum transfer coefficient td.cloud().UCoeff()[cellI] += np0*Spu; // Update sensible enthalpy transfer td.cloud().hsTrans()[cellI] += np0*dhsTrans; - - // Update sensible enthalpy coefficient td.cloud().hsCoeff()[cellI] += np0*Sph; } - - - // Remove the particle when mass falls below minimum threshold - // ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~ - - if (np0*mass1 < td.cloud().constProps().minParticleMass()) - { - td.keepParticle = false; - - if (td.cloud().solution().coupled()) - { - scalar dm = np0*mass1; - - // Absorb parcel into carrier phase - forAll(YGas_, i) - { - label gid = composition.localToGlobalCarrierId(GAS, i); - td.cloud().rhoTrans(gid)[cellI] += dm*YMix[GAS]*YGas_[i]; - } - forAll(YLiquid_, i) - { - label gid = composition.localToGlobalCarrierId(LIQ, i); - td.cloud().rhoTrans(gid)[cellI] += dm*YMix[LIQ]*YLiquid_[i]; - } -/* - // No mapping between solid components and carrier phase - forAll(YSolid_, i) - { - label gid = composition.localToGlobalCarrierId(SLD, i); - td.cloud().rhoTrans(gid)[cellI] += dm*YMix[SLD]*YSolid_[i]; - } -*/ - td.cloud().UTrans()[cellI] += dm*U1; - - td.cloud().hsTrans()[cellI] += dm*HsEff(td, pc, T1, idG, idL, idS); - - td.cloud().addToMassPhaseChange(dm); - } - } - - - // Set new particle properties - // ~~~~~~~~~~~~~~~~~~~~~~~~~~~ - - else - { - this->Cp_ = CpEff(td, pc, T1, idG, idL, idS); - this->T_ = T1; - this->U_ = U1; - - // Update particle density or diameter - if (td.cloud().constProps().constantVolume()) - { - this->rho_ = mass1/this->volume(); - } - else - { - this->d_ = cbrt(mass1/this->rho_*6.0/pi); - } - } } diff --git a/src/lagrangian/intermediate/parcels/Templates/ReactingParcel/ReactingParcel.C b/src/lagrangian/intermediate/parcels/Templates/ReactingParcel/ReactingParcel.C index 71330144640..21e7768aa77 100644 --- a/src/lagrangian/intermediate/parcels/Templates/ReactingParcel/ReactingParcel.C +++ b/src/lagrangian/intermediate/parcels/Templates/ReactingParcel/ReactingParcel.C @@ -267,21 +267,17 @@ void Foam::ReactingParcel<ParcelType>::calc // Define local properties at beginning of time step // ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~ + const scalar np0 = this->nParticle_; const scalar d0 = this->d_; const vector& U0 = this->U_; - const scalar rho0 = this->rho_; const scalar T0 = this->T_; - const scalar Cp0 = this->Cp_; const scalar mass0 = this->mass(); // Calc surface values - // ~~~~~~~~~~~~~~~~~~~ scalar Ts, rhos, mus, Prs, kappas; this->calcSurfaceValues(td, cellI, T0, Ts, rhos, mus, Prs, kappas); - - // Reynolds number scalar Res = this->Re(U0, d0, rhos, mus); @@ -291,16 +287,25 @@ void Foam::ReactingParcel<ParcelType>::calc // Explicit momentum source for particle vector Su = vector::zero; + // Linearised momentum source coefficient + scalar Spu = 0.0; + // Momentum transfer from the particle to the carrier phase vector dUTrans = vector::zero; // Explicit enthalpy source for particle scalar Sh = 0.0; + // Linearised enthalpy source coefficient + scalar Sph = 0.0; + // Sensible enthalpy transfer from the particle to the carrier phase scalar dhsTrans = 0.0; + // 1. Compute models that contribute to mass transfer - U, T held constant + // ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~ + // Phase change // ~~~~~~~~~~~~ @@ -338,22 +343,65 @@ void Foam::ReactingParcel<ParcelType>::calc Cs ); - // Correct surface values due to emitted species - correctSurfaceValues(td, cellI, Ts, Cs, rhos, mus, Prs, kappas); - Res = this->Re(U0, d0, rhos, mus); - // Update particle component mass and mass fractions - scalarField dMass(dMassPC); + // 2. Update the parcel properties due to change in mass + // ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~ + scalarField dMass(dMassPC); scalar mass1 = updateMassFraction(mass0, dMass, Y_); + this->Cp_ = composition.Cp(0, Y_, pc_, T0); + + // Update particle density or diameter + if (td.cloud().constProps().constantVolume()) + { + this->rho_ = mass1/this->volume(); + } + else + { + this->d_ = cbrt(mass1/this->rho_*6.0/pi); + } + + // Remove the particle when mass falls below minimum threshold + if (np0*mass1 < td.cloud().constProps().minParticleMass()) + { + td.keepParticle = false; + + if (td.cloud().solution().coupled()) + { + scalar dm = np0*mass1; + + // Absorb parcel into carrier phase + forAll(Y_, i) + { + scalar dmi = dm*Y_[i]; + label gid = composition.localToGlobalCarrierId(0, i); + scalar hs = composition.carrier().Hs(gid, T0); + + td.cloud().rhoTrans(gid)[cellI] += dmi; + td.cloud().hsTrans()[cellI] += dmi*hs; + } + td.cloud().UTrans()[cellI] += dm*U0; + + td.cloud().addToMassPhaseChange(dm); + } + + return; + } + + // Correct surface values due to emitted species + correctSurfaceValues(td, cellI, Ts, Cs, rhos, mus, Prs, kappas); + Res = this->Re(U0, this->d_, rhos, mus); + + + // 3. Compute heat- and momentum transfers + // ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~ // Heat transfer // ~~~~~~~~~~~~~ // Calculate new particle temperature - scalar Sph = 0.0; - scalar T1 = + this->T_ = this->calcHeatTransfer ( td, @@ -362,45 +410,29 @@ void Foam::ReactingParcel<ParcelType>::calc Res, Prs, kappas, - d0, - rho0, - T0, - Cp0, NCpW, Sh, dhsTrans, Sph ); + this->Cp_ = composition.Cp(0, Y_, pc_, T0); + // Motion // ~~~~~~ // Calculate new particle velocity - scalar Spu = 0.0; - vector U1 = - this->calcVelocity - ( - td, - dt, - cellI, - Res, - mus, - d0, - U0, - rho0, - mass0, - Su, - dUTrans, - Spu - ); + this->U_ = + this->calcVelocity(td, dt, cellI, Res, mus, mass1, Su, dUTrans, Spu); - // Accumulate carrier phase source terms - // ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~ + // 4. Accumulate carrier phase source terms + // ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~ + if (td.cloud().solution().coupled()) { - // Transfer mass lost to carrier mass and enthalpy sources + // Transfer mass lost to carrier mass, momentum and enthalpy sources forAll(dMass, i) { scalar dm = np0*dMass[i]; @@ -414,64 +446,12 @@ void Foam::ReactingParcel<ParcelType>::calc // Update momentum transfer td.cloud().UTrans()[cellI] += np0*dUTrans; - - // Update momentum transfer coefficient td.cloud().UCoeff()[cellI] += np0*Spu; // Update sensible enthalpy transfer td.cloud().hsTrans()[cellI] += np0*dhsTrans; - - // Update sensible enthalpy coefficient td.cloud().hsCoeff()[cellI] += np0*Sph; } - - - // Remove the particle when mass falls below minimum threshold - // ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~ - if (np0*mass1 < td.cloud().constProps().minParticleMass()) - { - td.keepParticle = false; - - if (td.cloud().solution().coupled()) - { - scalar dm = np0*mass1; - - // Absorb parcel into carrier phase - forAll(Y_, i) - { - scalar dmi = dm*Y_[i]; - label gid = composition.localToGlobalCarrierId(0, i); - scalar hs = composition.carrier().Hs(gid, T1); - - td.cloud().rhoTrans(gid)[cellI] += dmi; - td.cloud().hsTrans()[cellI] += dmi*hs; - } - td.cloud().UTrans()[cellI] += dm*U1; - - td.cloud().addToMassPhaseChange(dm); - } - } - - - // Set new particle properties - // ~~~~~~~~~~~~~~~~~~~~~~~~~~~ - - else - { - this->Cp_ = composition.Cp(0, Y_, pc_, T1); - this->T_ = T1; - this->U_ = U1; - - // Update particle density or diameter - if (td.cloud().constProps().constantVolume()) - { - this->rho_ = mass1/this->volume(); - } - else - { - this->d_ = cbrt(mass1/this->rho_*6.0/pi); - } - } } diff --git a/src/lagrangian/intermediate/parcels/Templates/ThermoParcel/ThermoParcel.C b/src/lagrangian/intermediate/parcels/Templates/ThermoParcel/ThermoParcel.C index dbc63e7c2e0..f7d37f3ad1a 100644 --- a/src/lagrangian/intermediate/parcels/Templates/ThermoParcel/ThermoParcel.C +++ b/src/lagrangian/intermediate/parcels/Templates/ThermoParcel/ThermoParcel.C @@ -172,21 +172,16 @@ void Foam::ThermoParcel<ParcelType>::calc // Define local properties at beginning of time step // ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~ const scalar np0 = this->nParticle_; - const scalar d0 = this->d_; - const vector U0 = this->U_; - const scalar rho0 = this->rho_; - const scalar T0 = this->T_; - const scalar Cp0 = this->Cp_; const scalar mass0 = this->mass(); // Calc surface values // ~~~~~~~~~~~~~~~~~~~ scalar Ts, rhos, mus, Pr, kappas; - calcSurfaceValues(td, cellI, T0, Ts, rhos, mus, Pr, kappas); + calcSurfaceValues(td, cellI, this->T_, Ts, rhos, mus, Pr, kappas); // Reynolds number - scalar Re = this->Re(U0, d0, rhos, mus); + scalar Re = this->Re(this->U_, this->d_, rhos, mus); // Sources @@ -195,12 +190,18 @@ void Foam::ThermoParcel<ParcelType>::calc // Explicit momentum source for particle vector Su = vector::zero; + // Linearised momentum source coefficient + scalar Spu = 0.0; + // Momentum transfer from the particle to the carrier phase vector dUTrans = vector::zero; // Explicit enthalpy source for particle scalar Sh = 0.0; + // Linearised enthalpy source coefficient + scalar Sph = 0.0; + // Sensible enthalpy transfer from the particle to the carrier phase scalar dhsTrans = 0.0; @@ -212,8 +213,7 @@ void Foam::ThermoParcel<ParcelType>::calc scalar NCpW = 0.0; // Calculate new particle temperature - scalar Sph = 0.0; - scalar T1 = + this->T_ = this->calcHeatTransfer ( td, @@ -222,10 +222,6 @@ void Foam::ThermoParcel<ParcelType>::calc Re, Pr, kappas, - d0, - rho0, - T0, - Cp0, NCpW, Sh, dhsTrans, @@ -237,23 +233,8 @@ void Foam::ThermoParcel<ParcelType>::calc // ~~~~~~ // Calculate new particle velocity - scalar Spu = 0.0; - vector U1 = - this->calcVelocity - ( - td, - dt, - cellI, - Re, - mus, - d0, - U0, - rho0, - mass0, - Su, - dUTrans, - Spu - ); + this->U_ = + this->calcVelocity(td, dt, cellI, Re, mus, mass0, Su, dUTrans, Spu); // Accumulate carrier phase source terms @@ -272,11 +253,6 @@ void Foam::ThermoParcel<ParcelType>::calc // Update sensible enthalpy coefficient td.cloud().hsCoeff()[cellI] += np0*Sph; } - - // Set new particle properties - // ~~~~~~~~~~~~~~~~~~~~~~~~~~~ - this->U_ = U1; - T_ = T1; } @@ -290,10 +266,6 @@ Foam::scalar Foam::ThermoParcel<ParcelType>::calcHeatTransfer const scalar Re, const scalar Pr, const scalar kappa, - const scalar d, - const scalar rho, - const scalar T, - const scalar Cp, const scalar NCpW, const scalar Sh, scalar& dhsTrans, @@ -302,9 +274,12 @@ Foam::scalar Foam::ThermoParcel<ParcelType>::calcHeatTransfer { if (!td.cloud().heatTransfer().active()) { - return T; + return T_; } + const scalar d = this->d(); + const scalar rho = this->rho(); + // Calc heat transfer coefficient scalar htc = td.cloud().heatTransfer().htc(d, Re, Pr, kappa, NCpW); @@ -313,7 +288,7 @@ Foam::scalar Foam::ThermoParcel<ParcelType>::calcHeatTransfer return max ( - T + dt*Sh/(this->volume(d)*rho*Cp), + T_ + dt*Sh/(this->volume(d)*rho*Cp_), td.cloud().constProps().TMin() ); } @@ -322,7 +297,7 @@ Foam::scalar Foam::ThermoParcel<ParcelType>::calcHeatTransfer const scalar As = this->areaS(d); scalar ap = Tc_ + Sh/As/htc; - scalar bp = 6.0*(Sh/As + htc*(Tc_ - T)); + scalar bp = 6.0*(Sh/As + htc*(Tc_ - T_)); if (td.cloud().radiation()) { tetIndices tetIs = this->currentTetIndices(); @@ -330,20 +305,20 @@ Foam::scalar Foam::ThermoParcel<ParcelType>::calcHeatTransfer const scalar sigma = physicoChemical::sigma.value(); const scalar epsilon = td.cloud().constProps().epsilon0(); - ap = (ap + epsilon*Gc/(4.0*htc))/(1.0 + epsilon*sigma*pow3(T)/htc); - bp += 6.0*(epsilon*(Gc/4.0 - sigma*pow4(T))); + ap = (ap + epsilon*Gc/(4.0*htc))/(1.0 + epsilon*sigma*pow3(T_)/htc); + bp += 6.0*(epsilon*(Gc/4.0 - sigma*pow4(T_))); } - bp /= rho*d*Cp*(ap - T) + ROOTVSMALL; + bp /= rho*d*Cp_*(ap - T_) + ROOTVSMALL; // Integrate to find the new parcel temperature IntegrationScheme<scalar>::integrationResult Tres = - td.cloud().TIntegrator().integrate(T, dt, ap*bp, bp); + td.cloud().TIntegrator().integrate(T_, dt, ap*bp, bp); scalar Tnew = max(Tres.value(), td.cloud().constProps().TMin()); Sph = dt*htc*As; - dhsTrans += Sph*(0.5*(T + Tnew) - Tc_); + dhsTrans += Sph*(0.5*(T_ + Tnew) - Tc_); return Tnew; } diff --git a/src/lagrangian/intermediate/parcels/Templates/ThermoParcel/ThermoParcel.H b/src/lagrangian/intermediate/parcels/Templates/ThermoParcel/ThermoParcel.H index c8fe14ad1a6..3b494e746e7 100644 --- a/src/lagrangian/intermediate/parcels/Templates/ThermoParcel/ThermoParcel.H +++ b/src/lagrangian/intermediate/parcels/Templates/ThermoParcel/ThermoParcel.H @@ -236,10 +236,6 @@ protected: const scalar Re, // Reynolds number const scalar Pr, // Prandtl number - surface const scalar kappa, // Thermal conductivity - surface - const scalar d, // diameter - const scalar rho, // density - const scalar T, // temperature - const scalar Cp, // specific heat capacity const scalar NCpW, // Sum of N*Cp*W of emission species const scalar Sh, // explicit particle enthalpy source scalar& dhsTrans, // sensible enthalpy transfer to carrier -- GitLab