33 template<
class BasePhaseSystem>
37 const phasePairKey& key
40 if (!iDmdt_.found(key))
42 return phaseSystem::dmdt(key);
47 return dmdtSign**iDmdt_[key];
51 template<
class BasePhaseSystem>
55 const phasePairKey& key
58 if (!wDmdt_.found(key))
60 return phaseSystem::dmdt(key);
65 return dmdtSign**wDmdt_[key];
71 template<
class BasePhaseSystem>
78 BasePhaseSystem(mesh),
79 volatile_(this->template lookupOrDefault<word>(
"volatile",
"none")),
82 saturationModel::
New(this->subDict(
"saturationModel"), mesh)
84 phaseChange_(this->
lookup(
"phaseChange"))
94 const phasePair& pair(phasePairIter());
110 this->
mesh().time().timeName(),
129 this->
mesh().time().timeName(),
148 this->
mesh().time().timeName(),
163 template<
class BasePhaseSystem>
171 template<
class BasePhaseSystem>
175 return saturationModel_();
179 template<
class BasePhaseSystem>
183 const phasePairKey& key
186 return BasePhaseSystem::dmdt(key) + this->iDmdt(key) + this->wDmdt(key);
190 template<
class BasePhaseSystem>
194 PtrList<volScalarField>
dmdts(BasePhaseSystem::dmdts());
198 const phasePair& pair = this->phasePairs_[iDmdtIter.key()];
201 this->addField(pair.phase1(),
"dmdt", iDmdt,
dmdts);
202 this->addField(pair.phase2(),
"dmdt", - iDmdt,
dmdts);
207 const phasePair& pair = this->phasePairs_[wDmdtIter.key()];
210 this->addField(pair.phase1(),
"dmdt", wDmdt,
dmdts);
211 this->addField(pair.phase2(),
"dmdt", - wDmdt,
dmdts);
218 template<
class BasePhaseSystem>
222 autoPtr<phaseSystem::heatTransferTable> eqnsPtr =
223 BasePhaseSystem::heatTransfer();
235 if (this->wMDotL_.found(phasePairIter.key()))
237 const phasePair& pair(phasePairIter());
244 const phaseModel& phase1 = pair.phase1();
245 const phaseModel& phase2 = pair.phase2();
247 *eqns[phase1.name()] +=
negPart(*this->wMDotL_[pair]);
248 *eqns[phase2.name()] -=
posPart(*this->wMDotL_[pair]);
256 template<
class BasePhaseSystem>
260 autoPtr<phaseSystem::massTransferTable> eqnsPtr =
272 const phasePair& pair(phasePairIter());
279 const phaseModel& phase = pair.phase1();
280 const phaseModel& otherPhase = pair.phase2();
282 const PtrList<volScalarField>& Yi = phase.Y();
286 if (Yi[i].member() != volatile_)
307 *eqns[otherName] -=
dmdt;
315 template<
class BasePhaseSystem>
319 typedef compressible::alphatPhaseChangeWallFunctionFvPatchScalarField
320 alphatPhaseChangeWallFunction;
324 typename BasePhaseSystem::heatTransferModelTable,
325 this->heatTransferModels_,
326 heatTransferModelIter
329 const phasePair& pair
331 this->phasePairs_[heatTransferModelIter.key()]
334 const phaseModel& phase1 = pair.phase1();
335 const phaseModel& phase2 = pair.phase2();
351 (
neg0(iDmdt)*hef2 +
pos(iDmdt)*he2)
352 - (
pos0(iDmdt)*hef1 +
neg(iDmdt)*he1)
362 Tf = saturationModel_->Tsat(phase1.thermo().p());
364 iDmdtNew = (H1*(Tf - T1) + H2*(Tf - T2))/L;
378 Tf = (H1*T1 + H2*T2 + iDmdtNew*L)/(H1 + H2);
380 Info<<
"Tf." << pair.name()
381 <<
": min = " <<
min(Tf.primitiveField())
382 <<
", mean = " <<
average(Tf.primitiveField())
383 <<
", max = " <<
max(Tf.primitiveField())
386 scalar iDmdtRelax(this->
mesh().fieldRelaxationFactor(
"iDmdt"));
387 iDmdt = (1 - iDmdtRelax)*iDmdt + iDmdtRelax*iDmdtNew;
391 Info<<
"iDmdt." << pair.name()
392 <<
": min = " <<
min(iDmdt.primitiveField())
393 <<
", mean = " <<
average(iDmdt.primitiveField())
394 <<
", max = " <<
max(iDmdt.primitiveField())
404 bool wallBoilingActive =
false;
408 const phaseModel& phase = iter();
409 const phaseModel& otherPhase = iter.otherPhase();
415 "alphat." + phase.name()
422 "alphat." + phase.name()
428 const fvPatch& currPatch = patches[
patchi];
432 isA<alphatPhaseChangeWallFunction>
434 alphat.boundaryField()[
patchi]
438 const alphatPhaseChangeWallFunction& PCpatch =
439 refCast<const alphatPhaseChangeWallFunction>
441 alphat.boundaryField()[
patchi]
444 phasePairKey key(phase.name(), otherPhase.name());
446 if (PCpatch.activePhasePair(key))
448 wallBoilingActive =
true;
462 const label faceCelli =
463 currPatch.faceCells()[facei];
464 wDmdt[faceCelli] -=
sign*patchDmdt[facei];
465 wMDotL[faceCelli] -=
sign*patchMDotL[facei];
473 if (wallBoilingActive)
475 Info<<
"wDmdt." << pair.name()
476 <<
": min = " <<
min(wDmdt.primitiveField())
477 <<
", mean = " <<
average(wDmdt.primitiveField())
478 <<
", max = " <<
max(wDmdt.primitiveField())
486 template<
class BasePhaseSystem>
dimensionedScalar sign(const dimensionedScalar &ds)
A simple container for copying or transferring objects of type <T>.
virtual Xfer< PtrList< volScalarField > > dmdts() const
Return the mass transfer rates for each phase.
#define forAll(list, i)
Loop across all elements in list.
intWM_LABEL_SIZE_t label
A label is an int32_t or int64_t as specified by the pre-processor macro WM_LABEL_SIZE.
HashPtrTable< fvScalarMatrix > massTransferTable
HashTable< autoPtr< phasePair >, phasePairKey, phasePairKey::hash > phasePairTable
dimensioned< Type > max(const dimensioned< Type > &, const dimensioned< Type > &)
ThermalPhaseChangePhaseSystem(const fvMesh &)
Construct from fvMesh.
static int compare(const Pair< word > &a, const Pair< word > &b)
Compare Pairs.
phaseSystem::massTransferTable & massTransfer(massTransferPtr())
Ostream & endl(Ostream &os)
Add newline and flush stream.
virtual tmp< volScalarField > dmdt(const phasePairKey &key) const
Return the mass transfer rate for a pair.
virtual autoPtr< phaseSystem::massTransferTable > massTransfer() const
Return the mass transfer matrices.
HashPtrTable< fvScalarMatrix > heatTransferTable
dimensioned< Type > domainIntegrate(const GeometricField< Type, fvPatchField, volMesh > &vf)
CGAL::Exact_predicates_exact_constructions_kernel K
dimensionedScalar posPart(const dimensionedScalar &ds)
dimensionedScalar neg(const dimensionedScalar &ds)
tmp< DimensionedField< TypeR, GeoMesh > > New(const tmp< DimensionedField< TypeR, GeoMesh >> &tdf1, const word &name, const dimensionSet &dimensions)
GeometricField< scalar, fvPatchField, volMesh > volScalarField
dimensionedScalar pos(const dimensionedScalar &ds)
const dimensionSet dimVolume(pow3(dimLength))
stressControl lookup("compactNormalStress") >> compactNormalStress
virtual void correctInterfaceThermo()
Correct the interface thermodynamics.
Field< scalar > scalarField
Specialisation of Field<T> for scalar.
static word groupName(Name name, const word &group)
dimensionedScalar neg0(const dimensionedScalar &ds)
virtual autoPtr< phaseSystem::heatTransferTable > heatTransfer() const
Return the heat transfer matrices.
Volume integrate volField creating a volField.
dimensioned< Type > average(const DimensionedField< Type, GeoMesh > &df)
forAllConstIter(PtrDictionary< phaseModel >, mixture.phases(), phase)
dimensionedScalar pos0(const dimensionedScalar &ds)
const Mesh & mesh() const
Return mesh.
virtual bool read()
Read base phaseProperties dictionary.
tmp< volScalarField > iDmdt(const phasePairKey &key) const
Return the interfacial mass transfer rate for a pair.
tmp< volScalarField > wDmdt(const phasePairKey &key) const
Return the boundary mass transfer rate for a pair.
dimensioned< Type > min(const dimensioned< Type > &, const dimensioned< Type > &)
word name(const complex &)
Return a string representation of a complex.
const dimensionSet dimEnergy
const dimensionSet dimDensity
const saturationModel & saturation() const
Return the saturationModel.
PtrList< fvPatch > fvPatchList
container classes for fvPatch
virtual ~ThermalPhaseChangePhaseSystem()
Destructor.
virtual Xfer< PtrList< volScalarField > > dmdts() const
Return the mass transfer rates for each phase.
dimensioned< scalar > dimensionedScalar
Dimensioned scalar obtained from generic dimensioned type.
const dimensionSet dimTime(0, 0, 1, 0, 0, 0, 0)
virtual tmp< volScalarField > dmdt(const phasePairKey &key) const
Return the mass transfer rate for a pair.
An auto-pointer similar to the STL auto_ptr but with automatic casting to a reference to the type and...
A class for managing temporary objects.
Calculate the matrix for implicit and explicit sources.
dimensionedScalar negPart(const dimensionedScalar &ds)