research/cantera
1
0
Fork 0
Code Issues Pull requests Projects Releases Packages Wiki Activity Actions

Initial attempt at exposing the phase stability routine to the

Browse source 
outside without actually doing an equilibrium solution.
This commit is contained in:
Harry Moffat 2010-12-16 22:44:29 +00:00
parent 32816b1fa7
commit fc7b7e8dc6
5 changed files with 52 additions and 26 deletions
Show stats Download patch file Download diff file Expand all files Collapse all files

28
Cantera/src/equil/vcs_MultiPhaseEquil.cpp
View file

@ -1559,18 +1559,14 @@ namespace VCSnonideal {
* states.
*/
m_mix->uploadMoleFractionsFromPhases();
int kGlob = 0;
for (int ip = 0; ip < m_vprob->NPhase; ip++) {
double phaseMole = 0.0;
Cantera::ThermoPhase &tref = m_mix->phase(ip);
int nspPhase = tref.nSpecies();
for (int k = 0; k < nspPhase; k++, kGlob++) {
phaseMole += m_vprob->w[kGlob];
}
//phaseMole *= 1.0E-3;
m_mix->setPhaseMoles(ip, phaseMole);
}
for (int i = 0; i < m_vprob->nspecies; i++) {
plogf("%d %15.3e\n", m_vprob->m_gibbsSpecies[i]);
}
m_mix->getChemPotentials(DATA_PTR(m_vprob->m_gibbsSpecies));
for (int i = 0; i < m_vprob->nspecies; i++) {
plogf("%d %15.3e\n", m_vprob->m_gibbsSpecies[i]);
}
double te = tickTock.secondsWC();
if (printLvl > 0) {
plogf("\n Results from vcs_PS:\n");
@ -1615,11 +1611,11 @@ namespace VCSnonideal {
if (m_vprob->w[i] <= 0.0) {
int iph = m_vprob->PhaseID[i];
vcs_VolPhase *VPhase = m_vprob->VPhaseList[iph];
if (VPhase->nSpecies() > 1) {
plogf(" -1.000e+300\n");
} else {
//if (VPhase->nSpecies() > 1) {
// plogf(" -1.000e+300\n");
//} else {
plogf("%15.3e\n", m_vprob->m_gibbsSpecies[i]);
}
//}
} else {
plogf("%15.3e\n", m_vprob->m_gibbsSpecies[i]);
}

2
Cantera/src/equil/vcs_solve.cpp
View file

@ -896,7 +896,7 @@ namespace VCSnonideal {
retn = VCS_PUB_BAD;
}
if (vPhase->PhaseName == pub_phase_ptr->PhaseName) {
if (vPhase->PhaseName != pub_phase_ptr->PhaseName) {
plogf("%sPhaseName value have changed:%s %s\n", yo.c_str(),
vPhase->PhaseName.c_str(),
pub_phase_ptr->PhaseName.c_str());

6
Cantera/src/equil/vcs_solve.h
View file

@ -1257,16 +1257,16 @@ private:
* the mole numbers of the component species. Therefore the following
* approximation is valid for a small component of an ideal phase:
*
* 0 = m_deltaGrxn_old(I) + log(molNum_new(I)/molNum_old(I))
* 0 = m_deltaGRxn_old(I) + log(molNum_new(I)/molNum_old(I))
*
* m_deltaGrxn_old contains the contribution from
* m_deltaGRxn_old contains the contribution from
*
* m_feSpecies_old(I) =
* m_SSfeSpecies(I) +
* log(ActCoeff[i] * molNum_old(I) / m_tPhaseMoles_old(iph))
* Thus,
*
* molNum_new(I)= molNum_old(I) * EXP(-m_deltaGrxn_old(I))
* molNum_new(I)= molNum_old(I) * EXP(-m_deltaGRxn_old(I))
*
* Most of this section is mainly restricting the update to reasonable
* values.

6
Cantera/src/equil/vcs_solve_TP.cpp
View file

@ -2015,16 +2015,16 @@ namespace VCSnonideal {
* the mole numbers of the component species. Therefore the following
* approximation is valid for a small component of an ideal phase:
*
* 0 = m_deltaGrxn_old(I) + log(molNum_new(I)/molNum_old(I))
* 0 = m_deltaGRxn_old(I) + log(molNum_new(I)/molNum_old(I))
*
* m_deltaGrxn_old contains the contribution from
* m_deltaGRxn_old contains the contribution from
*
* m_feSpecies_old(I) =
* m_SSfeSpecies(I) +
* log(ActCoeff[i] * molNum_old(I) / m_tPhaseMoles_old(iph))
* Thus,
*
* molNum_new(I)= molNum_old(I) * EXP(-m_deltaGrxn_old(I))
* molNum_new(I)= molNum_old(I) * EXP(-m_deltaGRxn_old(I))
*
* Most of this section is mainly restricting the update to reasonable
* values.

36
Cantera/src/equil/vcs_solve_phaseStability.cpp
View file

@ -41,7 +41,16 @@ namespace VCSnonideal {
int ifunc = 0;
int iStab = 0;
/*
* This function is called to create the private data
* using the public data.
*/
int nspecies0 = vprob->nspecies + 10;
int nelements0 = vprob->ne;
int nphase0 = vprob->NPhase;
vcs_initSizes(nspecies0, nelements0, nphase0);
if (ifunc < 0 || ifunc > 2) {
plogf("vcs: Unrecognized value of ifunc, %d: bailing!\n",
@ -131,8 +140,7 @@ namespace VCSnonideal {
* (all information concerning Temperature and Pressure has already
* been derived. The free energies are now in dimensionless form.)
*/
double funcVal;
iStab = vcs_solve_phaseStability(iphase, ifunc, funcVal, printLvl);
iStab = vcs_solve_phaseStability(iphase, ifunc, feStable, printLvl);
/*
@ -140,6 +148,23 @@ namespace VCSnonideal {
* the reverse of vcs_nondim to add back units.
*/
vcs_redim_TP();
/*
vcs_VolPhase *Vphase = m_VolPhaseList[iphase];
std::vector<double> mfPop = Vphase->moleFractions();
int nsp = Vphase->nSpecies();
vcs_VolPhase *VPphase = vprob->VPhaseList[iphase];
int kstart = Vphase->spGlobalIndexVCS(0);
for (int k = 0; k < nsp; k++) {
vprob->mf[kstart + k] = mfPop[k];
}
VPphase->setMoleFractionsState(Vphase->totalMoles(),
VCS_DATA_PTR(Vphase->moleFractions()),
VCS_STATECALC_TMP);
*/
vcs_prob_update(vprob);
/*
* Return the convergence success flag.
*/
@ -197,6 +222,11 @@ namespace VCSnonideal {
retn = vcs_basopt(FALSE, VCS_DATA_PTR(aw), VCS_DATA_PTR(sa),
VCS_DATA_PTR(sm), VCS_DATA_PTR(ss),
test, &usedZeroedSpecies);
vcs_dfe(VCS_STATECALC_OLD, 0, 0, m_numSpeciesRdc);
vcs_deltag(0, true, VCS_STATECALC_OLD);
phasePopPhaseIDs.clear();
iphasePop = vcs_popPhaseID(phasePopPhaseIDs);
funcVal = vcs_phaseStabilityTest(iph);