32 template<
class BasePhaseSystem>
35 const phaseSystem::dmdtfTable& dmdtfs,
36 phaseSystem::heatTransferTable& eqns
42 const phaseInterface interface(*
this, dmdtfIter.key());
48 const phaseModel& phase1 = interface.phase1();
49 const phaseModel& phase2 = interface.phase2();
50 const rhoThermo& thermo1 = phase1.thermo();
51 const rhoThermo& thermo2 = phase2.thermo();
60 *eqns[phase1.name()] +=
61 dmdtf*hs1 +
fvm::Sp(dmdtf12, he1) - dmdtf12*he1;
62 *eqns[phase2.name()] -=
63 dmdtf*hs2 +
fvm::Sp(dmdtf21, he2) - dmdtf21*he2;
66 *eqns[phase1.name()] += dmdtf21*(hs2 - hs1);
67 *eqns[phase2.name()] -= dmdtf12*(hs1 - hs2);
70 *eqns[phase1.name()] += dmdtf21*
K2 + dmdtf12*
K1;
71 *eqns[phase2.name()] -= dmdtf12*K1 + dmdtf21*
K2;
76 template<
class BasePhaseSystem>
79 const phaseSystem::dmidtfTable& dmidtfs,
80 phaseSystem::heatTransferTable& eqns
88 const phaseInterface interface(*
this, dmidtfIter.key());
90 const phaseModel& phase1 = interface.phase1();
91 const phaseModel& phase2 = interface.phase2();
92 const rhoThermo& thermo1 = phase1.thermo();
93 const rhoThermo& thermo2 = phase2.thermo();
94 const basicSpecieMixture* compositionPtr1 =
96 ? &refCast<const rhoReactionThermo>(thermo1).composition()
97 :
static_cast<const basicSpecieMixture*
>(
nullptr);
98 const basicSpecieMixture* compositionPtr2 =
100 ? &refCast<const rhoReactionThermo>(thermo2).composition()
101 :
static_cast<const basicSpecieMixture*
>(
nullptr);
110 forAllConstIter(HashPtrTable<volScalarField>, *dmidtfIter(), dmidtfJter)
112 const word& specie = dmidtfJter.key();
120 const label speciei1 =
121 compositionPtr1 ? compositionPtr1->species()[specie] : -1;
122 const label speciei2 =
123 compositionPtr2 ? compositionPtr2->species()[specie] : -1;
129 ? compositionPtr1->Hs(speciei1, thermo1.p(), thermo1.T())
130 : tmp<volScalarField>(hs1)
135 ? compositionPtr2->Hs(speciei2, thermo2.p(), thermo2.T())
136 : tmp<volScalarField>(hs2)
140 tmp<volScalarField> tYi1, tYi2;
145 ?
max(compositionPtr1->Y(speciei1), residualY_)
149 ?
max(compositionPtr2->Y(speciei2), residualY_)
154 *eqns[phase1.name()] += dmidtf*hsi1;
155 *eqns[phase2.name()] -= dmidtf*hsi2;
158 *eqns[phase1.name()] +=
159 fvm::Sp(dmidtf12/tYi1(), he1) - dmidtf12/tYi1()*he1;
160 *eqns[phase2.name()] -=
161 fvm::Sp(dmidtf21/tYi2(), he2) - dmidtf21/tYi2()*he2;
165 *eqns[phase1.name()] += dmidtf21*(hsi2 - hsi1);
166 *eqns[phase2.name()] -= dmidtf12*(hsi1 - hsi2);
169 *eqns[phase1.name()] += dmidtf21*
K2 + dmidtf12*
K1;
170 *eqns[phase2.name()] -= dmidtf12*K1 + dmidtf21*
K2;
176 template<
class BasePhaseSystem>
179 const phaseSystem::dmdtfTable& dmdtfs,
180 const phaseSystem::dmdtfTable& Tfs,
181 const latentHeatScheme scheme,
182 phaseSystem::heatTransferTable& eqns
188 const phaseInterface interface(*
this, dmdtfIter.key());
196 const phaseModel& phase1 = interface.phase1();
197 const phaseModel& phase2 = interface.phase2();
198 const rhoThermo& thermo1 = phase1.thermo();
199 const rhoThermo& thermo2 = phase2.thermo();
212 case latentHeatScheme::symmetric:
214 *eqns[phase1.name()] += dmdtf*hsf1;
215 *eqns[phase2.name()] -= dmdtf*hsf2;
219 case latentHeatScheme::upwind:
225 *eqns[phase1.name()] += dmdtf21*hsf1 + dmdtf12*hs1;
226 *eqns[phase2.name()] -= dmdtf12*hsf2 + dmdtf21*hs2;
231 *eqns[phase1.name()] +=
fvm::Sp(dmdtf12, he1) - dmdtf12*he1;
232 *eqns[phase2.name()] -=
fvm::Sp(dmdtf21, he2) - dmdtf21*he2;
235 *eqns[phase1.name()] += dmdtf21*
K2 + dmdtf12*
K1;
236 *eqns[phase2.name()] -= dmdtf12*K1 + dmdtf21*
K2;
241 template<
class BasePhaseSystem>
244 const phaseSystem::dmdtfTable& dmdtfs,
245 const phaseSystem::dmdtfTable& Tfs,
247 const latentHeatScheme scheme,
248 phaseSystem::heatTransferTable& eqns
254 const phaseInterface interface(*
this, dmdtfIter.key());
262 const phaseModel& phase1 = interface.phase1();
263 const phaseModel& phase2 = interface.phase2();
267 *eqns[phase1.name()] += ((1 - weight)*dmdtf12 + weight*dmdtf21)*L;
268 *eqns[phase2.name()] += ((1 - weight)*dmdtf21 + weight*dmdtf12)*L;
273 template<
class BasePhaseSystem>
276 const phaseSystem::dmdtfTable& dmdtfs,
277 const phaseSystem::dmdtfTable& Tfs,
279 const latentHeatScheme scheme,
280 phaseSystem::heatTransferTable& eqns
283 addDmdtHefsWithoutL(dmdtfs, Tfs, scheme, eqns);
284 addDmdtL(dmdtfs, Tfs, weight, scheme, eqns);
288 template<
class BasePhaseSystem>
291 const phaseSystem::dmidtfTable& dmidtfs,
292 const phaseSystem::dmdtfTable& Tfs,
293 const latentHeatScheme scheme,
294 phaseSystem::heatTransferTable& eqns
302 const phaseInterface interface(*
this, dmidtfIter.key());
306 const phaseModel& phase1 = interface.phase1();
307 const phaseModel& phase2 = interface.phase2();
308 const rhoThermo& thermo1 = phase1.thermo();
309 const rhoThermo& thermo2 = phase2.thermo();
310 const basicSpecieMixture* compositionPtr1 =
311 isA<rhoReactionThermo>(
thermo1)
312 ? &refCast<const rhoReactionThermo>(thermo1).composition()
313 :
static_cast<const basicSpecieMixture*
>(
nullptr);
314 const basicSpecieMixture* compositionPtr2 =
315 isA<rhoReactionThermo>(
thermo2)
316 ? &refCast<const rhoReactionThermo>(thermo2).composition()
317 :
static_cast<const basicSpecieMixture*
>(
nullptr);
328 forAllConstIter(HashPtrTable<volScalarField>, *dmidtfIter(), dmidtfJter)
330 const word& specie = dmidtfJter.key();
338 const label speciei1 =
339 compositionPtr1 ? compositionPtr1->species()[specie] : -1;
340 const label speciei2 =
341 compositionPtr2 ? compositionPtr2->species()[specie] : -1;
347 ? compositionPtr1->Hs(speciei1, thermo1.p(), Tf)
348 : tmp<volScalarField>(hsf1)
353 ? compositionPtr2->Hs(speciei2, thermo2.p(), Tf)
354 : tmp<volScalarField>(hsf2)
358 tmp<volScalarField> tYi1, tYi2;
359 if (this->residualY_ > 0)
363 ?
max(compositionPtr1->Y(speciei1), this->residualY_)
367 ?
max(compositionPtr2->Y(speciei2), this->residualY_)
374 case latentHeatScheme::symmetric:
376 *eqns[phase1.name()] += dmidtf*hsfi1;
377 *eqns[phase2.name()] -= dmidtf*hsfi2;
381 case latentHeatScheme::upwind:
387 ? compositionPtr1->Hs(speciei1, thermo1.p(), thermo1.T())
393 ? compositionPtr2->Hs(speciei2, thermo2.p(), thermo2.T())
397 *eqns[phase1.name()] += dmidtf21*hsfi1 + dmidtf12*hsi1;
398 *eqns[phase2.name()] -= dmidtf12*hsfi2 + dmidtf21*hsi2;
403 if (this->residualY_ > 0)
405 *eqns[phase1.name()] +=
406 fvm::Sp(dmidtf12/tYi1(), he1) - dmidtf12/tYi1()*he1;
408 if (this->residualY_ > 0)
410 *eqns[phase2.name()] -=
411 fvm::Sp(dmidtf21/tYi2(), he2) - dmidtf21/tYi2()*he2;
416 *eqns[phase1.name()] += dmidtf21*
K2 + dmidtf12*
K1;
417 *eqns[phase2.name()] -= dmidtf12*K1 + dmidtf21*
K2;
423 template<
class BasePhaseSystem>
426 const phaseSystem::dmidtfTable& dmidtfs,
427 const phaseSystem::dmdtfTable& Tfs,
429 const latentHeatScheme scheme,
430 phaseSystem::heatTransferTable& eqns
436 const phaseInterface interface(*
this, dmidtfIter.key());
440 const phaseModel& phase1 = interface.phase1();
441 const phaseModel& phase2 = interface.phase2();
444 forAllConstIter(HashPtrTable<volScalarField>, *dmidtfIter(), dmidtfJter)
446 const word& specie = dmidtfJter.key();
456 this->Li(interface, specie, dmidtf, Tf, scheme)
458 *eqns[phase1.name()] +=
459 ((1 - weight)*dmidtf12 + weight*dmidtf21)*Li;
460 *eqns[phase2.name()] +=
461 ((1 - weight)*dmidtf21 + weight*dmidtf12)*Li;
467 template<
class BasePhaseSystem>
470 const phaseSystem::dmidtfTable& dmidtfs,
471 const phaseSystem::dmdtfTable& Tfs,
473 const latentHeatScheme scheme,
474 phaseSystem::heatTransferTable& eqns
477 addDmidtHefsWithoutL(dmidtfs, Tfs, scheme, eqns);
478 addDmidtL(dmidtfs, Tfs, weight, scheme, eqns);
484 template<
class BasePhaseSystem>
490 heatTransferPhaseSystem(),
491 BasePhaseSystem(mesh),
492 residualY_(this->template lookupOrDefault<scalar>(
"residualY", -1))
498 template<
class BasePhaseSystem>
505 template<
class BasePhaseSystem>
509 const phaseInterface& interface,
512 const latentHeatScheme scheme
515 const rhoThermo& thermo1 = interface.phase1().thermo();
516 const rhoThermo& thermo2 = interface.phase2().thermo();
524 case latentHeatScheme::symmetric:
528 case latentHeatScheme::upwind:
535 neg0(dmdtf)*haf2 +
pos(dmdtf)*ha2
536 -
pos0(dmdtf)*haf1 -
neg(dmdtf)*ha1;
540 return tmp<volScalarField>(
nullptr);
544 template<
class BasePhaseSystem>
548 const phaseInterface& interface,
552 const latentHeatScheme scheme
555 const rhoThermo& thermo1 = interface.phase1().thermo();
556 const rhoThermo& thermo2 = interface.phase2().thermo();
564 case latentHeatScheme::symmetric:
568 case latentHeatScheme::upwind:
578 neg0(dmdtf)*haf2 +
pos(dmdtf)*ha2
579 -
pos0(dmdtf)*haf1 -
neg(dmdtf)*ha1;
583 return tmp<scalarField>(
nullptr);
587 template<
class BasePhaseSystem>
591 const phaseInterface& interface,
595 const latentHeatScheme scheme
598 const rhoThermo& thermo1 = interface.phase1().thermo();
599 const rhoThermo& thermo2 = interface.phase2().thermo();
600 const basicSpecieMixture* compositionPtr1 =
601 isA<rhoReactionThermo>(
thermo1)
602 ? &refCast<const rhoReactionThermo>(thermo1).composition()
603 :
static_cast<const basicSpecieMixture*
>(
nullptr);
604 const basicSpecieMixture* compositionPtr2 =
605 isA<rhoReactionThermo>(
thermo2)
606 ? &refCast<const rhoReactionThermo>(thermo2).composition()
607 :
static_cast<const basicSpecieMixture*
>(
nullptr);
608 const label speciei1 =
609 compositionPtr1 ? compositionPtr1->species()[specie] : -1;
610 const label speciei2 =
611 compositionPtr2 ? compositionPtr2->species()[specie] : -1;
617 ? compositionPtr1->Ha(speciei1, thermo1.p(), Tf)
618 : thermo1.ha(thermo1.p(), Tf)
623 ? compositionPtr2->Ha(speciei2, thermo2.p(), Tf)
624 : thermo2.ha(thermo1.p(), Tf)
629 case latentHeatScheme::symmetric:
631 return hafi2 - hafi1;
633 case latentHeatScheme::upwind:
639 ? compositionPtr1->Ha(speciei1, thermo1.p(), thermo1.T())
645 ? compositionPtr2->Ha(speciei2, thermo2.p(), thermo2.T())
650 neg0(dmdtf)*hafi2 +
pos(dmdtf)*hai2
651 -
pos0(dmdtf)*hafi1 -
neg(dmdtf)*hai1;
655 return tmp<volScalarField>(
nullptr);
659 template<
class BasePhaseSystem>
663 const phaseInterface& interface,
668 const latentHeatScheme scheme
671 const rhoThermo& thermo1 = interface.phase1().thermo();
672 const rhoThermo& thermo2 = interface.phase2().thermo();
673 const basicSpecieMixture* compositionPtr1 =
674 isA<rhoReactionThermo>(
thermo1)
675 ? &refCast<const rhoReactionThermo>(thermo1).composition()
676 :
static_cast<const basicSpecieMixture*
>(
nullptr);
677 const basicSpecieMixture* compositionPtr2 =
678 isA<rhoReactionThermo>(
thermo2)
679 ? &refCast<const rhoReactionThermo>(thermo2).composition()
680 :
static_cast<const basicSpecieMixture*
>(
nullptr);
681 const label speciei1 =
682 compositionPtr1 ? compositionPtr1->species()[specie] : -1;
683 const label speciei2 =
684 compositionPtr2 ? compositionPtr2->species()[specie] : -1;
693 ? compositionPtr1->Ha(speciei1, p1, Tf)
694 : thermo1.ha(Tf, cells)
699 ? compositionPtr2->Ha(speciei2, p2, Tf)
700 : thermo2.ha(Tf, cells)
705 case latentHeatScheme::symmetric:
707 return hafi2 - hafi1;
709 case latentHeatScheme::upwind:
718 ? compositionPtr1->Ha(speciei1, p1, T1)
719 : thermo1.ha(T1, cells)
724 ? compositionPtr2->Ha(speciei2, p2, T2)
725 : thermo2.ha(T2, cells)
729 neg0(dmdtf)*hafi2 +
pos(dmdtf)*hai2
730 -
pos0(dmdtf)*hafi1 -
neg(dmdtf)*hai1;
734 return tmp<scalarField>(
nullptr);
738 template<
class BasePhaseSystem>
virtual tmp< volScalarField > Li(const phaseInterface &interface, const word &member, const volScalarField &dmidtf, const volScalarField &Tf, const latentHeatScheme scheme) const
Return the latent heat for a given interface, specie, mass transfer.
layerAndWeight max(const layerAndWeight &a, const layerAndWeight &b)
FvWallInfoData< WallInfo, label > label
A label is an int32_t or int64_t as specified by the pre-processor macro WM_LABEL_SIZE.
autoPtr< CompressibleMomentumTransportModel > New(const volScalarField &rho, const volVectorField &U, const surfaceScalarField &phi, const viscosity &viscosity)
#define forAllConstIter(Container, container, iter)
Iterate across all elements in the container object of type.
void addDmdtHefsWithoutL(const phaseSystem::dmdtfTable &dmdtfs, const phaseSystem::dmdtfTable &Tfs, const latentHeatScheme scheme, phaseSystem::heatTransferTable &eqns) const
Add energy transfer terms which result from bulk phase changes,.
const dimensionSet dimless
mixture thermo1().correctRho(psi1 *(p_rgh - p_rgh_0))
dimensionedScalar posPart(const dimensionedScalar &ds)
dimensionedScalar neg(const dimensionedScalar &ds)
UList< label > labelUList
void addDmdtHefs(const phaseSystem::dmdtfTable &dmdtfs, phaseSystem::heatTransferTable &eqns) const
Add energy transfer terms which result from bulk mass transfers.
GeometricField< scalar, fvPatchField, volMesh > volScalarField
void addDmidtHefs(const phaseSystem::dmidtfTable &dmidtfs, phaseSystem::heatTransferTable &eqns) const
Add energy transfer terms which result from specie mass transfers.
dimensionedScalar pos(const dimensionedScalar &ds)
mixture thermo2().correctRho(psi2 *(p_rgh - p_rgh_0))
Field< scalar > scalarField
Specialisation of Field<T> for scalar.
virtual tmp< volScalarField > L(const phaseInterface &interface, const volScalarField &dmdtf, const volScalarField &Tf, const latentHeatScheme scheme) const
Return the latent heat for a given interface, mass transfer rate.
virtual ~HeatTransferPhaseSystem()
Destructor.
dimensionedScalar neg0(const dimensionedScalar &ds)
dimensionedScalar pos0(const dimensionedScalar &ds)
void addDmidtL(const phaseSystem::dmidtfTable &dmidtfs, const phaseSystem::dmdtfTable &Tfs, const scalar weight, const latentHeatScheme scheme, phaseSystem::heatTransferTable &eqns) const
Add latent heat terms which result from specie phase changes.
void addDmdtL(const phaseSystem::dmdtfTable &dmdtfs, const phaseSystem::dmdtfTable &Tfs, const scalar weight, const latentHeatScheme scheme, phaseSystem::heatTransferTable &eqns) const
Add latent heat terms which result from bulk phase changes.
const tmp< volScalarField::Internal > & Sp
void addDmidtHefsWithoutL(const phaseSystem::dmidtfTable &dmidtfs, const phaseSystem::dmdtfTable &Tfs, const latentHeatScheme scheme, phaseSystem::heatTransferTable &eqns) const
Add energy transfer terms which result from specie phase changes,.
dimensioned< scalar > dimensionedScalar
Dimensioned scalar obtained from generic dimensioned type.
A class for managing temporary objects.
HeatTransferPhaseSystem(const fvMesh &)
Construct from fvMesh.
virtual bool read()
Read base phaseProperties dictionary.
Calculate the matrix for implicit and explicit sources.
dimensionedScalar negPart(const dimensionedScalar &ds)