Go to the documentation of this file.
35 template<
class BasePhaseSystem>
42 if (!iDmdt_.found(key))
44 return phaseSystem::dmdt(key);
49 return dmdtSign**iDmdt_[key];
53 template<
class BasePhaseSystem>
60 if (!wDmdt_.found(key))
62 return phaseSystem::dmdt(key);
67 return dmdtSign**wDmdt_[key];
73 template<
class BasePhaseSystem>
80 BasePhaseSystem(
mesh),
81 volatile_(this->
template lookupOrDefault<word>(
"volatile",
"none")),
86 phaseChange_(this->
lookup(
"phaseChange"))
111 IOobject::groupName(
"iDmdt", pair.
name()),
112 this->
mesh().time().timeName(),
114 IOobject::READ_IF_PRESENT,
130 IOobject::groupName(
"wDmdt", pair.
name()),
131 this->
mesh().time().timeName(),
133 IOobject::READ_IF_PRESENT,
149 IOobject::groupName(
"wMDotL", pair.
name()),
150 this->
mesh().time().timeName(),
152 IOobject::READ_IF_PRESENT,
165 template<
class BasePhaseSystem>
173 template<
class BasePhaseSystem>
177 return saturationModel_();
181 template<
class BasePhaseSystem>
188 return BasePhaseSystem::dmdt(key) + this->iDmdt(key) + this->wDmdt(key);
192 template<
class BasePhaseSystem>
200 const phasePair& pair = this->phasePairs_[iDmdtIter.key()];
203 this->addField(pair.
phase1(),
"dmdt", iDmdt, dmdts);
204 this->addField(pair.
phase2(),
"dmdt", - iDmdt, dmdts);
209 const phasePair& pair = this->phasePairs_[wDmdtIter.key()];
212 this->addField(pair.
phase1(),
"dmdt", wDmdt, dmdts);
213 this->addField(pair.
phase2(),
"dmdt", - wDmdt, dmdts);
220 template<
class BasePhaseSystem>
225 BasePhaseSystem::heatTransfer();
237 if (this->wMDotL_.found(phasePairIter.key()))
254 phase1.thermo().he().member() ==
"e"
255 ||
phase2.thermo().he().member() ==
"e"
260 this->iDmdt(pair) + this->wDmdt(pair)
263 if (
phase1.thermo().he().member() ==
"e")
269 if (
phase2.thermo().he().member() ==
"e")
282 template<
class BasePhaseSystem>
312 if (Yi[i].member() != volatile_)
319 IOobject::groupName(volatile_,
phase.
name())
324 IOobject::groupName(volatile_, otherPhase.
name())
333 *eqns[otherName] -= dmdt;
341 template<
class BasePhaseSystem>
346 alphatPhaseChangeWallFunction;
350 typename BasePhaseSystem::heatTransferModelTable,
351 this->heatTransferModels_,
352 heatTransferModelIter
357 this->phasePairs_[heatTransferModelIter.key()]
405 (
neg0(iDmdt)*hf2 +
pos(iDmdt)*h2)
406 - (
pos0(iDmdt)*hf1 +
neg(iDmdt)*h1)
416 iDmdtNew = (H1*(Tsat - T1) + H2*(Tsat - T2))/
L;
430 Tf = (H1*T1 + H2*T2 + iDmdtNew*
L)/(H1 + H2);
432 Info<<
"Tf." << pair.name()
438 scalar iDmdtRelax(this->
mesh().fieldRelaxationFactor(
"iDmdt"));
439 iDmdt = (1 - iDmdtRelax)*iDmdt + iDmdtRelax*iDmdtNew;
443 Info<<
"iDmdt." << pair.name()
444 <<
": min = " <<
min(iDmdt.primitiveField())
445 <<
", mean = " <<
average(iDmdt.primitiveField())
446 <<
", max = " <<
max(iDmdt.primitiveField())
456 bool wallBoilingActive =
false;
484 isA<alphatPhaseChangeWallFunction>
490 const alphatPhaseChangeWallFunction& PCpatch =
491 refCast<const alphatPhaseChangeWallFunction>
498 if (PCpatch.activePhasePair(key))
500 wallBoilingActive =
true;
514 const label faceCelli =
516 wDmdt[faceCelli] -=
sign*patchDmdt[facei];
517 wMDotL[faceCelli] -=
sign*patchMDotL[facei];
525 if (wallBoilingActive)
527 Info<<
"wDmdt." << pair.name()
538 template<
class BasePhaseSystem>
Single incompressible phase derived from the phase-fraction. Used as part of the multiPhaseMixture fo...
Defines the attributes of an object for which implicit objectRegistry management is supported,...
const vector L(dict.get< vector >("L"))
Description for mass transfer between a pair of phases. The direction of the mass transfer is from th...
A class for handling words, derived from Foam::string.
virtual PtrList< volScalarField > dmdts() const
Return the mass transfer rates for each phase.
tmp< volScalarField > iDmdt(const phasePairKey &key) const
Return the interfacial mass transfer rate for a pair.
Single incompressible phase derived from the phase-fraction. Used as part of the multiPhaseMixture fo...
A class for managing temporary objects.
const dimensionSet dimEnergy
const dimensionSet dimDensity
dimensioned< Type > domainIntegrate(const GeometricField< Type, fvPatchField, volMesh > &vf)
dimensionedScalar posPart(const dimensionedScalar &ds)
const Internal::FieldType & primitiveField() const
Return a const-reference to the internal field.
dimensionedScalar neg0(const dimensionedScalar &ds)
const phaseModel & otherPhase() const
Return the other phase in this two-phase system.
Abstract base-class for all alphatWallFunctions supporting phase-change.
Ostream & endl(Ostream &os)
Add newline and flush stream.
dimensionedScalar pos0(const dimensionedScalar &ds)
dimensionedScalar sign(const dimensionedScalar &ds)
label min(const labelHashSet &set, label minValue=labelMax)
Find the min value in labelHashSet, optionally limited by second argument.
#define forAllConstIter(Container, container, iter)
Iterate across all elements in the container object.
#define forAll(list, i)
Loop across all elements in list.
virtual bool read()
Read base phaseProperties dictionary.
CGAL::Exact_predicates_exact_constructions_kernel K
const dimensionSet dimTime(0, 0, 1, 0, 0, 0, 0)
intWM_LABEL_SIZE_t label
A label is an int32_t or int64_t as specified by the pre-processor macro WM_LABEL_SIZE.
tmp< GeometricField< Type, faPatchField, areaMesh > > average(const GeometricField< Type, faePatchField, edgeMesh > &ssf)
Area-weighted average a edgeField creating a areaField.
virtual ~ThermalPhaseChangePhaseSystem()
Destructor.
messageStream Info
Information stream (uses stdout - output is on the master only)
ThermalPhaseChangePhaseSystem(const fvMesh &)
Construct from fvMesh.
word name(const complex &c)
Return string representation of complex.
const cellModel * lookup(const word &modelName)
Deprecated(2017-11) equivalent to cellModel::ptr static method.
A finiteVolume patch using a polyPatch and a fvBoundaryMesh.
virtual void correctInterfaceThermo()
Correct the interface thermodynamics.
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....
label max(const labelHashSet &set, label maxValue=labelMin)
Find the max value in labelHashSet, optionally limited by second argument.
Mesh data needed to do the Finite Volume discretisation.
bool ordered() const
Return the ordered flag.
Calculate the matrix for implicit and explicit sources.
const word & name() const
tmp< DimensionedField< TypeR, GeoMesh > > New(const tmp< DimensionedField< TypeR, GeoMesh >> &tdf1, const word &name, const dimensionSet &dimensions)
Global function forwards to reuseTmpDimensionedField::New.
const word & name() const
Return the name of this phase.
Pointer management similar to std::unique_ptr, with some additional methods and type checking.
dimensionedScalar negPart(const dimensionedScalar &ds)
An ordered pair of two objects of type <T> with first() and second() elements.
virtual word name() const
Pair name.
const saturationModel & saturation() const
Return the saturationModel.
Volume integrate volField creating a volField.
const polyBoundaryMesh & patches
const phaseModel & phase2() const
Return phase 2.
virtual tmp< volScalarField > dmdt(const phasePairKey &key) const
Return the mass transfer rate for a pair.
virtual const labelUList & faceCells() const
Return faceCells.
void max(const dimensioned< Type > &dt)
Use the maximum of the field and specified value.
const dimensionSet dimVolume(pow3(dimLength))
phaseSystem::massTransferTable & massTransfer(massTransferPtr())
const phaseModel & phase1() const
Return phase 1.
dimensionedScalar neg(const dimensionedScalar &ds)
tmp< volScalarField > wDmdt(const phasePairKey &key) const
Return the boundary mass transfer rate for a pair.
virtual autoPtr< phaseSystem::heatTransferTable > heatTransfer() const
Return the heat transfer matrices.
const Boundary & boundaryField() const
Return const-reference to the boundary field.
virtual autoPtr< phaseSystem::massTransferTable > massTransfer() const
Return the mass transfer matrices.
dimensionedScalar pos(const dimensionedScalar &ds)