32 template<
class BasePhaseSystem>
39 HeatAndMassTransferPhaseSystem<BasePhaseSystem>(mesh),
40 volatile_(this->
lookup(
"volatile")),
41 saturationModel_(saturationModel::
New(this->subDict(
"saturationModel"))),
42 massTransfer_(this->
lookup(
"massTransfer"))
52 const phasePair& pair(phasePairIter());
68 this->
mesh().time().timeName(),
83 template<
class BasePhaseSystem>
91 template<
class BasePhaseSystem>
95 return saturationModel_();
99 template<
class BasePhaseSystem>
103 typedef compressible::alphatPhaseChangeWallFunctionFvPatchScalarField
104 alphatPhaseChangeWallFunction;
106 autoPtr<phaseSystem::heatTransferTable> eqnsPtr =
119 const phasePair& pair(phasePairIter());
126 const phaseModel& phase = pair.phase1();
127 const phaseModel& otherPhase = pair.phase2();
134 phase.mesh().time().timeName(),
148 "alphat." + otherPhase.name()
155 "alphat." + otherPhase.name()
161 const fvPatch& currPatch = patches[
patchi];
165 isA<alphatPhaseChangeWallFunction>
167 alphat.boundaryField()[
patchi]
172 refCast<const alphatPhaseChangeWallFunction>
174 alphat.boundaryField()[
patchi]
179 label faceCelli = currPatch.faceCells()[facei];
180 mDotL[faceCelli] = patchMDotL[facei];
186 *eqns[otherPhase.name()] -= mDotL;
194 template<
class BasePhaseSystem>
199 autoPtr<phaseSystem::massTransferTable> eqnsPtr
208 const phaseModel& phase = this->phaseModels_[
phasei];
210 const PtrList<volScalarField>& Yi = phase.Y();
229 const phasePair& pair(phasePairIter());
235 const phaseModel& phase = pair.phase1();
236 const phaseModel& otherPhase = pair.phase2();
253 *eqns[otherName] += dmdt12 -
fvm::Sp(dmdt12, eqns[otherName]->
psi());
260 template<
class BasePhaseSystem>
263 typedef compressible::alphatPhaseChangeWallFunctionFvPatchScalarField
264 alphatPhaseChangeWallFunction;
266 BasePhaseSystem::correctThermo();
277 const phasePair& pair(phasePairIter());
284 const phaseModel& phase1 = pair.phase1();
285 const phaseModel& phase2 = pair.phase2();
287 Info<< phase1.name() <<
" min/max T " 288 <<
min(phase1.thermo().T()).value()
290 <<
max(phase1.thermo().T()).value()
293 Info<< phase2.name() <<
" min/max T " 294 <<
min(phase2.thermo().T()).value()
296 <<
max(phase2.thermo().T()).value()
317 this->heatTransferModels_[pair][pair.first()]->K(0)
322 this->heatTransferModels_[pair][pair.second()]->K(0)
325 Tf = saturationModel_->Tsat(phase1.thermo().p());
327 scalar iDmdtRelax(this->
mesh().fieldRelaxationFactor(
"iDmdt"));
330 (1 - iDmdtRelax)*iDmdt
331 + iDmdtRelax*(H1*(Tf - T1) + H2*(Tf - T2))
335 - (
neg(iDmdt)*he1 +
pos(iDmdt)*hef1),
339 Info<<
"iDmdt." << pair.name()
340 <<
": min = " <<
min(iDmdt.primitiveField())
341 <<
", mean = " <<
average(iDmdt.primitiveField())
342 <<
", max = " <<
max(iDmdt.primitiveField())
351 volScalarField H1(this->heatTransferModels_[pair][pair.first()]->K());
352 volScalarField H2(this->heatTransferModels_[pair][pair.second()]->K());
355 const scalar HLimit = 1
e-4;
356 H1.boundaryFieldRef() =
357 max(H1.boundaryField(), phase1.boundaryField()*HLimit);
358 H2.boundaryFieldRef() =
359 max(H2.boundaryField(), phase2.boundaryField()*HLimit);
366 - (
neg(iDmdt)*he1 +
pos(iDmdt)*hef1)
370 Tf = (H1*T1 + H2*T2 + mDotL)/(H1 + H2);
372 Info<<
"Tf." << pair.name()
373 <<
": min = " <<
min(Tf.primitiveField())
374 <<
", mean = " <<
average(Tf.primitiveField())
375 <<
", max = " <<
max(Tf.primitiveField())
384 this->
mesh().time().timeName(),
398 "alphat." + phase2.name()
405 "alphat." + phase2.name()
411 const fvPatch& currPatch = patches[
patchi];
415 isA<alphatPhaseChangeWallFunction>
417 alphat.boundaryField()[
patchi]
422 refCast<const alphatPhaseChangeWallFunction>
424 alphat.boundaryField()[
patchi]
429 label faceCelli = currPatch.faceCells()[facei];
430 wDmdt[faceCelli] += patchDmdt[facei];
435 Info<<
"wDmdt." << pair.name()
436 <<
": min = " <<
min(wDmdt.primitiveField())
437 <<
", mean = " <<
average(wDmdt.primitiveField())
438 <<
", max = " <<
max(wDmdt.primitiveField())
443 dmdt = wDmdt + iDmdt;
445 Info<<
"dmdt." << pair.name()
446 <<
": min = " <<
min(dmdt.primitiveField())
447 <<
", mean = " <<
average(dmdt.primitiveField())
448 <<
", max = " <<
max(dmdt.primitiveField())
455 template<
class BasePhaseSystem>
fvMatrix< scalar > fvScalarMatrix
#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.
const double e
Elementary charge.
HashTable< autoPtr< phasePair >, phasePairKey, phasePairKey::hash > phasePairTable
dimensioned< Type > max(const dimensioned< Type > &, const dimensioned< Type > &)
ThermalPhaseChangePhaseSystem(const fvMesh &)
Construct from fvMesh.
tmp< fvMatrix< Type > > Sp(const DimensionedField< scalar, volMesh > &, const GeometricField< Type, fvPatchField, volMesh > &)
virtual void correctThermo()
Correct the thermodynamics.
const saturationModel & saturation() const
Return the saturationModel.
Ostream & endl(Ostream &os)
Add newline and flush stream.
dimensioned< Type > domainIntegrate(const GeometricField< Type, fvPatchField, volMesh > &vf)
virtual autoPtr< phaseSystem::massTransferTable > massTransfer() const
Return the mass transfer matrices.
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
HashPtrTable< fvScalarMatrix, word, string::hash > heatTransferTable
dimensionedScalar pos(const dimensionedScalar &ds)
bool read(const char *, int32_t &)
stressControl lookup("compactNormalStress") >> compactNormalStress
Field< scalar > scalarField
Specialisation of Field<T> for scalar.
static word groupName(Name name, const word &group)
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)
virtual bool read()
Read base phaseProperties dictionary.
dimensioned< Type > min(const dimensioned< Type > &, const dimensioned< Type > &)
word name(const complex &)
Return a string representation of a complex.
const dimensionSet dimDensity
PtrList< fvPatch > fvPatchList
container classes for fvPatch
virtual ~ThermalPhaseChangePhaseSystem()
Destructor.
dimensioned< scalar > dimensionedScalar
Dimensioned scalar obtained from generic dimensioned type.
const Mesh & mesh() const
Return mesh.
const dimensionSet dimTime(0, 0, 1, 0, 0, 0, 0)
virtual autoPtr< phaseSystem::heatTransferTable > heatTransfer() const
Return the heat transfer matrices.
dimensioned< scalar > mag(const dimensioned< Type > &)
An auto-pointer similar to the STL auto_ptr but with automatic casting to a reference to the type and...
const volScalarField & psi
HashPtrTable< fvScalarMatrix, word, string::hash > massTransferTable
const dimensionSet dimMass(1, 0, 0, 0, 0, 0, 0)
dimensionedScalar negPart(const dimensionedScalar &ds)