36template<
class BasePhaseSystem>
43 if (!iDmdt_.found(key))
45 return phaseSystem::dmdt(key);
50 return dmdtSign**iDmdt_[key];
54template<
class BasePhaseSystem>
61 if (!wDmdt_.found(key))
63 return phaseSystem::dmdt(key);
68 return dmdtSign**wDmdt_[key];
74template<
class BasePhaseSystem>
81 BasePhaseSystem(
mesh),
82 volatile_(this->template getOrDefault<
word>(
"volatile",
"none")),
87 phaseChange_(this->
lookup(
"phaseChange"))
113 this->mesh().time().timeName(),
132 this->mesh().time().timeName(),
151 this->mesh().time().timeName(),
166template<
class BasePhaseSystem>
174template<
class BasePhaseSystem>
178 return saturationModel_();
182template<
class BasePhaseSystem>
189 return BasePhaseSystem::dmdt(key) + this->iDmdt(key) + this->wDmdt(key);
193template<
class BasePhaseSystem>
201 const phasePair& pair = this->phasePairs_[iDmdtIter.key()];
204 this->addField(pair.
phase1(),
"dmdt", iDmdt, dmdts);
205 this->addField(pair.
phase2(),
"dmdt", - iDmdt, dmdts);
210 const phasePair& pair = this->phasePairs_[wDmdtIter.key()];
213 this->addField(pair.
phase1(),
"dmdt", wDmdt, dmdts);
214 this->addField(pair.
phase2(),
"dmdt", - wDmdt, dmdts);
221template<
class BasePhaseSystem>
226 BasePhaseSystem::heatTransfer();
238 if (this->wMDotL_.found(phasePairIter.key()))
261 this->iDmdt(pair) + this->wDmdt(pair)
283template<
class BasePhaseSystem>
288 BasePhaseSystem::massTransfer();
313 if (Yi[i].member() != volatile_)
334 *eqns[otherName] -= dmdt;
342template<
class BasePhaseSystem>
347 alphatPhaseChangeWallFunction;
351 typename BasePhaseSystem::heatTransferModelTable,
352 this->heatTransferModels_,
353 heatTransferModelIter
358 this->phasePairs_[heatTransferModelIter.key()]
406 (
neg0(iDmdt)*hf2 +
pos(iDmdt)*h2)
407 - (
pos0(iDmdt)*hf1 +
neg(iDmdt)*h1)
417 iDmdtNew = (H1*(Tsat - T1) + H2*(Tsat - T2))/
L;
431 Tf = (H1*T1 + H2*T2 + iDmdtNew*
L)/(H1 + H2);
439 scalar iDmdtRelax(this->
mesh().fieldRelaxationFactor(
"iDmdt"));
440 iDmdt = (1 - iDmdtRelax)*iDmdt + iDmdtRelax*iDmdtNew;
457 bool wallBoilingActive =
false;
485 isA<alphatPhaseChangeWallFunction>
491 const alphatPhaseChangeWallFunction& PCpatch =
492 refCast<const alphatPhaseChangeWallFunction>
499 if (PCpatch.activePhasePair(key))
501 wallBoilingActive =
true;
515 const label faceCelli =
517 wDmdt[faceCelli] -=
sign*patchDmdt[facei];
518 wMDotL[faceCelli] -=
sign*patchMDotL[facei];
526 if (wallBoilingActive)
539template<
class BasePhaseSystem>
542 if (BasePhaseSystem::read())
CGAL::Exact_predicates_exact_constructions_kernel K
tmp< volScalarField > dmdt() const
Return the blended mass transfer rate.
const dimensionSet & dimensions() const
Return dimensions.
const Mesh & mesh() const
Return mesh.
const Internal::FieldType & primitiveField() const
Return a const-reference to the internal field.
void max(const dimensioned< Type > &dt)
Use the maximum of the field and specified value.
const Boundary & boundaryField() const
Return const-reference to the boundary field.
bool set(const Key &key, T *ptr)
Assign a new entry, overwriting existing entries.
Defines the attributes of an object for which implicit objectRegistry management is supported,...
static word member(const word &name)
Return member (name without the extension)
static word groupName(StringType base, const word &group)
Create dot-delimited name.group string.
An ordered pair of two objects of type <T> with first() and second() elements.
A list of pointers to objects of type <T>, with allocation/deallocation management of the pointers....
Class to provide interfacial heat and mass transfer between a number of phases according the interfac...
virtual ~ThermalPhaseChangePhaseSystem()
Destructor.
virtual autoPtr< phaseSystem::heatTransferTable > heatTransfer() const
Return the heat transfer matrices.
virtual void correctInterfaceThermo()
Correct the interface thermodynamics.
tmp< volScalarField > wDmdt(const phasePairKey &key) const
Return the boundary mass transfer rate for a pair.
wMDotLTable wMDotL_
Boundary thermal energy transfer rate.
virtual autoPtr< phaseSystem::massTransferTable > massTransfer() const
Return the mass transfer matrices.
iDmdtTable iDmdt_
Interfacial Mass transfer rate.
tmp< volScalarField > iDmdt(const phasePairKey &key) const
Return the interfacial mass transfer rate for a pair.
virtual PtrList< volScalarField > dmdts() const
Return the mass transfer rates for each phase.
wDmdtTable wDmdt_
Boundary Mass transfer rate.
virtual bool read()
Read base phaseProperties dictionary.
const saturationModel & saturation() const
Return the saturationModel.
Pointer management similar to std::unique_ptr, with some additional methods and type checking.
virtual const volScalarField & T() const
Temperature [K].
virtual volScalarField & p()
Pressure [Pa].
virtual volScalarField & he()=0
Enthalpy/Internal energy [J/kg].
Abstract base-class for all alphatWallFunctions supporting phase-change.
Mesh data needed to do the Finite Volume discretisation.
const fvBoundaryMesh & boundary() const
Return reference to boundary mesh.
A finiteVolume patch using a polyPatch and a fvBoundaryMesh.
virtual const labelUList & faceCells() const
Return faceCells.
Single incompressible phase derived from the phase-fraction. Used as part of the multiPhaseMixture fo...
const dimensionedScalar & rho() const
const phaseModel & otherPhase() const
Return the other phase in this two-phase system.
const word & name() const
virtual const rhoThermo & thermo() const =0
Return the thermophysical model.
An ordered or unorder pair of phase names. Typically specified as follows.
bool ordered() const noexcept
Return the ordered flag.
Description for mass transfer between a pair of phases. The direction of the mass transfer is from th...
virtual word name() const
Pair name.
const multiphaseInter::phaseModel & phase1() const
const multiphaseInter::phaseModel & phase2() const
Single incompressible phase derived from the phase-fraction. Used as part of the multiPhaseMixture fo...
const word & name() const
Lookup type of boundary radiation properties.
virtual tmp< volScalarField > rho() const
Density [kg/m^3].
A class for managing temporary objects.
A class for handling words, derived from Foam::string.
const polyBoundaryMesh & patches
Volume integrate volField creating a volField.
Calculate the finiteVolume matrix for implicit and explicit sources.
dimensioned< Type > domainIntegrate(const GeometricField< Type, fvPatchField, volMesh > &vf)
dimensioned< scalar > dimensionedScalar
Dimensioned scalar obtained from generic dimensioned type.
dimensionedScalar pos(const dimensionedScalar &ds)
dimensionedScalar pos0(const dimensionedScalar &ds)
dimensionedScalar sign(const dimensionedScalar &ds)
const dimensionSet dimEnergy
const dimensionSet dimTime(0, 0, 1, 0, 0, 0, 0)
messageStream Info
Information stream (stdout output on master, null elsewhere)
Ostream & endl(Ostream &os)
Add newline and flush stream.
dimensionedScalar negPart(const dimensionedScalar &ds)
dimensionedScalar neg(const dimensionedScalar &ds)
Type gAverage(const FieldField< Field, Type > &f)
tmp< DimensionedField< TypeR, GeoMesh > > New(const tmp< DimensionedField< TypeR, GeoMesh > > &tdf1, const word &name, const dimensionSet &dimensions)
Global function forwards to reuseTmpDimensionedField::New.
dimensionedScalar neg0(const dimensionedScalar &ds)
Type gMin(const FieldField< Field, Type > &f)
const dimensionSet dimVolume(pow3(dimLength))
word name(const expressions::valueTypeCode typeCode)
A word representation of a valueTypeCode. Empty for INVALID.
const dimensionSet dimDensity
Type gMax(const FieldField< Field, Type > &f)
dimensionedScalar posPart(const dimensionedScalar &ds)
#define forAll(list, i)
Loop across all elements in list.
#define forAllConstIter(Container, container, iter)
Iterate across all elements in the container object.
const vector L(dict.get< vector >("L"))