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

Variable name changing. no change in algorithm

Browse source
This commit is contained in:
Harry Moffat 2008-05-18 15:37:20 +00:00
parent 7a99ca115f
commit 61a3d2b2c4
11 changed files with 201 additions and 177 deletions
Show stats Download patch file Download diff file Expand all files Collapse all files

14
Cantera/src/equil/vcs_elem.cpp
View file

@ -209,7 +209,7 @@ namespace VCSnonideal {
double l2before = 0.0, l2after = 0.0;
std::vector<double> ga_save(m_numElemConstraints, 0.0);
vcs_dcopy(VCS_DATA_PTR(ga_save), VCS_DATA_PTR(m_elemAbundances), m_numElemConstraints);
if (vcs_debug_print_lvl >= 2) {
if (m_debug_print_lvl >= 2) {
plogf(" --- vcsc_elcorr: Element abundances correction routine");
if (m_numElemConstraints != m_numComponents) {
plogf(" (m_numComponents != m_numElemConstraints)");
@ -313,8 +313,9 @@ namespace VCSnonideal {
if (m_molNumSpecies_old[kspec] > maxPermissible) {
#ifdef DEBUG_MODE
if (vcs_debug_print_lvl >= 3) {
plogf(" --- vcs_elcorr: Reduced species %s from %g to %g due to %s max bounds constraint\n",
if (m_debug_print_lvl >= 3) {
plogf(" --- vcs_elcorr: Reduced species %s from %g to %g "
"due to %s max bounds constraint\n",
m_speciesName[kspec].c_str(), m_molNumSpecies_old[kspec],
maxPermissible, m_elementName[i].c_str());
}
@ -329,8 +330,9 @@ namespace VCSnonideal {
m_rxnStatus[kspec] = VCS_SPECIES_ZEROEDMS;
}
#ifdef DEBUG_MODE
if (vcs_debug_print_lvl >= 2) {
plogf(" --- vcs_elcorr: Zeroed species %s and changed status to %d due to max bounds constraint\n",
if (m_debug_print_lvl >= 2) {
plogf(" --- vcs_elcorr: Zeroed species %s and changed "
"status to %d due to max bounds constraint\n",
m_speciesName[kspec].c_str(), m_rxnStatus[irxn]);
}
#endif
@ -581,7 +583,7 @@ namespace VCSnonideal {
l2after += SQUARE(m_elemAbundances[i] - m_elemAbundancesGoal[i]);
}
l2after = sqrt(l2after/m_numElemConstraints);
if (vcs_debug_print_lvl >= 2) {
if (m_debug_print_lvl >= 2) {
plogf(" --- Elem_Abund: Correct Initial "
" Final\n");
for (i = 0; i < m_numElemConstraints; ++i) {

4
Cantera/src/equil/vcs_elem_rearrange.cpp
View file

@ -66,7 +66,7 @@ namespace VCSnonideal {
int ncomponents = m_numComponents;
double test = -1.0E10;
#ifdef DEBUG_MODE
if (vcs_debug_print_lvl >= 2) {
if (m_debug_print_lvl >= 2) {
plogf(" "); for(i=0; i<77; i++) plogf("-"); plogf("\n");
plogf(" --- Subroutine elem_rearrange() called to ");
plogf("check stoich. coefficent matrix\n");
@ -187,7 +187,7 @@ namespace VCSnonideal {
/* ****************************************** */
if (jr != k) {
#ifdef DEBUG_MODE
if (vcs_debug_print_lvl >= 2) {
if (m_debug_print_lvl >= 2) {
plogf(" --- "); plogf("%-2.2s", (m_elementName[k]).c_str());
plogf("(%9.2g) replaces ", m_elemAbundancesGoal[k]);
plogf("%-2.2s", (m_elementName[jr]).c_str());

32
Cantera/src/equil/vcs_inest.cpp
View file

@ -53,8 +53,8 @@ namespace VCSnonideal {
double *molNum = VCS_DATA_PTR(m_molNumSpecies_old);
double TMolesMultiphase;
double *xtphMax = VCS_DATA_PTR(TmpPhase);
double *xtphMin = VCS_DATA_PTR(TmpPhase2);
double *xtphMax = VCS_DATA_PTR(m_TmpPhase);
double *xtphMin = VCS_DATA_PTR(m_TmpPhase2);
ikl = 0;
lt = 0;
@ -98,7 +98,7 @@ namespace VCSnonideal {
#endif
#ifdef DEBUG_MODE
if (vcs_debug_print_lvl >= 2) {
if (m_debug_print_lvl >= 2) {
plogf("%s Mole Numbers returned from linear programming (vcs_inest initial guess):\n",
pprefix);
plogf("%s SPECIES MOLE_NUMBER -SS_ChemPotential\n", pprefix);
@ -209,7 +209,7 @@ namespace VCSnonideal {
}
vcs_deltag(0, true);
#ifdef DEBUG_MODE
if (vcs_debug_print_lvl >= 2) {
if (m_debug_print_lvl >= 2) {
for (kspec = 0; kspec < nspecies; ++kspec) {
plogf("%s", pprefix); plogf("%-12.12s", m_speciesName[kspec].c_str());
if (kspec < m_numComponents)
@ -263,7 +263,7 @@ namespace VCSnonideal {
}
}
#ifdef DEBUG_MODE
if (vcs_debug_print_lvl >= 2) {
if (m_debug_print_lvl >= 2) {
for (kspec = 0; kspec < nspecies; ++kspec) {
if (m_speciesUnknownType[kspec] != VCS_SPECIES_TYPE_INTERFACIALVOLTAGE) {
plogf("%sdirection (", pprefix); plogf("%-12.12s", m_speciesName[kspec].c_str());
@ -371,7 +371,7 @@ namespace VCSnonideal {
finished:
;
#ifdef DEBUG_MODE
if (vcs_debug_print_lvl >= 2) {
if (m_debug_print_lvl >= 2) {
plogf("%s Final Mole Numbers produced by inest:\n",
pprefix);
plogf("%s SPECIES MOLE_NUMBER\n", pprefix);
@ -406,18 +406,20 @@ namespace VCSnonideal {
vcs_elab();
if (vcs_elabcheck(0)) {
#ifdef DEBUG_MODE
if (vcs_debug_print_lvl >= 2) {
if (m_debug_print_lvl >= 2) {
plogf("%s Initial guess passed element abundances on input\n", pprefix);
plogf("%s m_doEstimateEquil = 1 so will use the input mole numbers as estimates", pprefix);
plogf("%s m_doEstimateEquil = 1 so will use the input mole "
"numbers as estimates", pprefix);
plogendl();
}
#endif
return retn;
#ifdef DEBUG_MODE
} else {
if (vcs_debug_print_lvl >= 2) {
if (m_debug_print_lvl >= 2) {
plogf("%s Initial guess failed element abundances on input\n", pprefix);
plogf("%s m_doEstimateEquil = 1 so will discard input mole numbers and find our own estimate", pprefix);
plogf("%s m_doEstimateEquil = 1 so will discard input "
"mole numbers and find our own estimate", pprefix);
plogendl();
}
#endif
@ -440,7 +442,7 @@ namespace VCSnonideal {
* Go get the estimate of the solution
*/
#ifdef DEBUG_MODE
if (vcs_debug_print_lvl >= 2) {
if (m_debug_print_lvl >= 2) {
plogf("%sGo find an initial estimate for the equilibrium problem",
pprefix);
plogendl();
@ -466,7 +468,7 @@ namespace VCSnonideal {
int rangeCheck = vcs_elabcheck(1);
if (!vcs_elabcheck(0)) {
#ifdef DEBUG_MODE
if (vcs_debug_print_lvl >= 2) {
if (m_debug_print_lvl >= 2) {
plogf("%sInitial guess failed element abundances\n", pprefix);
plogf("%sCall vcs_elcorr to attempt fix", pprefix);
plogendl();
@ -483,7 +485,7 @@ namespace VCSnonideal {
retn = -1;
} else {
#ifdef DEBUG_MODE
if (vcs_debug_print_lvl >= 2) {
if (m_debug_print_lvl >= 2) {
if (rangeCheck) {
plogf("%sInitial guess now satisfies element abundances", pprefix);
plogendl();
@ -500,7 +502,7 @@ namespace VCSnonideal {
}
else {
#ifdef DEBUG_MODE
if (vcs_debug_print_lvl >= 2) {
if (m_debug_print_lvl >= 2) {
if (rangeCheck) {
plogf("%sInitial guess satisfies element abundances", pprefix);
plogendl();
@ -516,7 +518,7 @@ namespace VCSnonideal {
}
#ifdef DEBUG_MODE
if (vcs_debug_print_lvl >= 2) {
if (m_debug_print_lvl >= 2) {
plogf("%sTotal Dimensionless Gibbs Free Energy = %15.7E", pprefix,
vcs_Total_Gibbs(VCS_DATA_PTR(m_molNumSpecies_old), VCS_DATA_PTR(m_feSpecies_curr),
VCS_DATA_PTR(m_tPhaseMoles_old)));

12
Cantera/src/equil/vcs_nondim.cpp
View file

@ -106,8 +106,8 @@ namespace VCSnonideal {
void VCS_SOLVE::vcs_nondim_TP() {
int i;
double tf;
if (UnitsState == VCS_DIMENSIONAL_G) {
UnitsState = VCS_NONDIMENSIONAL_G;
if (m_unitsState == VCS_DIMENSIONAL_G) {
m_unitsState = VCS_NONDIMENSIONAL_G;
tf = 1.0 / vcs_nondimMult_TP(m_VCS_UnitsFormat, m_temperature);
for (i = 0; i < m_numSpeciesTot; ++i) {
/*
@ -122,7 +122,7 @@ namespace VCSnonideal {
m_feSpecies_old[i] *= tf;
}
Faraday_dim = vcs_nondim_Farad(m_VCS_UnitsFormat, m_temperature);
m_Faraday_dim = vcs_nondim_Farad(m_VCS_UnitsFormat, m_temperature);
if (m_VCS_UnitsFormat == VCS_UNITS_MKS) {
for (i = 0; i < m_numSpeciesTot; ++i) {
if (m_speciesUnknownType[i] != VCS_SPECIES_TYPE_INTERFACIALVOLTAGE) {
@ -152,8 +152,8 @@ namespace VCSnonideal {
{
int i;
double tf;
if (UnitsState != VCS_DIMENSIONAL_G) {
UnitsState = VCS_DIMENSIONAL_G;
if (m_unitsState != VCS_DIMENSIONAL_G) {
m_unitsState = VCS_DIMENSIONAL_G;
tf = vcs_nondimMult_TP(m_VCS_UnitsFormat, m_temperature);
for (i = 0; i < m_numSpeciesTot; ++i) {
/*
@ -166,7 +166,7 @@ namespace VCSnonideal {
m_deltaGRxn_old[i] *= tf;
m_feSpecies_old[i] *= tf;
}
Faraday_dim *= tf;
m_Faraday_dim *= tf;
}
if (m_VCS_UnitsFormat == VCS_UNITS_MKS) {
for (i = 0; i < m_numSpeciesTot; ++i) {

5
Cantera/src/equil/vcs_prep.cpp
View file

@ -119,9 +119,8 @@ namespace VCSnonideal {
double *aw, *sa, *sm, *ss;
bool modifiedSoln = false;
#ifdef DEBUG_MODE
vcs_debug_print_lvl = printLvl;
#endif
m_debug_print_lvl = printLvl;
/*
* Calculate the Single Species status of phases

16
Cantera/src/equil/vcs_report.cpp
View file

@ -51,7 +51,7 @@ namespace VCSnonideal {
int nspecies = m_numSpeciesTot;
double g;
char originalUnitsState = UnitsState;
char originalUnitsState = m_unitsState;
std::vector<int> sortindex(nspecies,0);
@ -83,7 +83,7 @@ namespace VCSnonideal {
* -> For the printouts from this routine, we will use nondimensional
* representations. This may be expanded in the future.
*/
if (UnitsState == VCS_DIMENSIONAL_G) {
if (m_unitsState == VCS_DIMENSIONAL_G) {
vcs_nondim_TP();
}
@ -321,7 +321,7 @@ namespace VCSnonideal {
plogf("%14.7E ", log(m_actCoeffSpecies_old[l]));
double tpmoles = m_tPhaseMoles_old[pid];
double phi = m_phasePhi[pid];
double eContrib = phi * m_chargeSpecies[l] * Faraday_dim;
double eContrib = phi * m_chargeSpecies[l] * m_Faraday_dim;
double lx = 0.0;
if (m_speciesUnknownType[l] == VCS_SPECIES_TYPE_INTERFACIALVOLTAGE) {
lx = 0.0;
@ -330,20 +330,20 @@ namespace VCSnonideal {
lx = log(m_molNumSpecies_old[l]) - log(tpmoles);
} else {
lx = m_feSpecies_curr[l] - m_SSfeSpecies[l]
- log(m_actCoeffSpecies_old[l]) + SpecLnMnaught[l];
- log(m_actCoeffSpecies_old[l]) + m_lnMnaughtSpecies[l];
}
}
plogf("%14.7E |", lx);
plogf("%14.7E | ", eContrib);
double tmp = m_SSfeSpecies[l] + log(m_actCoeffSpecies_old[l])
+ lx - SpecLnMnaught[l] + eContrib;
+ lx - m_lnMnaughtSpecies[l] + eContrib;
if (fabs(m_feSpecies_curr[l] - tmp) > 1.0E-8) {
plogf("\n\t\twe have a problem - doesn't add up\n");
exit(-1);
}
plogf(" %12.4E |", m_feSpecies_curr[l]);
if (SpecLnMnaught[l] != 0.0) {
plogf(" (%14.7E)", - SpecLnMnaught[l]);
if (m_lnMnaughtSpecies[l] != 0.0) {
plogf(" (%14.7E)", - m_lnMnaughtSpecies[l]);
}
plogf("\n");
}
@ -372,7 +372,7 @@ namespace VCSnonideal {
* Set the Units state of the system back to where it was when we
* entered the program.
*/
if (originalUnitsState != UnitsState) {
if (originalUnitsState != m_unitsState) {
if (originalUnitsState == VCS_DIMENSIONAL_G ) vcs_redim_TP();
else vcs_nondim_TP();
}

4
Cantera/src/equil/vcs_rxnadj.cpp
View file

@ -432,7 +432,7 @@ namespace VCSnonideal {
if (dx_orig > 0.0) {
dx = 0.0;
#ifdef DEBUG_MODE
if (vcs_debug_print_lvl >= 2) {
if (m_debug_print_lvl >= 2) {
//plogf(" --- %s :Warning possible error dx>0 dg > 0\n", SpName[kspec]);
}
sprintf(ANOTE,"Rxn reduced to zero step size in line search: dx>0 dg > 0");
@ -443,7 +443,7 @@ namespace VCSnonideal {
if (dx_orig < 0.0) {
dx = 0.0;
#ifdef DEBUG_MODE
if (vcs_debug_print_lvl >= 2) {
if (m_debug_print_lvl >= 2) {
//plogf(" --- %s :Warning possible error dx<0 dg < 0\n", SpName[kspec]);
}
sprintf(ANOTE,"Rxn reduced to zero step size in line search: dx<0 dg < 0");

14
Cantera/src/equil/vcs_setMolesLinProg.cpp
View file

@ -62,7 +62,7 @@ int VCS_SOLVE::vcs_setMolesLinProg() {
#ifdef DEBUG_MODE
std::string pprefix(" --- seMolesLinProg ");
if (vcs_debug_print_lvl >= 2) {
if (m_debug_print_lvl >= 2) {
plogf(" --- call setInitialMoles\n");
}
#endif
@ -97,7 +97,7 @@ int VCS_SOLVE::vcs_setMolesLinProg() {
}
#ifdef DEBUG_MODE
if (vcs_debug_print_lvl >= 2) {
if (m_debug_print_lvl >= 2) {
printProgress(m_speciesName, m_molNumSpecies_old, m_SSfeSpecies);
}
#endif
@ -106,7 +106,7 @@ int VCS_SOLVE::vcs_setMolesLinProg() {
if (!vcs_elabcheck(0)) {
#ifdef DEBUG_MODE
if (vcs_debug_print_lvl >= 2) {
if (m_debug_print_lvl >= 2) {
plogf("%s Mole numbers failing element abundances\n", pprefix.c_str());
plogf("%sCall vcs_elcorr to attempt fix\n", pprefix.c_str());
}
@ -131,7 +131,7 @@ int VCS_SOLVE::vcs_setMolesLinProg() {
if (retn != VCS_SUCCESS) return retn;
#ifdef DEBUG_MODE
if (vcs_debug_print_lvl >= 2) {
if (m_debug_print_lvl >= 2) {
plogf("iteration %d\n", iter);
}
#endif
@ -170,7 +170,7 @@ int VCS_SOLVE::vcs_setMolesLinProg() {
if (!redo) {
if (delta_xi < 1.0e-10 && (m_molNumSpecies_old[ik] >= 1.0E-10)) {
#ifdef DEBUG_MODE
if (vcs_debug_print_lvl >= 2) {
if (m_debug_print_lvl >= 2) {
plogf(" --- Component too small: %s\n", m_speciesName[jcomp].c_str());
}
#endif
@ -207,14 +207,14 @@ int VCS_SOLVE::vcs_setMolesLinProg() {
// Update the phase objects with the contents of the m_molNumSpecies_old vector
// vcs_updateVP(0);
#ifdef DEBUG_MODE
if (vcs_debug_print_lvl >= 2) {
if (m_debug_print_lvl >= 2) {
printProgress(m_speciesName, m_molNumSpecies_old, m_SSfeSpecies);
}
#endif
}
#ifdef DEBUG_MODE
if (vcs_debug_print_lvl == 1) {
if (m_debug_print_lvl == 1) {
printProgress(m_speciesName, m_molNumSpecies_old, m_SSfeSpecies);
plogf(" --- setInitialMoles end\n");
}

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

@ -49,12 +49,12 @@ namespace VCSnonideal {
m_tolmin(0.0),
m_tolmaj2(0.0),
m_tolmin2(0.0),
UnitsState(VCS_DIMENSIONAL_G),
m_unitsState(VCS_DIMENSIONAL_G),
m_useActCoeffJac(0),
m_totalVol(0.0),
Faraday_dim(1.602e-19 * 6.022136736e26),
m_Faraday_dim(1.602e-19 * 6.022136736e26),
m_VCount(0),
vcs_debug_print_lvl(0),
m_debug_print_lvl(0),
m_timing_print_lvl(1),
m_VCS_UnitsFormat(VCS_UNITS_UNITLESS)
{
@ -134,8 +134,8 @@ namespace VCSnonideal {
m_tPhaseMoles_old.resize(nphase0, 0.0);
m_tPhaseMoles_new.resize(nphase0, 0.0);
m_deltaPhaseMoles.resize(nphase0, 0.0);
TmpPhase.resize(nphase0, 0.0);
TmpPhase2.resize(nphase0, 0.0);
m_TmpPhase.resize(nphase0, 0.0);
m_TmpPhase2.resize(nphase0, 0.0);
m_formulaMatrix.resize(nelements, nspecies0);
@ -178,12 +178,12 @@ namespace VCSnonideal {
* Malloc space for activity coefficients for all species
* -> Set it equal to one.
*/
SpecActConvention.resize(nspecies0, 0);
PhaseActConvention.resize(nphase0, 0);
SpecLnMnaught.resize(nspecies0, 0.0);
m_actConventionSpecies.resize(nspecies0, 0);
m_phaseActConvention.resize(nphase0, 0);
m_lnMnaughtSpecies.resize(nspecies0, 0.0);
m_actCoeffSpecies_new.resize(nspecies0, 1.0);
m_actCoeffSpecies_old.resize(nspecies0, 1.0);
CurrPhAC.resize(nphase0, 0);
m_phaseACAreCurrent.resize(nphase0, 0);
m_wtSpecies.resize(nspecies0, 0.0);
m_chargeSpecies.resize(nspecies0, 0.0);
m_speciesThermoList.resize(nspecies0, (VCS_SPECIES_THERMO *)0);
@ -489,7 +489,9 @@ namespace VCSnonideal {
m_numPhases = nph;
#ifdef DEBUG_MODE
vcs_debug_print_lvl = pub->vcs_debug_print_lvl;
m_debug_print_lvl = pub->vcs_debug_print_lvl;
#else
m_debug_print_lvl = MIN(2, pub->vcs_debug_print_lvl);
#endif
/*
@ -714,7 +716,7 @@ namespace VCSnonideal {
*/
for (iph = 0; iph < nph; iph++) {
Vphase = m_VolPhaseList[iph];
PhaseActConvention[iph] = Vphase->ActivityConvention;
m_phaseActConvention[iph] = Vphase->ActivityConvention;
if (Vphase->ActivityConvention != 0) {
/*
* We assume here that species 0 is the solvent.
@ -729,8 +731,8 @@ namespace VCSnonideal {
double mnaught = m_wtSpecies[iSolvent] / 1000.;
for (int k = 1; k < Vphase->NVolSpecies; k++) {
int kspec = Vphase->IndSpecies[k];
SpecActConvention[kspec] = Vphase->ActivityConvention;
SpecLnMnaught[kspec] = log(mnaught);
m_actConventionSpecies[kspec] = Vphase->ActivityConvention;
m_lnMnaughtSpecies[kspec] = log(mnaught);
}
}
}

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

@ -1092,10 +1092,10 @@ public:
std::vector<double> m_tPhaseMoles_new;
//! Temporary vector of length NPhase
std::vector<double> TmpPhase;
std::vector<double> m_TmpPhase;
//! Temporary vector of length NPhase
std::vector<double> TmpPhase2;
std::vector<double> m_TmpPhase2;
//! Change in the total moles in each phase
/*!
@ -1272,7 +1272,7 @@ public:
/*!
*. The default is to have this unitless
*/
char UnitsState;
char m_unitsState;
//! specifies the activity convention of the phase containing the species
/*!
@ -1281,7 +1281,7 @@ public:
* 1 = molality based
* length = number of species
*/
std::vector<int> SpecActConvention;
std::vector<int> m_actConventionSpecies;
//! specifies the activity convention of the phase.
/*!
@ -1289,7 +1289,7 @@ public:
* 1 = molality based
* length = number of phases
*/
std::vector<int> PhaseActConvention;
std::vector<int> m_phaseActConvention;
//! specifies the ln(Mnaught) used to calculate the chemical potentials
/*!
@ -1297,7 +1297,7 @@ public:
* this will be equal to 0.0
* length = number of species
*/
std::vector<double> SpecLnMnaught;
std::vector<double> m_lnMnaughtSpecies;
//! Molar-based Activity Coefficients for Species
/*!
@ -1327,7 +1327,7 @@ public:
//! This boolean indicates whether the activity coefficients for a phase
//! are current.
std::vector<int> CurrPhAC;
std::vector<int> m_phaseACAreCurrent;
//! Molecular weight of each species
/*!
@ -1377,12 +1377,27 @@ public:
/*!
* F / RT (1/volt)
*/
double Faraday_dim;
double m_Faraday_dim;
//! Timing and iteration counters for the vcs object
VCS_COUNTERS *m_VCount;
int vcs_debug_print_lvl;
//! Debug printing lvl
/*!
* Levels correspond to the following guidlines
* - 0 No printing at all
* - 1 Serious warnings or fatal errors get one line
* - 2 one line per eacdh successful vcs package call
* - 3 one line per every successful solve_TP calculation
* - 4 one line for every successful operation -> solve_TP gets a summary report
* - 5 each iteration in solve_TP gets a report with one line per species
* - 6 Each decision in solve_TP gets a line in addition to 4
* - 10 Additionally Hessian matrix is printed out
*
* Levels of printing above 4 are only accessible when DEBUG_MODE is turned on
*/
int m_debug_print_lvl;
//! printing level of timing information
/*!

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

@ -128,7 +128,7 @@ namespace VCSnonideal {
/*
* Set the debug print lvl to the same as the print lvl.
*/
vcs_debug_print_lvl = printDetails;
m_debug_print_lvl = printDetails;
#endif
if (printDetails > 0 && print_lvl == 0) {
print_lvl = 1;
@ -320,7 +320,7 @@ namespace VCSnonideal {
if (m_rxnStatus[irxn] == VCS_SPECIES_MINOR) {
m_rxnStatus[irxn] = VCS_SPECIES_MAJOR;
#ifdef DEBUG_MODE
if (vcs_debug_print_lvl >= 2) {
if (m_debug_print_lvl >= 2) {
plogf(" --- Minor species changed to major: ");
plogf("%-12s", m_speciesName[kspec].c_str());
plogendl();
@ -335,7 +335,7 @@ namespace VCSnonideal {
lec = FALSE;
if (! vcs_elabcheck(0)) {
#ifdef DEBUG_MODE
if (vcs_debug_print_lvl >= 2) {
if (m_debug_print_lvl >= 2) {
plogf(" --- Element Abundance check failed");
plogendl();
}
@ -345,7 +345,7 @@ namespace VCSnonideal {
}
#ifdef DEBUG_MODE
else {
if (vcs_debug_print_lvl >= 2) {
if (m_debug_print_lvl >= 2) {
plogf(" --- Element Abundance check passed");
plogendl();
}
@ -439,7 +439,7 @@ namespace VCSnonideal {
}
} else {
#ifdef DEBUG_MODE
if (vcs_debug_print_lvl >= 2) {
if (m_debug_print_lvl >= 2) {
if (im) {
plogf(" --- vcs_RxnStepSizes not called because all"
"species are minors\n");
@ -490,7 +490,7 @@ namespace VCSnonideal {
*
*/
#ifdef DEBUG_MODE
if (vcs_debug_print_lvl >= 2) {
if (m_debug_print_lvl >= 2) {
plogf(" --- Main Loop Treatment of each non-component species ");
if (iti == 0) plogf("- Full Calculation:\n");
else plogf("- Major Components Calculation:\n");
@ -525,7 +525,7 @@ namespace VCSnonideal {
/********************************************************************/
bool resurrect = true;
#ifdef DEBUG_MODE
if (vcs_debug_print_lvl >= 3) {
if (m_debug_print_lvl >= 3) {
plogf(" --- %s currently zeroed (SpStatus=%-2d):",
m_speciesName[kspec].c_str(), m_rxnStatus[irxn]);
plogf("%3d DG = %11.4E WT = %11.4E W = %11.4E DS = %11.4E\n",
@ -551,7 +551,7 @@ namespace VCSnonideal {
sprintf(ANOTE, "Species stays zeroed even though dg neg:DG = %11.4E, ds zeroed ",
m_deltaGRxn_new[irxn]);
}
//if (vcs_debug_print_lvl >= 2) {
//if (m_debug_print_lvl >= 2) {
//plogf(" --- "); plogf("%-12s", m_speciesName[kspec]);
//plogf("%3d%11.4E%11.4E%11.4E | %s\n",
// m_rxnStatus[irxn], w[kspec], wt[kspec],
@ -584,7 +584,7 @@ namespace VCSnonideal {
if (Vphase->Existence == 0) Vphase->Existence = 1;
--m_numRxnMinorZeroed;
#ifdef DEBUG_MODE
if (vcs_debug_print_lvl >= 2) {
if (m_debug_print_lvl >= 2) {
plogf(" --- Zeroed species changed to major: ");
plogf("%-12s\n", m_speciesName[kspec].c_str());
}
@ -623,7 +623,7 @@ namespace VCSnonideal {
dx = 0.0;
#ifdef DEBUG_MODE
sprintf(ANOTE,"minor species not considered");
if (vcs_debug_print_lvl >= 2) {
if (m_debug_print_lvl >= 2) {
plogf(" --- "); plogf("%-12s", m_speciesName[kspec].c_str());
plogf("%3d%11.4E%11.4E%11.4E | %s",
m_rxnStatus[irxn], m_molNumSpecies_old[kspec], m_molNumSpecies_new[kspec],
@ -659,7 +659,7 @@ namespace VCSnonideal {
/***** MOLE NUMBER */
/*******************************************************************/
#ifdef DEBUG_MODE
if (vcs_debug_print_lvl >= 2) {
if (m_debug_print_lvl >= 2) {
plogf(" --- Delete minor species in multispec phase: %-12s",
m_speciesName[kspec].c_str());
plogendl();
@ -707,7 +707,7 @@ namespace VCSnonideal {
dx = 0.0;
#ifdef DEBUG_MODE
sprintf(ANOTE, "major species is converged");
if (vcs_debug_print_lvl >= 2) {
if (m_debug_print_lvl >= 2) {
plogf(" --- "); plogf("%-12s", m_speciesName[kspec].c_str());
plogf("%3d%11.4E%11.4E%11.4E | %s",
m_rxnStatus[irxn], m_molNumSpecies_old[kspec], m_molNumSpecies_new[kspec],
@ -785,7 +785,7 @@ namespace VCSnonideal {
iph = m_phaseID[kspec];
if (m_molNumSpecies_new[kspec] < (m_tPhaseMoles_old[iph] * 0.01)) {
#ifdef DEBUG_MODE
if (vcs_debug_print_lvl >= 2) {
if (m_debug_print_lvl >= 2) {
plogf(" --- Major species changed to minor: ");
plogf("%-12s", m_speciesName[kspec].c_str());
plogendl();
@ -879,7 +879,7 @@ namespace VCSnonideal {
m_deltaGRxn_old[ll] = m_deltaGRxn_new[ll];
}
#ifdef DEBUG_MODE
if (vcs_debug_print_lvl >= 2) {
if (m_debug_print_lvl >= 2) {
if (m_rxnStatus[irxn] >= 0) {
plogf(" --- SS species changed to zeroedss: ");
plogf("%-12s", m_speciesName[kspec].c_str());
@ -958,7 +958,7 @@ namespace VCSnonideal {
L_MAIN_LOOP_END: ;
#endif
#ifdef DEBUG_MODE
if (vcs_debug_print_lvl >= 2) {
if (m_debug_print_lvl >= 2) {
m_molNumSpecies_new[kspec] = m_molNumSpecies_old[kspec] + m_deltaMolNumSpecies[kspec];
plogf(" --- "); plogf("%-12.12s", m_speciesName[kspec].c_str());
plogf("%3d%11.4E%11.4E%11.4E | %s",
@ -972,7 +972,7 @@ namespace VCSnonideal {
} /**************** END OF MAIN LOOP OVER FORMATION REACTIONS ************/
#ifdef DEBUG_MODE
if (vcs_debug_print_lvl >= 2) {
if (m_debug_print_lvl >= 2) {
for (k = 0; k < m_numComponents; k++) {
plogf(" --- "); plogf("%-12.12s", m_speciesName[k].c_str());
plogf(" c%11.4E%11.4E%11.4E |\n",
@ -1018,7 +1018,7 @@ namespace VCSnonideal {
/* Reduce the size of the step by the multiplicative factor, par */
par *= 0.99;
#ifdef DEBUG_MODE
if (vcs_debug_print_lvl >= 2) {
if (m_debug_print_lvl >= 2) {
plogf(" --- Reduction in step size due to component ");
plogf("%s", m_speciesName[ll].c_str());
plogf(" going negative = %11.3E", par);
@ -1229,7 +1229,7 @@ namespace VCSnonideal {
++(m_VCount->Its);
++it1;
#ifdef DEBUG_MODE
if (vcs_debug_print_lvl >= 2) {
if (m_debug_print_lvl >= 2) {
plogf(" --- Increment counter increased, step is accepted: %4d",
m_VCount->Its);
plogendl();
@ -1274,7 +1274,7 @@ namespace VCSnonideal {
}
if (soldel) {
#ifdef DEBUG_MODE
if (vcs_debug_print_lvl >= 1) {
if (m_debug_print_lvl >= 1) {
plogf(" --- Setting microscopic phase %d to zero", iph);
plogendl();
}
@ -1317,14 +1317,14 @@ namespace VCSnonideal {
/***************** CHECK FOR ELEMENT ABUNDANCE****************************/
/*************************************************************************/
#ifdef DEBUG_MODE
if (vcs_debug_print_lvl >= 2) {
if (m_debug_print_lvl >= 2) {
plogf(" --- Normal element abundance check");
}
#endif
vcs_elab();
if (! vcs_elabcheck(0)) {
#ifdef DEBUG_MODE
if (vcs_debug_print_lvl >= 2) {
if (m_debug_print_lvl >= 2) {
plogf(" - failed -> redoing element abundances.");
plogendl();
}
@ -1336,7 +1336,7 @@ namespace VCSnonideal {
}
#ifdef DEBUG_MODE
else {
if (vcs_debug_print_lvl >= 2) {
if (m_debug_print_lvl >= 2) {
plogf(" - passed");
plogendl();
}
@ -1389,7 +1389,7 @@ namespace VCSnonideal {
if (doSwap) {
if (m_stoichCoeffRxnMatrix[i][j] != 0.0) {
#ifdef DEBUG_MODE
if (vcs_debug_print_lvl >= 2) {
if (m_debug_print_lvl >= 2) {
plogf(" --- Get a new basis because %s", m_speciesName[l].c_str());
plogf(" is better than comp %s", m_speciesName[j].c_str());
plogf(" and share nonzero stoic: %-9.1f",
@ -1408,7 +1408,7 @@ namespace VCSnonideal {
if (m_stoichCoeffRxnMatrix[i][j] != 0.0) {
if (dg[i] < 0.0) {
#ifdef DEBUG_MODE
if (vcs_debug_print_lvl >= 2) {
if (m_debug_print_lvl >= 2) {
plogf(" --- Get a new basis because %s", m_speciesName[l].c_str());
plogf(" has dg < 0.0 and comp %s has zero mole num", m_speciesName[j].c_str());
plogf(" and share nonzero stoic: %-9.1f",
@ -1449,7 +1449,7 @@ namespace VCSnonideal {
if (doSwap) {
if (m_stoichCoeffRxnMatrix[i][j] != 0.0) {
#ifdef DEBUG_MODE
if (vcs_debug_print_lvl >= 2) {
if (m_debug_print_lvl >= 2) {
plogf(" --- Get a new basis because ");
plogf("%s", m_speciesName[l].c_str());
plogf(" is better than comp ");
@ -1468,7 +1468,7 @@ namespace VCSnonideal {
if (m_stoichCoeffRxnMatrix[i][j] != 0.0) {
if (dg[i] < 0.0) {
#ifdef DEBUG_MODE
if (vcs_debug_print_lvl >= 2) {
if (m_debug_print_lvl >= 2) {
plogf(" --- Get a new basis because %s", m_speciesName[l].c_str());
plogf(" has dg < 0.0 and comp %s has zero mole num", m_speciesName[j].c_str());
plogf(" and share nonzero stoic: %-9.1f",
@ -1486,7 +1486,7 @@ namespace VCSnonideal {
}
}
#ifdef DEBUG_MODE
if (vcs_debug_print_lvl >= 2) {
if (m_debug_print_lvl >= 2) {
plogf(" --- Check for an optimum basis passed");
plogendl();
}
@ -1500,7 +1500,7 @@ namespace VCSnonideal {
*/
if (iti == 0) {
#ifdef DEBUG_MODE
if (vcs_debug_print_lvl >= 2) {
if (m_debug_print_lvl >= 2) {
plogf(" --- Reevaluate major-minor status of noncomponents:\n");
}
#endif
@ -1511,7 +1511,7 @@ namespace VCSnonideal {
int speciesType = vcs_species_type(kspec);
if (speciesType < VCS_SPECIES_MINOR) {
#ifdef DEBUG_MODE
if (vcs_debug_print_lvl >= 2) {
if (m_debug_print_lvl >= 2) {
if (m_rxnStatus[irxn] >= VCS_SPECIES_MINOR) {
plogf(" --- major/minor species is now zeroed out: %s\n",
m_speciesName[kspec].c_str());
@ -1521,7 +1521,7 @@ namespace VCSnonideal {
++m_numRxnMinorZeroed;
} else if (speciesType == VCS_SPECIES_MINOR) {
#ifdef DEBUG_MODE
if (vcs_debug_print_lvl >= 2) {
if (m_debug_print_lvl >= 2) {
if (m_rxnStatus[irxn] != VCS_SPECIES_MINOR) {
if (m_rxnStatus[irxn] == VCS_SPECIES_MAJOR) {
plogf(" --- Noncomponent turned from major to minor: ");
@ -1539,7 +1539,7 @@ namespace VCSnonideal {
} else if (speciesType == VCS_SPECIES_MAJOR) {
if (m_rxnStatus[irxn] != VCS_SPECIES_MAJOR) {
#ifdef DEBUG_MODE
if (vcs_debug_print_lvl >= 2) {
if (m_debug_print_lvl >= 2) {
if (m_rxnStatus[irxn] == VCS_SPECIES_MINOR) {
plogf(" --- Noncomponent turned from minor to major: ");
} else if (m_rxnStatus[irxn] == VCS_SPECIES_COMPONENT) {
@ -1574,7 +1574,7 @@ namespace VCSnonideal {
L_EQUILIB_CHECK: ;
if (! im) {
#ifdef DEBUG_MODE
if (vcs_debug_print_lvl >= 2) {
if (m_debug_print_lvl >= 2) {
plogf(" --- Equilibrium check for major species: ");
}
#endif
@ -1589,7 +1589,7 @@ namespace VCSnonideal {
goto L_RETURN_BLOCK;
} else {
#ifdef DEBUG_MODE
if (vcs_debug_print_lvl >= 2) {
if (m_debug_print_lvl >= 2) {
plogf("%s failed\n", m_speciesName[m_indexRxnToSpecies[irxn]].c_str());
}
#endif
@ -1607,7 +1607,7 @@ namespace VCSnonideal {
}
}
#ifdef DEBUG_MODE
if (vcs_debug_print_lvl >= 2) {
if (m_debug_print_lvl >= 2) {
plogf(" MAJOR SPECIES CONVERGENCE achieved");
plogendl();
}
@ -1615,7 +1615,7 @@ namespace VCSnonideal {
}
#ifdef DEBUG_MODE
else {
if (vcs_debug_print_lvl >= 2) {
if (m_debug_print_lvl >= 2) {
plogf(" MAJOR SPECIES CONVERGENCE achieved "
"(because there are no major species)");
plogendl();
@ -1641,7 +1641,7 @@ namespace VCSnonideal {
uptodate_minors = TRUE;
}
#ifdef DEBUG_MODE
if (vcs_debug_print_lvl >= 2) {
if (m_debug_print_lvl >= 2) {
plogf(" --- Equilibrium check for minor species: ");
}
#endif
@ -1656,7 +1656,7 @@ namespace VCSnonideal {
goto L_RETURN_BLOCK;
}
#ifdef DEBUG_MODE
if (vcs_debug_print_lvl >= 2) {
if (m_debug_print_lvl >= 2) {
plogf("%s failed\n", m_speciesName[m_indexRxnToSpecies[irxn]].c_str());
}
#endif
@ -1669,7 +1669,7 @@ namespace VCSnonideal {
}
}
#ifdef DEBUG_MODE
if (vcs_debug_print_lvl >= 2) {
if (m_debug_print_lvl >= 2) {
plogf(" CONVERGENCE achieved\n");
}
#endif
@ -1687,7 +1687,7 @@ namespace VCSnonideal {
if (lec) {
if (!giveUpOnElemAbund) {
#ifdef DEBUG_MODE
if (vcs_debug_print_lvl >= 2) {
if (m_debug_print_lvl >= 2) {
plogf(" --- Check the Full Element Abundances: ");
}
#endif
@ -1698,7 +1698,7 @@ namespace VCSnonideal {
*/
if (! vcs_elabcheck(1) ) {
#ifdef DEBUG_MODE
if (vcs_debug_print_lvl >= 2) {
if (m_debug_print_lvl >= 2) {
if (! vcs_elabcheck(0)) {
plogf(" failed\n");
} else {
@ -1710,7 +1710,7 @@ namespace VCSnonideal {
goto L_ELEM_ABUND_CHECK;
}
#ifdef DEBUG_MODE
if (vcs_debug_print_lvl >= 2) {
if (m_debug_print_lvl >= 2) {
plogf(" passed\n");
}
#endif
@ -1790,7 +1790,7 @@ namespace VCSnonideal {
* Probably an unrecoverable range error
*/
#ifdef DEBUG_MODE
if (vcs_debug_print_lvl >= 2) {
if (m_debug_print_lvl >= 2) {
plogf(" --- vcs_solve_tp: RANGE SPACE ERROR ENCOUNTERED\n");
plogf(" --- vcs_solve_tp: - Giving up on NE Element Abundance satisfaction \n");
plogf(" --- vcs_solve_tp: - However, NC Element Abundance criteria is satisfied \n");
@ -1872,7 +1872,7 @@ namespace VCSnonideal {
MajorSpeciesHaveConverged = true;
iti = 0;
#ifdef DEBUG_MODE
if (vcs_debug_print_lvl >= 1) {
if (m_debug_print_lvl >= 1) {
plogf(" --- add_all_deleted(): some rxns not converged. RETURNING TO LOOP!");
plogendl();
}
@ -2088,7 +2088,7 @@ namespace VCSnonideal {
* Voltage calculation
* HKM -> Need to check the sign
*/
dx = m_deltaGRxn_new[irxn]/ Faraday_dim;
dx = m_deltaGRxn_new[irxn]/ m_Faraday_dim;
#ifdef DEBUG_MODE
sprintf(ANOTE,"voltage species alternative calc");
#endif
@ -2318,7 +2318,7 @@ namespace VCSnonideal {
int *phaseID = VCS_DATA_PTR(m_phaseID);
double dx;
#ifdef DEBUG_MODE
if (vcs_debug_print_lvl >= 2) {
if (m_debug_print_lvl >= 2) {
plogf(" --- Add back a deleted species: %-12s\n", m_speciesName[kspec].c_str());
}
#endif
@ -2398,7 +2398,7 @@ namespace VCSnonideal {
*/
Vphase->Existence = 0;
#ifdef DEBUG_MODE
if (vcs_debug_print_lvl >= 2) {
if (m_debug_print_lvl >= 2) {
plogf(" --- delete_multiphase %d, %s\n", iph, Vphase->PhaseName.c_str());
}
#endif
@ -2466,7 +2466,7 @@ namespace VCSnonideal {
++(m_numRxnRdc);
++(m_numSpeciesRdc);
#ifdef DEBUG_MODE
if (vcs_debug_print_lvl >= 2) {
if (m_debug_print_lvl >= 2) {
plogf(" --- Make %s", m_speciesName[kspec].c_str());
plogf(" an active but zeroed species because its phase "
"was zeroed\n");
@ -2498,9 +2498,9 @@ namespace VCSnonideal {
int VCS_SOLVE::recheck_deleted(void)
{
int iph, kspec, irxn, npb;
double *xtcutoff = VCS_DATA_PTR(TmpPhase);
double *xtcutoff = VCS_DATA_PTR(m_TmpPhase);
#ifdef DEBUG_MODE
if (vcs_debug_print_lvl >= 2) {
if (m_debug_print_lvl >= 2) {
plogf(" --- Start rechecking deleted species in multispec phases\n");
}
#endif
@ -2594,8 +2594,8 @@ namespace VCSnonideal {
for (kspec = m_numSpeciesRdc; kspec < m_numSpeciesTot; ++kspec) {
iph = m_phaseID[kspec];
m_feSpecies_curr[kspec] = (m_SSfeSpecies[kspec] + log(m_actCoeffSpecies_old[kspec])
- SpecLnMnaught[kspec]
+ m_chargeSpecies[kspec] * Faraday_dim * m_phasePhi[iph]);
- m_lnMnaughtSpecies[kspec]
+ m_chargeSpecies[kspec] * m_Faraday_dim * m_phasePhi[iph]);
}
/*
* Recalculate the DeltaG's of the formation reactions for the
@ -2613,7 +2613,7 @@ namespace VCSnonideal {
retn = delta_species(kspec, &dx);
if (retn == 0) {
#ifdef DEBUG_MODE
if (vcs_debug_print_lvl) {
if (m_debug_print_lvl) {
plogf(" --- add_deleted(): delta_species() failed for species %s (%d) with mol number %g\n",
m_speciesName[kspec].c_str(), kspec, dx);
}
@ -2623,7 +2623,7 @@ namespace VCSnonideal {
retn = delta_species(kspec, &dx);
#ifdef DEBUG_MODE
if (retn == 0) {
if (vcs_debug_print_lvl) {
if (m_debug_print_lvl) {
plogf(" --- add_deleted(): delta_species() failed for species %s (%d) with mol number %g\n",
m_speciesName[kspec].c_str(), kspec, dx);
}
@ -2632,7 +2632,7 @@ namespace VCSnonideal {
}
}
#ifdef DEBUG_MODE
if (vcs_debug_print_lvl >= 2) {
if (m_debug_print_lvl >= 2) {
if (retn != 0) {
plogf(" --- add_deleted(): species %s added back in with mol number %g",
m_speciesName[kspec].c_str(), dx);
@ -2657,7 +2657,7 @@ namespace VCSnonideal {
if (fabs(m_deltaGRxn_old[irxn]) > m_tolmin) {
retn++;
#ifdef DEBUG_MODE
if (vcs_debug_print_lvl >= 2) {
if (m_debug_print_lvl >= 2) {
plogf(" --- add_deleted(): species %s with mol number %g not converged: DG = %g",
m_speciesName[kspec].c_str(), m_molNumSpecies_old[kspec],
m_deltaGRxn_old[irxn]);
@ -2719,7 +2719,7 @@ namespace VCSnonideal {
}
#ifdef DEBUG_MODE
if (vcs_debug_print_lvl >= 2) {
if (m_debug_print_lvl >= 2) {
plogf(" --- subroutine FORCE: Beginning Slope = %g\n", s1);
plogf(" --- subroutine FORCE: End Slope = %g\n", s2);
}
@ -2727,7 +2727,7 @@ namespace VCSnonideal {
if (s1 > 0.0) {
#ifdef DEBUG_MODE
if (vcs_debug_print_lvl >= 2) {
if (m_debug_print_lvl >= 2) {
plogf(" --- subroutine FORCE produced no adjustments,");
if (s1 < 1.0E-40) {
plogf(" s1 positive but really small");
@ -2742,7 +2742,7 @@ namespace VCSnonideal {
if (s2 <= 0.0) {
#ifdef DEBUG_MODE
if (vcs_debug_print_lvl >= 2) {
if (m_debug_print_lvl >= 2) {
plogf(" --- subroutine FORCE produced no adjustments, s2 < 0");
plogendl();
}
@ -2759,14 +2759,14 @@ namespace VCSnonideal {
}
if (al >= 0.95 || al < 0.0) {
#ifdef DEBUG_MODE
if (vcs_debug_print_lvl >= 2) {
if (m_debug_print_lvl >= 2) {
plogf(" --- subroutine FORCE produced no adjustments (al = %g)\n", al);
}
#endif
return FALSE;
}
#ifdef DEBUG_MODE
if (vcs_debug_print_lvl >= 2) {
if (m_debug_print_lvl >= 2) {
plogf(" --- subroutine FORCE produced a damping factor = %g\n", al);
}
#endif
@ -2775,7 +2775,7 @@ namespace VCSnonideal {
/* **** ADJUST MOLE NUMBERS, CHEM. POT *************** */
/* *************************************************** */
#ifdef DEBUG_MODE
if (vcs_debug_print_lvl >= 2) {
if (m_debug_print_lvl >= 2) {
vcs_dcopy(VCS_DATA_PTR(m_deltaGRxn_tmp), VCS_DATA_PTR(m_deltaGRxn_new),
m_numRxnRdc);
}
@ -2791,7 +2791,7 @@ namespace VCSnonideal {
vcs_updateVP(1);
#ifdef DEBUG_MODE
if (vcs_debug_print_lvl >= 2) {
if (m_debug_print_lvl >= 2) {
plogf(" --- subroutine FORCE adjusted the mole "
"numbers, AL = %10.3f\n", al);
}
@ -2820,7 +2820,7 @@ namespace VCSnonideal {
#ifdef DEBUG_MODE
if (vcs_debug_print_lvl >= 2) {
if (m_debug_print_lvl >= 2) {
plogf(" --- subroutine FORCE: Adj End Slope = %g", s2);
plogendl();
}
@ -2856,7 +2856,7 @@ namespace VCSnonideal {
double *dnPhase_irxn;
#ifdef DEBUG_MODE
char ANOTE[128];
if (vcs_debug_print_lvl >= 2) {
if (m_debug_print_lvl >= 2) {
plogf(" "); for (j = 0; j < 82; j++) plogf("-"); plogf("\n");
plogf(" --- Subroutine vcs_RxnStepSizes called - Details:\n");
plogf(" "); for (j = 0; j < 82; j++) plogf("-"); plogf("\n");
@ -2941,7 +2941,7 @@ namespace VCSnonideal {
if (fabs(m_deltaGRxn_new[irxn]) <= m_tolmaj2) {
#ifdef DEBUG_MODE
sprintf(ANOTE,"Skipped: superconverged DG = %11.3E", m_deltaGRxn_new[irxn]);
if (vcs_debug_print_lvl >= 2) {
if (m_debug_print_lvl >= 2) {
plogf(" --- %-12.12s", m_speciesName[kspec].c_str());
plogf(" %12.4E %12.4E %12.4E | %s\n",
m_molNumSpecies_old[kspec], m_deltaMolNumSpecies[kspec],
@ -2958,7 +2958,7 @@ namespace VCSnonideal {
#ifdef DEBUG_MODE
sprintf(ANOTE,"Skipped: IC = %3d and DG >0: %11.3E",
m_rxnStatus[irxn], m_deltaGRxn_new[irxn]);
if (vcs_debug_print_lvl >= 2) {
if (m_debug_print_lvl >= 2) {
plogf(" --- %-12.12s", m_speciesName[kspec].c_str());
plogf(" %12.4E %12.4E %12.4E | %s\n",
m_molNumSpecies_old[kspec], m_deltaMolNumSpecies[kspec],
@ -3103,7 +3103,7 @@ namespace VCSnonideal {
Vphase->Existence = 0;
m_tPhaseMoles_old[iph] = 0.0;
#ifdef DEBUG_MODE
if (vcs_debug_print_lvl >= 2) {
if (m_debug_print_lvl >= 2) {
plogf(" --- vcs_RxnStepSizes Special section to delete %s\n",
m_speciesName[k].c_str());
plogf(" --- Immediate return - Restart iteration\n");
@ -3121,7 +3121,7 @@ namespace VCSnonideal {
}
} /* End of regular processing */
#ifdef DEBUG_MODE
if (vcs_debug_print_lvl >= 2) {
if (m_debug_print_lvl >= 2) {
plogf(" --- %-12.12s", m_speciesName[kspec].c_str());
plogf(" %12.4E %12.4E %12.4E | %s\n",
m_molNumSpecies_old[kspec], m_deltaMolNumSpecies[kspec],
@ -3131,7 +3131,7 @@ namespace VCSnonideal {
} /* End of loop over m_speciesUnknownType */
} /* End of loop over non-component stoichiometric formation reactions */
#ifdef DEBUG_MODE
if (vcs_debug_print_lvl >= 2) {
if (m_debug_print_lvl >= 2) {
plogf(" "); vcs_print_line("-", 82);
}
#endif
@ -3177,7 +3177,7 @@ namespace VCSnonideal {
}
#ifdef DEBUG_MODE
if (vcs_debug_print_lvl >= 2) {
if (m_debug_print_lvl >= 2) {
plogf(" --- Subroutine vcs_deltag called for ");
if (l < 0) {
plogf("major noncomponents\n");
@ -3391,12 +3391,12 @@ namespace VCSnonideal {
double *dptr, *scrxn_ptr;
Cantera::clockWC tickTock;
#ifdef DEBUG_MODE
if (vcs_debug_print_lvl >= 2) {
if (m_debug_print_lvl >= 2) {
plogf(" "); for(i=0; i<77; i++) plogf("-"); plogf("\n");
plogf(" --- Subroutine BASOPT called to ");
if (ifirst) plogf("calculate the number of components\n");
else plogf("reevaluate the components\n");
if (vcs_debug_print_lvl >= 2) {
if (m_debug_print_lvl >= 2) {
plogf("\n");
plogf(" --- Formula Matrix used in BASOPT calculation\n");
plogf(" --- Active | ");
@ -3581,7 +3581,7 @@ namespace VCSnonideal {
m_indexRxnToSpecies[i] = ncTrial + i;
}
#ifdef DEBUG_MODE
if (vcs_debug_print_lvl >= 2) {
if (m_debug_print_lvl >= 2) {
plogf(" --- Total number of components found = %3d (ne = %d)\n ",
ncTrial, m_numElemConstraints);
}
@ -3647,7 +3647,7 @@ namespace VCSnonideal {
/* ****************************************** */
if (jr != k) {
#ifdef DEBUG_MODE
if (vcs_debug_print_lvl >= 2) {
if (m_debug_print_lvl >= 2) {
plogf(" --- %-12.12s", (m_speciesName[k]).c_str());
plogf("(%9.2g) replaces %-12.12s", m_molNumSpecies_old[k], m_speciesName[jr].c_str());
plogf("(%9.2g) as component %3d\n", m_molNumSpecies_old[jr], jr);
@ -3658,7 +3658,7 @@ namespace VCSnonideal {
}
#ifdef DEBUG_MODE
else {
if (vcs_debug_print_lvl >= 2) {
if (m_debug_print_lvl >= 2) {
plogf(" --- %-12.12s", m_speciesName[k].c_str());
plogf("(%9.2g) remains ", m_molNumSpecies_old[k]);
plogf(" as component %3d\n", jr);
@ -3798,7 +3798,7 @@ namespace VCSnonideal {
#ifdef DEBUG_MODE
if (vcs_debug_print_lvl >= 2) {
if (m_debug_print_lvl >= 2) {
plogf(" --- Components:");
for (j = 0; j < ncTrial; j++) {
plogf(" %3d ", j);
@ -3962,7 +3962,7 @@ namespace VCSnonideal {
if (negChangeComp > 0.0) {
if (m_molNumSpecies_old[j] < 1.0E-60) {
#ifdef DEBUG_MODE
if (vcs_debug_print_lvl >= 2) {
if (m_debug_print_lvl >= 2) {
plogf(" --- %s would have popped back into existance but"
" needed component %s is zero\n",
m_speciesName[kspec].c_str(), m_speciesName[j].c_str());
@ -3986,7 +3986,7 @@ namespace VCSnonideal {
double maxPermissible = m_elemAbundancesGoal[j] / atomComp;
if (maxPermissible < VCS_DELETE_MINORSPECIES_CUTOFF) {
#ifdef DEBUG_MODE
if (vcs_debug_print_lvl >= 2) {
if (m_debug_print_lvl >= 2) {
plogf(" --- %s would have popped back into existance but"
" needed element %s is zero\n",
m_speciesName[kspec].c_str(), (m_elementName[j]).c_str());
@ -4105,7 +4105,7 @@ namespace VCSnonideal {
int k, kspec;
#ifdef DEBUG_MODE
//if (vcs_debug_print_lvl >= 2) {
//if (m_debug_print_lvl >= 2) {
// plogf(" --- Subroutine vcs_chemPotPhase called for phase %d\n",
// iph);
//}
@ -4124,7 +4124,7 @@ namespace VCSnonideal {
Vphase->sendToVCSActCoeff(ac);
double phi = Vphase->electricPotential();
double Faraday_phi = Faraday_dim * phi;
double Faraday_phi = m_Faraday_dim * phi;
for (k = 0; k < nkk; k++) {
kspec = Vphase->IndSpecies[k];
@ -4152,10 +4152,10 @@ namespace VCSnonideal {
mu_i[kspec] = m_SSfeSpecies[kspec] + m_chargeSpecies[kspec] * Faraday_phi;
} else if (molNum[kspec] <= VCS_DELETE_MINORSPECIES_CUTOFF) {
mu_i[kspec] = m_SSfeSpecies[kspec] + log(ac[kspec] * VCS_DELETE_MINORSPECIES_CUTOFF)
- tlogMoles - SpecLnMnaught[kspec] + m_chargeSpecies[kspec] * Faraday_phi;
- tlogMoles - m_lnMnaughtSpecies[kspec] + m_chargeSpecies[kspec] * Faraday_phi;
} else {
mu_i[kspec] = m_SSfeSpecies[kspec] + log(ac[kspec] * molNum[kspec])
- tlogMoles - SpecLnMnaught[kspec] + m_chargeSpecies[kspec] * Faraday_phi;
- tlogMoles - m_lnMnaughtSpecies[kspec] + m_chargeSpecies[kspec] * Faraday_phi;
}
}
}
@ -4304,7 +4304,7 @@ namespace VCSnonideal {
#endif
#ifdef DEBUG_MODE
if (vcs_debug_print_lvl >= 2) {
if (m_debug_print_lvl >= 2) {
if (ll == 0) {
if (lbot != 0) {
plogf(" --- Subroutine vcs_dfe called for one species: ");
@ -4328,7 +4328,7 @@ namespace VCSnonideal {
tPhMoles_ptr = VCS_DATA_PTR(m_tPhaseMoles_new);
actCoeff_ptr = VCS_DATA_PTR(m_actCoeffSpecies_new);
}
tlogMoles = VCS_DATA_PTR(TmpPhase);
tlogMoles = VCS_DATA_PTR(m_TmpPhase);
/*
* Might as well recalculate the phase mole vector
* and compare to the storred one. They should be correct.
@ -4363,7 +4363,7 @@ namespace VCSnonideal {
* Zero the indicator that that tells us the activity coefficients
* are current
*/
vcs_izero(VCS_DATA_PTR(CurrPhAC), m_numPhases);
vcs_izero(VCS_DATA_PTR(m_phaseACAreCurrent), m_numPhases);
if (ll != 0) {
l1 = lbot;
@ -4378,14 +4378,14 @@ namespace VCSnonideal {
* not current
*/
for (iphase = 0; iphase < m_numPhases; iphase++) {
if (!CurrPhAC[iphase]) {
if (!m_phaseACAreCurrent[iphase]) {
Vphase = m_VolPhaseList[iphase];
if (!Vphase->SingleSpecies) {
Vphase->setMolesFromVCS(z);
Vphase->sendToVCSActCoeff(VCS_DATA_PTR(actCoeff_ptr));
}
m_phasePhi[iphase] = Vphase->electricPotential();
CurrPhAC[iphase] = 1;
m_phaseACAreCurrent[iphase] = 1;
}
}
/* ************************************************************** */
@ -4411,7 +4411,7 @@ namespace VCSnonideal {
}
#endif
m_feSpecies_curr[kspec] =
m_SSfeSpecies[kspec] + m_chargeSpecies[kspec] * Faraday_dim * m_phasePhi[iphase];
m_SSfeSpecies[kspec] + m_chargeSpecies[kspec] * m_Faraday_dim * m_phasePhi[iphase];
} else {
if (m_SSPhase[kspec]) {
m_feSpecies_curr[kspec] = m_SSfeSpecies[kspec];
@ -4421,15 +4421,16 @@ namespace VCSnonideal {
if (tPhMoles_ptr[iph] > 0.0) {
m_feSpecies_curr[kspec] = m_SSfeSpecies[kspec]
+ log(actCoeff_ptr[kspec] * VCS_DELETE_MINORSPECIES_CUTOFF)
- tlogMoles[m_phaseID[kspec]] - SpecLnMnaught[kspec]
+ m_chargeSpecies[kspec] * Faraday_dim * m_phasePhi[iphase];
- tlogMoles[m_phaseID[kspec]] - m_lnMnaughtSpecies[kspec]
+ m_chargeSpecies[kspec] * m_Faraday_dim * m_phasePhi[iphase];
} else {
m_feSpecies_curr[kspec] = m_SSfeSpecies[kspec];
}
} else {
m_feSpecies_curr[kspec] = m_SSfeSpecies[kspec] + log(actCoeff_ptr[kspec] * z[kspec])
- tlogMoles[m_phaseID[kspec]] - SpecLnMnaught[kspec]
+ m_chargeSpecies[kspec] * Faraday_dim * m_phasePhi[iphase];
m_feSpecies_curr[kspec] = m_SSfeSpecies[kspec]
+ log(actCoeff_ptr[kspec] * z[kspec])
- tlogMoles[m_phaseID[kspec]] - m_lnMnaughtSpecies[kspec]
+ m_chargeSpecies[kspec] * m_Faraday_dim * m_phasePhi[iphase];
}
}
}
@ -4454,7 +4455,8 @@ namespace VCSnonideal {
}
#endif
m_feSpecies_curr[kspec] =
m_SSfeSpecies[kspec] + m_chargeSpecies[kspec] * Faraday_dim * m_phasePhi[iphase];
m_SSfeSpecies[kspec]
+ m_chargeSpecies[kspec] * m_Faraday_dim * m_phasePhi[iphase];
} else {
if (m_SSPhase[kspec]) {
m_feSpecies_curr[kspec] = m_SSfeSpecies[kspec];
@ -4464,15 +4466,15 @@ namespace VCSnonideal {
if (tPhMoles_ptr[iph] > 0.0) {
m_feSpecies_curr[kspec] = m_SSfeSpecies[kspec]
+ log(actCoeff_ptr[kspec] * VCS_DELETE_MINORSPECIES_CUTOFF)
- tlogMoles[m_phaseID[kspec]] - SpecLnMnaught[kspec]
+ m_chargeSpecies[kspec] * Faraday_dim * m_phasePhi[iphase]; ;
- tlogMoles[m_phaseID[kspec]] - m_lnMnaughtSpecies[kspec]
+ m_chargeSpecies[kspec] * m_Faraday_dim * m_phasePhi[iphase]; ;
} else {
m_feSpecies_curr[kspec] = m_SSfeSpecies[kspec];
}
} else {
m_feSpecies_curr[kspec] = m_SSfeSpecies[kspec] + log(actCoeff_ptr[kspec] * z[kspec])
- tlogMoles[m_phaseID[kspec]] - SpecLnMnaught[kspec]
+ m_chargeSpecies[kspec] * Faraday_dim * m_phasePhi[iphase];
- tlogMoles[m_phaseID[kspec]] - m_lnMnaughtSpecies[kspec]
+ m_chargeSpecies[kspec] * m_Faraday_dim * m_phasePhi[iphase];
}
}
}
@ -4498,7 +4500,8 @@ namespace VCSnonideal {
}
#endif
m_feSpecies_curr[kspec] =
m_SSfeSpecies[kspec] + m_chargeSpecies[kspec] * Faraday_dim * m_phasePhi[iphase]; ;
m_SSfeSpecies[kspec] +
m_chargeSpecies[kspec] * m_Faraday_dim * m_phasePhi[iphase]; ;
} else {
if (m_SSPhase[kspec]) {
m_feSpecies_curr[kspec] = m_SSfeSpecies[kspec];
@ -4508,14 +4511,15 @@ namespace VCSnonideal {
if (tPhMoles_ptr[iph] > 0.0) {
m_feSpecies_curr[kspec] = m_SSfeSpecies[kspec]
+ log(actCoeff_ptr[kspec] * VCS_DELETE_MINORSPECIES_CUTOFF)
- tlogMoles[m_phaseID[kspec]] - SpecLnMnaught[kspec];
- tlogMoles[m_phaseID[kspec]] - m_lnMnaughtSpecies[kspec];
} else {
m_feSpecies_curr[kspec] = m_SSfeSpecies[kspec];
}
} else {
st_ptr = m_speciesThermoList[kspec];
m_feSpecies_curr[kspec] = m_SSfeSpecies[kspec] + log(actCoeff_ptr[kspec] * z[kspec])
- tlogMoles[m_phaseID[kspec]] - SpecLnMnaught[kspec];
m_feSpecies_curr[kspec] = m_SSfeSpecies[kspec]
+ log(actCoeff_ptr[kspec] * z[kspec])
- tlogMoles[m_phaseID[kspec]] - m_lnMnaughtSpecies[kspec];
}
}
}
@ -4741,8 +4745,8 @@ namespace VCSnonideal {
SWAP(m_phaseID[k1], m_phaseID[k2], j);
SWAP(m_speciesMapIndex[k1], m_speciesMapIndex[k2], j);
SWAP(m_speciesLocalPhaseIndex[k1], m_speciesLocalPhaseIndex[k2], j);
SWAP(SpecActConvention[k1], SpecActConvention[k2], j);
SWAP(SpecLnMnaught[k1], SpecLnMnaught[k2], t1);
SWAP(m_actConventionSpecies[k1], m_actConventionSpecies[k2], j);
SWAP(m_lnMnaughtSpecies[k1], m_lnMnaughtSpecies[k2], t1);
SWAP(m_actCoeffSpecies_new[k1], m_actCoeffSpecies_new[k2], t1);
SWAP(m_actCoeffSpecies_old[k1], m_actCoeffSpecies_old[k2], t1);
SWAP(m_wtSpecies[k1], m_wtSpecies[k2], t1);
@ -4823,7 +4827,7 @@ namespace VCSnonideal {
vcs_VolPhase *vPhase = m_VolPhaseList[iphase];
#ifdef DEBUG_MODE
if (vcs_debug_print_lvl >= 2) {
if (m_debug_print_lvl >= 2) {
plogf(" --- Subroutine vcs_deltag_Phase called for phase %d\n",
iphase);
}