31 template<
class CloudType>
35 return *cloudCopyPtr_;
39 template<
class CloudType>
40 inline const typename CloudType::particleType::constantProperties&
47 template<
class CloudType>
48 inline typename CloudType::particleType::constantProperties&
55 template<
class CloudType>
59 return *compositionModel_;
63 template<
class CloudType>
67 return *phaseChangeModel_;
71 template<
class CloudType>
75 return *phaseChangeModel_;
79 template<
class CloudType>
83 const typename parcelType::trackingData&
td 88 const label celli =
p.cell();
90 const scalar m =
p.nParticle()*
p.mass();
92 this->rhokTrans()[celli] += m;
94 this->UTrans()[celli] += m*
p.U();
96 const scalar pc =
td.pc();
97 const scalar
T =
p.T();
98 const auto&
Y =
p.Y();
102 const scalar dm = m*
p.Y[i];
103 const label gid = comp.localToCarrierId(0, i);
104 const scalar hs = comp.carrier().Hs(gid, pc,
T);
106 this->rhoTrans(gid)[celli] += dm;
107 this->hsTrans()[celli] += dm*hs;
112 template<
class CloudType>
120 template<
class CloudType>
129 template<
class CloudType>
137 template<
class CloudType>
146 if (this->
solution().semiImplicit(
"Yi"))
150 IOobject::scopedName(this->
name(),
"rhoTrans"),
151 IOobject::NO_REGISTER,
155 auto& sourceField = trhoTrans.ref();
157 sourceField.primitiveFieldRef() =
158 rhoTrans_[i]/(this->db().time().deltaTValue()*this->
mesh().V());
163 fvm::Sp(
neg(sourceField)*sourceField/(Yi + YiSMALL), Yi)
164 +
pos0(sourceField)*sourceField;
169 auto& fvm = tfvm.ref();
171 fvm.source() = -rhoTrans_[i]/this->db().time().deltaTValue();
181 template<
class CloudType>
187 IOobject::scopedName(this->
name(),
"rhoTrans"),
188 IOobject::NO_REGISTER,
199 rhoi = rhoTrans_[i]/(this->db().time().deltaTValue()*this->
mesh().V());
206 template<
class CloudType>
212 IOobject::scopedName(this->
name(),
"rhoTrans"),
213 IOobject::NO_REGISTER,
226 sourceField += rhoTrans_[i];
229 sourceField /= this->db().time().deltaTValue()*this->
mesh().V();
236 template<
class CloudType>
244 IOobject::scopedName(this->
name(),
"rhoTrans"),
245 IOobject::NO_REGISTER,
251 if (this->
solution().semiImplicit(
"rho"))
255 sourceField += rhoTrans_[i];
257 sourceField /= this->db().time().deltaTValue()*this->
mesh().V();
264 auto& fvm = tfvm.ref();
268 sourceField += rhoTrans_[i];
271 fvm.source() = -trhoTrans()/this->db().time().deltaT();
Templated phase change model class.
basicSpecieMixture & composition
wallPoints::trackData td(isBlockedFace, regionToBlockSize)
const PtrList< volScalarField::Internal > & rhoTrans() const
Return const access to mass source fields.
tmp< DimensionedField< TypeR, GeoMesh > > New(const tmp< DimensionedField< TypeR, GeoMesh >> &tf1, const word &name, const dimensionSet &dimensions, const bool initCopy=false)
Global function forwards to reuseTmpDimensionedField::New.
Generic dimensioned Type class.
void transferToCarrier(const parcelType &p, const typename parcelType::trackingData &td)
Transfer the effect of parcel to the carrier phase.
const dimensionSet dimless
Dimensionless.
zeroField SuSp(const Foam::zero, const GeometricField< Type, faPatchField, areaMesh > &)
A no-op source.
tmp< volScalarField::Internal > Srho() const
Return tmp total mass source for carrier phase.
dimensionedScalar neg(const dimensionedScalar &ds)
#define forAll(list, i)
Loop across all elements in list.
const dimensionSet dimVolume(pow3(dimLength))
const ReactingCloud & cloudCopy() const
Return a reference to the cloud copy.
word name(const expressions::valueTypeCode typeCode)
A word representation of a valueTypeCode. Empty for expressions::valueTypeCode::INVALID.
const parcelType::constantProperties & constProps() const
Return the constant properties.
dimensionedScalar pos0(const dimensionedScalar &ds)
Templated base class for reacting cloud.
dimensioned< scalar > dimensionedScalar
Dimensioned scalar obtained from generic dimensioned type.
A list of pointers to objects of type <T>, with allocation/deallocation management of the pointers...
PtrList< volScalarField > & Y
Field with dimensions and associated with geometry type GeoMesh which is used to size the field and a...
const dimensionSet dimTime(0, 0, 1, 0, 0, 0, 0)
const CompositionModel< ReactingCloud< CloudType > > & composition() const
Return const access to reacting composition model.
Selector class for relaxation factors, solver type and solution.
const dimensionSet dimMass(1, 0, 0, 0, 0, 0, 0)
A class for managing temporary objects.
Templated reacting parcel composition model class Consists of carrier species (via thermo package)...
tmp< fvScalarMatrix > SYi(const label i, volScalarField &Yi) const
Return mass source term for specie i - specie eqn.
const PhaseChangeModel< ReactingCloud< CloudType > > & phaseChange() const
Return const access to reacting phase change model.
static constexpr const zero Zero
Global zero (0)