33 template<
class CloudType>
36 heatTransferModel_.reset
40 this->subModelProperties(),
45 compositionModel_.reset
49 this->subModelProperties(),
59 this->
solution().integrationSchemes()
65 this->subModelProperties().lookup(
"radiation") >> radiation_;
76 this->
name() +
":radAreaP",
77 this->db().time().
name(),
79 IOobject::READ_IF_PRESENT,
93 this->
name() +
":radT4",
94 this->db().time().
name(),
96 IOobject::READ_IF_PRESENT,
110 this->
name() +
":radAreaPT4",
111 this->db().time().
name(),
113 IOobject::READ_IF_PRESENT,
124 template<
class CloudType>
127 CloudType::cloudReset(
c);
129 heatTransferModel_.reset(
c.heatTransferModel_.ptr());
130 compositionModel_.reset(
c.compositionModel_.ptr());
132 TIntegrator_.reset(
c.TIntegrator_.ptr());
134 radiation_ =
c.radiation_;
140 template<
class CloudType>
152 cloudCopyPtr_(nullptr),
153 constProps_(this->particleProperties()),
154 carrierThermo_(carrierThermo),
155 thermo_(carrierThermo_),
156 T_(carrierThermo.
T()),
157 p_(carrierThermo.
p()),
158 heatTransferModel_(nullptr),
159 compositionModel_(nullptr),
160 TIntegrator_(nullptr),
164 radAreaPT4_(nullptr),
171 this->
name() +
":hsTrans",
172 this->db().time().
name(),
187 this->
name() +
":hsCoeff",
188 this->db().time().
name(),
206 if (this->
solution().resetSourcesOnStartup())
213 template<
class CloudType>
221 cloudCopyPtr_(nullptr),
222 constProps_(
c.constProps_),
223 carrierThermo_(
c.carrierThermo_),
227 heatTransferModel_(
c.heatTransferModel_->
clone()),
228 compositionModel_(
c.compositionModel_->
clone()),
229 TIntegrator_(
c.TIntegrator_->
clone()),
230 radiation_(
c.radiation_),
233 radAreaPT4_(nullptr),
240 this->
name() +
":hsTrans",
241 this->db().time().
name(),
256 this->
name() +
":hsCoeff",
257 this->db().time().
name(),
275 this->
name() +
":radAreaP",
292 this->
name() +
":radT4",
309 this->
name() +
":radAreaPT4",
323 template<
class CloudType>
332 cloudCopyPtr_(nullptr),
334 carrierThermo_(
c.carrierThermo_),
338 heatTransferModel_(nullptr),
339 compositionModel_(
c.compositionModel_->
clone()),
340 TIntegrator_(nullptr),
344 radAreaPT4_(nullptr),
352 template<
class CloudType>
359 template<
class CloudType>
365 CloudType::setParcelThermoProperties(parcel);
367 parcel.T() = constProps_.T0();
368 parcel.Cp() = constProps_.Cp0();
372 template<
class CloudType>
376 const label injectori
379 CloudType::checkParcelProperties(parcel, injectori);
383 template<
class CloudType>
396 template<
class CloudType>
399 cloudReset(cloudCopyPtr_());
400 cloudCopyPtr_.clear();
404 template<
class CloudType>
407 CloudType::resetSourceTerms();
408 hsTrans_->primitiveFieldRef() = 0.0;
409 hsCoeff_->primitiveFieldRef() = 0.0;
413 radAreaP_->primitiveFieldRef() = 0.0;
414 radT4_->primitiveFieldRef() = 0.0;
415 radAreaPT4_->primitiveFieldRef() = 0.0;
420 template<
class CloudType>
426 CloudType::relaxSources(cloudOldTime);
433 this->
relax(radAreaP_(), cloudOldTime.
radAreaP(),
"radiation");
434 this->
relax(radT4_(), cloudOldTime.
radT4(),
"radiation");
440 template<
class CloudType>
443 CloudType::scaleSources();
445 this->scale(hsTrans_(),
"h");
446 this->scale(hsCoeff_(),
"h");
450 this->scale(radAreaP_(),
"radiation");
451 this->scale(radT4_(),
"radiation");
452 this->scale(radAreaPT4_(),
"radiation");
457 template<
class CloudType>
460 CloudType::preEvolve();
462 this->pAmbient() = carrierThermo_.p().average().value();
466 template<
class CloudType>
471 typename parcelType::trackingData td(*
this);
473 this->
solve(*
this, td);
478 template<
class CloudType>
483 Info<<
" Temperature min/max = " << Tmin() <<
", " << Tmax()
488 template<
class CloudType>
491 if (compositionModel_.valid())
493 CloudType::particleType::writeFields(*
this, this->composition());
497 CloudType::particleType::writeFields(*
this);
void deleteLostParticles()
Remove lost particles from cloud and delete.
Templated reacting parcel composition model class Consists of carrier species (via thermo package),...
Templated base class for dsmc cloud.
ParcelType parcelType
Type of parcel the cloud was instantiated for.
void info() const
Print cloud information.
Field with dimensions and associated with geometry type GeoMesh which is used to size the field and a...
Generic GeometricField class.
Templated heat transfer model class.
IOobject defines the attributes of an object for which implicit objectRegistry management is supporte...
const objectRegistry & db() const
Return the local objectRegistry.
const word & name() const
Return name.
Templated base class for thermodynamic cloud.
volScalarField::Internal & radAreaP()
Radiation sum of parcel projected areas [m^2].
void setModels()
Set cloud sub-models.
volScalarField::Internal & radAreaPT4()
Radiation sum of parcel projected area*temperature^4 [m2K4].
void storeState()
Store the current cloud state.
volScalarField::Internal & radT4()
Radiation sum of parcel temperature^4 [K4].
void relaxSources(const ThermoCloud< CloudType > &cloudOldTime)
Apply relaxation to (steady state) cloud sources.
autoPtr< volScalarField::Internal > hsCoeff_
Coefficient for carrier phase hs equation [W/K].
autoPtr< volScalarField::Internal > radT4_
Radiation sum of parcel temperature^4.
void scaleSources()
Apply scaling to (transient) cloud sources.
ThermoCloud(const word &cloudName, const volScalarField &rho, const volVectorField &U, const dimensionedVector &g, const fluidThermo &carrierThermo, const bool readFields=true)
Construct given carrier fields and thermo.
Switch radiation_
Include radiation.
virtual void writeFields() const
Write the field data for the cloud.
void cloudReset(ThermoCloud< CloudType > &c)
Reset state of cloud.
void evolve()
Evolve the cloud.
autoPtr< volScalarField::Internal > radAreaPT4_
Radiation sum of parcel projected areas * temperature^4.
void info()
Print cloud information.
void restoreState()
Reset the current cloud to the previously stored state.
void checkParcelProperties(parcelType &parcel, const label injectori)
Check parcel properties.
void preEvolve()
Pre-evolve.
void resetSourceTerms()
Reset the cloud source terms.
virtual ~ThermoCloud()
Destructor.
autoPtr< volScalarField::Internal > hsTrans_
Sensible enthalpy transfer [J/kg].
autoPtr< volScalarField::Internal > radAreaP_
Radiation sum of parcel projected areas.
void setParcelThermoProperties(parcelType &parcel)
Set parcel thermo properties.
Base-class for fluid thermodynamic properties.
Mesh data needed to do the Finite Volume discretisation.
const Time & time() const
Return time.
Selector class for relaxation factors, solver type and solution.
A class for handling words, derived from string.
autoPtr< CompressibleMomentumTransportModel > New(const volScalarField &rho, const volVectorField &U, const surfaceScalarField &phi, const viscosity &viscosity)
const dimensionedScalar c
Speed of light in a vacuum.
intWM_LABEL_SIZE_t label
A label is an int32_t or int64_t as specified by the pre-processor macro WM_LABEL_SIZE.
const dimensionSet dimEnergy
dimensionedSymmTensor sqr(const dimensionedVector &dv)
Ostream & endl(Ostream &os)
Add newline and flush stream.
word name(const bool)
Return a word representation of a bool.
void readFields(const typename GeoFieldType::Mesh &mesh, const IOobjectList &objects, const HashSet< word > &selectedFields, LIFOStack< regIOobject * > &storedObjects)
Read the selected GeometricFields of the specified type.
const dimensionSet dimLength
const dimensionSet dimTemperature
dimensionedScalar pow4(const dimensionedScalar &ds)
const dimensionSet dimArea
void T(FieldField< Field, Type > &f1, const FieldField< Field, Type > &f2)
dimensioned< scalar > dimensionedScalar
Dimensioned scalar obtained from generic dimensioned type.
SolverPerformance< Type > solve(fvMatrix< Type > &, const word &)
Solve returning the solution statistics given convergence tolerance.
fluidMulticomponentThermo & thermo
const word cloudName(propsDict.lookup("cloudName"))