Go to the documentation of this file.
34 template<
class CloudType>
37 atomizationModel_.reset
41 this->subModelProperties(),
50 this->subModelProperties(),
57 template<
class CloudType>
63 CloudType::cloudReset(
c);
65 atomizationModel_.reset(
c.atomizationModel_.ptr());
66 breakupModel_.reset(
c.breakupModel_.ptr());
72 template<
class CloudType>
85 cloudCopyPtr_(
nullptr),
86 averageParcelMass_(0.0),
87 atomizationModel_(
nullptr),
88 breakupModel_(
nullptr)
94 averageParcelMass_ = this->injectors().averageParcelMass();
99 this->deleteLostParticles();
102 Info <<
"Average parcel mass: " << averageParcelMass_ <<
endl;
105 if (this->
solution().resetSourcesOnStartup())
107 CloudType::resetSourceTerms();
112 template<
class CloudType>
121 cloudCopyPtr_(
nullptr),
122 averageParcelMass_(
c.averageParcelMass_),
123 atomizationModel_(
c.atomizationModel_->clone()),
124 breakupModel_(
c.breakupModel_->clone())
128 template<
class CloudType>
138 cloudCopyPtr_(
nullptr),
139 averageParcelMass_(0.0),
140 atomizationModel_(
nullptr),
141 breakupModel_(
nullptr)
147 template<
class CloudType>
154 template<
class CloudType>
158 const scalar lagrangianDt
161 CloudType::setParcelThermoProperties(parcel, lagrangianDt);
167 const scalar pc = this->
p()[parcel.cell()];
170 parcel.Cp() = liqMix.
Cp(pc, parcel.T(), X);
171 parcel.rho() = liqMix.
rho(pc, parcel.T(), X);
172 parcel.sigma() = liqMix.
sigma(pc, parcel.T(), X);
173 parcel.mu() = liqMix.
mu(pc, parcel.T(), X);
177 template<
class CloudType>
181 const scalar lagrangianDt,
182 const bool fullyDescribed
185 CloudType::checkParcelProperties(parcel, lagrangianDt, fullyDescribed);
188 parcel.position0() = parcel.position();
189 parcel.d0() = parcel.d();
191 parcel.y() = breakup().y0();
192 parcel.yDot() = breakup().yDot0();
194 parcel.liquidCore() = atomization().initLiquidCore();
198 template<
class CloudType>
205 clone(this->
name() +
"Copy").ptr()
211 template<
class CloudType>
214 cloudReset(cloudCopyPtr_());
215 cloudCopyPtr_.clear();
219 template<
class CloudType>
224 typename parcelType::trackingData td(*
this);
226 this->
solve(*
this, td);
231 template<
class CloudType>
235 scalar d32 = 1.0e+6*this->Dij(3, 2);
236 scalar d10 = 1.0e+6*this->Dij(1, 0);
237 scalar dMax = 1.0e+6*this->Dmax();
238 scalar pen = this->penetration(0.95);
240 Info <<
" D10, D32, Dmax (mu) = " << d10 <<
", " << d32
241 <<
", " << dMax <<
nl
242 <<
" Liquid penetration 95% mass (m) = " << pen <<
endl;
Templated base class for spray cloud.
Selector class for relaxation factors, solver type and solution.
const word cloudName(propsDict.get< word >("cloud"))
A class for handling words, derived from Foam::string.
Thermo package for (S)olids (L)iquids and (G)ases Takes reference to thermo package,...
void storeState()
Store the current cloud state.
Basic thermodynamics type based on the use of fitting functions for cp, h, s obtained from the templa...
Templated atomization model class.
scalar mu(const scalar p, const scalar T, const scalarField &X) const
Calculate the mixture viscosity [Pa s].
Ostream & endl(Ostream &os)
Add newline and flush stream.
Foam::DSMCCloud ::particleType parcelType
Type of parcel the cloud was instantiated for.
void restoreState()
Reset the current cloud to the previously stored state.
basicSpecieMixture & composition
scalar rho(const scalar p, const scalar T, const scalarField &X) const
Calculate the mixture density [kg/m^3].
Virtual abstract base class for templated SprayCloud.
void info()
Print cloud information.
scalarField X(const scalarField &Y) const
Returns the mole fractions corresponding to the given mass fractions.
void cloudReset(SprayCloud< CloudType > &c)
Reset state of cloud.
scalar sigma(const scalar p, const scalar T, const scalarField &X) const
Estimate mixture surface tension [N/m].
messageStream Info
Information stream (stdout output on master, null elsewhere)
void checkParcelProperties(parcelType &parcel, const scalar lagrangianDt, const bool fullyDescribed)
Check parcel properties.
DSMCCloud< dsmcParcel > CloudType
Mesh data needed to do the Finite Volume discretisation.
const uniformDimensionedVectorField & g
PtrList< volScalarField > & Y
Templated break-up model class.
void readFields(const typename GeoFieldType::Mesh &mesh, const IOobjectList &objects, const wordHashSet &selectedFields, LIFOStack< regIOobject * > &storedObjects)
Read the selected GeometricFields of the templated type.
tmp< DimensionedField< TypeR, GeoMesh > > New(const tmp< DimensionedField< TypeR, GeoMesh >> &tdf1, const word &name, const dimensionSet &dimensions)
Global function forwards to reuseTmpDimensionedField::New.
virtual ~SprayCloud()
Destructor.
scalar Cp(const scalar p, const scalar T, const scalarField &X) const
Calculate the mixture heat capacity [J/(kg K)].
void setParcelThermoProperties(parcelType &parcel, const scalar lagrangianDt)
Set parcel thermo properties.
void evolve()
Evolve the spray (inject, move)
const dimensionedScalar c
Speed of light in a vacuum.
word name(const expressions::valueTypeCode typeCode)
A word representation of a valueTypeCode. Empty for INVALID.
void setModels()
Set cloud sub-models.