32 template<
class BasePhaseSystem>
35 const phaseSystem::dmdtfTable& dmdtfs,
36 phaseSystem::heatTransferTable& eqns
42 const phasePairKey& key = dmdtfIter.key();
43 const phasePair& pair(this->phasePairs_[key]);
49 const phaseModel& phase1 = pair.phase1();
50 const phaseModel& phase2 = pair.phase2();
51 const rhoThermo& thermo1 = phase1.thermo();
52 const rhoThermo& thermo2 = phase2.thermo();
61 *eqns[phase1.name()] +=
62 dmdtf*hs1 +
fvm::Sp(dmdtf12, he1) - dmdtf12*he1;
63 *eqns[phase2.name()] -=
64 dmdtf*hs2 +
fvm::Sp(dmdtf21, he2) - dmdtf21*he2;
67 *eqns[phase1.name()] += dmdtf21*(hs2 - hs1);
68 *eqns[phase2.name()] -= dmdtf12*(hs1 - hs2);
71 *eqns[phase1.name()] += dmdtf21*
K2 + dmdtf12*
K1;
72 *eqns[phase2.name()] -= dmdtf12*K1 + dmdtf21*
K2;
77 template<
class BasePhaseSystem>
80 const phaseSystem::dmidtfTable& dmidtfs,
81 phaseSystem::heatTransferTable& eqns
89 const phasePairKey& key = dmidtfIter.key();
90 const phasePair& pair(this->phasePairs_[key]);
92 const phaseModel& phase1 = pair.phase1();
93 const phaseModel& phase2 = pair.phase2();
94 const rhoThermo& thermo1 = phase1.thermo();
95 const rhoThermo& thermo2 = phase2.thermo();
96 const basicSpecieMixture* compositionPtr1 =
98 ? &refCast<const rhoReactionThermo>(thermo1).composition()
99 :
static_cast<const basicSpecieMixture*
>(
nullptr);
100 const basicSpecieMixture* compositionPtr2 =
101 isA<rhoReactionThermo>(
thermo2)
102 ? &refCast<const rhoReactionThermo>(thermo2).composition()
103 :
static_cast<const basicSpecieMixture*
>(
nullptr);
112 forAllConstIter(HashPtrTable<volScalarField>, *dmidtfIter(), dmidtfJter)
114 const word& specie = dmidtfJter.key();
125 const label speciei1 =
126 compositionPtr1 ? compositionPtr1->species()[specie] : -1;
127 const label speciei2 =
128 compositionPtr2 ? compositionPtr2->species()[specie] : -1;
134 ? compositionPtr1->Hs(speciei1, thermo1.p(), thermo1.T())
135 : tmp<volScalarField>(hs1)
140 ? compositionPtr2->Hs(speciei2, thermo2.p(), thermo2.T())
141 : tmp<volScalarField>(hs2)
145 tmp<volScalarField> tYi1, tYi2;
150 ?
max(compositionPtr1->Y(speciei1), residualY_)
154 ?
max(compositionPtr2->Y(speciei2), residualY_)
159 *eqns[phase1.name()] += dmidtf*hsi1;
160 *eqns[phase2.name()] -= dmidtf*hsi2;
163 *eqns[phase1.name()] +=
164 fvm::Sp(dmidtf12/tYi1(), he1) - dmidtf12/tYi1()*he1;
165 *eqns[phase2.name()] -=
166 fvm::Sp(dmidtf21/tYi2(), he2) - dmidtf21/tYi2()*he2;
170 *eqns[phase1.name()] += dmidtf21*(hsi2 - hsi1);
171 *eqns[phase2.name()] -= dmidtf12*(hsi1 - hsi2);
174 *eqns[phase1.name()] += dmidtf21*
K2 + dmidtf12*
K1;
175 *eqns[phase2.name()] -= dmidtf12*K1 + dmidtf21*
K2;
181 template<
class BasePhaseSystem>
184 const phaseSystem::dmdtfTable& dmdtfs,
185 const phaseSystem::dmdtfTable& Tfs,
186 const latentHeatScheme scheme,
187 phaseSystem::heatTransferTable& eqns
193 const phasePairKey& key = dmdtfIter.key();
194 const phasePair& pair(this->phasePairs_[key]);
202 const phaseModel& phase1 = pair.phase1();
203 const phaseModel& phase2 = pair.phase2();
204 const rhoThermo& thermo1 = phase1.thermo();
205 const rhoThermo& thermo2 = phase2.thermo();
218 case latentHeatScheme::symmetric:
220 *eqns[phase1.name()] += dmdtf*hsf1;
221 *eqns[phase2.name()] -= dmdtf*hsf2;
225 case latentHeatScheme::upwind:
231 *eqns[phase1.name()] += dmdtf21*hsf1 + dmdtf12*hs1;
232 *eqns[phase2.name()] -= dmdtf12*hsf2 + dmdtf21*hs2;
237 *eqns[phase1.name()] +=
fvm::Sp(dmdtf12, he1) - dmdtf12*he1;
238 *eqns[phase2.name()] -=
fvm::Sp(dmdtf21, he2) - dmdtf21*he2;
241 *eqns[phase1.name()] += dmdtf21*
K2 + dmdtf12*
K1;
242 *eqns[phase2.name()] -= dmdtf12*K1 + dmdtf21*
K2;
247 template<
class BasePhaseSystem>
250 const phaseSystem::dmdtfTable& dmdtfs,
251 const phaseSystem::dmdtfTable& Tfs,
253 const latentHeatScheme scheme,
254 phaseSystem::heatTransferTable& eqns
260 const phasePairKey& key = dmdtfIter.key();
261 const phasePair& pair(this->phasePairs_[key]);
269 const phaseModel& phase1 = pair.phase1();
270 const phaseModel& phase2 = pair.phase2();
274 *eqns[phase1.name()] += ((1 - weight)*dmdtf12 + weight*dmdtf21)*L;
275 *eqns[phase2.name()] += ((1 - weight)*dmdtf21 + weight*dmdtf12)*L;
280 template<
class BasePhaseSystem>
283 const phaseSystem::dmdtfTable& dmdtfs,
284 const phaseSystem::dmdtfTable& Tfs,
286 const latentHeatScheme scheme,
287 phaseSystem::heatTransferTable& eqns
290 addDmdtHefsWithoutL(dmdtfs, Tfs, scheme, eqns);
291 addDmdtL(dmdtfs, Tfs, weight, scheme, eqns);
295 template<
class BasePhaseSystem>
298 const phaseSystem::dmidtfTable& dmidtfs,
299 const phaseSystem::dmdtfTable& Tfs,
300 const latentHeatScheme scheme,
301 phaseSystem::heatTransferTable& eqns
309 const phasePairKey& key = dmidtfIter.key();
310 const phasePair& pair(this->phasePairs_[key]);
314 const phaseModel& phase1 = pair.phase1();
315 const phaseModel& phase2 = pair.phase2();
316 const rhoThermo& thermo1 = phase1.thermo();
317 const rhoThermo& thermo2 = phase2.thermo();
318 const basicSpecieMixture* compositionPtr1 =
319 isA<rhoReactionThermo>(
thermo1)
320 ? &refCast<const rhoReactionThermo>(thermo1).composition()
321 :
static_cast<const basicSpecieMixture*
>(
nullptr);
322 const basicSpecieMixture* compositionPtr2 =
323 isA<rhoReactionThermo>(
thermo2)
324 ? &refCast<const rhoReactionThermo>(thermo2).composition()
325 :
static_cast<const basicSpecieMixture*
>(
nullptr);
336 forAllConstIter(HashPtrTable<volScalarField>, *dmidtfIter(), dmidtfJter)
338 const word& specie = dmidtfJter.key();
349 const label speciei1 =
350 compositionPtr1 ? compositionPtr1->species()[specie] : -1;
351 const label speciei2 =
352 compositionPtr2 ? compositionPtr2->species()[specie] : -1;
358 ? compositionPtr1->Hs(speciei1, thermo1.p(), Tf)
359 : tmp<volScalarField>(hsf1)
364 ? compositionPtr2->Hs(speciei2, thermo2.p(), Tf)
365 : tmp<volScalarField>(hsf2)
369 tmp<volScalarField> tYi1, tYi2;
370 if (this->residualY_ > 0)
374 ?
max(compositionPtr1->Y(speciei1), this->residualY_)
378 ?
max(compositionPtr2->Y(speciei2), this->residualY_)
385 case latentHeatScheme::symmetric:
387 *eqns[phase1.name()] += dmidtf*hsfi1;
388 *eqns[phase2.name()] -= dmidtf*hsfi2;
392 case latentHeatScheme::upwind:
398 ? compositionPtr1->Hs(speciei1, thermo1.p(), thermo1.T())
404 ? compositionPtr2->Hs(speciei2, thermo2.p(), thermo2.T())
408 *eqns[phase1.name()] += dmidtf21*hsfi1 + dmidtf12*hsi1;
409 *eqns[phase2.name()] -= dmidtf12*hsfi2 + dmidtf21*hsi2;
414 if (this->residualY_ > 0)
416 *eqns[phase1.name()] +=
417 fvm::Sp(dmidtf12/tYi1(), he1) - dmidtf12/tYi1()*he1;
419 if (this->residualY_ > 0)
421 *eqns[phase2.name()] -=
422 fvm::Sp(dmidtf21/tYi2(), he2) - dmidtf21/tYi2()*he2;
427 *eqns[phase1.name()] += dmidtf21*
K2 + dmidtf12*
K1;
428 *eqns[phase2.name()] -= dmidtf12*K1 + dmidtf21*
K2;
434 template<
class BasePhaseSystem>
437 const phaseSystem::dmidtfTable& dmidtfs,
438 const phaseSystem::dmdtfTable& Tfs,
440 const latentHeatScheme scheme,
441 phaseSystem::heatTransferTable& eqns
447 const phasePairKey& key = dmidtfIter.key();
448 const phasePair& pair(this->phasePairs_[key]);
452 const phaseModel& phase1 = pair.phase1();
453 const phaseModel& phase2 = pair.phase2();
456 forAllConstIter(HashPtrTable<volScalarField>, *dmidtfIter(), dmidtfJter)
458 const word& specie = dmidtfJter.key();
469 const volScalarField Li(this->Li(pair, specie, dmidtf, Tf, scheme));
470 *eqns[phase1.name()] +=
471 ((1 - weight)*dmidtf12 + weight*dmidtf21)*Li;
472 *eqns[phase2.name()] +=
473 ((1 - weight)*dmidtf21 + weight*dmidtf12)*Li;
479 template<
class BasePhaseSystem>
482 const phaseSystem::dmidtfTable& dmidtfs,
483 const phaseSystem::dmdtfTable& Tfs,
485 const latentHeatScheme scheme,
486 phaseSystem::heatTransferTable& eqns
489 addDmidtHefsWithoutL(dmidtfs, Tfs, scheme, eqns);
490 addDmidtL(dmidtfs, Tfs, weight, scheme, eqns);
496 template<
class BasePhaseSystem>
502 heatTransferPhaseSystem(),
503 BasePhaseSystem(mesh),
504 residualY_(this->template lookupOrDefault<scalar>(
"residualY", -1))
510 template<
class BasePhaseSystem>
517 template<
class BasePhaseSystem>
521 const phasePair& pair,
524 const latentHeatScheme scheme
527 const rhoThermo& thermo1 = pair.phase1().thermo();
528 const rhoThermo& thermo2 = pair.phase2().thermo();
536 case latentHeatScheme::symmetric:
540 case latentHeatScheme::upwind:
547 neg0(dmdtf)*haf2 +
pos(dmdtf)*ha2
548 -
pos0(dmdtf)*haf1 -
neg(dmdtf)*ha1;
552 return tmp<volScalarField>(
nullptr);
556 template<
class BasePhaseSystem>
560 const phasePair& pair,
564 const latentHeatScheme scheme
567 const rhoThermo& thermo1 = pair.phase1().thermo();
568 const rhoThermo& thermo2 = pair.phase2().thermo();
576 case latentHeatScheme::symmetric:
580 case latentHeatScheme::upwind:
590 neg0(dmdtf)*haf2 +
pos(dmdtf)*ha2
591 -
pos0(dmdtf)*haf1 -
neg(dmdtf)*ha1;
595 return tmp<scalarField>(
nullptr);
599 template<
class BasePhaseSystem>
603 const phasePair& pair,
607 const latentHeatScheme scheme
610 const rhoThermo& thermo1 = pair.phase1().thermo();
611 const rhoThermo& thermo2 = pair.phase2().thermo();
612 const basicSpecieMixture* compositionPtr1 =
613 isA<rhoReactionThermo>(
thermo1)
614 ? &refCast<const rhoReactionThermo>(thermo1).composition()
615 :
static_cast<const basicSpecieMixture*
>(
nullptr);
616 const basicSpecieMixture* compositionPtr2 =
617 isA<rhoReactionThermo>(
thermo2)
618 ? &refCast<const rhoReactionThermo>(thermo2).composition()
619 :
static_cast<const basicSpecieMixture*
>(
nullptr);
620 const label speciei1 =
621 compositionPtr1 ? compositionPtr1->species()[specie] : -1;
622 const label speciei2 =
623 compositionPtr2 ? compositionPtr2->species()[specie] : -1;
629 ? compositionPtr1->Ha(speciei1, thermo1.p(), Tf)
630 : thermo1.ha(thermo1.p(), Tf)
635 ? compositionPtr2->Ha(speciei2, thermo2.p(), Tf)
636 : thermo2.ha(thermo1.p(), Tf)
641 case latentHeatScheme::symmetric:
643 return hafi2 - hafi1;
645 case latentHeatScheme::upwind:
651 ? compositionPtr1->Ha(speciei1, thermo1.p(), thermo1.T())
657 ? compositionPtr2->Ha(speciei2, thermo2.p(), thermo2.T())
662 neg0(dmdtf)*hafi2 +
pos(dmdtf)*hai2
663 -
pos0(dmdtf)*hafi1 -
neg(dmdtf)*hai1;
667 return tmp<volScalarField>(
nullptr);
671 template<
class BasePhaseSystem>
675 const phasePair& pair,
680 const latentHeatScheme scheme
683 const rhoThermo& thermo1 = pair.phase1().thermo();
684 const rhoThermo& thermo2 = pair.phase2().thermo();
685 const basicSpecieMixture* compositionPtr1 =
686 isA<rhoReactionThermo>(
thermo1)
687 ? &refCast<const rhoReactionThermo>(thermo1).composition()
688 :
static_cast<const basicSpecieMixture*
>(
nullptr);
689 const basicSpecieMixture* compositionPtr2 =
690 isA<rhoReactionThermo>(
thermo2)
691 ? &refCast<const rhoReactionThermo>(thermo2).composition()
692 :
static_cast<const basicSpecieMixture*
>(
nullptr);
693 const label speciei1 =
694 compositionPtr1 ? compositionPtr1->species()[specie] : -1;
695 const label speciei2 =
696 compositionPtr2 ? compositionPtr2->species()[specie] : -1;
705 ? compositionPtr1->Ha(speciei1, p1, Tf)
706 : thermo1.ha(Tf, cells)
711 ? compositionPtr2->Ha(speciei2, p2, Tf)
712 : thermo2.ha(Tf, cells)
717 case latentHeatScheme::symmetric:
719 return hafi2 - hafi1;
721 case latentHeatScheme::upwind:
730 ? compositionPtr1->Ha(speciei1, p1, T1)
731 : thermo1.ha(T1, cells)
736 ? compositionPtr2->Ha(speciei2, p2, T2)
737 : thermo2.ha(T2, cells)
741 neg0(dmdtf)*hafi2 +
pos(dmdtf)*hai2
742 -
pos0(dmdtf)*hafi1 -
neg(dmdtf)*hai1;
746 return tmp<scalarField>(
nullptr);
750 template<
class BasePhaseSystem>
virtual tmp< volScalarField > L(const phasePair &pair, const volScalarField &dmdtf, const volScalarField &Tf, const latentHeatScheme scheme) const
Return the latent heat for a given pair, mass transfer rate (used.
intWM_LABEL_SIZE_t label
A label is an int32_t or int64_t as specified by the pre-processor macro WM_LABEL_SIZE.
dimensioned< Type > max(const dimensioned< Type > &, const dimensioned< 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,.
static int compare(const Pair< word > &a, const Pair< word > &b)
Compare Pairs.
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))
autoPtr< BasicCompressibleMomentumTransportModel > New(const volScalarField &rho, const volVectorField &U, const surfaceScalarField &phi, const typename BasicCompressibleMomentumTransportModel::transportModel &transport)
Field< scalar > scalarField
Specialisation of Field<T> for scalar.
virtual ~HeatTransferPhaseSystem()
Destructor.
dimensionedScalar neg0(const dimensionedScalar &ds)
forAllConstIter(PtrDictionary< phaseModel >, mixture.phases(), phase)
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,.
virtual tmp< volScalarField > Li(const phasePair &pair, const word &member, const volScalarField &dmidtf, const volScalarField &Tf, const latentHeatScheme scheme) const
Return the latent heat for a given pair, specie, mass transfer rate.
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)