33 template<
class ParcelType>
34 template<
class TrackData>
42 ParcelType::setCellValues(td, dt, celli);
46 Cpc_ = td.CpInterp().interpolate(this->position(), tetIs);
48 Tc_ = td.TInterp().interpolate(this->position(), tetIs);
50 if (Tc_ < td.cloud().constProps().TMin())
55 <<
"Limiting observed temperature in cell " << celli <<
" to " 56 << td.cloud().constProps().TMin() <<
nl <<
endl;
59 Tc_ = td.cloud().constProps().TMin();
64 template<
class ParcelType>
65 template<
class TrackData>
73 this->Uc_ += td.cloud().UTrans()[celli]/this->massCell(celli);
75 const scalar CpMean = td.CpInterp().psi()[celli];
76 Tc_ += td.cloud().hsTrans()[celli]/(CpMean*this->massCell(celli));
78 if (Tc_ < td.cloud().constProps().TMin())
83 <<
"Limiting observed temperature in cell " << celli <<
" to " 84 << td.cloud().constProps().TMin() <<
nl <<
endl;
87 Tc_ = td.cloud().constProps().TMin();
92 template<
class ParcelType>
93 template<
class TrackData>
107 Ts = (2.0*T + Tc_)/3.0;
109 if (Ts < td.cloud().constProps().TMin())
114 <<
"Limiting parcel surface temperature to " 115 << td.cloud().constProps().TMin() <<
nl <<
endl;
118 Ts = td.cloud().constProps().TMin();
122 const scalar TRatio = Tc_/Ts;
124 rhos = this->rhoc_*TRatio;
127 mus = td.muInterp().interpolate(this->position(), tetIs)/TRatio;
128 kappas = td.kappaInterp().interpolate(this->position(), tetIs)/TRatio;
130 Pr = Cpc_*mus/kappas;
131 Pr =
max(ROOTVSMALL, Pr);
135 template<
class ParcelType>
136 template<
class TrackData>
146 const scalar np0 = this->nParticle_;
147 const scalar mass0 = this->mass();
150 const scalar T0 = this->T_;
155 scalar Ts, rhos, mus,
Pr, kappas;
156 calcSurfaceValues(td, celli, this->T_, Ts, rhos, mus, Pr, kappas);
159 scalar Re = this->
Re(this->U_, this->d_, rhos, mus);
181 scalar dhsTrans = 0.0;
192 this->calcHeatTransfer
212 this->calcVelocity(td, dt, celli, Re, mus, mass0, Su, dUTrans, Spu);
217 if (td.cloud().solution().coupled())
220 td.cloud().UTrans()[celli] += np0*dUTrans;
223 td.cloud().UCoeff()[celli] += np0*Spu;
226 td.cloud().hsTrans()[celli] += np0*dhsTrans;
229 td.cloud().hsCoeff()[celli] += np0*Sph;
232 if (td.cloud().radiation())
234 const scalar ap = this->areaP();
235 const scalar T4 =
pow4(T0);
236 td.cloud().radAreaP()[celli] += dt*np0*ap;
237 td.cloud().radT4()[celli] += dt*np0*T4;
238 td.cloud().radAreaPT4()[celli] += dt*np0*ap*T4;
244 template<
class ParcelType>
245 template<
class TrackData>
260 if (!td.cloud().heatTransfer().active())
265 const scalar d = this->d();
266 const scalar
rho = this->
rho();
269 scalar htc = td.cloud().heatTransfer().htc(d, Re, Pr, kappa, NCpW);
271 if (
mag(htc) < ROOTVSMALL && !td.cloud().radiation())
276 T_ + dt*Sh/(this->volume(d)*rho*Cp_),
277 td.cloud().constProps().TMin()
281 htc =
max(htc, ROOTVSMALL);
282 const scalar As = this->areaS(d);
284 scalar ap = Tc_ + Sh/(As*htc);
285 scalar bp = 6.0*(Sh/As + htc*(Tc_ - T_));
286 if (td.cloud().radiation())
289 const scalar Gc = td.GInterp().interpolate(this->position(), tetIs);
291 const scalar
epsilon = td.cloud().constProps().epsilon0();
294 scalar
s = epsilon*(Gc/4.0 - sigma*
pow4(T_));
299 bp /= rho*d*Cp_*(ap - T_) + ROOTVSMALL;
303 td.cloud().TIntegrator().integrate(T_, dt, ap*bp, bp);
311 td.cloud().constProps().TMin()
313 td.cloud().constProps().TMax()
318 dhsTrans += Sph*(Tres.average() - Tc_);
326 template<
class ParcelType>
340 template<
class ParcelType>
dimensionedScalar Pr("Pr", dimless, laminarTransport)
void calcSurfaceValues(TrackData &td, const label celli, const scalar T, scalar &Ts, scalar &rhos, scalar &mus, scalar &Pr, scalar &kappas) const
Calculate surface thermo properties.
intWM_LABEL_SIZE_t label
A label is an int32_t or int64_t as specified by the pre-processor macro WM_LABEL_SIZE.
void cellValueSourceCorrection(TrackData &td, const scalar dt, const label celli)
Correct cell values using latest transfer information.
dimensioned< Type > max(const dimensioned< Type > &, const dimensioned< Type > &)
scalar Cp_
Specific heat capacity [J/(kg.K)].
scalar calcHeatTransfer(TrackData &td, const scalar dt, const label celli, const scalar Re, const scalar Pr, const scalar kappa, const scalar NCpW, const scalar Sh, scalar &dhsTrans, scalar &Sph)
Calculate new particle temperature.
void calc(TrackData &td, const scalar dt, const label celli)
Update parcel properties over the time interval.
Helper class to supply results of integration.
Ostream & endl(Ostream &os)
Add newline and flush stream.
const dimensionedScalar sigma
Stefan-Boltzmann constant: default SI units: [W/m2/K4].
const Type & value() const
Return const reference to value.
void setCellValues(TrackData &td, const scalar dt, const label celli)
Set cell values.
scalar T_
Temperature [K].
gmvFile<< "tracers "<< particles.size()<< nl;forAllConstIter(Cloud< passiveParticle >, particles, iter){gmvFile<< iter().position().x()<< " ";}gmvFile<< nl;forAllConstIter(Cloud< passiveParticle >, particles, iter){gmvFile<< iter().position().y()<< " ";}gmvFile<< nl;forAllConstIter(Cloud< passiveParticle >, particles, iter){gmvFile<< iter().position().z()<< " ";}gmvFile<< nl;forAll(lagrangianScalarNames, i){word name=lagrangianScalarNames[i];IOField< scalar > s(IOobject(name, runTime.timeName(), cloud::prefix, mesh, IOobject::MUST_READ, IOobject::NO_WRITE))
ThermoParcel(const polyMesh &mesh, const vector &position, const label celli, const label tetFacei, const label tetPtI)
Construct from owner, position, and cloud owner.
Storage and named access for the indices of a tet which is part of the decomposition of a cell...
scalar Cpc_
Specific heat capacity [J/(kg.K)].
dimensioned< Type > min(const dimensioned< Type > &, const dimensioned< Type > &)
#define WarningInFunction
Report a warning using Foam::Warning.
dimensionedScalar pow4(const dimensionedScalar &ds)
dimensioned< scalar > mag(const dimensioned< Type > &)
Mesh consisting of general polyhedral cells.
tmp< GeometricField< Type, fvPatchField, volMesh > > Su(const GeometricField< Type, fvPatchField, volMesh > &su, const GeometricField< Type, fvPatchField, volMesh > &vf)
scalarField Re(const UList< complex > &cf)
Thermodynamic parcel class with one/two-way coupling with the continuous phase. Includes Kinematic pa...
scalar Tc_
Temperature [K].