28 #include "heatTransferModel.H" 34 template<
class BasePhaseSystem>
37 const phaseSystem::dmdtfTable& dmdtfs,
38 phaseSystem::heatTransferTable& eqns
43 const phasePairKey& key = dmdtfIter.key();
44 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();
59 *eqns[phase1.name()] += dmdtf21*he2 +
fvm::Sp(dmdtf12, he1);
60 *eqns[phase2.name()] -= dmdtf12*he1 +
fvm::Sp(dmdtf21, he2);
63 *eqns[phase1.name()] += dmdtf21*
K2 + dmdtf12*
K1;
64 *eqns[phase2.name()] -= dmdtf12*K1 + dmdtf21*
K2;
69 template<
class BasePhaseSystem>
72 const phaseSystem::dmidtfTable& dmidtfs,
73 phaseSystem::heatTransferTable& eqns
78 const phasePairKey& key = dmidtfIter.key();
79 const phasePair& pair(this->phasePairs_[key]);
81 const phaseModel& phase1 = pair.phase1();
82 const phaseModel& phase2 = pair.phase2();
83 const rhoThermo& thermo1 = phase1.thermo();
84 const rhoThermo& thermo2 = phase2.thermo();
90 forAllConstIter(HashPtrTable<volScalarField>, *dmidtfIter(), dmidtfJter)
92 const word& member = dmidtfJter.key();
103 if (isA<rhoReactionThermo>(thermo1))
105 const basicSpecieMixture& composition1 =
110 composition1.species()[member],
116 if (isA<rhoReactionThermo>(thermo2))
118 const basicSpecieMixture& composition2 =
123 composition2.species()[member],
130 *eqns[phase1.name()] += dmidtf21*hei2 +
fvm::Sp(dmidtf12, he1);
131 *eqns[phase2.name()] -= dmidtf12*hei1 +
fvm::Sp(dmidtf21, he2);
134 *eqns[phase1.name()] += dmidtf21*
K2 + dmidtf12*
K1;
135 *eqns[phase2.name()] -= dmidtf12*K1 + dmidtf21*
K2;
143 template<
class BasePhaseSystem>
150 BasePhaseSystem(mesh)
152 this->generatePairsAndSubModels
163 template<
class BasePhaseSystem>
171 template<
class BasePhaseSystem>
176 autoPtr<phaseSystem::heatTransferTable> eqnsPtr
185 const phaseModel& phase = this->phaseModels_[
phasei];
197 heatTransferModelTable,
199 heatTransferModelIter
204 const phasePair& pair(this->phasePairs_[heatTransferModelIter.key()]);
208 const phaseModel& phase = iter();
209 const phaseModel& otherPhase = iter.otherPhase();
214 *eqns[phase.name()] +=
215 K*(otherPhase.thermo().T() - phase.thermo().T() +
he/Cpv)
224 template<
class BasePhaseSystem>
OneResistanceHeatTransferPhaseSystem(const fvMesh &)
Construct from fvMesh.
void addDmidtHef(const phaseSystem::dmidtfTable &dmidtfs, phaseSystem::heatTransferTable &eqns) const
Add energy transfer terms which result from specie mass transfers.
fvMatrix< scalar > fvScalarMatrix
#define forAll(list, i)
Loop across all elements in list.
basicSpecieMixture & composition
static int compare(const Pair< word > &a, const Pair< word > &b)
Compare Pairs.
HashPtrTable< fvScalarMatrix > heatTransferTable
tmp< fvMatrix< Type > > Sp(const volScalarField::Internal &, const GeometricField< Type, fvPatchField, volMesh > &)
CGAL::Exact_predicates_exact_constructions_kernel K
dimensionedScalar posPart(const dimensionedScalar &ds)
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
mixture thermo2().correctRho(psi2 *(p_rgh - p_rgh_0))
mixture thermo1().correctRho(psi1 *(p_rgh - p_rgh_0))
forAllConstIter(PtrDictionary< phaseModel >, mixture.phases(), phase)
const dimensionSet dimEnergy
virtual autoPtr< phaseSystem::heatTransferTable > heatTransfer() const
Return the heat transfer matrices.
const dimensionSet dimTime(0, 0, 1, 0, 0, 0, 0)
An auto-pointer similar to the STL auto_ptr but with automatic casting to a reference to the type and...
Calculate the matrix for implicit and explicit sources.
virtual bool read()
Read base phaseProperties dictionary.
dimensionedScalar negPart(const dimensionedScalar &ds)
virtual ~OneResistanceHeatTransferPhaseSystem()
Destructor.