35 template<
class ParcelType>
36 template<
class TrackCloudType>
39 TrackCloudType&
cloud,
43 ParcelType::setCellValues(
cloud, td);
51 if (td.
Tc() <
cloud.constProps().TMin())
56 <<
"Limiting observed temperature in cell " << this->
cell()
57 <<
" to " <<
cloud.constProps().TMin() <<
nl <<
endl;
60 td.
Tc() =
cloud.constProps().TMin();
65 template<
class ParcelType>
66 template<
class TrackCloudType>
69 TrackCloudType& cloud,
74 const label celli = this->cell();
75 const scalar massCell = this->massCell(td);
77 td.Uc() += cloud.UTrans()[celli]/massCell;
82 const scalar CpMean = td.CpInterp().psi()[celli];
83 td.Tc() += cloud.hsTrans()[celli]/(CpMean*massCell);
85 if (td.Tc() < cloud.constProps().TMin())
90 <<
"Limiting observed temperature in cell " << celli
91 <<
" to " << cloud.constProps().TMin() <<
nl <<
endl;
94 td.Tc() = cloud.constProps().TMin();
99 template<
class ParcelType>
100 template<
class TrackCloudType>
103 TrackCloudType& cloud,
114 Ts = (2.0*
T + td.Tc())/3.0;
116 if (Ts < cloud.constProps().TMin())
121 <<
"Limiting parcel surface temperature to " 122 << cloud.constProps().TMin() <<
nl <<
endl;
125 Ts = cloud.constProps().TMin();
129 const scalar TRatio = td.Tc()/Ts;
131 rhos = td.rhoc()*TRatio;
133 tetIndices tetIs = this->currentTetIndices();
134 mus = td.muInterp().interpolate(this->
coordinates(), tetIs)/TRatio;
135 kappas = td.kappaInterp().interpolate(this->
coordinates(), tetIs)/TRatio;
137 Pr = td.Cpc()*mus/kappas;
142 template<
class ParcelType>
143 template<
class TrackCloudType>
146 TrackCloudType&
cloud,
153 const scalar np0 = this->nParticle_;
154 const scalar mass0 = this->mass();
157 const scalar
T0 = this->T_;
162 scalar Ts, rhos, mus,
Pr, kappas;
163 calcSurfaceValues(
cloud, td, this->T_, Ts, rhos, mus,
Pr, kappas);
166 scalar
Re = this->
Re(rhos, this->U_, td.Uc(), this->d_, mus);
188 scalar dhsTrans = 0.0;
199 this->calcHeatTransfer
219 this->calcVelocity(
cloud, td, dt,
Re, mus, mass0,
Su, dUTrans, Spu);
224 if (
cloud.solution().coupled())
227 cloud.UTrans()[this->
cell()] += np0*dUTrans;
233 cloud.hsTrans()[this->
cell()] += np0*dhsTrans;
236 cloud.hsCoeff()[this->
cell()] += np0*Sph;
239 if (
cloud.radiation())
241 const scalar ap = this->areaP();
242 const scalar T4 =
pow4(
T0);
243 cloud.radAreaP()[this->
cell()] += dt*np0*ap;
244 cloud.radT4()[this->
cell()] += dt*np0*T4;
245 cloud.radAreaPT4()[this->
cell()] += dt*np0*ap*T4;
251 template<
class ParcelType>
252 template<
class TrackCloudType>
255 TrackCloudType& cloud,
267 if (!cloud.heatTransfer().active())
272 const scalar d = this->d();
273 const scalar
rho = this->
rho();
274 const scalar As = this->areaS(d);
275 const scalar V = this->
volume(d);
276 const scalar m =
rho*V;
279 scalar htc = cloud.heatTransfer().htc(d,
Re,
Pr,
kappa, NCpW);
282 const scalar bcp = htc*As/(m*Cp_);
283 const scalar acp = bcp*td.Tc();
285 if (cloud.radiation())
287 const tetIndices tetIs = this->currentTetIndices();
288 const scalar Gc = td.GInterp().interpolate(this->
coordinates(), tetIs);
290 const scalar
epsilon = cloud.constProps().epsilon0();
297 const scalar deltaT = cloud.TIntegrator().delta(T_, dt, acp + ancp, bcp);
298 const scalar deltaTncp = ancp*dt;
299 const scalar deltaTcp = deltaT - deltaTncp;
302 scalar Tnew = T_ + deltaT;
303 Tnew =
min(
max(Tnew, cloud.constProps().TMin()), cloud.constProps().TMax());
305 dhsTrans -= m*Cp_*deltaTcp;
315 template<
class ParcelType>
318 const ThermoParcel<ParcelType>&
p 327 template<
class ParcelType>
330 const ThermoParcel<ParcelType>&
p,
Different types of constants.
dimensionedScalar Pr("Pr", dimless, laminarTransport)
const Type & value() const noexcept
Return const reference to value.
label max(const labelHashSet &set, label maxValue=labelMin)
Find the max value in labelHashSet, optionally limited by second argument.
void calc(TrackCloudType &cloud, trackingData &td, const scalar dt)
Update parcel properties over the time interval.
constexpr char nl
The newline '\n' character (0x0a)
Ostream & endl(Ostream &os)
Add newline and flush stream.
const dimensionedScalar sigma
Stefan-Boltzmann constant: default SI units: [W/m2/K4].
ThermoParcel(const polyMesh &mesh, const barycentric &coordinates, const label celli, const label tetFacei, const label tetPti)
Construct from mesh, coordinates and topology.
const dimensionedScalar kappa
Coulomb constant: default SI units: [N.m2/C2].
dimensionedScalar sigma("sigma", dimMass/sqr(dimTime), transportProperties)
const interpolation< scalar > & TInterp() const
Return const access to the interpolator for continuous.
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.
A cloud is a registry collection of lagrangian particles.
void calcSurfaceValues(TrackCloudType &cloud, trackingData &td, const scalar T, scalar &Ts, scalar &rhos, scalar &mus, scalar &Pr, scalar &kappas) const
Calculate surface thermo properties.
void cellValueSourceCorrection(TrackCloudType &cloud, trackingData &td, const scalar dt)
Correct cell values using latest transfer information.
scalarField Re(const UList< complex > &cf)
Extract real component.
Storage and named access for the indices of a tet which is part of the decomposition of a cell...
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...
int debug
Static debugging option.
scalar Tc() const
Return the continuous phase temperature.
void T(FieldField< Field, Type > &f1, const FieldField< Field, Type > &f2)
const interpolation< scalar > & CpInterp() const
Return const access to the interpolator for continuous.
PtrList< coordinateSystem > coordinates(solidRegions.size())
scalar calcHeatTransfer(TrackCloudType &cloud, trackingData &td, const scalar dt, const scalar Re, const scalar Pr, const scalar kappa, const scalar NCpW, const scalar Sh, scalar &dhsTrans, scalar &Sph)
Calculate new particle temperature.
#define WarningInFunction
Report a warning using Foam::Warning.
scalar Cpc() const
Return the continuous phase specific heat capacity.
A cell is defined as a list of faces with extra functionality.
dimensionedScalar pow4(const dimensionedScalar &ds)
void setCellValues(TrackCloudType &cloud, trackingData &td)
Set cell values.
static constexpr const zero Zero
Global zero (0)