diff --git a/src/equil/vcs_MultiPhaseEquil.cpp b/src/equil/vcs_MultiPhaseEquil.cpp index ba2f69f4c..3de7d7017 100644 --- a/src/equil/vcs_MultiPhaseEquil.cpp +++ b/src/equil/vcs_MultiPhaseEquil.cpp @@ -60,19 +60,15 @@ int vcs_MultiPhaseEquil::equilibrate_TV(int XY, doublereal xtarget, { // doublereal dt = 1.0e3; doublereal Vtarget = m_mix->volume(); - doublereal dVdP; if ((XY != TV) && (XY != HV) && (XY != UV) && (XY != SV)) { throw CanteraError("vcs_MultiPhaseEquil::equilibrate_TV", "Wrong XY flag:" + int2str(XY)); } int maxiter = 100; int iSuccess = 0; - int innerXY; - double Pnow; if (XY == TV) { m_mix->setTemperature(xtarget); } - double Pnew; int strt = estimateEquil; double P1 = 0.0; double V1 = 0.0; @@ -83,24 +79,21 @@ int vcs_MultiPhaseEquil::equilibrate_TV(int XY, doublereal xtarget, doublereal Vnow, Verr; int printLvlSub = std::max(0, printLvl - 1); for (int n = 0; n < maxiter; n++) { - Pnow = m_mix->pressure(); + double Pnow = m_mix->pressure(); switch (XY) { case TV: iSuccess = equilibrate_TP(strt, printLvlSub, err, maxsteps, loglevel); break; case HV: - innerXY = HP; - iSuccess = equilibrate_HP(xtarget, innerXY, Tlow, Thigh, strt, + iSuccess = equilibrate_HP(xtarget, HP, Tlow, Thigh, strt, printLvlSub, err, maxsteps, loglevel); break; case UV: - innerXY = UP; - iSuccess = equilibrate_HP(xtarget, innerXY, Tlow, Thigh, strt, + iSuccess = equilibrate_HP(xtarget, UP, Tlow, Thigh, strt, printLvlSub, err, maxsteps, loglevel); break; case SV: - innerXY = SP; iSuccess = equilibrate_SP(xtarget, Tlow, Thigh, strt, printLvlSub, err, maxsteps, loglevel); break; @@ -127,9 +120,10 @@ int vcs_MultiPhaseEquil::equilibrate_TV(int XY, doublereal xtarget, if (Verr < err) { goto done; } + double Pnew; // find dV/dP if (n > 1) { - dVdP = (V2 - V1) / (P2 - P1); + double dVdP = (V2 - V1) / (P2 - P1); if (dVdP == 0.0) { throw CanteraError("vcs_MultiPhase::equilibrate_TV", "dVdP == 0.0"); @@ -145,7 +139,7 @@ int vcs_MultiPhaseEquil::equilibrate_TV(int XY, doublereal xtarget, } else { m_mix->setPressure(Pnow*1.01); - dVdP = (m_mix->volume() - Vnow)/(0.01*Pnow); + double dVdP = (m_mix->volume() - Vnow)/(0.01*Pnow); Pnew = Pnow + 0.5*(Vtarget - Vnow)/dVdP; if (Pnew < 0.5* Pnow) { Pnew = 0.5 * Pnow; @@ -188,11 +182,9 @@ int vcs_MultiPhaseEquil::equilibrate_HP(doublereal Htarget, Thigh = 2.0 * m_mix->maxTemp(); } - doublereal cpb = 1.0, dT, dTa, dTmax, Tnew; - doublereal Hnow; + doublereal cpb = 1.0, Tnew; doublereal Hlow = Undef; doublereal Hhigh = Undef; - doublereal Herr, HConvErr; doublereal Tnow = m_mix->temperature(); int printLvlSub = std::max(printLvl - 1, 0); @@ -203,11 +195,7 @@ int vcs_MultiPhaseEquil::equilibrate_HP(doublereal Htarget, Tnow = m_mix->temperature(); iSuccess = equilibrate_TP(strt, printLvlSub, err, maxsteps, loglevel); strt = 0; - if (XY == UP) { - Hnow = m_mix->IntEnergy(); - } else { - Hnow = m_mix->enthalpy(); - } + double Hnow = (XY == UP) ? m_mix->IntEnergy() : m_mix->enthalpy(); double pmoles[10]; pmoles[0] = m_mix->phaseMoles(0); double Tmoles = pmoles[0]; @@ -237,6 +225,7 @@ int vcs_MultiPhaseEquil::equilibrate_HP(doublereal Htarget, Hhigh = Hnow; } } + double dT, dTa, dTmax, Tnew; if (Hlow != Undef && Hhigh != Undef) { cpb = (Hhigh - Hlow)/(Thigh - Tlow); dT = (Htarget - Hnow)/cpb; @@ -257,8 +246,8 @@ int vcs_MultiPhaseEquil::equilibrate_HP(doublereal Htarget, } double acpb = std::max(fabs(cpb), 1.0E-6); double denom = std::max(fabs(Htarget), acpb); - Herr = Htarget - Hnow; - HConvErr = fabs((Herr)/denom); + double Herr = Htarget - Hnow; + double HConvErr = fabs((Herr)/denom); if (printLvl > 0) { plogf(" equilibrate_HP: It = %d, Tcurr = %g Hcurr = %g, Htarget = %g\n", n, Tnow, Hnow, Htarget); @@ -320,10 +309,8 @@ int vcs_MultiPhaseEquil::equilibrate_SP(doublereal Starget, } doublereal cpb = 1.0, dT, dTa, dTmax, Tnew; - doublereal Snow; doublereal Slow = Undef; doublereal Shigh = Undef; - doublereal Serr, SConvErr; doublereal Tnow = m_mix->temperature(); if (Tnow < Tlow) { Tlow = Tnow; @@ -341,7 +328,7 @@ int vcs_MultiPhaseEquil::equilibrate_SP(doublereal Starget, Tnow = m_mix->temperature(); int iSuccess = equilibrate_TP(strt, printLvlSub, err, maxsteps, loglevel); strt = 0; - Snow = m_mix->entropy(); + double Snow = m_mix->entropy(); double pmoles[10]; pmoles[0] = m_mix->phaseMoles(0); double Tmoles = pmoles[0]; @@ -396,8 +383,8 @@ int vcs_MultiPhaseEquil::equilibrate_SP(doublereal Starget, double acpb = std::max(fabs(cpb), 1.0E-6); double denom = std::max(fabs(Starget), acpb); - Serr = Starget - Snow; - SConvErr = fabs((Serr)/denom); + double Serr = Starget - Snow; + double SConvErr = fabs((Serr)/denom); if (printLvl > 0) { plogf(" equilibrate_SP: It = %d, Tcurr = %g Scurr = %g, Starget = %g\n", n, Tnow, Snow, Starget); @@ -441,10 +428,9 @@ int vcs_MultiPhaseEquil::equilibrate(int XY, int estimateEquil, int printLvl, doublereal err, int maxsteps, int loglevel) { - int iSuccess; doublereal xtarget; if (XY == TP) { - iSuccess = equilibrate_TP(estimateEquil, printLvl, err, maxsteps, loglevel); + return equilibrate_TP(estimateEquil, printLvl, err, maxsteps, loglevel); } else if (XY == HP || XY == UP) { if (XY == HP) { xtarget = m_mix->enthalpy(); @@ -453,36 +439,34 @@ int vcs_MultiPhaseEquil::equilibrate(int XY, int estimateEquil, } double Tlow = 0.5 * m_mix->minTemp(); double Thigh = 2.0 * m_mix->maxTemp(); - iSuccess = equilibrate_HP(xtarget, XY, Tlow, Thigh, - estimateEquil, printLvl, err, maxsteps, loglevel); + return equilibrate_HP(xtarget, XY, Tlow, Thigh, + estimateEquil, printLvl, err, maxsteps, loglevel); } else if (XY == SP) { xtarget = m_mix->entropy(); double Tlow = 0.5 * m_mix->minTemp(); double Thigh = 2.0 * m_mix->maxTemp(); - iSuccess = equilibrate_SP(xtarget, Tlow, Thigh, - estimateEquil, printLvl, err, maxsteps, loglevel); - + return equilibrate_SP(xtarget, Tlow, Thigh, + estimateEquil, printLvl, err, maxsteps, loglevel); } else if (XY == TV) { xtarget = m_mix->temperature(); - iSuccess = equilibrate_TV(XY, xtarget, - estimateEquil, printLvl, err, maxsteps, loglevel); + return equilibrate_TV(XY, xtarget, + estimateEquil, printLvl, err, maxsteps, loglevel); } else if (XY == HV) { xtarget = m_mix->enthalpy(); - iSuccess = equilibrate_TV(XY, xtarget, - estimateEquil, printLvl, err, maxsteps, loglevel); + return equilibrate_TV(XY, xtarget, + estimateEquil, printLvl, err, maxsteps, loglevel); } else if (XY == UV) { xtarget = m_mix->IntEnergy(); - iSuccess = equilibrate_TV(XY, xtarget, - estimateEquil, printLvl, err, maxsteps, loglevel); + return equilibrate_TV(XY, xtarget, + estimateEquil, printLvl, err, maxsteps, loglevel); } else if (XY == SV) { xtarget = m_mix->entropy(); - iSuccess = equilibrate_TV(XY, xtarget, estimateEquil, + return equilibrate_TV(XY, xtarget, estimateEquil, printLvl, err, maxsteps, loglevel); } else { throw CanteraError(" vcs_MultiPhaseEquil::equilibrate", "Unsupported Option"); } - return iSuccess; } int vcs_MultiPhaseEquil::equilibrate_TP(int estimateEquil, @@ -631,12 +615,7 @@ int vcs_MultiPhaseEquil::equilibrate_TP(int estimateEquil, void vcs_MultiPhaseEquil::reportCSV(const std::string& reportFile) { - size_t k; - size_t istart; - size_t nSpecies; - double vol = 0.0; - string sName; size_t nphase = m_vprob.NPhase; FILE* FP = fopen(reportFile.c_str(), "w"); @@ -657,19 +636,18 @@ void vcs_MultiPhaseEquil::reportCSV(const std::string& reportFile) std::vector mu0; std::vector molalities; - vol = 0.0; for (size_t iphase = 0; iphase < nphase; iphase++) { - istart = m_mix->speciesIndex(0, iphase); + size_t istart = m_mix->speciesIndex(0, iphase); Cantera::ThermoPhase& tref = m_mix->phase(iphase); - nSpecies = tref.nSpecies(); + size_t nSpecies = tref.nSpecies(); VolPM.resize(nSpecies, 0.0); tref.getPartialMolarVolumes(VCS_DATA_PTR(VolPM)); vcs_VolPhase* volP = m_vprob.VPhaseList[iphase]; double TMolesPhase = volP->totalMoles(); double VolPhaseVolumes = 0.0; - for (k = 0; k < nSpecies; k++) { + for (size_t k = 0; k < nSpecies; k++) { VolPhaseVolumes += VolPM[k] * mf[istart + k]; } VolPhaseVolumes *= TMolesPhase; @@ -685,14 +663,13 @@ void vcs_MultiPhaseEquil::reportCSV(const std::string& reportFile) fprintf(FP,"Number VCS iterations = %d\n", m_vprob.m_Iterations); for (size_t iphase = 0; iphase < nphase; iphase++) { - istart = m_mix->speciesIndex(0, iphase); + size_t istart = m_mix->speciesIndex(0, iphase); Cantera::ThermoPhase& tref = m_mix->phase(iphase); - Cantera::ThermoPhase* tp = &tref; string phaseName = tref.name(); vcs_VolPhase* volP = m_vprob.VPhaseList[iphase]; double TMolesPhase = volP->totalMoles(); //AssertTrace(TMolesPhase == m_mix->phaseMoles(iphase)); - nSpecies = tref.nSpecies(); + size_t nSpecies = tref.nSpecies(); activity.resize(nSpecies, 0.0); ac.resize(nSpecies, 0.0); @@ -701,15 +678,15 @@ void vcs_MultiPhaseEquil::reportCSV(const std::string& reportFile) VolPM.resize(nSpecies, 0.0); molalities.resize(nSpecies, 0.0); - int actConvention = tp->activityConvention(); - tp->getActivities(VCS_DATA_PTR(activity)); - tp->getActivityCoefficients(VCS_DATA_PTR(ac)); - tp->getStandardChemPotentials(VCS_DATA_PTR(mu0)); + int actConvention = tref.activityConvention(); + tref.getActivities(VCS_DATA_PTR(activity)); + tref.getActivityCoefficients(VCS_DATA_PTR(ac)); + tref.getStandardChemPotentials(VCS_DATA_PTR(mu0)); - tp->getPartialMolarVolumes(VCS_DATA_PTR(VolPM)); - tp->getChemPotentials(VCS_DATA_PTR(mu)); + tref.getPartialMolarVolumes(VCS_DATA_PTR(VolPM)); + tref.getChemPotentials(VCS_DATA_PTR(mu)); double VolPhaseVolumes = 0.0; - for (k = 0; k < nSpecies; k++) { + for (size_t k = 0; k < nSpecies; k++) { VolPhaseVolumes += VolPM[k] * mf[istart + k]; } VolPhaseVolumes *= TMolesPhase; @@ -717,9 +694,9 @@ void vcs_MultiPhaseEquil::reportCSV(const std::string& reportFile) if (actConvention == 1) { - MolalityVPSSTP* mTP = static_cast(tp); + MolalityVPSSTP* mTP = static_cast(&tref); mTP->getMolalities(VCS_DATA_PTR(molalities)); - tp->getChemPotentials(VCS_DATA_PTR(mu)); + tref.getChemPotentials(VCS_DATA_PTR(mu)); if (iphase == 0) { fprintf(FP," Name, Phase, PhaseMoles, Mole_Fract, " @@ -730,8 +707,8 @@ void vcs_MultiPhaseEquil::reportCSV(const std::string& reportFile) " , , ," " (J/kmol), (J/kmol), (kmol), (m**3/kmol), (m**3)\n"); } - for (k = 0; k < nSpecies; k++) { - sName = tp->speciesName(k); + for (size_t k = 0; k < nSpecies; k++) { + std::string sName = tref.speciesName(k); fprintf(FP,"%12s, %11s, %11.3e, %11.3e, %11.3e, %11.3e, %11.3e," "%11.3e, %11.3e, %11.3e, %11.3e, %11.3e\n", sName.c_str(), @@ -752,11 +729,11 @@ void vcs_MultiPhaseEquil::reportCSV(const std::string& reportFile) " , , ," " (J/kmol), (J/kmol), (kmol), (m**3/kmol), (m**3)\n"); } - for (k = 0; k < nSpecies; k++) { + for (size_t k = 0; k < nSpecies; k++) { molalities[k] = 0.0; } - for (k = 0; k < nSpecies; k++) { - sName = tp->speciesName(k); + for (size_t k = 0; k < nSpecies; k++) { + std::string sName = tref.speciesName(k); fprintf(FP,"%12s, %11s, %11.3e, %11.3e, %11.3e, %11.3e, %11.3e, " "%11.3e, %11.3e,% 11.3e, %11.3e, %11.3e\n", sName.c_str(), @@ -772,8 +749,8 @@ void vcs_MultiPhaseEquil::reportCSV(const std::string& reportFile) /* * Check consistency: These should be equal */ - tp->getChemPotentials(fe+istart); - for (k = 0; k < nSpecies; k++) { + tref.getChemPotentials(fe+istart); + for (size_t k = 0; k < nSpecies; k++) { if (!vcs_doubleEqual(fe[istart+k], mu[k])) { fprintf(FP,"ERROR: incompatibility!\n"); fclose(FP); @@ -828,9 +805,6 @@ int vcs_Cantera_to_vprob(Cantera::MultiPhase* mphase, vprob->Vol = mphase->volume(); vprob->Title = "MultiPhase Object"; - Cantera::ThermoPhase* tPhase = 0; - - bool gasPhase; int printLvl = vprob->m_printLvl; /* @@ -842,7 +816,7 @@ int vcs_Cantera_to_vprob(Cantera::MultiPhase* mphase, /* * Get the thermophase object - assume volume phase */ - tPhase = &(mphase->phase(iphase)); + Cantera::ThermoPhase* tPhase = &(mphase->phase(iphase)); size_t nelem = tPhase->nElements(); /* @@ -850,11 +824,7 @@ int vcs_Cantera_to_vprob(Cantera::MultiPhase* mphase, * current phase. */ int eos = tPhase->eosType(); - if (eos == cIdealGas) { - gasPhase = true; - } else { - gasPhase = false; - } + bool gasPhase = (eos == cIdealGas); /* * Find out the number of species in the phase @@ -1226,10 +1196,9 @@ int vcs_Cantera_update_vprob(Cantera::MultiPhase* mphase, vprob->T = mphase->temperature(); vprob->PresPA = mphase->pressure(); vprob->Vol = mphase->volume(); - Cantera::ThermoPhase* tPhase = 0; for (size_t iphase = 0; iphase < totNumPhases; iphase++) { - tPhase = &(mphase->phase(iphase)); + Cantera::ThermoPhase* tPhase = &(mphase->phase(iphase)); vcs_VolPhase* volPhase = vprob->VPhaseList[iphase]; /* * Set the electric potential of the volume phase from the diff --git a/src/equil/vcs_TP.cpp b/src/equil/vcs_TP.cpp index 6c2c13835..f9e125917 100644 --- a/src/equil/vcs_TP.cpp +++ b/src/equil/vcs_TP.cpp @@ -8,7 +8,6 @@ namespace VCSnonideal { int VCS_SOLVE::vcs_TP(int ipr, int ip1, int maxit, double T_arg, double pres_arg) { - int retn, iconv; /* * Store the temperature and pressure in the private global variables */ @@ -18,7 +17,7 @@ int VCS_SOLVE::vcs_TP(int ipr, int ip1, int maxit, double T_arg, double pres_arg * Evaluate the standard state free energies * at the current temperatures and pressures. */ - iconv = vcs_evalSS_TP(ipr, ip1, m_temperature, pres_arg); + int iconv = vcs_evalSS_TP(ipr, ip1, m_temperature, pres_arg); /* * Prepare the problem data: @@ -35,7 +34,7 @@ int VCS_SOLVE::vcs_TP(int ipr, int ip1, int maxit, double T_arg, double pres_arg * If so, go get one. If not, then */ if (m_doEstimateEquil) { - retn = vcs_inest_TP(); + int retn = vcs_inest_TP(); if (retn != VCS_SUCCESS) { plogf("vcs_inest_TP returned a failure flag\n"); } diff --git a/src/equil/vcs_VolPhase.cpp b/src/equil/vcs_VolPhase.cpp index d1391940c..95c657017 100644 --- a/src/equil/vcs_VolPhase.cpp +++ b/src/equil/vcs_VolPhase.cpp @@ -481,8 +481,6 @@ void vcs_VolPhase::setMoleFractionsState(const double totalMoles, void vcs_VolPhase::setMolesFromVCS(const int stateCalc, const double* molesSpeciesVCS) { - size_t kglob; - double tmp; v_totalMoles = m_totalMolesInert; if (molesSpeciesVCS == 0) { @@ -524,15 +522,15 @@ void vcs_VolPhase::setMolesFromVCS(const int stateCalc, for (size_t k = 0; k < m_numSpecies; k++) { if (m_speciesUnknownType[k] != VCS_SPECIES_TYPE_INTERFACIALVOLTAGE) { - kglob = IndSpecies[k]; + size_t kglob = IndSpecies[k]; v_totalMoles += std::max(0.0, molesSpeciesVCS[kglob]); } } if (v_totalMoles > 0.0) { for (size_t k = 0; k < m_numSpecies; k++) { if (m_speciesUnknownType[k] != VCS_SPECIES_TYPE_INTERFACIALVOLTAGE) { - kglob = IndSpecies[k]; - tmp = std::max(0.0, molesSpeciesVCS[kglob]); + size_t kglob = IndSpecies[k]; + double tmp = std::max(0.0, molesSpeciesVCS[kglob]); Xmol_[k] = tmp / v_totalMoles; } } @@ -551,7 +549,7 @@ void vcs_VolPhase::setMolesFromVCS(const int stateCalc, * in the equation system */ if (m_phiVarIndex != npos) { - kglob = IndSpecies[m_phiVarIndex]; + size_t kglob = IndSpecies[m_phiVarIndex]; if (m_numSpecies == 1) { Xmol_[m_phiVarIndex] = 1.0; } else { @@ -1200,9 +1198,6 @@ static bool chargeNeutralityElement(const Cantera::ThermoPhase* const tPhase) size_t vcs_VolPhase::transferElementsFM(const Cantera::ThermoPhase* const tPhase) { - size_t e, k, eT; - std::string ename; - size_t eFound = npos; size_t nebase = tPhase->nElements(); size_t ne = nebase; size_t ns = tPhase->nSpecies(); @@ -1227,6 +1222,7 @@ size_t vcs_VolPhase::transferElementsFM(const Cantera::ThermoPhase* const tPhase m_elementType[ChargeNeutralityElement] = VCS_ELEM_TYPE_CHARGENEUTRALITY; } + size_t eFound = npos; if (hasChargedSpecies(tPhase)) { if (cne) { /* @@ -1238,18 +1234,16 @@ size_t vcs_VolPhase::transferElementsFM(const Cantera::ThermoPhase* const tPhase * by toggling the ElActive variable. If we find we need it * later, we will retoggle ElActive to true. */ - for (eT = 0; eT < nebase; eT++) { - ename = tPhase->elementName(eT); - if (ename == "E") { + for (size_t eT = 0; eT < nebase; eT++) { + if (tPhase->elementName(eT) == "E") { eFound = eT; m_elementActive[eT] = 0; m_elementType[eT] = VCS_ELEM_TYPE_ELECTRONCHARGE; } } } else { - for (eT = 0; eT < nebase; eT++) { - ename = tPhase->elementName(eT); - if (ename == "E") { + for (size_t eT = 0; eT < nebase; eT++) { + if (tPhase->elementName(eT) == "E") { eFound = eT; m_elementType[eT] = VCS_ELEM_TYPE_ELECTRONCHARGE; } @@ -1273,10 +1267,9 @@ size_t vcs_VolPhase::transferElementsFM(const Cantera::ThermoPhase* const tPhase elemResize(ne); - e = 0; - for (eT = 0; eT < nebase; eT++) { - ename = tPhase->elementName(eT); - m_elementNames[e] = ename; + size_t e = 0; + for (size_t eT = 0; eT < nebase; eT++) { + m_elementNames[e] = tPhase->elementName(eT); m_elementType[e] = tPhase->elementType(eT); e++; } @@ -1288,15 +1281,14 @@ size_t vcs_VolPhase::transferElementsFM(const Cantera::ThermoPhase* const tPhase sss << "phase" << VP_ID_; pname = sss.str(); } - ename = "cn_" + pname; e = ChargeNeutralityElement; - m_elementNames[e] = ename; + m_elementNames[e] = "cn_" + pname; } double* const* const fm = m_formulaMatrix.baseDataAddr(); - for (k = 0; k < ns; k++) { + for (size_t k = 0; k < ns; k++) { e = 0; - for (eT = 0; eT < nebase; eT++) { + for (size_t eT = 0; eT < nebase; eT++) { fm[e][k] = tPhase->nAtoms(k, eT); e++; } @@ -1306,7 +1298,7 @@ size_t vcs_VolPhase::transferElementsFM(const Cantera::ThermoPhase* const tPhase } if (cne) { - for (k = 0; k < ns; k++) { + for (size_t k = 0; k < ns; k++) { fm[ChargeNeutralityElement][k] = tPhase->charge(k); } } diff --git a/src/equil/vcs_elem.cpp b/src/equil/vcs_elem.cpp index 961454b5a..2269fcf30 100644 --- a/src/equil/vcs_elem.cpp +++ b/src/equil/vcs_elem.cpp @@ -30,9 +30,6 @@ void VCS_SOLVE::vcs_elab() bool VCS_SOLVE::vcs_elabcheck(int ibound) { size_t top = m_numComponents; - double eval, scale; - int numNonZero; - bool multisign = false; if (ibound) { top = m_numElemConstraints; } @@ -51,23 +48,21 @@ bool VCS_SOLVE::vcs_elabcheck(int ibound) "Problem with charge neutrality condition"); } if (m_elemAbundancesGoal[i] == 0.0 || (m_elType[i] == VCS_ELEM_TYPE_ELECTRONCHARGE)) { - scale = VCS_DELETE_MINORSPECIES_CUTOFF; + double scale = VCS_DELETE_MINORSPECIES_CUTOFF; /* * Find out if the constraint is a multisign constraint. * If it is, then we have to worry about roundoff error * in the addition of terms. We are limited to 13 * digits of finite arithmetic accuracy. */ - numNonZero = 0; - multisign = false; + bool multisign = false; for (size_t kspec = 0; kspec < m_numSpeciesTot; kspec++) { - eval = m_formulaMatrix[i][kspec]; + double eval = m_formulaMatrix[i][kspec]; if (eval < 0.0) { multisign = true; } if (eval != 0.0) { scale = std::max(scale, fabs(eval * m_molNumSpecies_old[kspec])); - numNonZero++; } } if (multisign) { @@ -112,12 +107,9 @@ void VCS_SOLVE::vcs_elabPhase(size_t iphase, double* const elemAbundPhase) int VCS_SOLVE::vcs_elcorr(double aa[], double x[]) { - int retn = 0, its; - bool goodSpec; - double xx, par, saveDir, dir; + int retn = 0; #ifdef DEBUG_MODE - double l2before = 0.0, l2after = 0.0; std::vector ga_save(m_numElemConstraints, 0.0); vcs_dcopy(VCS_DATA_PTR(ga_save), VCS_DATA_PTR(m_elemAbundances), m_numElemConstraints); if (m_debug_print_lvl >= 2) { @@ -131,7 +123,7 @@ int VCS_SOLVE::vcs_elcorr(double aa[], double x[]) for (size_t i = 0; i < m_numElemConstraints; ++i) { x[i] = m_elemAbundances[i] - m_elemAbundancesGoal[i]; } - l2before = 0.0; + double l2before = 0.0; for (size_t i = 0; i < m_numElemConstraints; ++i) { l2before += x[i] * x[i]; } @@ -144,12 +136,10 @@ int VCS_SOLVE::vcs_elcorr(double aa[], double x[]) * formula matrix, and no other species have zero values either. * */ - int numNonZero = 0; bool changed = false; - bool multisign = false; for (size_t i = 0; i < m_numElemConstraints; ++i) { - numNonZero = 0; - multisign = false; + int numNonZero = 0; + bool multisign = false; for (size_t kspec = 0; kspec < m_numSpeciesTot; kspec++) { if (m_speciesUnknownType[kspec] != VCS_SPECIES_TYPE_INTERFACIALVOLTAGE) { double eval = m_formulaMatrix[i][kspec]; @@ -286,10 +276,10 @@ int VCS_SOLVE::vcs_elcorr(double aa[], double x[]) /* * Now apply the new direction without creating negative species. */ - par = 0.5; + double par = 0.5; for (size_t i = 0; i < m_numComponents; ++i) { if (m_molNumSpecies_old[i] > 0.0) { - xx = -x[i] / m_molNumSpecies_old[i]; + double xx = -x[i] / m_molNumSpecies_old[i]; if (par < xx) { par = xx; } @@ -350,10 +340,10 @@ int VCS_SOLVE::vcs_elcorr(double aa[], double x[]) if (m_speciesUnknownType[kspec] == VCS_SPECIES_TYPE_INTERFACIALVOLTAGE) { continue; } - saveDir = 0.0; - goodSpec = true; + double saveDir = 0.0; + bool goodSpec = true; for (size_t i = 0; i < m_numComponents; ++i) { - dir = m_formulaMatrix[i][kspec] * (m_elemAbundancesGoal[i] - m_elemAbundances[i]); + double dir = m_formulaMatrix[i][kspec] * (m_elemAbundancesGoal[i] - m_elemAbundances[i]); if (fabs(dir) > 1.0E-10) { if (dir > 0.0) { if (saveDir < 0.0) { @@ -375,8 +365,8 @@ int VCS_SOLVE::vcs_elcorr(double aa[], double x[]) } } if (goodSpec) { - its = 0; - xx = 0.0; + int its = 0; + double xx = 0.0; for (size_t i = 0; i < m_numComponents; ++i) { if (m_formulaMatrix[i][kspec] != 0.0) { xx += (m_elemAbundancesGoal[i] - m_elemAbundances[i]) / m_formulaMatrix[i][kspec]; @@ -500,15 +490,15 @@ L_CLEANUP: ; vcs_tmoles(); #ifdef DEBUG_MODE - l2after = 0.0; - for (int i = 0; i < m_numElemConstraints; ++i) { + double l2after = 0.0; + for (size_t i = 0; i < m_numElemConstraints; ++i) { l2after += SQUARE(m_elemAbundances[i] - m_elemAbundancesGoal[i]); } l2after = sqrt(l2after/m_numElemConstraints); if (m_debug_print_lvl >= 2) { plogf(" --- Elem_Abund: Correct Initial " " Final\n"); - for (int i = 0; i < m_numElemConstraints; ++i) { + for (size_t i = 0; i < m_numElemConstraints; ++i) { plogf(" --- "); plogf("%-2.2s", m_elementName[i].c_str()); plogf(" %20.12E %20.12E %20.12E\n", m_elemAbundancesGoal[i], ga_save[i], m_elemAbundances[i]); diff --git a/src/equil/vcs_elem_rearrange.cpp b/src/equil/vcs_elem_rearrange.cpp index a13dd9fd0..3ae62e4f6 100644 --- a/src/equil/vcs_elem_rearrange.cpp +++ b/src/equil/vcs_elem_rearrange.cpp @@ -19,14 +19,11 @@ namespace VCSnonideal int VCS_SOLVE::vcs_elem_rearrange(double* const aw, double* const sa, double* const sm, double* const ss) { - size_t j, k, l, i, jl, ml, jr, ielem; - bool lindep; size_t ncomponents = m_numComponents; - double test = -1.0E10; #ifdef DEBUG_MODE if (m_debug_print_lvl >= 2) { plogf(" "); - for (i=0; i<77; i++) { + for (size_t i=0; i<77; i++) { plogf("-"); } plogf("\n"); @@ -41,10 +38,11 @@ int VCS_SOLVE::vcs_elem_rearrange(double* const aw, double* const sa, * Use a temporary work array for the element numbers * Also make sure the value of test is unique. */ - lindep = false; + bool lindep = false; + double test = -1.0E10; do { lindep = false; - for (i = 0; i < m_numElemConstraints; ++i) { + for (size_t i = 0; i < m_numElemConstraints; ++i) { test -= 1.0; aw[i] = m_elemAbundancesGoal[i]; if (test == aw[i]) { @@ -57,9 +55,10 @@ int VCS_SOLVE::vcs_elem_rearrange(double* const aw, double* const sa, * Top of a loop of some sort based on the index JR. JR is the * current number independent elements found. */ - jr = npos; + size_t jr = npos; do { ++jr; + size_t k; /* * Top of another loop point based on finding a linearly * independent species @@ -70,7 +69,7 @@ int VCS_SOLVE::vcs_elem_rearrange(double* const aw, double* const sa, * for the largest remaining species. Return its identity in K. */ k = m_numElemConstraints; - for (ielem = jr; ielem < m_numElemConstraints; ielem++) { + for (size_t ielem = jr; ielem < m_numElemConstraints; ielem++) { if (m_elementActive[ielem]) { if (aw[ielem] != test) { k = ielem; @@ -98,13 +97,13 @@ int VCS_SOLVE::vcs_elem_rearrange(double* const aw, double* const sa, * Modified Gram-Schmidt Method, p. 202 Dalquist * QR factorization of a matrix without row pivoting. */ - jl = jr; + size_t jl = jr; /* * Fill in the row for the current element, k, under consideration * The row will contain the Formula matrix value for that element * from the current component. */ - for (j = 0; j < ncomponents; ++j) { + for (size_t j = 0; j < ncomponents; ++j) { sm[j + jr*ncomponents] = m_formulaMatrix[k][j]; } if (jl > 0) { @@ -114,9 +113,9 @@ int VCS_SOLVE::vcs_elem_rearrange(double* const aw, double* const sa, * (this is slightly different than Dalquist) * R_JA_JA = 1 */ - for (j = 0; j < jl; ++j) { + for (size_t j = 0; j < jl; ++j) { ss[j] = 0.0; - for (i = 0; i < ncomponents; ++i) { + for (size_t i = 0; i < ncomponents; ++i) { ss[j] += sm[i + jr*ncomponents] * sm[i + j*ncomponents]; } ss[j] /= sa[j]; @@ -125,8 +124,8 @@ int VCS_SOLVE::vcs_elem_rearrange(double* const aw, double* const sa, * Now make the new column, (*,JR), orthogonal to the * previous columns */ - for (j = 0; j < jl; ++j) { - for (l = 0; l < ncomponents; ++l) { + for (size_t j = 0; j < jl; ++j) { + for (size_t l = 0; l < ncomponents; ++l) { sm[l + jr*ncomponents] -= ss[j] * sm[l + j*ncomponents]; } } @@ -137,7 +136,7 @@ int VCS_SOLVE::vcs_elem_rearrange(double* const aw, double* const sa, * It will be used in the denominator in future row calcs. */ sa[jr] = 0.0; - for (ml = 0; ml < ncomponents; ++ml) { + for (size_t ml = 0; ml < ncomponents; ++ml) { sa[jr] += SQUARE(sm[ml + jr*ncomponents]); } /* **************************************************** */ @@ -181,8 +180,6 @@ void VCS_SOLVE::vcs_switch_elem_pos(size_t ipos, size_t jpos) if (ipos == jpos) { return; } - size_t j; - vcs_VolPhase* volPhase; #ifdef DEBUG_MODE if (ipos > (m_numElemConstraints - 1) || jpos > (m_numElemConstraints - 1)) { @@ -197,7 +194,7 @@ void VCS_SOLVE::vcs_switch_elem_pos(size_t ipos, size_t jpos) * to reflect the switch in the element positions. */ for (size_t iph = 0; iph < m_numPhases; iph++) { - volPhase = m_VolPhaseList[iph]; + vcs_VolPhase* volPhase = m_VolPhaseList[iph]; for (size_t e = 0; e < volPhase->nElemConstraints(); e++) { if (volPhase->elemGlobalIndex(e) == ipos) { volPhase->setElemGlobalIndex(e, jpos); @@ -212,7 +209,7 @@ void VCS_SOLVE::vcs_switch_elem_pos(size_t ipos, size_t jpos) std::swap(m_elementMapIndex[ipos], m_elementMapIndex[jpos]); std::swap(m_elType[ipos], m_elType[jpos]); std::swap(m_elementActive[ipos], m_elementActive[jpos]); - for (j = 0; j < m_numSpeciesTot; ++j) { + for (size_t j = 0; j < m_numSpeciesTot; ++j) { std::swap(m_formulaMatrix[ipos][j], m_formulaMatrix[jpos][j]); } std::swap(m_elementName[ipos], m_elementName[jpos]); diff --git a/src/equil/vcs_equilibrate.cpp b/src/equil/vcs_equilibrate.cpp index 44a980dc5..25085ac8f 100644 --- a/src/equil/vcs_equilibrate.cpp +++ b/src/equil/vcs_equilibrate.cpp @@ -36,7 +36,6 @@ int vcs_equilibrate(thermo_t& s, const char* XY, { MultiPhase* m = 0; int retn = 1; - int retnSub = 0; if (solver == 2) { m = new MultiPhase; @@ -78,7 +77,7 @@ int vcs_equilibrate(thermo_t& s, const char* XY, if (estimateEquil == 0) { useThermoPhaseElementPotentials = true; } - retnSub = e->equilibrate(s, XY, + int retnSub = e->equilibrate(s, XY, useThermoPhaseElementPotentials, loglevel-1); if (retnSub < 0) { delete e; diff --git a/src/equil/vcs_inest.cpp b/src/equil/vcs_inest.cpp index b210e6803..f2d8ee5c9 100644 --- a/src/equil/vcs_inest.cpp +++ b/src/equil/vcs_inest.cpp @@ -22,22 +22,10 @@ static char pprefix[20] = " --- vcs_inest: "; void VCS_SOLVE::vcs_inest(double* const aw, double* const sa, double* const sm, double* const ss, double test) { - size_t lt, ikl, kspec, iph, irxn; - double s; - double s1 = 0.0; - double xl, par; - bool finished, conv; size_t nspecies = m_numSpeciesTot; size_t nrxn = m_numRxnTot; // double *molNum = VCS_DATA_PTR(m_molNumSpecies_old); - double TMolesMultiphase; - double* xtphMax = VCS_DATA_PTR(m_TmpPhase); - double* xtphMin = VCS_DATA_PTR(m_TmpPhase2); - - ikl = 0; - lt = 0; - /* * CALL ROUTINE TO SOLVE MAX(CC*molNum) SUCH THAT AX*molNum = BB @@ -52,7 +40,7 @@ void VCS_SOLVE::vcs_inest(double* const aw, double* const sa, double* const sm, plogf("%s Mole Numbers returned from linear programming (vcs_inest initial guess):\n", pprefix); plogf("%s SPECIES MOLE_NUMBER -SS_ChemPotential\n", pprefix); - for (kspec = 0; kspec < nspecies; ++kspec) { + for (size_t kspec = 0; kspec < nspecies; ++kspec) { plogf("%s ", pprefix); plogf("%-12.12s", m_speciesName[kspec].c_str()); plogf(" %15.5g %12.3g\n", m_molNumSpecies_old[kspec], -m_SSfeSpecies[kspec]); @@ -62,31 +50,27 @@ void VCS_SOLVE::vcs_inest(double* const aw, double* const sa, double* const sm, pprefix); plogendl(); plogf("%s Element Goal Actual\n", pprefix); - int jj = 0; for (size_t j = 0; j < m_numElemConstraints; j++) { if (m_elementActive[j]) { double tmp = 0.0; - for (kspec = 0; kspec < nspecies; ++kspec) { + for (size_t kspec = 0; kspec < nspecies; ++kspec) { tmp += m_formulaMatrix[j][kspec] * m_molNumSpecies_old[kspec]; } plogf("%s ", pprefix); plogf(" %-9.9s", (m_elementName[j]).c_str()); plogf(" %12.3g %12.3g\n", m_elemAbundancesGoal[j], tmp); - jj++; } } plogendl(); } #endif - /* * Make sure all species have positive definite mole numbers * Set voltages to zero for now, until we figure out what to do */ vcs_dzero(VCS_DATA_PTR(m_deltaMolNumSpecies), nspecies); - for (kspec = 0; kspec < nspecies; ++kspec) { - iph = m_phaseID[kspec]; + for (size_t kspec = 0; kspec < nspecies; ++kspec) { if (m_speciesUnknownType[kspec] != VCS_SPECIES_TYPE_INTERFACIALVOLTAGE) { if (m_molNumSpecies_old[kspec] <= 0.0) { /* @@ -105,6 +89,7 @@ void VCS_SOLVE::vcs_inest(double* const aw, double* const sa, double* const sm, * Now find the optimized basis that spans the stoichiometric * coefficient matrix */ + bool conv; (void) vcs_basopt(false, aw, sa, sm, ss, test, &conv); /* ***************************************************************** */ @@ -118,22 +103,22 @@ void VCS_SOLVE::vcs_inest(double* const aw, double* const sa, double* const sm, /* * m_tPhaseMoles_new[] will consist of just the component moles */ - for (iph = 0; iph < m_numPhases; iph++) { + for (size_t iph = 0; iph < m_numPhases; iph++) { m_tPhaseMoles_new[iph] = TPhInertMoles[iph] + 1.0E-20; } - for (kspec = 0; kspec < m_numComponents; ++kspec) { + for (size_t kspec = 0; kspec < m_numComponents; ++kspec) { if (m_speciesUnknownType[kspec] == VCS_SPECIES_TYPE_MOLNUM) { m_tPhaseMoles_new[m_phaseID[kspec]] += m_molNumSpecies_old[kspec]; } } - TMolesMultiphase = 0.0; - for (iph = 0; iph < m_numPhases; iph++) { + double TMolesMultiphase = 0.0; + for (size_t iph = 0; iph < m_numPhases; iph++) { if (! m_VolPhaseList[iph]->m_singleSpecies) { TMolesMultiphase += m_tPhaseMoles_new[iph]; } } vcs_dcopy(VCS_DATA_PTR(m_molNumSpecies_new), VCS_DATA_PTR(m_molNumSpecies_old), nspecies); - for (kspec = 0; kspec < m_numComponents; ++kspec) { + for (size_t kspec = 0; kspec < m_numComponents; ++kspec) { if (m_speciesUnknownType[kspec] != VCS_SPECIES_TYPE_MOLNUM) { m_molNumSpecies_new[kspec] = 0.0; } @@ -141,10 +126,10 @@ void VCS_SOLVE::vcs_inest(double* const aw, double* const sa, double* const sm, vcs_dcopy(VCS_DATA_PTR(m_feSpecies_new), VCS_DATA_PTR(m_SSfeSpecies), nspecies); - for (kspec = 0; kspec < m_numComponents; ++kspec) { + for (size_t kspec = 0; kspec < m_numComponents; ++kspec) { if (m_speciesUnknownType[kspec] == VCS_SPECIES_TYPE_MOLNUM) { if (! m_SSPhase[kspec]) { - iph = m_phaseID[kspec]; + size_t iph = m_phaseID[kspec]; m_feSpecies_new[kspec] += log(m_molNumSpecies_new[kspec] / m_tPhaseMoles_old[iph]); } } else { @@ -154,7 +139,7 @@ void VCS_SOLVE::vcs_inest(double* const aw, double* const sa, double* const sm, vcs_deltag(0, true, VCS_STATECALC_NEW); #ifdef DEBUG_MODE if (m_debug_print_lvl >= 2) { - for (kspec = 0; kspec < nspecies; ++kspec) { + for (size_t kspec = 0; kspec < nspecies; ++kspec) { plogf("%s", pprefix); plogf("%-12.12s", m_speciesName[kspec].c_str()); if (kspec < m_numComponents) @@ -169,13 +154,15 @@ void VCS_SOLVE::vcs_inest(double* const aw, double* const sa, double* const sm, /* ********************************************************** */ /* **** ESTIMATE REACTION ADJUSTMENTS *********************** */ /* ********************************************************** */ + double* xtphMax = VCS_DATA_PTR(m_TmpPhase); + double* xtphMin = VCS_DATA_PTR(m_TmpPhase2); vcs_dzero(VCS_DATA_PTR(m_deltaPhaseMoles), m_numPhases); - for (iph = 0; iph < m_numPhases; iph++) { + for (size_t iph = 0; iph < m_numPhases; iph++) { xtphMax[iph] = log(m_tPhaseMoles_new[iph] * 1.0E32); xtphMin[iph] = log(m_tPhaseMoles_new[iph] * 1.0E-32); } - for (irxn = 0; irxn < nrxn; ++irxn) { - kspec = m_indexRxnToSpecies[irxn]; + for (size_t irxn = 0; irxn < nrxn; ++irxn) { + size_t kspec = m_indexRxnToSpecies[irxn]; /* * For single species phases, we will not estimate the * mole numbers. If the phase exists, it stays. If it @@ -183,7 +170,7 @@ void VCS_SOLVE::vcs_inest(double* const aw, double* const sa, double* const sm, * existence here. */ if (! m_SSPhase[kspec]) { - iph = m_phaseID[kspec]; + size_t iph = m_phaseID[kspec]; if (m_deltaGRxn_new[irxn] > xtphMax[iph]) { m_deltaGRxn_new[irxn] = 0.8 * xtphMax[iph]; } @@ -213,7 +200,7 @@ void VCS_SOLVE::vcs_inest(double* const aw, double* const sa, double* const sm, } #ifdef DEBUG_MODE if (m_debug_print_lvl >= 2) { - for (kspec = 0; kspec < nspecies; ++kspec) { + for (size_t kspec = 0; kspec < nspecies; ++kspec) { if (m_speciesUnknownType[kspec] != VCS_SPECIES_TYPE_INTERFACIALVOLTAGE) { plogf("%sdirection (", pprefix); plogf("%-12.12s", m_speciesName[kspec].c_str()); @@ -233,8 +220,8 @@ void VCS_SOLVE::vcs_inest(double* const aw, double* const sa, double* const sm, /* *********************************************************** */ /* **** KEEP COMPONENT SPECIES POSITIVE ********************** */ /* *********************************************************** */ - par = 0.5; - for (kspec = 0; kspec < m_numComponents; ++kspec) { + double par = 0.5; + for (size_t kspec = 0; kspec < m_numComponents; ++kspec) { if (m_speciesUnknownType[kspec] != VCS_SPECIES_TYPE_INTERFACIALVOLTAGE) { if (par < -m_deltaMolNumSpecies[kspec] / m_molNumSpecies_new[kspec]) { par = -m_deltaMolNumSpecies[kspec] / m_molNumSpecies_new[kspec]; @@ -250,16 +237,19 @@ void VCS_SOLVE::vcs_inest(double* const aw, double* const sa, double* const sm, /* ******************************************** */ /* **** CALCULATE NEW MOLE NUMBERS ************ */ /* ******************************************** */ - finished = false; + bool finished = false; + size_t lt = 0; + size_t ikl = 0; + double s1 = 0.0; do { - for (kspec = 0; kspec < m_numComponents; ++kspec) { + for (size_t kspec = 0; kspec < m_numComponents; ++kspec) { if (m_speciesUnknownType[kspec] != VCS_SPECIES_TYPE_INTERFACIALVOLTAGE) { m_molNumSpecies_old[kspec] = m_molNumSpecies_new[kspec] + par * m_deltaMolNumSpecies[kspec]; } else { m_deltaMolNumSpecies[kspec] = 0.0; } } - for (kspec = m_numComponents; kspec < nspecies; ++kspec) { + for (size_t kspec = m_numComponents; kspec < nspecies; ++kspec) { if (m_speciesUnknownType[kspec] != VCS_SPECIES_TYPE_INTERFACIALVOLTAGE) { if (m_deltaMolNumSpecies[kspec] != 0.0) { m_molNumSpecies_old[kspec] = m_deltaMolNumSpecies[kspec] * par; @@ -279,7 +269,8 @@ void VCS_SOLVE::vcs_inest(double* const aw, double* const sa, double* const sm, /* ******************************************* */ vcs_setFlagsVolPhases(false, VCS_STATECALC_OLD); vcs_dfe(VCS_STATECALC_OLD, 0, 0, nspecies); - for (kspec = 0, s = 0.0; kspec < nspecies; ++kspec) { + double s = 0.0; + for (size_t kspec = 0; kspec < nspecies; ++kspec) { s += m_deltaMolNumSpecies[kspec] * m_feSpecies_old[kspec]; } if (s == 0.0) { @@ -304,7 +295,7 @@ void VCS_SOLVE::vcs_inest(double* const aw, double* const sa, double* const sm, /* **************************************************** */ /* **** FIT PARABOLA THROUGH HALF AND FULL STEPS ****** */ /* **************************************************** */ - xl = (1.0 - s / (s1 - s)) * 0.5; + double xl = (1.0 - s / (s1 - s)) * 0.5; if (xl < 0.0) { /* *************************************************** */ /* *** POOR DIRECTION, REDUCE STEP SIZE TO 0.2 ******* */ @@ -332,7 +323,7 @@ finished: plogf("%s Final Mole Numbers produced by inest:\n", pprefix); plogf("%s SPECIES MOLE_NUMBER\n", pprefix); - for (kspec = 0; kspec < nspecies; ++kspec) { + for (size_t kspec = 0; kspec < nspecies; ++kspec) { plogf("%s ", pprefix); plogf("%-12.12s", m_speciesName[kspec].c_str()); plogf(" %g", m_molNumSpecies_old[kspec]); @@ -345,9 +336,7 @@ finished: int VCS_SOLVE::vcs_inest_TP() { int retn = 0; - double test; Cantera::clockWC tickTock; - test = -1.0E20; if (m_doEstimateEquil > 0) { /* @@ -398,6 +387,7 @@ int VCS_SOLVE::vcs_inest_TP() plogendl(); } #endif + double test = -1.0E20; vcs_inest(VCS_DATA_PTR(aw), VCS_DATA_PTR(sa), VCS_DATA_PTR(sm), VCS_DATA_PTR(ss), test); /* @@ -472,8 +462,7 @@ int VCS_SOLVE::vcs_inest_TP() /* * Record time */ - double tsecond = tickTock.secondsWC(); - m_VCount->T_Time_inest += tsecond; + m_VCount->T_Time_inest += tickTock.secondsWC(); (m_VCount->T_Calls_Inest)++; return retn; } diff --git a/src/equil/vcs_nondim.cpp b/src/equil/vcs_nondim.cpp index fde4c7405..7761a9ff9 100644 --- a/src/equil/vcs_nondim.cpp +++ b/src/equil/vcs_nondim.cpp @@ -18,7 +18,6 @@ namespace VCSnonideal { double VCS_SOLVE::vcs_nondim_Farad(int mu_units, double TKelvin) const { - double Farad; if (TKelvin <= 0.0) { TKelvin = 293.15; } @@ -26,59 +25,47 @@ double VCS_SOLVE::vcs_nondim_Farad(int mu_units, double TKelvin) const case VCS_UNITS_MKS: case VCS_UNITS_KJMOL: case VCS_UNITS_KCALMOL: - Farad = Cantera::ElectronCharge * Cantera::Avogadro / + return Cantera::ElectronCharge * Cantera::Avogadro / (TKelvin * Cantera::GasConstant); - break; case VCS_UNITS_UNITLESS: - Farad = Cantera::ElectronCharge * Cantera::Avogadro; - break; + return Cantera::ElectronCharge * Cantera::Avogadro; case VCS_UNITS_KELVIN: - Farad = Cantera::ElectronCharge * Cantera::Avogadro/ TKelvin; - break; + return Cantera::ElectronCharge * Cantera::Avogadro/ TKelvin; default: plogf("vcs_nondim_Farad error: unknown units: %d\n", mu_units); plogendl(); exit(EXIT_FAILURE); } - return Farad; } double VCS_SOLVE::vcs_nondimMult_TP(int mu_units, double TKelvin) const { - double rt; if (TKelvin <= 0.0) { TKelvin = 293.15; } switch (mu_units) { case VCS_UNITS_KCALMOL: - rt = TKelvin * Cantera::GasConst_cal_mol_K * 1e-3; - break; + return TKelvin * Cantera::GasConst_cal_mol_K * 1e-3; case VCS_UNITS_UNITLESS: - rt = 1.0; - break; + return 1.0; case VCS_UNITS_KJMOL: - rt = TKelvin * Cantera::GasConstant * 1e-6; - break; + return TKelvin * Cantera::GasConstant * 1e-6; case VCS_UNITS_KELVIN: - rt = TKelvin; - break; + return TKelvin; case VCS_UNITS_MKS: - rt = TKelvin * Cantera::GasConstant; - break; + return TKelvin * Cantera::GasConstant; default: plogf("vcs_nondimMult_TP error: unknown units: %d\n", mu_units); plogendl(); exit(EXIT_FAILURE); } - return rt; } void VCS_SOLVE::vcs_nondim_TP() { - double tf; if (m_unitsState == VCS_DIMENSIONAL_G) { m_unitsState = VCS_NONDIMENSIONAL_G; - tf = 1.0 / vcs_nondimMult_TP(m_VCS_UnitsFormat, m_temperature); + double tf = 1.0 / vcs_nondimMult_TP(m_VCS_UnitsFormat, m_temperature); for (size_t i = 0; i < m_numSpeciesTot; ++i) { /* * Modify the standard state and total chemical potential data, @@ -166,10 +153,9 @@ void VCS_SOLVE::vcs_nondim_TP() void VCS_SOLVE::vcs_redim_TP(void) { - double tf; if (m_unitsState != VCS_DIMENSIONAL_G) { m_unitsState = VCS_DIMENSIONAL_G; - tf = vcs_nondimMult_TP(m_VCS_UnitsFormat, m_temperature); + double tf = vcs_nondimMult_TP(m_VCS_UnitsFormat, m_temperature); for (size_t i = 0; i < m_numSpeciesTot; ++i) { /* * Modify the standard state and total chemical potential data, diff --git a/src/equil/vcs_phaseStability.cpp b/src/equil/vcs_phaseStability.cpp index 1889c415b..8e7c3de46 100644 --- a/src/equil/vcs_phaseStability.cpp +++ b/src/equil/vcs_phaseStability.cpp @@ -139,10 +139,6 @@ bool VCS_SOLVE::vcs_popPhasePossible(const size_t iphasePop) const int VCS_SOLVE::vcs_phasePopDeterminePossibleList() { int nfound = 0; - vcs_VolPhase* Vphase = 0; - double stoicC; - double molComp; - std::vector linkedPhases; phasePopProblemLists_.clear(); /* @@ -161,19 +157,17 @@ int VCS_SOLVE::vcs_phasePopDeterminePossibleList() */ for (size_t j = 0; j < m_numComponents; j++) { if (m_elType[j] == VCS_ELEM_TYPE_ABSPOS) { - molComp = m_molNumSpecies_old[j]; - if (molComp <= 0.0) { + if (m_molNumSpecies_old[j] <= 0.0) { std::vector &jList = zeroedComponentLinkedPhasePops[j]; size_t iph = m_phaseID[j]; jList.push_back(iph); for (size_t irxn = 0; irxn < m_numRxnTot; irxn++) { size_t kspec = irxn + m_numComponents; iph = m_phaseID[kspec]; - Vphase = m_VolPhaseList[iph]; + vcs_VolPhase* Vphase = m_VolPhaseList[iph]; int existence = Vphase->exists(); if (existence < 0) { - stoicC = m_stoichCoeffRxnMatrix[irxn][j]; - if (stoicC > 0.0) { + if (m_stoichCoeffRxnMatrix[irxn][j] > 0.0) { if (std::find(jList.begin(), jList.end(), iph) != jList.end()) { jList.push_back(iph); } @@ -190,15 +184,15 @@ int VCS_SOLVE::vcs_phasePopDeterminePossibleList() * Cut out components which have a pos stoichiometric value with another species in the phase. */ std::vector< std::vector > zeroedPhaseLinkedZeroComponents(m_numPhases); + std::vector linkedPhases; /* * The logic below calculates zeroedPhaseLinkedZeroComponents */ for (size_t iph = 0; iph < m_numPhases; iph++) { std::vector &iphList = zeroedPhaseLinkedZeroComponents[iph]; iphList.clear(); - Vphase = m_VolPhaseList[iph]; - int existence = Vphase->exists(); - if (existence < 0) { + vcs_VolPhase* Vphase = m_VolPhaseList[iph]; + if (Vphase->exists() < 0) { linkedPhases.clear(); size_t nsp = Vphase->nSpecies(); @@ -208,10 +202,8 @@ int VCS_SOLVE::vcs_phasePopDeterminePossibleList() for (size_t j = 0; j < m_numComponents; j++) { if (m_elType[j] == VCS_ELEM_TYPE_ABSPOS) { - molComp = m_molNumSpecies_old[j]; - if (molComp <= 0.0) { - stoicC = m_stoichCoeffRxnMatrix[irxn][j]; - if (stoicC < 0.0) { + if (m_molNumSpecies_old[j] <= 0.0) { + if (m_stoichCoeffRxnMatrix[irxn][j] < 0.0) { bool foundPos = false; for (size_t kk = 0; kk < nsp; kk++) { size_t kkspec = Vphase->spGlobalIndexVCS(kk); @@ -240,9 +232,8 @@ int VCS_SOLVE::vcs_phasePopDeterminePossibleList() * */ for (size_t iph = 0; iph < m_numPhases; iph++) { - Vphase = m_VolPhaseList[iph]; - int existence = Vphase->exists(); - if (existence < 0) { + vcs_VolPhase* Vphase = m_VolPhaseList[iph]; + if (Vphase->exists() < 0) { std::vector &iphList = zeroedPhaseLinkedZeroComponents[iph]; std::vector popProblem(0); popProblem.push_back(iph); @@ -266,11 +257,8 @@ int VCS_SOLVE::vcs_phasePopDeterminePossibleList() size_t VCS_SOLVE::vcs_popPhaseID(std::vector & phasePopPhaseIDs) { size_t iphasePop = npos; - size_t irxn, kspec; doublereal FephaseMax = -1.0E30; doublereal Fephase = -1.0E30; - vcs_VolPhase* Vphase = 0; - #ifdef DEBUG_MODE char anote[128]; @@ -281,7 +269,7 @@ size_t VCS_SOLVE::vcs_popPhaseID(std::vector & phasePopPhaseIDs) } #endif for (size_t iph = 0; iph < m_numPhases; iph++) { - Vphase = m_VolPhaseList[iph]; + vcs_VolPhase* Vphase = m_VolPhaseList[iph]; int existence = Vphase->exists(); #ifdef DEBUG_MODE strcpy(anote, ""); @@ -303,8 +291,8 @@ size_t VCS_SOLVE::vcs_popPhaseID(std::vector & phasePopPhaseIDs) * Single Phase Stability Resolution * ***********************************************************************/ - kspec = Vphase->spGlobalIndexVCS(0); - irxn = kspec - m_numComponents; + size_t kspec = Vphase->spGlobalIndexVCS(0); + size_t irxn = kspec - m_numComponents; doublereal deltaGRxn = m_deltaGRxn_old[irxn]; Fephase = exp(-deltaGRxn) - 1.0; if (Fephase > 0.0) { @@ -403,16 +391,13 @@ int VCS_SOLVE::vcs_popPhaseRxnStepSizes(const size_t iphasePop) size_t irxn = kspec - m_numComponents; std::vector creationGlobalRxnNumbers; - doublereal s; // Calculate the initial moles of the phase being born. // Here we set it to 10x of the value which would cause the phase to be // zeroed out within the algorithm. We may later adjust the value. doublereal tPhaseMoles = 10. * m_totalMolNum * VCS_DELETE_PHASE_CUTOFF; - #ifdef DEBUG_MODE - int existence = Vphase->exists(); - if (existence > 0) { + if (Vphase->exists() > 0) { printf("ERROR vcs_popPhaseRxnStepSizes called for a phase that exists!"); exit(-1); } @@ -424,7 +409,7 @@ int VCS_SOLVE::vcs_popPhaseRxnStepSizes(const size_t iphasePop) #endif // Section for a single-species phase if (Vphase->m_singleSpecies) { - s = 0.0; + double s = 0.0; double* dnPhase_irxn = m_deltaMolNumPhase[irxn]; for (size_t j = 0; j < m_numComponents; ++j) { if (!m_SSPhase[j]) { @@ -595,10 +580,8 @@ double VCS_SOLVE::vcs_phaseStabilityTest(const size_t iph) /* * We will use the _new state calc here */ - size_t kspec, irxn, k, i, kc, kc_spec; vcs_VolPhase* Vphase = m_VolPhaseList[iph]; const size_t nsp = Vphase->nSpecies(); - doublereal deltaGRxn; int minNumberIterations = 3; if (nsp <= 1) { minNumberIterations = 1; @@ -635,14 +618,14 @@ double VCS_SOLVE::vcs_phaseStabilityTest(const size_t iph) std::vector componentList; - for (k = 0; k < nsp; k++) { - kspec = Vphase->spGlobalIndexVCS(k); + for (size_t k = 0; k < nsp; k++) { + size_t kspec = Vphase->spGlobalIndexVCS(k); if (kspec < m_numComponents) { componentList.push_back(k); } } - for (k = 0; k < m_numComponents; k++) { + for (size_t k = 0; k < m_numComponents; k++) { m_feSpecies_Deficient[k] = m_feSpecies_old[k]; } normUpdate = 0.1 * vcs_l2norm(fracDelta_new); @@ -663,7 +646,7 @@ double VCS_SOLVE::vcs_phaseStabilityTest(const size_t iph) } #endif - for (k = 0; k < nsp; k++) { + for (size_t k = 0; k < nsp; k++) { if (fracDelta_new[k] < 1.0E-13) { fracDelta_new[k] = 1.0E-13; } @@ -678,14 +661,14 @@ double VCS_SOLVE::vcs_phaseStabilityTest(const size_t iph) // Given a set of fracDelta's, we calculate the fracDelta's // for the component species, if any - for (i = 0; i < componentList.size(); i++) { - kc = componentList[i]; - kc_spec = Vphase->spGlobalIndexVCS(kc); + for (size_t i = 0; i < componentList.size(); i++) { + size_t kc = componentList[i]; + size_t kc_spec = Vphase->spGlobalIndexVCS(kc); fracDelta_old[kc] = 0.0; - for (k = 0; k < nsp; k++) { - kspec = Vphase->spGlobalIndexVCS(k); + for (size_t k = 0; k < nsp; k++) { + size_t kspec = Vphase->spGlobalIndexVCS(k); if (kspec >= m_numComponents) { - irxn = kspec - m_numComponents; + size_t irxn = kspec - m_numComponents; fracDelta_old[kc] += m_stoichCoeffRxnMatrix[irxn][kc_spec] * fracDelta_old[k]; } } @@ -693,7 +676,7 @@ double VCS_SOLVE::vcs_phaseStabilityTest(const size_t iph) // Now, calculate the predicted mole fractions, X_est[k] double sumFrac = 0.0; - for (k = 0; k < nsp; k++) { + for (size_t k = 0; k < nsp; k++) { sumFrac += fracDelta_old[k]; } // Necessary because this can be identically zero. -> we need to fix this algorithm! @@ -701,10 +684,9 @@ double VCS_SOLVE::vcs_phaseStabilityTest(const size_t iph) sumFrac = 1.0; } double sum_Xcomp = 0.0; - for (k = 0; k < nsp; k++) { + for (size_t k = 0; k < nsp; k++) { X_est[k] = fracDelta_old[k] / sumFrac; - kc_spec = Vphase->spGlobalIndexVCS(k); - if (kc_spec < m_numComponents) { + if (Vphase->spGlobalIndexVCS(k) < m_numComponents) { sum_Xcomp += X_est[k]; } } @@ -725,9 +707,9 @@ double VCS_SOLVE::vcs_phaseStabilityTest(const size_t iph) * First calculate altered chemical potentials for component species * belonging to this phase. */ - for (i = 0; i < componentList.size(); i++) { - kc = componentList[i]; - kc_spec = Vphase->spGlobalIndexVCS(kc); + for (size_t i = 0; i < componentList.size(); i++) { + size_t kc = componentList[i]; + size_t kc_spec = Vphase->spGlobalIndexVCS(kc); if (X_est[kc] > VCS_DELETE_MINORSPECIES_CUTOFF) { m_feSpecies_Deficient[kc_spec] = m_feSpecies_old[kc_spec] + log(m_actCoeffSpecies_new[kc_spec] * X_est[kc]); @@ -737,14 +719,12 @@ double VCS_SOLVE::vcs_phaseStabilityTest(const size_t iph) } } - for (i = 0; i < componentList.size(); i++) { - kc = componentList[i]; - kc_spec = Vphase->spGlobalIndexVCS(kc); - - for (k = 0; k < Vphase->nSpecies(); k++) { - kspec = Vphase->spGlobalIndexVCS(k); + for (size_t i = 0; i < componentList.size(); i++) { + size_t kc_spec = Vphase->spGlobalIndexVCS(componentList[i]); + for (size_t k = 0; k < Vphase->nSpecies(); k++) { + size_t kspec = Vphase->spGlobalIndexVCS(k); if (kspec >= m_numComponents) { - irxn = kspec - m_numComponents; + size_t irxn = kspec - m_numComponents; if (i == 0) { m_deltaGRxn_Deficient[irxn] = m_deltaGRxn_old[irxn]; } @@ -763,11 +743,11 @@ double VCS_SOLVE::vcs_phaseStabilityTest(const size_t iph) */ sum = 0.0; funcPhaseStability = sum_Xcomp - 1.0; - for (k = 0; k < nsp; k++) { - kspec = Vphase->spGlobalIndexVCS(k); + for (size_t k = 0; k < nsp; k++) { + size_t kspec = Vphase->spGlobalIndexVCS(k); if (kspec >= m_numComponents) { - irxn = kspec - m_numComponents; - deltaGRxn = m_deltaGRxn_Deficient[irxn]; + size_t irxn = kspec - m_numComponents; + double deltaGRxn = m_deltaGRxn_Deficient[irxn]; if (deltaGRxn > 50.0) { deltaGRxn = 50.0; } @@ -785,11 +765,10 @@ double VCS_SOLVE::vcs_phaseStabilityTest(const size_t iph) /* * Calculate the raw estimate of the new fracs */ - for (k = 0; k < nsp; k++) { - kspec = Vphase->spGlobalIndexVCS(k); + for (size_t k = 0; k < nsp; k++) { + size_t kspec = Vphase->spGlobalIndexVCS(k); double b = E_phi[k] / sum * (1.0 - sum_Xcomp); if (kspec >= m_numComponents) { - irxn = kspec - m_numComponents; fracDelta_raw[k] = b; } } @@ -797,14 +776,14 @@ double VCS_SOLVE::vcs_phaseStabilityTest(const size_t iph) // Given a set of fracDelta's, we calculate the fracDelta's // for the component species, if any - for (i = 0; i < componentList.size(); i++) { - kc = componentList[i]; - kc_spec = Vphase->spGlobalIndexVCS(kc); + for (size_t i = 0; i < componentList.size(); i++) { + size_t kc = componentList[i]; + size_t kc_spec = Vphase->spGlobalIndexVCS(kc); fracDelta_raw[kc] = 0.0; - for (k = 0; k < nsp; k++) { - kspec = Vphase->spGlobalIndexVCS(k); + for (size_t k = 0; k < nsp; k++) { + size_t kspec = Vphase->spGlobalIndexVCS(k); if (kspec >= m_numComponents) { - irxn = kspec - m_numComponents; + size_t irxn = kspec - m_numComponents; fracDelta_raw[kc] += m_stoichCoeffRxnMatrix[irxn][kc_spec] * fracDelta_raw[k]; } } @@ -816,14 +795,14 @@ double VCS_SOLVE::vcs_phaseStabilityTest(const size_t iph) * Now possibly dampen the estimate. */ doublereal sumADel = 0.0; - for (k = 0; k < nsp; k++) { + for (size_t k = 0; k < nsp; k++) { delFrac[k] = fracDelta_raw[k] - fracDelta_old[k]; sumADel += fabs(delFrac[k]); } normUpdate = vcs_l2norm(delFrac); dirProd = 0.0; - for (k = 0; k < nsp; k++) { + for (size_t k = 0; k < nsp; k++) { dirProd += fracDelta_old[k] * delFrac[k]; } bool crossedSign = false; @@ -850,7 +829,7 @@ double VCS_SOLVE::vcs_phaseStabilityTest(const size_t iph) } } - for (k = 0; k < nsp; k++) { + for (size_t k = 0; k < nsp; k++) { if (fabs(damp * delFrac[k]) > 0.3*fabs(fracDelta_old[k])) { damp = std::max(0.3*fabs(fracDelta_old[k]) / fabs(delFrac[k]), 1.0E-8/fabs(delFrac[k])); @@ -870,7 +849,7 @@ double VCS_SOLVE::vcs_phaseStabilityTest(const size_t iph) damp = 0.000001; } - for (k = 0; k < nsp; k++) { + for (size_t k = 0; k < nsp; k++) { fracDelta_new[k] = fracDelta_old[k] + damp * (delFrac[k]); } diff --git a/src/equil/vcs_prep.cpp b/src/equil/vcs_prep.cpp index 89b61a8ea..dd6493c17 100644 --- a/src/equil/vcs_prep.cpp +++ b/src/equil/vcs_prep.cpp @@ -19,9 +19,6 @@ namespace VCSnonideal { void VCS_SOLVE::vcs_SSPhase() { - size_t iph; - vcs_VolPhase* Vphase; - std::vector numPhSpecies(m_numPhases, 0); for (size_t kspec = 0; kspec < m_numSpeciesTot; ++kspec) { @@ -32,8 +29,8 @@ void VCS_SOLVE::vcs_SSPhase() * has been earmarked as a multispecies phase. * Treat that species as a single-species phase */ - for (iph = 0; iph < m_numPhases; iph++) { - Vphase = m_VolPhaseList[iph]; + for (size_t iph = 0; iph < m_numPhases; iph++) { + vcs_VolPhase* Vphase = m_VolPhaseList[iph]; Vphase->m_singleSpecies = false; if (TPhInertMoles[iph] > 0.0) { Vphase->setExistence(2); @@ -52,8 +49,8 @@ void VCS_SOLVE::vcs_SSPhase() * single species phase or not. */ for (size_t kspec = 0; kspec < m_numSpeciesTot; kspec++) { - iph = m_phaseID[kspec]; - Vphase = m_VolPhaseList[iph]; + size_t iph = m_phaseID[kspec]; + vcs_VolPhase* Vphase = m_VolPhaseList[iph]; if (Vphase->m_singleSpecies) { m_SSPhase[kspec] = true; } else { @@ -64,13 +61,7 @@ void VCS_SOLVE::vcs_SSPhase() int VCS_SOLVE::vcs_prep_oneTime(int printLvl) { - size_t kspec, i; int retn = VCS_SUCCESS; - double pres, test; - double* aw, *sa, *sm, *ss; - bool modifiedSoln = false; - bool conv; - m_debug_print_lvl = printLvl; /* @@ -90,11 +81,11 @@ int VCS_SOLVE::vcs_prep_oneTime(int printLvl) } m_numRxnRdc = m_numRxnTot; m_numSpeciesRdc = m_numSpeciesTot; - for (i = 0; i < m_numRxnRdc; ++i) { + for (size_t i = 0; i < m_numRxnRdc; ++i) { m_indexRxnToSpecies[i] = m_numElemConstraints + i; } - for (kspec = 0; kspec < m_numSpeciesTot; ++kspec) { + for (size_t kspec = 0; kspec < m_numSpeciesTot; ++kspec) { size_t pID = m_phaseID[kspec]; size_t spPhIndex = m_speciesLocalPhaseIndex[kspec]; vcs_VolPhase* vPhase = m_VolPhaseList[pID]; @@ -130,23 +121,20 @@ int VCS_SOLVE::vcs_prep_oneTime(int printLvl) * * For voltage unknowns, set these to zero for the moment. */ - test = -1.0e-10; + double test = -1.0e-10; + bool modifiedSoln = false; if (m_doEstimateEquil < 0) { double sum = 0.0; - for (kspec = 0; kspec < m_numSpeciesTot; ++kspec) { + for (size_t kspec = 0; kspec < m_numSpeciesTot; ++kspec) { if (m_speciesUnknownType[kspec] == VCS_SPECIES_TYPE_MOLNUM) { sum += fabs(m_molNumSpecies_old[kspec]); } } if (fabs(sum) < 1.0E-6) { modifiedSoln = true; - if (m_pressurePA <= 0.0) { - pres = 1.01325E5; - } else { - pres = m_pressurePA; - } + double pres = (m_pressurePA <= 0.0) ? 1.01325E5 : m_pressurePA; retn = vcs_evalSS_TP(0, 0, m_temperature, pres); - for (kspec = 0; kspec < m_numSpeciesTot; ++kspec) { + for (size_t kspec = 0; kspec < m_numSpeciesTot; ++kspec) { if (m_speciesUnknownType[kspec] == VCS_SPECIES_TYPE_MOLNUM) { m_molNumSpecies_old[kspec] = - m_SSfeSpecies[kspec]; } else { @@ -163,14 +151,15 @@ int VCS_SOLVE::vcs_prep_oneTime(int printLvl) * reaction matrix. */ std::vector awSpace(m_numSpeciesTot + (m_numElemConstraints + 2)*(m_numElemConstraints), 0.0); - aw = VCS_DATA_PTR(awSpace); + double* aw = VCS_DATA_PTR(awSpace); if (aw == NULL) { plogf("vcs_prep_oneTime: failed to get memory: global bailout\n"); return VCS_NOMEMORY; } - sa = aw + m_numSpeciesTot; - sm = sa + m_numElemConstraints; - ss = sm + (m_numElemConstraints)*(m_numElemConstraints); + double* sa = aw + m_numSpeciesTot; + double* sm = sa + m_numElemConstraints; + double* ss = sm + (m_numElemConstraints)*(m_numElemConstraints); + bool conv; retn = vcs_basopt(true, aw, sa, sm, ss, test, &conv); if (retn != VCS_SUCCESS) { plogf("vcs_prep_oneTime:"); @@ -182,7 +171,7 @@ int VCS_SOLVE::vcs_prep_oneTime(int printLvl) if (m_numSpeciesTot >= m_numComponents) { m_numRxnTot = m_numRxnRdc = m_numSpeciesTot - m_numComponents; - for (i = 0; i < m_numRxnRdc; ++i) { + for (size_t i = 0; i < m_numRxnRdc; ++i) { m_indexRxnToSpecies[i] = m_numComponents + i; } } else { @@ -209,7 +198,7 @@ int VCS_SOLVE::vcs_prep_oneTime(int printLvl) // If we mucked up the solution unknowns because they were all // zero to start with, set them back to zero here if (modifiedSoln) { - for (kspec = 0; kspec < m_numSpeciesTot; ++kspec) { + for (size_t kspec = 0; kspec < m_numSpeciesTot; ++kspec) { m_molNumSpecies_old[kspec] = 0.0; } } diff --git a/src/equil/vcs_prob.cpp b/src/equil/vcs_prob.cpp index 005a04381..520b200ac 100644 --- a/src/equil/vcs_prob.cpp +++ b/src/equil/vcs_prob.cpp @@ -184,7 +184,6 @@ static void print_char(const char letter, const int num) void VCS_PROB::prob_report(int print_lvl) { m_printLvl = print_lvl; - vcs_VolPhase* Vphase = 0; /* * Printout the species information: PhaseID's and mole nums */ @@ -215,7 +214,7 @@ void VCS_PROB::prob_report(int print_lvl) plogf(" species phaseID phaseName "); plogf(" Initial_Estimated_Moles Species_Type\n"); for (size_t i = 0; i < nspecies; i++) { - Vphase = VPhaseList[PhaseID[i]]; + vcs_VolPhase* Vphase = VPhaseList[PhaseID[i]]; plogf("%16s %5d %16s", SpName[i].c_str(), PhaseID[i], Vphase->PhaseName.c_str()); if (iest >= 0) { @@ -245,7 +244,7 @@ void VCS_PROB::prob_report(int print_lvl) plogf(" TMolesInert TKmoles\n"); for (size_t iphase = 0; iphase < NPhase; iphase++) { - Vphase = VPhaseList[iphase]; + vcs_VolPhase* Vphase = VPhaseList[iphase]; std::string EOS_cstr = string16_EOSType(Vphase->m_eqnState); plogf("%16s %5d %5d %8d ", Vphase->PhaseName.c_str(), Vphase->VP_ID_, Vphase->m_singleSpecies, Vphase->m_gasPhase); @@ -288,7 +287,7 @@ void VCS_PROB::prob_report(int print_lvl) plogf(" Species (phase) " " SS0ChemPot StarChemPot\n"); for (size_t iphase = 0; iphase < NPhase; iphase++) { - Vphase = VPhaseList[iphase]; + vcs_VolPhase* Vphase = VPhaseList[iphase]; Vphase->setState_TP(T, PresPA); for (size_t kindex = 0; kindex < Vphase->nSpecies(); kindex++) { size_t kglob = Vphase->spGlobalIndexVCS(kindex); @@ -317,24 +316,20 @@ void VCS_PROB::prob_report(int print_lvl) void VCS_PROB::addPhaseElements(vcs_VolPhase* volPhase) { - size_t e, eVP; - size_t foundPos = npos; size_t neVP = volPhase->nElemConstraints(); - std::string en; - std::string enVP; /* * Loop through the elements in the vol phase object */ - for (eVP = 0; eVP < neVP; eVP++) { - foundPos = npos; - enVP = volPhase->elementName(eVP); + for (size_t eVP = 0; eVP < neVP; eVP++) { + size_t foundPos = npos; + std::string enVP = volPhase->elementName(eVP); /* * Search for matches with the existing elements. * If found, then fill in the entry in the global * mapping array. */ - for (e = 0; e < ne; e++) { - en = ElName[e]; + for (size_t e = 0; e < ne; e++) { + std::string en = ElName[e]; if (!strcmp(enVP.c_str(), en.c_str())) { volPhase->setElemGlobalIndex(eVP, e); foundPos = e; @@ -343,7 +338,7 @@ void VCS_PROB::addPhaseElements(vcs_VolPhase* volPhase) if (foundPos == npos) { int elType = volPhase->elementType(eVP); int elactive = volPhase->elementActive(eVP); - e = addElement(enVP.c_str(), elType, elactive); + size_t e = addElement(enVP.c_str(), elType, elactive); volPhase->setElemGlobalIndex(eVP, e); } } @@ -366,7 +361,6 @@ size_t VCS_PROB::addElement(const char* elNameNew, int elType, int elactive) size_t VCS_PROB::addOnePhaseSpecies(vcs_VolPhase* volPhase, size_t k, size_t kT) { - size_t e, eVP; if (kT > nspecies) { /* * Need to expand the number of species here @@ -375,8 +369,8 @@ size_t VCS_PROB::addOnePhaseSpecies(vcs_VolPhase* volPhase, size_t k, size_t kT) exit(EXIT_FAILURE); } double const* const* const fm = volPhase->getFormulaMatrix(); - for (eVP = 0; eVP < volPhase->nElemConstraints(); eVP++) { - e = volPhase->elemGlobalIndex(eVP); + for (size_t eVP = 0; eVP < volPhase->nElemConstraints(); eVP++) { + size_t e = volPhase->elemGlobalIndex(eVP); #ifdef DEBUG_MODE if (e == npos) { exit(EXIT_FAILURE); @@ -394,18 +388,11 @@ size_t VCS_PROB::addOnePhaseSpecies(vcs_VolPhase* volPhase, size_t k, size_t kT) void VCS_PROB::reportCSV(const std::string& reportFile) { - size_t k; - size_t istart; - - double vol = 0.0; - string sName; - FILE* FP = fopen(reportFile.c_str(), "w"); if (!FP) { plogf("Failure to open file\n"); exit(EXIT_FAILURE); } - double Temp = T; std::vector volPM(nspecies, 0.0); std::vector activity(nspecies, 0.0); @@ -414,10 +401,10 @@ void VCS_PROB::reportCSV(const std::string& reportFile) std::vector mu0(nspecies, 0.0); std::vector molalities(nspecies, 0.0); - vol = 0.0; + double vol = 0.0; size_t iK = 0; for (size_t iphase = 0; iphase < NPhase; iphase++) { - istart = iK; + size_t istart = iK; vcs_VolPhase* volP = VPhaseList[iphase]; //const Cantera::ThermoPhase *tptr = volP->ptrThermoPhase(); size_t nSpeciesPhase = volP->nSpecies(); @@ -426,7 +413,7 @@ void VCS_PROB::reportCSV(const std::string& reportFile) double TMolesPhase = volP->totalMoles(); double VolPhaseVolumes = 0.0; - for (k = 0; k < nSpeciesPhase; k++) { + for (size_t k = 0; k < nSpeciesPhase; k++) { iK++; VolPhaseVolumes += volPM[istart + k] * mf[istart + k]; } @@ -436,7 +423,7 @@ void VCS_PROB::reportCSV(const std::string& reportFile) fprintf(FP,"--------------------- VCS_MULTIPHASE_EQUIL FINAL REPORT" " -----------------------------\n"); - fprintf(FP,"Temperature = %11.5g kelvin\n", Temp); + fprintf(FP,"Temperature = %11.5g kelvin\n", T); fprintf(FP,"Pressure = %11.5g Pascal\n", PresPA); fprintf(FP,"Total Volume = %11.5g m**3\n", vol); fprintf(FP,"Number Basis optimizations = %d\n", m_NumBasisOptimizations); @@ -444,7 +431,7 @@ void VCS_PROB::reportCSV(const std::string& reportFile) iK = 0; for (size_t iphase = 0; iphase < NPhase; iphase++) { - istart = iK; + size_t istart = iK; vcs_VolPhase* volP = VPhaseList[iphase]; const Cantera::ThermoPhase* tp = volP->ptrThermoPhase(); @@ -469,13 +456,12 @@ void VCS_PROB::reportCSV(const std::string& reportFile) tp->getPartialMolarVolumes(VCS_DATA_PTR(volPM)); tp->getChemPotentials(VCS_DATA_PTR(mu)); double VolPhaseVolumes = 0.0; - for (k = 0; k < nSpeciesPhase; k++) { + for (size_t k = 0; k < nSpeciesPhase; k++) { VolPhaseVolumes += volPM[k] * mf[istart + k]; } VolPhaseVolumes *= TMolesPhase; vol += VolPhaseVolumes; - if (actConvention == 1) { const Cantera::MolalityVPSSTP* mTP = static_cast(tp); tp->getChemPotentials(VCS_DATA_PTR(mu)); @@ -491,8 +477,8 @@ void VCS_PROB::reportCSV(const std::string& reportFile) " , , ," " (J/kmol), (J/kmol), (kmol), (m**3/kmol), (m**3)\n"); } - for (k = 0; k < nSpeciesPhase; k++) { - sName = tp->speciesName(k); + for (size_t k = 0; k < nSpeciesPhase; k++) { + std::string sName = tp->speciesName(k); fprintf(FP,"%12s, %11s, %11.3e, %11.3e, %11.3e, %11.3e, %11.3e," "%11.3e, %11.3e, %11.3e, %11.3e, %11.3e\n", sName.c_str(), @@ -513,11 +499,11 @@ void VCS_PROB::reportCSV(const std::string& reportFile) " , , ," " (J/kmol), (J/kmol), (kmol), (m**3/kmol), (m**3)\n"); } - for (k = 0; k < nSpeciesPhase; k++) { + for (size_t k = 0; k < nSpeciesPhase; k++) { molalities[k] = 0.0; } - for (k = 0; k < nSpeciesPhase; k++) { - sName = tp->speciesName(k); + for (size_t k = 0; k < nSpeciesPhase; k++) { + std::string sName = tp->speciesName(k); fprintf(FP,"%12s, %11s, %11.3e, %11.3e, %11.3e, %11.3e, %11.3e, " "%11.3e, %11.3e,% 11.3e, %11.3e, %11.3e\n", sName.c_str(), @@ -534,7 +520,7 @@ void VCS_PROB::reportCSV(const std::string& reportFile) * Check consistency: These should be equal */ tp->getChemPotentials(VCS_DATA_PTR(m_gibbsSpecies)+istart); - for (k = 0; k < nSpeciesPhase; k++) { + for (size_t k = 0; k < nSpeciesPhase; k++) { if (!vcs_doubleEqual(m_gibbsSpecies[istart+k], mu[k])) { fprintf(FP,"ERROR: incompatibility!\n"); fclose(FP); diff --git a/src/equil/vcs_rearrange.cpp b/src/equil/vcs_rearrange.cpp index ad07e4cb9..beeeb7080 100644 --- a/src/equil/vcs_rearrange.cpp +++ b/src/equil/vcs_rearrange.cpp @@ -15,17 +15,16 @@ namespace VCSnonideal { int VCS_SOLVE::vcs_rearrange() { - size_t i, l, j; size_t k1 = 0; /* - Loop over all of the species */ - for (i = 0; i < m_numSpeciesTot; ++i) { + for (size_t i = 0; i < m_numSpeciesTot; ++i) { /* * Find the index of I in the index vector m_speciesIndexVector[]. * Call it k1 and continue. */ - for (j = 0; j < m_numSpeciesTot; ++j) { - l = m_speciesMapIndex[j]; + for (size_t j = 0; j < m_numSpeciesTot; ++j) { + size_t l = m_speciesMapIndex[j]; k1 = j; if (l == i) { break; diff --git a/src/equil/vcs_report.cpp b/src/equil/vcs_report.cpp index d43495d5f..9e089b4bd 100644 --- a/src/equil/vcs_report.cpp +++ b/src/equil/vcs_report.cpp @@ -29,13 +29,10 @@ static void print_line(const std::string& schar, size_t num) int VCS_SOLVE::vcs_report(int iconv) { bool printActualMoles = true, inertYes = false; - size_t i, j, l, k, kspec; size_t nspecies = m_numSpeciesTot; - double g; char originalUnitsState = m_unitsState; - std::vector sortindex(nspecies,0); std::vector xy(nspecies,0.0); @@ -43,7 +40,7 @@ int VCS_SOLVE::vcs_report(int iconv) /* **** SORT DEPENDENT SPECIES IN DECREASING ORDER OF MOLES ***** */ /* ************************************************************** */ - for (i = 0; i < nspecies; ++i) { + for (size_t i = 0; i < nspecies; ++i) { sortindex[i] = i; xy[i] = m_molNumSpecies_old[i]; } @@ -52,8 +49,8 @@ int VCS_SOLVE::vcs_report(int iconv) * and the sort index vector, sortindex, according to * the magnitude of the mole fraction vector. */ - for (l = m_numComponents; l < m_numSpeciesRdc; ++l) { - k = vcs_optMax(VCS_DATA_PTR(xy), 0, l, m_numSpeciesRdc); + for (size_t l = m_numComponents; l < m_numSpeciesRdc; ++l) { + size_t k = vcs_optMax(VCS_DATA_PTR(xy), 0, l, m_numSpeciesRdc); if (k != l) { std::swap(xy[k], xy[l]); std::swap(sortindex[k], sortindex[l]); @@ -113,7 +110,7 @@ int VCS_SOLVE::vcs_report(int iconv) plogf(" Species Equilibrium kmoles "); plogf("Mole Fraction ChemPot/RT SpecUnkType\n"); print_line("-", 80); - for (i = 0; i < m_numComponents; ++i) { + for (size_t i = 0; i < m_numComponents; ++i) { plogf(" %-12.12s", m_speciesName[i].c_str()); print_space(13); plogf("%14.7E %14.7E %12.4E", m_molNumSpecies_old[i] * molScale, @@ -121,8 +118,8 @@ int VCS_SOLVE::vcs_report(int iconv) plogf(" %3d", m_speciesUnknownType[i]); plogf("\n"); } - for (i = m_numComponents; i < m_numSpeciesRdc; ++i) { - l = sortindex[i]; + for (size_t i = m_numComponents; i < m_numSpeciesRdc; ++i) { + size_t l = sortindex[i]; plogf(" %-12.12s", m_speciesName[l].c_str()); print_space(13); @@ -139,7 +136,7 @@ int VCS_SOLVE::vcs_report(int iconv) } plogf("\n"); } - for (i = 0; i < m_numPhases; i++) { + for (size_t i = 0; i < m_numPhases; i++) { if (TPhInertMoles[i] > 0.0) { inertYes = true; if (i == 0) { @@ -154,15 +151,15 @@ int VCS_SOLVE::vcs_report(int iconv) } if (m_numSpeciesRdc != nspecies) { plogf("\n SPECIES WITH LESS THAN 1.0E-32 KMOLES:\n\n"); - for (kspec = m_numSpeciesRdc; kspec < nspecies; ++kspec) { + for (size_t kspec = m_numSpeciesRdc; kspec < nspecies; ++kspec) { plogf(" %-12.12s", m_speciesName[kspec].c_str()); // Note m_deltaGRxn_new[] stores in kspec slot not irxn slot, after solve plogf(" %14.7E %14.7E %12.4E", m_molNumSpecies_old[kspec]*molScale, m_molNumSpecies_new[kspec]*molScale, m_deltaGRxn_new[kspec]); - if (m_speciesUnknownType[i] == VCS_SPECIES_TYPE_MOLNUM) { + if (m_speciesUnknownType[kspec] == VCS_SPECIES_TYPE_MOLNUM) { plogf(" KMol_Num"); - } else if (m_speciesUnknownType[i] == VCS_SPECIES_TYPE_INTERFACIALVOLTAGE) { + } else if (m_speciesUnknownType[kspec] == VCS_SPECIES_TYPE_INTERFACIALVOLTAGE) { plogf(" Voltage"); } else { plogf(" Unknown"); @@ -180,17 +177,17 @@ int VCS_SOLVE::vcs_report(int iconv) plogf("\n"); print_line("-", m_numComponents*10 + 45); plogf(" |ComponentID|"); - for (j = 0; j < m_numComponents; j++) { + for (size_t j = 0; j < m_numComponents; j++) { plogf(" %3d", j); } plogf(" | |\n"); plogf(" | Components|"); - for (j = 0; j < m_numComponents; j++) { + for (size_t j = 0; j < m_numComponents; j++) { plogf(" %10.10s", m_speciesName[j].c_str()); } plogf(" | |\n"); plogf(" NonComponent | Moles |"); - for (j = 0; j < m_numComponents; j++) { + for (size_t j = 0; j < m_numComponents; j++) { plogf(" %10.3g", m_molNumSpecies_old[j] * molScale); } plogf(" | DG/RT Rxn |\n"); @@ -200,7 +197,7 @@ int VCS_SOLVE::vcs_report(int iconv) plogf(" %3d ", kspec); plogf("%-10.10s", m_speciesName[kspec].c_str()); plogf("|%10.3g |", m_molNumSpecies_old[kspec]*molScale); - for (j = 0; j < m_numComponents; j++) { + for (size_t j = 0; j < m_numComponents; j++) { plogf(" %6.2f", m_stoichCoeffRxnMatrix[irxn][j]); } plogf(" |%10.3g |", m_deltaGRxn_new[irxn]); @@ -221,17 +218,17 @@ int VCS_SOLVE::vcs_report(int iconv) plogf("\n"); print_line("-", m_numElemConstraints*10 + 58); plogf(" | ElementID |"); - for (j = 0; j < m_numElemConstraints; j++) { + for (size_t j = 0; j < m_numElemConstraints; j++) { plogf(" %3d", j); } plogf(" | |\n"); plogf(" | Element |"); - for (j = 0; j < m_numElemConstraints; j++) { + for (size_t j = 0; j < m_numElemConstraints; j++) { plogf(" %10.10s", (m_elementName[j]).c_str()); } plogf(" | |\n"); plogf(" PhaseName |KMolTarget |"); - for (j = 0; j < m_numElemConstraints; j++) { + for (size_t j = 0; j < m_numElemConstraints; j++) { plogf(" %10.3g", m_elemAbundancesGoal[j]); } plogf(" | Gibbs Total |\n"); @@ -249,7 +246,7 @@ int VCS_SOLVE::vcs_report(int iconv) } } vcs_elabPhase(iphase, VCS_DATA_PTR(gaPhase)); - for (j = 0; j < m_numElemConstraints; j++) { + for (size_t j = 0; j < m_numElemConstraints; j++) { plogf(" %10.3g", gaPhase[j]); gaTPhase[j] += gaPhase[j]; } @@ -260,7 +257,7 @@ int VCS_SOLVE::vcs_report(int iconv) } print_line("-", m_numElemConstraints*10 + 58); plogf(" TOTAL |%10.3e |", totalMoles); - for (j = 0; j < m_numElemConstraints; j++) { + for (size_t j = 0; j < m_numElemConstraints; j++) { plogf(" %10.3g", gaTPhase[j]); } plogf(" | %18.11E |\n", gibbsTotal); @@ -278,8 +275,8 @@ int VCS_SOLVE::vcs_report(int iconv) * energy of zero */ - g = vcs_Total_Gibbs(VCS_DATA_PTR(m_molNumSpecies_old), VCS_DATA_PTR(m_feSpecies_old), - VCS_DATA_PTR(m_tPhaseMoles_old)); + double g = vcs_Total_Gibbs(VCS_DATA_PTR(m_molNumSpecies_old), VCS_DATA_PTR(m_feSpecies_old), + VCS_DATA_PTR(m_tPhaseMoles_old)); plogf("\n\tTotal Dimensionless Gibbs Free Energy = G/RT = %15.7E\n", g); if (inertYes) { plogf("\t\t(Inert species have standard free energy of zero)\n"); @@ -287,7 +284,7 @@ int VCS_SOLVE::vcs_report(int iconv) plogf("\nElemental Abundances (kmol): "); plogf(" Actual Target Type ElActive\n"); - for (i = 0; i < m_numElemConstraints; ++i) { + for (size_t i = 0; i < m_numElemConstraints; ++i) { print_space(26); plogf("%-2.2s", (m_elementName[i]).c_str()); plogf("%20.12E %20.12E", m_elemAbundances[i]*molScale, m_elemAbundancesGoal[i]*molScale); @@ -311,8 +308,8 @@ int VCS_SOLVE::vcs_report(int iconv) plogf("| (MolNum ChemPot)|"); plogf("\n"); print_line("-", 147); - for (i = 0; i < nspecies; ++i) { - l = sortindex[i]; + for (size_t i = 0; i < nspecies; ++i) { + size_t l = sortindex[i]; size_t pid = m_phaseID[l]; plogf(" %-12.12s", m_speciesName[l].c_str()); plogf(" %14.7E ", m_molNumSpecies_old[l]*molScale); @@ -351,7 +348,7 @@ int VCS_SOLVE::vcs_report(int iconv) plogf("| %20.9E |", m_feSpecies_old[l] * m_molNumSpecies_old[l] * molScale); plogf("\n"); } - for (i = 0; i < 125; i++) { + for (size_t i = 0; i < 125; i++) { plogf(" "); } plogf(" %20.9E\n", g); diff --git a/src/equil/vcs_rxnadj.cpp b/src/equil/vcs_rxnadj.cpp index c7f0d0cc9..5ac1e65cb 100644 --- a/src/equil/vcs_rxnadj.cpp +++ b/src/equil/vcs_rxnadj.cpp @@ -19,12 +19,8 @@ namespace VCSnonideal size_t VCS_SOLVE::vcs_RxnStepSizes(int& forceComponentCalc, size_t& kSpecial) { - size_t kspec, iph; size_t iphDel = npos; - double s, xx, dss; size_t k = 0; - vcs_VolPhase* Vphase = 0; - double* dnPhase_irxn; #ifdef DEBUG_MODE char ANOTE[128]; if (m_debug_print_lvl >= 2) { @@ -59,7 +55,7 @@ size_t VCS_SOLVE::vcs_RxnStepSizes(int& forceComponentCalc, size_t& kSpecial) sprintf(ANOTE, "Normal Calc"); #endif - kspec = m_indexRxnToSpecies[irxn]; + size_t kspec = m_indexRxnToSpecies[irxn]; if (m_speciesStatus[kspec] == VCS_SPECIES_ZEROEDPHASE) { m_deltaMolNumSpecies[kspec] = 0.0; @@ -68,7 +64,7 @@ size_t VCS_SOLVE::vcs_RxnStepSizes(int& forceComponentCalc, size_t& kSpecial) #endif } else if (m_speciesUnknownType[kspec] != VCS_SPECIES_TYPE_INTERFACIALVOLTAGE) { - dnPhase_irxn = m_deltaMolNumPhase[irxn]; + double* dnPhase_irxn = m_deltaMolNumPhase[irxn]; if (m_molNumSpecies_old[kspec] == 0.0 && (!m_SSPhase[kspec])) { /********************************************************************/ @@ -84,10 +80,10 @@ size_t VCS_SOLVE::vcs_RxnStepSizes(int& forceComponentCalc, size_t& kSpecial) * First decide if this species is part of a multiphase that * is nontrivial in size. */ - iph = m_phaseID[kspec]; + size_t iph = m_phaseID[kspec]; double tphmoles = m_tPhaseMoles_old[iph]; double trphmoles = tphmoles / m_totalMolNum; - Vphase = m_VolPhaseList[iph]; + vcs_VolPhase* Vphase = m_VolPhaseList[iph]; if (Vphase->exists() && (trphmoles > VCS_DELETE_PHASE_CUTOFF)) { m_deltaMolNumSpecies[kspec] = m_totalMolNum * VCS_SMALL_MULTIPHASE_SPECIES; if (m_speciesStatus[kspec] == VCS_SPECIES_STOICHZERO) { @@ -169,6 +165,7 @@ size_t VCS_SOLVE::vcs_RxnStepSizes(int& forceComponentCalc, size_t& kSpecial) /* * Start of the regular processing */ + double s; if (m_SSPhase[kspec]) { s = 0.0; } else { @@ -182,7 +179,7 @@ size_t VCS_SOLVE::vcs_RxnStepSizes(int& forceComponentCalc, size_t& kSpecial) } } for (size_t j = 0; j < m_numPhases; j++) { - Vphase = m_VolPhaseList[j]; + vcs_VolPhase* Vphase = m_VolPhaseList[j]; if (!Vphase->m_singleSpecies) { if (m_tPhaseMoles_old[j] > 0.0) { s -= SQUARE(dnPhase_irxn[j]) / m_tPhaseMoles_old[j]; @@ -254,12 +251,13 @@ size_t VCS_SOLVE::vcs_RxnStepSizes(int& forceComponentCalc, size_t& kSpecial) * will zero out first. * -> The species to be zeroed out will be "k". */ + double dss; if (m_deltaGRxn_new[irxn] > 0.0) { dss = m_molNumSpecies_old[kspec]; k = kspec; for (size_t j = 0; j < m_numComponents; ++j) { if (m_stoichCoeffRxnMatrix[irxn][j] > 0.0) { - xx = m_molNumSpecies_old[j] / m_stoichCoeffRxnMatrix[irxn][j]; + double xx = m_molNumSpecies_old[j] / m_stoichCoeffRxnMatrix[irxn][j]; if (xx < dss) { dss = xx; k = j; @@ -271,7 +269,7 @@ size_t VCS_SOLVE::vcs_RxnStepSizes(int& forceComponentCalc, size_t& kSpecial) dss = 1.0e10; for (size_t j = 0; j < m_numComponents; ++j) { if (m_stoichCoeffRxnMatrix[irxn][j] < 0.0) { - xx = -m_molNumSpecies_old[j] / m_stoichCoeffRxnMatrix[irxn][j]; + double xx = -m_molNumSpecies_old[j] / m_stoichCoeffRxnMatrix[irxn][j]; if (xx < dss) { dss = xx; k = j; @@ -382,16 +380,11 @@ size_t VCS_SOLVE::vcs_RxnStepSizes(int& forceComponentCalc, size_t& kSpecial) int VCS_SOLVE::vcs_rxn_adj_cg() { - size_t irxn, j; - size_t k = 0; - size_t kspec; int soldel = 0; - double s, xx, dss; - double* dnPhase_irxn; #ifdef DEBUG_MODE char ANOTE[128]; plogf(" "); - for (j = 0; j < 77; j++) { + for (size_t j = 0; j < 77; j++) { plogf("-"); } plogf("\n --- Subroutine rxn_adj_cg() called\n"); @@ -404,13 +397,13 @@ int VCS_SOLVE::vcs_rxn_adj_cg() * We also evaluate whether the matrix is appropriate for * this algorithm. If not, we bail out. */ - for (irxn = 0; irxn < m_numRxnRdc; ++irxn) { + for (size_t irxn = 0; irxn < m_numRxnRdc; ++irxn) { #ifdef DEBUG_MODE sprintf(ANOTE, "Normal Calc"); #endif - kspec = m_indexRxnToSpecies[irxn]; - dnPhase_irxn = m_deltaMolNumPhase[irxn]; + size_t kspec = m_indexRxnToSpecies[irxn]; + double* dnPhase_irxn = m_deltaMolNumPhase[irxn]; if (m_molNumSpecies_old[kspec] == 0.0 && (!m_SSPhase[kspec])) { /* *******************************************************************/ @@ -471,17 +464,18 @@ int VCS_SOLVE::vcs_rxn_adj_cg() /* * Start of the regular processing */ + double s; if (m_SSPhase[kspec]) { s = 0.0; } else { s = 1.0 / m_molNumSpecies_old[kspec]; } - for (j = 0; j < m_numComponents; ++j) { + for (size_t j = 0; j < m_numComponents; ++j) { if (!m_SSPhase[j]) { s += SQUARE(m_stoichCoeffRxnMatrix[irxn][j]) / m_molNumSpecies_old[j]; } } - for (j = 0; j < m_numPhases; j++) { + for (size_t j = 0; j < m_numPhases; j++) { if (!(m_VolPhaseList[j])->m_singleSpecies) { if (m_tPhaseMoles_old[j] > 0.0) { s -= SQUARE(dnPhase_irxn[j]) / m_tPhaseMoles_old[j]; @@ -502,12 +496,14 @@ int VCS_SOLVE::vcs_rxn_adj_cg() * Then, we need to follow the reaction to see which species * will zero out first. */ + size_t k; + double dss; if (m_deltaGRxn_new[irxn] > 0.0) { dss = m_molNumSpecies_old[kspec]; k = kspec; - for (j = 0; j < m_numComponents; ++j) { + for (size_t j = 0; j < m_numComponents; ++j) { if (m_stoichCoeffRxnMatrix[irxn][j] > 0.0) { - xx = m_molNumSpecies_old[j] / m_stoichCoeffRxnMatrix[irxn][j]; + double xx = m_molNumSpecies_old[j] / m_stoichCoeffRxnMatrix[irxn][j]; if (xx < dss) { dss = xx; k = j; @@ -517,9 +513,9 @@ int VCS_SOLVE::vcs_rxn_adj_cg() dss = -dss; } else { dss = 1.0e10; - for (j = 0; j < m_numComponents; ++j) { + for (size_t j = 0; j < m_numComponents; ++j) { if (m_stoichCoeffRxnMatrix[irxn][j] < 0.0) { - xx = -m_molNumSpecies_old[j] / m_stoichCoeffRxnMatrix[irxn][j]; + double xx = -m_molNumSpecies_old[j] / m_stoichCoeffRxnMatrix[irxn][j]; if (xx < dss) { dss = xx; k = j; @@ -538,7 +534,7 @@ int VCS_SOLVE::vcs_rxn_adj_cg() if (dss != 0.0) { m_molNumSpecies_old[kspec] += dss; m_tPhaseMoles_old[m_phaseID[kspec]] += dss; - for (j = 0; j < m_numComponents; ++j) { + for (size_t j = 0; j < m_numComponents; ++j) { m_molNumSpecies_old[j] += dss * m_stoichCoeffRxnMatrix[irxn][j]; m_tPhaseMoles_old[m_phaseID[j]] += dss * m_stoichCoeffRxnMatrix[irxn][j]; } @@ -582,7 +578,7 @@ int VCS_SOLVE::vcs_rxn_adj_cg() #ifdef DEBUG_MODE plogf(" "); - for (j = 0; j < 77; j++) { + for (size_t j = 0; j < 77; j++) { plogf("-"); } plogf("\n"); @@ -610,27 +606,24 @@ double VCS_SOLVE::vcs_Hessian_diag_adj(size_t irxn, double hessianDiag_Ideal) double VCS_SOLVE::vcs_Hessian_actCoeff_diag(size_t irxn) { - size_t kspec, k, l, kph; - double s; - double* sc_irxn; - kspec = m_indexRxnToSpecies[irxn]; - kph = m_phaseID[kspec]; + size_t kspec = m_indexRxnToSpecies[irxn]; + size_t kph = m_phaseID[kspec]; double np_kspec = m_tPhaseMoles_old[kph]; if (np_kspec < 1.0E-13) { np_kspec = 1.0E-13; } - sc_irxn = m_stoichCoeffRxnMatrix[irxn]; + double* sc_irxn = m_stoichCoeffRxnMatrix[irxn]; /* * First the diagonal term of the Jacobian */ - s = m_np_dLnActCoeffdMolNum[kspec][kspec] / np_kspec; + double s = m_np_dLnActCoeffdMolNum[kspec][kspec] / np_kspec; /* * Next, the other terms. Note this only a loop over the components * So, it's not too expensive to calculate. */ - for (l = 0; l < m_numComponents; l++) { + for (size_t l = 0; l < m_numComponents; l++) { if (!m_SSPhase[l]) { - for (k = 0; k < m_numComponents; ++k) { + for (size_t k = 0; k < m_numComponents; ++k) { if (m_phaseID[k] == m_phaseID[l]) { double np = m_tPhaseMoles_old[m_phaseID[k]]; if (np > 0.0) { @@ -697,7 +690,6 @@ double VCS_SOLVE::vcs_line_search(const size_t irxn, const double dx_orig) #endif { int its = 0; - size_t k; size_t kspec = m_indexRxnToSpecies[irxn]; const int MAXITS = 10; double dx = dx_orig; @@ -705,8 +697,6 @@ double VCS_SOLVE::vcs_line_search(const size_t irxn, const double dx_orig) double* molNumBase = VCS_DATA_PTR(m_molNumSpecies_old); double* acBase = VCS_DATA_PTR(m_actCoeffSpecies_old); double* ac = VCS_DATA_PTR(m_actCoeffSpecies_new); - double molSum = 0.0; - double slope; /* * Calculate the deltaG value at the dx = 0.0 point */ @@ -743,9 +733,9 @@ double VCS_SOLVE::vcs_line_search(const size_t irxn, const double dx_orig) } vcs_dcopy(VCS_DATA_PTR(m_molNumSpecies_new), molNumBase, m_numSpeciesRdc); - molSum = molNumBase[kspec]; + double molSum = molNumBase[kspec]; m_molNumSpecies_new[kspec] = molNumBase[kspec] + dx_orig; - for (k = 0; k < m_numComponents; k++) { + for (size_t k = 0; k < m_numComponents; k++) { m_molNumSpecies_new[k] = molNumBase[k] + sc_irxn[k] * dx_orig; molSum += molNumBase[k]; } @@ -769,7 +759,7 @@ double VCS_SOLVE::vcs_line_search(const size_t irxn, const double dx_orig) */ if (fabs(deltaG1) < 0.8 * forig) { if (deltaG1 * deltaGOrig < 0.0) { - slope = (deltaG1 - deltaGOrig) / dx_orig; + double slope = (deltaG1 - deltaGOrig) / dx_orig; dx = -deltaGOrig / slope; } else { dx = dx_orig; @@ -786,7 +776,7 @@ double VCS_SOLVE::vcs_line_search(const size_t irxn, const double dx_orig) */ dx *= 0.5; m_molNumSpecies_new[kspec] = molNumBase[kspec] + dx; - for (k = 0; k < m_numComponents; k++) { + for (size_t k = 0; k < m_numComponents; k++) { m_molNumSpecies_new[k] = molNumBase[k] + sc_irxn[k] * dx; } vcs_setFlagsVolPhases(false, VCS_STATECALC_NEW); @@ -806,7 +796,7 @@ double VCS_SOLVE::vcs_line_search(const size_t irxn, const double dx_orig) */ if (fabs(deltaG) / forig < (1.0 - 0.1 * dx / dx_orig)) { if (deltaG * deltaGOrig < 0.0) { - slope = (deltaG - deltaGOrig) / dx; + double slope = (deltaG - deltaGOrig) / dx; dx = -deltaGOrig / slope; } goto finalize; diff --git a/src/equil/vcs_solve.cpp b/src/equil/vcs_solve.cpp index 9304e026f..c911734b2 100644 --- a/src/equil/vcs_solve.cpp +++ b/src/equil/vcs_solve.cpp @@ -222,15 +222,14 @@ VCS_SOLVE::~VCS_SOLVE() void VCS_SOLVE::vcs_delete_memory() { - size_t j; size_t nspecies = m_numSpeciesTot; - for (j = 0; j < m_numPhases; j++) { + for (size_t j = 0; j < m_numPhases; j++) { delete m_VolPhaseList[j]; m_VolPhaseList[j] = 0; } - for (j = 0; j < nspecies; j++) { + for (size_t j = 0; j < nspecies; j++) { delete m_speciesThermoList[j]; m_speciesThermoList[j] = 0; } @@ -247,7 +246,6 @@ void VCS_SOLVE::vcs_delete_memory() int VCS_SOLVE::vcs(VCS_PROB* vprob, int ifunc, int ipr, int ip1, int maxit) { int retn = 0, iconv = 0; - size_t nspecies0, nelements0, nphase0; Cantera::clockWC tickTock; int iprintTime = std::max(ipr, ip1); @@ -266,9 +264,9 @@ int VCS_SOLVE::vcs(VCS_PROB* vprob, int ifunc, int ipr, int ip1, int maxit) * This function is called to create the private data * using the public data. */ - nspecies0 = vprob->nspecies + 10; - nelements0 = vprob->ne; - nphase0 = vprob->NPhase; + size_t nspecies0 = vprob->nspecies + 10; + size_t nelements0 = vprob->ne; + size_t nphase0 = vprob->NPhase; vcs_initSizes(nspecies0, nelements0, nphase0); @@ -380,7 +378,6 @@ int VCS_SOLVE::vcs(VCS_PROB* vprob, int ifunc, int ipr, int ip1, int maxit) int VCS_SOLVE::vcs_prob_specifyFully(const VCS_PROB* pub) { - vcs_VolPhase* Vphase = 0; const char* ser = "vcs_pub_to_priv ERROR :ill defined interface -> bailout:\n\t"; @@ -566,7 +563,7 @@ int VCS_SOLVE::vcs_prob_specifyFully(const VCS_PROB* pub) * TPhInertMoles[] -> must be copied over here */ for (size_t iph = 0; iph < nph; iph++) { - Vphase = pub->VPhaseList[iph]; + vcs_VolPhase* Vphase = pub->VPhaseList[iph]; TPhInertMoles[iph] = Vphase->totalMolesInert(); } @@ -626,7 +623,7 @@ int VCS_SOLVE::vcs_prob_specifyFully(const VCS_PROB* pub) numPhSp[iph]++; } for (size_t iph = 0; iph < nph; iph++) { - Vphase = pub->VPhaseList[iph]; + vcs_VolPhase* Vphase = pub->VPhaseList[iph]; if (numPhSp[iph] != Vphase->nSpecies()) { plogf("%sNumber of species in phase %d, %s, doesn't match\n", ser, iph, Vphase->PhaseName.c_str()); @@ -703,7 +700,7 @@ int VCS_SOLVE::vcs_prob_specifyFully(const VCS_PROB* pub) * It should point to the species thermo pointer in the private * data space. */ - Vphase = m_VolPhaseList[iph]; + vcs_VolPhase* Vphase = m_VolPhaseList[iph]; for (size_t k = 0; k < Vphase->nSpecies(); k++) { vcs_SpeciesProperties* sProp = Vphase->speciesProperty(k); size_t kT = Vphase->spGlobalIndexVCS(k); @@ -715,7 +712,7 @@ int VCS_SOLVE::vcs_prob_specifyFully(const VCS_PROB* pub) * Specify the Activity Convention information */ for (size_t iph = 0; iph < nph; iph++) { - Vphase = m_VolPhaseList[iph]; + vcs_VolPhase* Vphase = m_VolPhaseList[iph]; m_phaseActConvention[iph] = Vphase->p_activityConvention; if (Vphase->p_activityConvention != 0) { /* @@ -762,10 +759,8 @@ int VCS_SOLVE::vcs_prob_specifyFully(const VCS_PROB* pub) int VCS_SOLVE::vcs_prob_specify(const VCS_PROB* pub) { - size_t kspec, k, i, j, iph; string yo("vcs_prob_specify ERROR: "); int retn = VCS_SUCCESS; - bool status_change = false; m_temperature = pub->T; m_pressurePA = pub->PresPA; @@ -779,8 +774,8 @@ int VCS_SOLVE::vcs_prob_specify(const VCS_PROB* pub) m_tolmaj2 = 0.01 * m_tolmaj; m_tolmin2 = 0.01 * m_tolmin; - for (kspec = 0; kspec < m_numSpeciesTot; ++kspec) { - k = m_speciesMapIndex[kspec]; + for (size_t kspec = 0; kspec < m_numSpeciesTot; ++kspec) { + size_t k = m_speciesMapIndex[kspec]; m_molNumSpecies_old[kspec] = pub->w[k]; m_molNumSpecies_new[kspec] = pub->mf[k]; m_feSpecies_old[kspec] = pub->m_gibbsSpecies[k]; @@ -789,8 +784,8 @@ int VCS_SOLVE::vcs_prob_specify(const VCS_PROB* pub) /* * Transfer the element abundance goals to the solve object */ - for (i = 0; i < m_numElemConstraints; i++) { - j = m_elementMapIndex[i]; + for (size_t i = 0; i < m_numElemConstraints; i++) { + size_t j = m_elementMapIndex[i]; m_elemAbundancesGoal[i] = pub->gai[j]; } @@ -813,7 +808,8 @@ int VCS_SOLVE::vcs_prob_specify(const VCS_PROB* pub) * condition. */ - for (iph = 0; iph < m_numPhases; iph++) { + bool status_change = false; + for (size_t iph = 0; iph < m_numPhases; iph++) { vcs_VolPhase* vPhase = m_VolPhaseList[iph]; vcs_VolPhase* pub_phase_ptr = pub->VPhaseList[iph]; diff --git a/src/equil/vcs_solve_TP.cpp b/src/equil/vcs_solve_TP.cpp index 1eeefc776..62f97c765 100644 --- a/src/equil/vcs_solve_TP.cpp +++ b/src/equil/vcs_solve_TP.cpp @@ -61,13 +61,10 @@ void VCS_SOLVE::checkDelta1(double* const dsLocal, int VCS_SOLVE::vcs_solve_TP(int print_lvl, int printDetails, int maxit) { int retn = VCS_SUCCESS, soldel, solveFail; - double test, RT; - size_t j, k, l, l1, kspec, irxn, i; - bool conv = false, allMinorZeroedSpecies = false, forced, lec; - size_t iph; - double dx, xx, par; + size_t kspec; + bool conv = false, allMinorZeroedSpecies = false; size_t it1 = 0; - size_t npb, iti, lnospec; + size_t npb, iti; bool dofast; int rangeErrorFound = 0; bool giveUpOnElemAbund = false; @@ -76,10 +73,6 @@ int VCS_SOLVE::vcs_solve_TP(int print_lvl, int printDetails, int maxit) bool justDeletedMultiPhase = false; bool usedZeroedSpecies; /* return flag from basopt indicating that one of the components had a zero concentration */ - vcs_VolPhase* Vphase; - double* sc_irxn = NULL; /* Stoichiometric coefficients for cur rxn */ - double* dnPhase_irxn; - double atomComp; size_t iphasePop; int forceComponentCalc = 1; size_t iphaseDelete; /* integer that determines which phase is being deleted */ @@ -145,7 +138,7 @@ int VCS_SOLVE::vcs_solve_TP(int print_lvl, int printDetails, int maxit) plogf(" PRESSURE%22.8g %3s\n", m_pressurePA, "Pa "); plogf(" TEMPERATURE%19.3f K\n", m_temperature); - Vphase = m_VolPhaseList[0]; + vcs_VolPhase* Vphase = m_VolPhaseList[0]; if (Vphase->nSpecies() > 0) { plogf(" PHASE1 INERTS%17.3f\n", TPhInertMoles[0]); } @@ -154,7 +147,7 @@ int VCS_SOLVE::vcs_solve_TP(int print_lvl, int printDetails, int maxit) } plogf("\n ELEMENTAL ABUNDANCES CORRECT"); plogf(" FROM ESTIMATE Type\n\n"); - for (i = 0; i < m_numElemConstraints; ++i) { + for (size_t i = 0; i < m_numElemConstraints; ++i) { print_space(26); plogf("%-2.2s", (m_elementName[i]).c_str()); plogf("%20.12E%20.12E %3d\n", m_elemAbundancesGoal[i], m_elemAbundances[i], @@ -187,14 +180,14 @@ int VCS_SOLVE::vcs_solve_TP(int print_lvl, int printDetails, int maxit) print_space(41); plogf(" STAN_CHEM_POT EQUILIBRIUM_EST. Species_Type\n\n"); print_space(20); - for (i = 0; i < m_numElemConstraints; ++i) { + for (size_t i = 0; i < m_numElemConstraints; ++i) { plogf("%-4.4s ", m_elementName[i].c_str()); } plogf(" PhaseID\n"); - RT = vcs_nondimMult_TP(m_VCS_UnitsFormat, m_temperature); - for (i = 0; i < m_numSpeciesTot; ++i) { + double RT = vcs_nondimMult_TP(m_VCS_UnitsFormat, m_temperature); + for (size_t i = 0; i < m_numSpeciesTot; ++i) { plogf(" %-18.18s", m_speciesName[i].c_str()); - for (j = 0; j < m_numElemConstraints; ++j) { + for (size_t j = 0; j < m_numElemConstraints; ++j) { plogf("% -7.3g ", m_formulaMatrix[j][i]); } plogf(" %3d ", m_phaseID[i]); @@ -211,13 +204,13 @@ int VCS_SOLVE::vcs_solve_TP(int print_lvl, int printDetails, int maxit) } } - for (i = 0; i < m_numSpeciesTot; ++i) { + for (size_t i = 0; i < m_numSpeciesTot; ++i) { if (m_molNumSpecies_old[i] < 0.0) { plogf("On Input species %-12s has a " "negative MF, setting it small", m_speciesName[i].c_str()); plogendl(); - iph = m_phaseID[i]; + size_t iph = m_phaseID[i]; double tmp = m_tPhaseMoles_old[iph] * VCS_RELDELETE_SPECIES_CUTOFF * 10; if (VCS_DELETE_MINORSPECIES_CUTOFF*10. > tmp) { tmp = VCS_DELETE_MINORSPECIES_CUTOFF*10.; @@ -246,7 +239,7 @@ int VCS_SOLVE::vcs_solve_TP(int print_lvl, int printDetails, int maxit) */ L_COMPONENT_CALC: ; - test = -1.0e-10; + double test = -1.0e-10; retn = vcs_basopt(false, VCS_DATA_PTR(aw), VCS_DATA_PTR(sa), VCS_DATA_PTR(sm), VCS_DATA_PTR(ss), test, &usedZeroedSpecies); @@ -260,6 +253,7 @@ L_COMPONENT_CALC: // Turn off the force componentCalc flag forceComponentCalc = 0; + bool lec = false; if (conv) { goto L_RETURN_BLOCK; } @@ -269,7 +263,6 @@ L_COMPONENT_CALC: /************** EVALUATE INITIAL SPECIES STATUS VECTOR *******************/ /*************************************************************************/ allMinorZeroedSpecies = vcs_evaluate_speciesType(); - lec = false; /*************************************************************************/ /************** EVALUATE THE ELELEMT ABUNDANCE CHECK ******************/ @@ -468,13 +461,13 @@ L_MAINLOOP_ALL_SPECIES: } #endif - for (k = 0; k < m_numSpeciesTot; k++) { + for (size_t k = 0; k < m_numSpeciesTot; k++) { m_molNumSpecies_new[k] = m_molNumSpecies_old[k] + m_deltaMolNumSpecies[k]; - iph = m_phaseID[k]; + size_t iph = m_phaseID[k]; m_tPhaseMoles_new[iph] += m_deltaMolNumSpecies[k]; } if (kspec >= m_numComponents) { - if (m_molNumSpecies_new[k] != 0.0) { + if (m_molNumSpecies_new[m_numSpeciesTot] != 0.0) { printf("vcs_solve_tp:: we shouldn't be here!\n"); exit(EXIT_FAILURE); } @@ -518,14 +511,15 @@ L_MAINLOOP_ALL_SPECIES: plogf(" KMoles Tent_KMoles Rxn_Adj | Comment \n"); } #endif - for (irxn = 0; irxn < m_numRxnRdc; irxn++) { + for (size_t irxn = 0; irxn < m_numRxnRdc; irxn++) { kspec = m_indexRxnToSpecies[irxn]; - sc_irxn = m_stoichCoeffRxnMatrix[irxn]; - iph = m_phaseID[kspec]; - Vphase = m_VolPhaseList[iph]; + double* sc_irxn = m_stoichCoeffRxnMatrix[irxn]; + size_t iph = m_phaseID[kspec]; + vcs_VolPhase* Vphase = m_VolPhaseList[iph]; #ifdef DEBUG_MODE ANOTE[0] = '\0'; #endif + double dx; if (iphasePop != npos) { if (iph == iphasePop) { dx = m_deltaMolNumSpecies[kspec]; @@ -586,7 +580,7 @@ L_MAINLOOP_ALL_SPECIES: for (size_t j = 0; j < m_numElemConstraints; ++j) { int elType = m_elType[j]; if (elType == VCS_ELEM_TYPE_ABSPOS) { - atomComp = m_formulaMatrix[j][kspec]; + double atomComp = m_formulaMatrix[j][kspec]; if (atomComp > 0.0) { double maxPermissible = m_elemAbundancesGoal[j] / atomComp; if (maxPermissible < VCS_DELETE_MINORSPECIES_CUTOFF) { @@ -714,7 +708,7 @@ L_MAINLOOP_ALL_SPECIES: * branch to the code where we reevaluate the deletion * of all species. */ - lnospec = vcs_delete_species(kspec); + size_t lnospec = vcs_delete_species(kspec); if (lnospec) { goto L_RECHECK_DELETED; } @@ -830,7 +824,7 @@ L_MAINLOOP_ALL_SPECIES: * a little more so that the component values can only be * reduced by two 99%, */ - for (j = 0; j < m_numComponents; ++j) { + for (size_t j = 0; j < m_numComponents; ++j) { if (sc_irxn[j] != 0.0) { wx[j] = m_molNumSpecies_old[j] + sc_irxn[j] * dx; if (wx[j] <= m_molNumSpecies_old[j] * 0.01 - 1.0E-150) { @@ -889,7 +883,7 @@ L_MAINLOOP_ALL_SPECIES: m_deltaMolNumSpecies[kspec] = dx; m_molNumSpecies_new[kspec] = 0.0; - for (k = 0; k < m_numComponents; ++k) { + for (size_t k = 0; k < m_numComponents; ++k) { m_deltaMolNumSpecies[k] = 0.0; } for (iph = 0; iph < m_numPhases; iph++) { @@ -945,14 +939,14 @@ L_MAINLOOP_ALL_SPECIES: exit(EXIT_FAILURE); } #endif - for (k = 0; k < m_numComponents; ++k) { + for (size_t k = 0; k < m_numComponents; ++k) { m_deltaMolNumSpecies[k] += sc_irxn[k] * dx; } /* * Calculate the tentative change in the total number of * moles in all of the phases */ - dnPhase_irxn = m_deltaMolNumPhase[irxn]; + double* dnPhase_irxn = m_deltaMolNumPhase[irxn]; for (iph = 0; iph < m_numPhases; iph++) { m_deltaPhaseMoles[iph] += dx * dnPhase_irxn[iph]; } @@ -998,7 +992,7 @@ L_MAIN_LOOP_END_NO_PRINT: #ifdef DEBUG_MODE if (m_debug_print_lvl >= 2) { - for (k = 0; k < m_numComponents; k++) { + for (size_t k = 0; k < m_numComponents; k++) { plogf(" --- "); plogf("%-12.12s", m_speciesName[k].c_str()); plogf(" c%11.4E%11.4E%11.4E |\n", @@ -1021,10 +1015,10 @@ L_MAIN_LOOP_END_NO_PRINT: * * */ - par = 0.5; - for (k = 0; k < m_numComponents; ++k) { + double par = 0.5; + for (size_t k = 0; k < m_numComponents; ++k) { if (m_molNumSpecies_old[k] > 0.0) { - xx = -m_deltaMolNumSpecies[k] / m_molNumSpecies_old[k]; + double xx = -m_deltaMolNumSpecies[k] / m_molNumSpecies_old[k]; if (par < xx) { par = xx; #ifdef DEBUG_MODE @@ -1037,7 +1031,7 @@ L_MAIN_LOOP_END_NO_PRINT: * If we are here, we then do a step which violates element * conservation. */ - iph = m_phaseID[k]; + size_t iph = m_phaseID[k]; m_deltaPhaseMoles[iph] -= m_deltaMolNumSpecies[k]; m_deltaMolNumSpecies[k] = 0.0; } @@ -1055,10 +1049,10 @@ L_MAIN_LOOP_END_NO_PRINT: plogendl(); } #endif - for (i = 0; i < m_numSpeciesTot; ++i) { + for (size_t i = 0; i < m_numSpeciesTot; ++i) { m_deltaMolNumSpecies[i] *= par; } - for (iph = 0; iph < m_numPhases; iph++) { + for (size_t iph = 0; iph < m_numPhases; iph++) { m_deltaPhaseMoles[iph] *= par; } } else { @@ -1089,7 +1083,7 @@ L_MAIN_LOOP_END_NO_PRINT: /* * Calculate the tentative total mole numbers for each phase */ - for (iph = 0; iph < m_numPhases; iph++) { + for (size_t iph = 0; iph < m_numPhases; iph++) { m_tPhaseMoles_new[iph] = m_tPhaseMoles_old[iph] + m_deltaPhaseMoles[iph]; } @@ -1125,7 +1119,7 @@ L_MAIN_LOOP_END_NO_PRINT: plogendl(); } - forced = vcs_globStepDamp(); + bool forced = vcs_globStepDamp(); /* * Print out the changes to the solution that FORCER produced @@ -1136,13 +1130,13 @@ L_MAIN_LOOP_END_NO_PRINT: plogf(" --- FORCER SUBROUTINE changed the solution:\n"); plogf(" --- SPECIES Status INIT MOLES TENT_MOLES"); plogf(" FINAL KMOLES INIT_DEL_G/RT TENT_DEL_G/RT FINAL_DELTA_G/RT\n"); - for (i = 0; i < m_numComponents; ++i) { + for (size_t i = 0; i < m_numComponents; ++i) { plogf(" --- %-12.12s", m_speciesName[i].c_str()); plogf(" %14.6E %14.6E %14.6E\n", m_molNumSpecies_old[i], m_molNumSpecies_old[i] + m_deltaMolNumSpecies[i], m_molNumSpecies_new[i]); } for (kspec = m_numComponents; kspec < m_numSpeciesRdc; ++kspec) { - irxn = kspec - m_numComponents; + size_t irxn = kspec - m_numComponents; plogf(" --- %-12.12s", m_speciesName[kspec].c_str()); plogf(" %2d %14.6E%14.6E%14.6E%14.6E%14.6E%14.6E\n", m_speciesStatus[kspec], m_molNumSpecies_old[kspec], @@ -1186,14 +1180,14 @@ L_MAIN_LOOP_END_NO_PRINT: plogf("\n"); plogf(" --- Species Status Initial_KMoles Final_KMoles Initial_Mu/RT"); plogf(" Mu/RT Init_Del_G/RT Delta_G/RT\n"); - for (i = 0; i < m_numComponents; ++i) { + for (size_t i = 0; i < m_numComponents; ++i) { plogf(" --- %-12.12s", m_speciesName[i].c_str()); plogf(" "); plogf("%14.6E%14.6E%14.6E%14.6E\n", m_molNumSpecies_old[i], m_molNumSpecies_new[i], m_feSpecies_old[i], m_feSpecies_new[i]); } - for (i = m_numComponents; i < m_numSpeciesRdc; ++i) { - l1 = i - m_numComponents; + for (size_t i = m_numComponents; i < m_numSpeciesRdc; ++i) { + size_t l1 = i - m_numComponents; plogf(" --- %-12.12s", m_speciesName[i].c_str()); plogf(" %2d %14.6E%14.6E%14.6E%14.6E%14.6E%14.6E\n", m_speciesStatus[i], m_molNumSpecies_old[i], @@ -1201,7 +1195,7 @@ L_MAIN_LOOP_END_NO_PRINT: m_deltaGRxn_old[l1], m_deltaGRxn_new[l1]); } for (kspec = m_numSpeciesRdc; kspec < m_numSpeciesTot; ++kspec) { - l1 = kspec - m_numComponents; + size_t l1 = kspec - m_numComponents; plogf(" --- %-12.12s", m_speciesName[kspec].c_str()); plogf(" %2d %14.6E%14.6E%14.6E%14.6E%14.6E%14.6E\n", m_speciesStatus[kspec], m_molNumSpecies_old[kspec], @@ -1218,8 +1212,8 @@ L_MAIN_LOOP_END_NO_PRINT: plogf(" --- Phase_Name KMoles(after update)\n"); plogf(" --- "); vcs_print_line("-", 50); - for (iph = 0; iph < m_numPhases; iph++) { - Vphase = m_VolPhaseList[iph]; + for (size_t iph = 0; iph < m_numPhases; iph++) { + vcs_VolPhase* Vphase = m_VolPhaseList[iph]; plogf(" --- %18s = %15.7E\n", Vphase->PhaseName.c_str(), m_tPhaseMoles_new[iph]); } plogf(" "); @@ -1287,18 +1281,11 @@ L_MAIN_LOOP_END_NO_PRINT: * absolute zero. */ justDeletedMultiPhase = false; - for (iph = 0; iph < m_numPhases; iph++) { - Vphase = m_VolPhaseList[iph]; - if (!(Vphase->m_singleSpecies)) { + for (size_t iph = 0; iph < m_numPhases; iph++) { + if (!(m_VolPhaseList[iph]->m_singleSpecies)) { if (m_tPhaseMoles_old[iph] != 0.0 && m_tPhaseMoles_old[iph]/m_totalMolNum <= VCS_DELETE_PHASE_CUTOFF) { soldel = 1; - - for (kspec = 0; kspec < m_numSpeciesRdc; kspec++) { - if (m_phaseID[kspec] == iph && m_molNumSpecies_old[kspec] > 0.0) { - irxn = kspec - m_numComponents; - } - } if (soldel) { #ifdef DEBUG_MODE if (m_debug_print_lvl >= 1) { @@ -1382,7 +1369,7 @@ L_MAIN_LOOP_END_NO_PRINT: * the code logic. */ dofast = (m_numComponents != 1); - for (i = 1; i < m_numComponents; ++i) { + for (size_t i = 1; i < m_numComponents; ++i) { if ((m_molNumSpecies_old[i - 1] * m_spSize[i-1]) < (m_molNumSpecies_old[i] * m_spSize[i])) { dofast = false; break; @@ -1390,10 +1377,10 @@ L_MAIN_LOOP_END_NO_PRINT: } dofast = false; if (dofast) { - for (i = 0; i < m_numRxnRdc; ++i) { - l = m_indexRxnToSpecies[i]; + for (size_t i = 0; i < m_numRxnRdc; ++i) { + size_t l = m_indexRxnToSpecies[i]; if (m_speciesUnknownType[l] != VCS_SPECIES_TYPE_INTERFACIALVOLTAGE) { - for (j = m_numComponents - 1; j != npos; j--) { + for (size_t j = m_numComponents - 1; j != npos; j--) { bool doSwap = false; if (m_SSPhase[j]) { doSwap = (m_molNumSpecies_old[l] * m_spSize[l]) > @@ -1455,10 +1442,10 @@ L_MAIN_LOOP_END_NO_PRINT: } } } else { - for (i = 0; i < m_numRxnRdc; ++i) { - l = m_indexRxnToSpecies[i]; + for (size_t i = 0; i < m_numRxnRdc; ++i) { + size_t l = m_indexRxnToSpecies[i]; if (m_speciesUnknownType[l] != VCS_SPECIES_TYPE_INTERFACIALVOLTAGE) { - for (j = 0; j < m_numComponents; ++j) { + for (size_t j = 0; j < m_numComponents; ++j) { bool doSwap = false; if (m_SSPhase[j]) { doSwap = (m_molNumSpecies_old[l] * m_spSize[l]) > @@ -1541,7 +1528,7 @@ L_MAIN_LOOP_END_NO_PRINT: } #endif m_numRxnMinorZeroed = 0; - for (irxn = 0; irxn < m_numRxnRdc; irxn++) { + for (size_t irxn = 0; irxn < m_numRxnRdc; irxn++) { kspec = m_indexRxnToSpecies[irxn]; int speciesType = vcs_species_type(kspec); @@ -1616,7 +1603,7 @@ L_EQUILIB_CHECK: plogf(" --- Equilibrium check for major species: "); } #endif - for (irxn = 0; irxn < m_numRxnRdc; ++irxn) { + for (size_t irxn = 0; irxn < m_numRxnRdc; ++irxn) { kspec = irxn + m_numComponents; if (m_speciesStatus[kspec] == VCS_SPECIES_MAJOR && (fabs(m_deltaGRxn_new[irxn]) > m_tolmaj)) { if (m_VCount->Its >= maxit) { @@ -1680,7 +1667,7 @@ L_EQUILIB_CHECK: plogf(" --- Equilibrium check for minor species: "); } #endif - for (irxn = 0; irxn < m_numRxnRdc; ++irxn) { + for (size_t irxn = 0; irxn < m_numRxnRdc; ++irxn) { kspec = irxn + m_numComponents; if (m_speciesStatus[kspec] == VCS_SPECIES_MINOR && (fabs(m_deltaGRxn_new[irxn]) > m_tolmin)) { if (m_VCount->Its >= maxit) { @@ -1947,7 +1934,7 @@ L_RETURN_BLOCK_B: if (m_SSPhase[kspec]) { m_molNumSpecies_new[kspec] = 1.0; } else { - iph = m_phaseID[kspec]; + size_t iph = m_phaseID[kspec]; if (m_tPhaseMoles_old[iph] != 0.0) { m_molNumSpecies_new[kspec] = m_molNumSpecies_old[kspec] / m_tPhaseMoles_old[iph]; } else { @@ -1957,9 +1944,8 @@ L_RETURN_BLOCK_B: * This contains the mole fraction that would be true if * the phase just pops into existence. */ - i = m_speciesLocalPhaseIndex[kspec]; - Vphase = m_VolPhaseList[iph]; - m_molNumSpecies_new[kspec] = Vphase->molefraction(i); + size_t i = m_speciesLocalPhaseIndex[kspec]; + m_molNumSpecies_new[kspec] = m_VolPhaseList[iph]->molefraction(i); } } } @@ -2110,7 +2096,6 @@ int VCS_SOLVE::delta_species(const size_t kspec, double* const delta_ptr) { size_t irxn = kspec - m_numComponents; int retn = 1; - double tmp; double delta = *delta_ptr; #ifdef DEBUG_MODE if (kspec < m_numComponents) { @@ -2128,7 +2113,7 @@ int VCS_SOLVE::delta_species(const size_t kspec, double* const delta_ptr) double* sc_irxn = m_stoichCoeffRxnMatrix[irxn]; for (size_t j = 0; j < m_numComponents; ++j) { if (m_molNumSpecies_old[j] > 0.0) { - tmp = sc_irxn[j] * dx; + double tmp = sc_irxn[j] * dx; if (-tmp > m_molNumSpecies_old[j]) { retn = 0; dx = std::min(dx, - m_molNumSpecies_old[j] / sc_irxn[j]); @@ -2155,7 +2140,7 @@ int VCS_SOLVE::delta_species(const size_t kspec, double* const delta_ptr) vcs_setFlagsVolPhase(iph, false, VCS_STATECALC_OLD); for (size_t j = 0; j < m_numComponents; ++j) { - tmp = sc_irxn[j] * dx; + double tmp = sc_irxn[j] * dx; if (tmp != 0.0) { iph = m_phaseID[j]; m_molNumSpecies_old[j] += tmp; @@ -2275,10 +2260,8 @@ int VCS_SOLVE::vcs_delete_species(const size_t kspec) void VCS_SOLVE::vcs_reinsert_deleted(size_t kspec) { - size_t k; // int irxn = kspec - m_numComponents; size_t iph = m_phaseID[kspec]; - double dx; #ifdef DEBUG_MODE if (m_debug_print_lvl >= 2) { plogf(" --- Add back a deleted species: %-12s\n", m_speciesName[kspec].c_str()); @@ -2289,7 +2272,7 @@ void VCS_SOLVE::vcs_reinsert_deleted(size_t kspec) * this adjusts m_molNumSpecies_old[] and m_tPhaseMoles_old[] * HKM -> make this a relative mole number! */ - dx = m_tPhaseMoles_old[iph] * VCS_RELDELETE_SPECIES_CUTOFF * 10.; + double dx = m_tPhaseMoles_old[iph] * VCS_RELDELETE_SPECIES_CUTOFF * 10.; delta_species(kspec, &dx); m_speciesStatus[kspec] = VCS_SPECIES_MINOR; @@ -2314,7 +2297,7 @@ void VCS_SOLVE::vcs_reinsert_deleted(size_t kspec) if (! m_SSPhase[kspec]) { if (Vphase->exists() == VCS_PHASE_EXIST_NO) { Vphase->setExistence(VCS_PHASE_EXIST_YES); - for (k = 0; k < m_numSpeciesTot; k++) { + for (size_t k = 0; k < m_numSpeciesTot; k++) { if (m_phaseID[k] == iph) { if (m_speciesStatus[k] != VCS_SPECIES_DELETED) { m_speciesStatus[k] = VCS_SPECIES_MINOR; @@ -2340,8 +2323,6 @@ void VCS_SOLVE::vcs_reinsert_deleted(size_t kspec) bool VCS_SOLVE::vcs_delete_multiphase(const size_t iph) { - size_t kspec, irxn; - double dx; vcs_VolPhase* Vphase = m_VolPhaseList[iph]; bool successful = true; /* @@ -2357,13 +2338,13 @@ bool VCS_SOLVE::vcs_delete_multiphase(const size_t iph) /* * Loop over all of the species in the phase. */ - for (kspec = m_numComponents; kspec < m_numSpeciesRdc; ++kspec) { + for (size_t kspec = m_numComponents; kspec < m_numSpeciesRdc; ++kspec) { if (m_phaseID[kspec] == iph) { if (m_speciesUnknownType[kspec] != VCS_SPECIES_TYPE_INTERFACIALVOLTAGE) { /* * calculate an extent of rxn, dx, that zeroes out the species. */ - dx = - (m_molNumSpecies_old[kspec]); + double dx = - (m_molNumSpecies_old[kspec]); double dxTent = dx; int retn = delta_species(kspec, &dxTent); @@ -2409,8 +2390,8 @@ bool VCS_SOLVE::vcs_delete_multiphase(const size_t iph) } #endif if (m_molNumSpecies_old[kcomp] != 0.0) { - for (kspec = m_numComponents; kspec < m_numSpeciesRdc; ++kspec) { - irxn = kspec - m_numComponents; + for (size_t kspec = m_numComponents; kspec < m_numSpeciesRdc; ++kspec) { + size_t irxn = kspec - m_numComponents; if (m_phaseID[kspec] != iph) { if (m_stoichCoeffRxnMatrix[irxn][kcomp] != 0.0) { dxWant = -m_molNumSpecies_old[kcomp] / m_stoichCoeffRxnMatrix[irxn][kcomp]; @@ -2464,9 +2445,8 @@ bool VCS_SOLVE::vcs_delete_multiphase(const size_t iph) * existence or not is checked as part of the main iteration * loop. */ - for (kspec = m_numSpeciesRdc; kspec < m_numSpeciesTot; ++kspec) { + for (size_t kspec = m_numSpeciesRdc; kspec < m_numSpeciesTot; ++kspec) { if (m_phaseID[kspec] == iph) { - irxn = kspec - m_numComponents; m_molNumSpecies_old[kspec] = 0.0; m_molNumSpecies_new[kspec] = 0.0; m_deltaMolNumSpecies[kspec] = 0.0; @@ -2506,8 +2486,6 @@ bool VCS_SOLVE::vcs_delete_multiphase(const size_t iph) int VCS_SOLVE::vcs_recheck_deleted() { - int npb; - size_t iph, kspec, irxn; double* xtcutoff = VCS_DATA_PTR(m_TmpPhase); #ifdef DEBUG_MODE if (m_debug_print_lvl >= 2) { @@ -2523,8 +2501,8 @@ int VCS_SOLVE::vcs_recheck_deleted() * for formation reactions. * Note: fe[] here includes everything except for the ln(x[i]) term */ - for (kspec = m_numSpeciesRdc; kspec < m_numSpeciesTot; ++kspec) { - iph = m_phaseID[kspec]; + for (size_t kspec = m_numSpeciesRdc; kspec < m_numSpeciesTot; ++kspec) { + size_t iph = m_phaseID[kspec]; m_feSpecies_new[kspec] = (m_SSfeSpecies[kspec] + log(m_actCoeffSpecies_old[kspec]) - m_lnMnaughtSpecies[kspec] + m_chargeSpecies[kspec] * m_Faraday_dim * m_phasePhi[iph]); @@ -2536,7 +2514,7 @@ int VCS_SOLVE::vcs_recheck_deleted() */ vcs_deltag(0, true, VCS_STATECALC_NEW); - for (iph = 0; iph < m_numPhases; iph++) { + for (size_t iph = 0; iph < m_numPhases; iph++) { if (m_tPhaseMoles_old[iph] > 0.0) { xtcutoff[iph] = log(1.0 / VCS_RELDELETE_SPECIES_CUTOFF); } else { @@ -2570,10 +2548,10 @@ int VCS_SOLVE::vcs_recheck_deleted() * Thus, we need to amend the code. Also nonideal solutions will tend to * complicate matters severely also. */ - npb = 0; - for (irxn = m_numRxnRdc; irxn < m_numRxnTot; ++irxn) { - kspec = m_indexRxnToSpecies[irxn]; - iph = m_phaseID[kspec]; + int npb = 0; + for (size_t irxn = m_numRxnRdc; irxn < m_numRxnTot; ++irxn) { + size_t kspec = m_indexRxnToSpecies[irxn]; + size_t iph = m_phaseID[kspec]; if (m_tPhaseMoles_old[iph] == 0.0) { if (m_deltaGRxn_new[irxn] < 0.0) { vcs_reinsert_deleted(kspec); @@ -2602,10 +2580,9 @@ bool VCS_SOLVE::recheck_deleted_phase(const int iphase) if (Vphase->exists() == VCS_PHASE_EXIST_ZEROEDPHASE) { return false; } - size_t irxn, kspec; if (Vphase->m_singleSpecies) { - kspec = Vphase->spGlobalIndexVCS(0); - irxn = kspec + m_numComponents; + size_t kspec = Vphase->spGlobalIndexVCS(0); + size_t irxn = kspec + m_numComponents; if (m_deltaGRxn_old[irxn] < 0.0) { return true; } @@ -2614,8 +2591,8 @@ bool VCS_SOLVE::recheck_deleted_phase(const int iphase) double phaseDG = 1.0; for (size_t kk = 0; kk < Vphase->nSpecies(); kk++) { - kspec = Vphase->spGlobalIndexVCS(kk); - irxn = kspec + m_numComponents; + size_t kspec = Vphase->spGlobalIndexVCS(kk); + size_t irxn = kspec + m_numComponents; if (m_deltaGRxn_old[irxn] > 50.0) { m_deltaGRxn_old[irxn] = 50.0; } @@ -2633,7 +2610,7 @@ bool VCS_SOLVE::recheck_deleted_phase(const int iphase) size_t VCS_SOLVE::vcs_add_all_deleted() { - size_t iph, kspec, retn; + size_t retn; if (m_numSpeciesRdc == m_numSpeciesTot) { return 0; } @@ -2646,8 +2623,8 @@ size_t VCS_SOLVE::vcs_add_all_deleted() vcs_dcopy(VCS_DATA_PTR(m_molNumSpecies_new), VCS_DATA_PTR(m_molNumSpecies_old), m_numSpeciesTot); for (int cits = 0; cits < 3; cits++) { - for (kspec = m_numSpeciesRdc; kspec < m_numSpeciesTot; ++kspec) { - iph = m_phaseID[kspec]; + for (size_t kspec = m_numSpeciesRdc; kspec < m_numSpeciesTot; ++kspec) { + size_t iph = m_phaseID[kspec]; vcs_VolPhase* Vphase = m_VolPhaseList[iph]; if (m_molNumSpecies_new[kspec] == 0.0) { m_molNumSpecies_new[kspec] = VCS_DELETE_MINORSPECIES_CUTOFF * 1.0E-10; @@ -2663,8 +2640,8 @@ size_t VCS_SOLVE::vcs_add_all_deleted() */ vcs_deltag(0, true, VCS_STATECALC_NEW); for (size_t irxn = m_numRxnRdc; irxn < m_numRxnTot; ++irxn) { - kspec = m_indexRxnToSpecies[irxn]; - iph = m_phaseID[kspec]; + size_t kspec = m_indexRxnToSpecies[irxn]; + size_t iph = m_phaseID[kspec]; if (m_tPhaseMoles_old[iph] > 0.0) { double maxDG = std::min(m_deltaGRxn_new[irxn], 690.0); double dx = m_tPhaseMoles_old[iph] * exp(- maxDG); @@ -2673,8 +2650,8 @@ size_t VCS_SOLVE::vcs_add_all_deleted() } } for (size_t irxn = m_numRxnRdc; irxn < m_numRxnTot; ++irxn) { - kspec = m_indexRxnToSpecies[irxn]; - iph = m_phaseID[kspec]; + size_t kspec = m_indexRxnToSpecies[irxn]; + size_t iph = m_phaseID[kspec]; if (m_tPhaseMoles_old[iph] > 0.0) { double dx = m_molNumSpecies_new[kspec]; retn = delta_species(kspec, &dx); @@ -2718,8 +2695,8 @@ size_t VCS_SOLVE::vcs_add_all_deleted() retn = 0; for (size_t irxn = m_numRxnRdc; irxn < m_numRxnTot; ++irxn) { - kspec = m_indexRxnToSpecies[irxn]; - iph = m_phaseID[kspec]; + size_t kspec = m_indexRxnToSpecies[irxn]; + size_t iph = m_phaseID[kspec]; if (m_tPhaseMoles_old[iph] > 0.0) { if (fabs(m_deltaGRxn_old[irxn]) > m_tolmin) { if (((m_molNumSpecies_old[kspec] * exp(-m_deltaGRxn_old[irxn])) > @@ -2743,16 +2720,14 @@ size_t VCS_SOLVE::vcs_add_all_deleted() bool VCS_SOLVE::vcs_globStepDamp() { - double s1, s2, al; - size_t irxn, kspec, iph; double* dptr = VCS_DATA_PTR(m_deltaGRxn_new); /* *************************************************** */ /* **** CALCULATE SLOPE AT END OF THE STEP ********** */ /* *************************************************** */ - s2 = 0.0; - for (irxn = 0; irxn < m_numRxnRdc; ++irxn) { - kspec = irxn + m_numComponents; + double s2 = 0.0; + for (size_t irxn = 0; irxn < m_numRxnRdc; ++irxn) { + size_t kspec = irxn + m_numComponents; if (m_speciesUnknownType[kspec] != VCS_SPECIES_TYPE_INTERFACIALVOLTAGE) { s2 += dptr[irxn] * m_deltaMolNumSpecies[kspec]; } @@ -2762,10 +2737,10 @@ bool VCS_SOLVE::vcs_globStepDamp() /* *************************************************** */ /* **** CALCULATE ORIGINAL SLOPE ********************* */ /* ************************************************** */ - s1 = 0.0; + double s1 = 0.0; dptr = VCS_DATA_PTR(m_deltaGRxn_old); - for (irxn = 0; irxn < m_numRxnRdc; ++irxn) { - kspec = irxn + m_numComponents; + for (size_t irxn = 0; irxn < m_numRxnRdc; ++irxn) { + size_t kspec = irxn + m_numComponents; if (m_speciesUnknownType[kspec] != VCS_SPECIES_TYPE_INTERFACIALVOLTAGE) { s1 += dptr[irxn] * m_deltaMolNumSpecies[kspec]; } @@ -2806,7 +2781,7 @@ bool VCS_SOLVE::vcs_globStepDamp() /* *************************************************** */ /* **** FIT PCJ2822ARABOLA ********************************* */ /* *************************************************** */ - al = 1.0; + double al = 1.0; if (fabs(s1 -s2) > 1.0E-200) { al = s1 / (s1 - s2); } @@ -2835,11 +2810,11 @@ bool VCS_SOLVE::vcs_globStepDamp() #endif dptr = VCS_DATA_PTR(m_molNumSpecies_new); - for (kspec = 0; kspec < m_numSpeciesRdc; ++kspec) { + for (size_t kspec = 0; kspec < m_numSpeciesRdc; ++kspec) { m_molNumSpecies_new[kspec] = m_molNumSpecies_old[kspec] + al * m_deltaMolNumSpecies[kspec]; } - for (iph = 0; iph < m_numPhases; iph++) { + for (size_t iph = 0; iph < m_numPhases; iph++) { m_tPhaseMoles_new[iph] = m_tPhaseMoles_old[iph] + al * m_deltaPhaseMoles[iph]; } vcs_updateVP(VCS_STATECALC_NEW); @@ -2867,8 +2842,8 @@ bool VCS_SOLVE::vcs_globStepDamp() dptr = VCS_DATA_PTR(m_deltaGRxn_new); s2 = 0.0; - for (irxn = 0; irxn < m_numRxnRdc; ++irxn) { - kspec = irxn + m_numComponents; + for (size_t irxn = 0; irxn < m_numRxnRdc; ++irxn) { + size_t kspec = irxn + m_numComponents; if (m_speciesUnknownType[kspec] != VCS_SPECIES_TYPE_INTERFACIALVOLTAGE) { s2 += dptr[irxn] * m_deltaMolNumSpecies[kspec]; } @@ -2887,9 +2862,8 @@ bool VCS_SOLVE::vcs_globStepDamp() int VCS_SOLVE::vcs_basopt(const bool doJustComponents, double aw[], double sa[], double sm[], double ss[], double test, bool* const usedZeroedSpecies) { - size_t j, k, l, i, jl, ml, jr, irxn, kspec; + size_t k; bool lindep; - size_t ncTrial; size_t juse = npos; size_t jlose = npos; double* dptr, *scrxn_ptr; @@ -2897,7 +2871,7 @@ int VCS_SOLVE::vcs_basopt(const bool doJustComponents, double aw[], double sa[], #ifdef DEBUG_MODE if (m_debug_print_lvl >= 2) { plogf(" "); - for (i=0; i<77; i++) { + for (size_t i=0; i<77; i++) { plogf("-"); } plogf("\n"); @@ -2911,12 +2885,12 @@ int VCS_SOLVE::vcs_basopt(const bool doJustComponents, double aw[], double sa[], plogf("\n"); plogf(" --- Formula Matrix used in BASOPT calculation\n"); plogf(" --- Active | "); - for (j = 0; j < m_numElemConstraints; j++) { + for (size_t j = 0; j < m_numElemConstraints; j++) { plogf(" %1d ", m_elementActive[j]); } plogf("\n"); plogf(" --- Species | "); - for (j = 0; j < m_numElemConstraints; j++) { + for (size_t j = 0; j < m_numElemConstraints; j++) { plogf(" "); vcs_print_stringTrunc(m_elementName[j].c_str(), 8, 1); } @@ -2925,7 +2899,7 @@ int VCS_SOLVE::vcs_basopt(const bool doJustComponents, double aw[], double sa[], plogf(" --- "); vcs_print_stringTrunc(m_speciesName[k].c_str(), 11, 1); plogf(" | "); - for (j = 0; j < m_numElemConstraints; j++) { + for (size_t j = 0; j < m_numElemConstraints; j++) { plogf(" %8.2g", m_formulaMatrix[j][k]); } plogf("\n"); @@ -2940,7 +2914,7 @@ int VCS_SOLVE::vcs_basopt(const bool doJustComponents, double aw[], double sa[], * It's equal to the minimum of the number of elements and the * number of total species. */ - ncTrial = std::min(m_numElemConstraints, m_numSpeciesTot); + size_t ncTrial = std::min(m_numElemConstraints, m_numSpeciesTot); m_numComponents = ncTrial; *usedZeroedSpecies = false; vector_int ipiv(ncTrial); @@ -2959,7 +2933,7 @@ int VCS_SOLVE::vcs_basopt(const bool doJustComponents, double aw[], double sa[], } } - jr = npos; + size_t jr = npos; /* * Top of a loop of some sort based on the index JR. JR is the * current number of component species found. @@ -3024,7 +2998,7 @@ int VCS_SOLVE::vcs_basopt(const bool doJustComponents, double aw[], double sa[], size_t kfound = npos; int minNonZeroes = 100000; int nonZeroesKspec = 0; - for (kspec = ncTrial; kspec < m_numSpeciesTot; kspec++) { + for (size_t kspec = ncTrial; kspec < m_numSpeciesTot; kspec++) { if (aw[kspec] >= 0.0) { if (m_speciesUnknownType[kspec] != VCS_SPECIES_TYPE_INTERFACIALVOLTAGE) { maxConcPossKspec = 1.0E10; @@ -3064,9 +3038,9 @@ int VCS_SOLVE::vcs_basopt(const bool doJustComponents, double aw[], double sa[], if (kfound == npos) { double gmin = 0.0; kfound = k; - for (kspec = ncTrial; kspec < m_numSpeciesTot; kspec++) { + for (size_t kspec = ncTrial; kspec < m_numSpeciesTot; kspec++) { if (aw[kspec] >= 0.0) { - irxn = kspec - ncTrial; + size_t irxn = kspec - ncTrial; if (m_deltaGRxn_new[irxn] < gmin) { gmin = m_deltaGRxn_new[irxn]; kfound = kspec; @@ -3089,7 +3063,7 @@ int VCS_SOLVE::vcs_basopt(const bool doJustComponents, double aw[], double sa[], m_numRxnTot = m_numSpeciesTot - ncTrial; m_numRxnRdc = m_numRxnTot - numPreDeleted; m_numSpeciesRdc = m_numSpeciesTot - numPreDeleted; - for (i = 0; i < m_numSpeciesTot; ++i) { + for (size_t i = 0; i < m_numSpeciesTot; ++i) { m_indexRxnToSpecies[i] = ncTrial + i; } #ifdef DEBUG_MODE @@ -3112,8 +3086,8 @@ int VCS_SOLVE::vcs_basopt(const bool doJustComponents, double aw[], double sa[], * Modified Gram-Schmidt Method, p. 202 Dalquist * QR factorization of a matrix without row pivoting. */ - jl = jr; - for (j = 0; j < m_numElemConstraints; ++j) { + size_t jl = jr; + for (size_t j = 0; j < m_numElemConstraints; ++j) { sm[j + jr*m_numElemConstraints] = m_formulaMatrix[j][k]; } if (jl > 0) { @@ -3123,9 +3097,9 @@ int VCS_SOLVE::vcs_basopt(const bool doJustComponents, double aw[], double sa[], * (this is slightly different than Dalquist) * R_JA_JA = 1 */ - for (j = 0; j < jl; ++j) { + for (size_t j = 0; j < jl; ++j) { ss[j] = 0.0; - for (i = 0; i < m_numElemConstraints; ++i) { + for (size_t i = 0; i < m_numElemConstraints; ++i) { ss[j] += sm[i + jr*m_numElemConstraints] * sm[i + j*m_numElemConstraints]; } ss[j] /= sa[j]; @@ -3134,8 +3108,8 @@ int VCS_SOLVE::vcs_basopt(const bool doJustComponents, double aw[], double sa[], * Now make the new column, (*,JR), orthogonal to the * previous columns */ - for (j = 0; j < jl; ++j) { - for (l = 0; l < m_numElemConstraints; ++l) { + for (size_t j = 0; j < jl; ++j) { + for (size_t l = 0; l < m_numElemConstraints; ++l) { sm[l + jr*m_numElemConstraints] -= ss[j] * sm[l + j*m_numElemConstraints]; } } @@ -3145,7 +3119,7 @@ int VCS_SOLVE::vcs_basopt(const bool doJustComponents, double aw[], double sa[], * It will be used in the denominator in future row calcs. */ sa[jr] = 0.0; - for (ml = 0; ml < m_numElemConstraints; ++ml) { + for (size_t ml = 0; ml < m_numElemConstraints; ++ml) { sa[jr] += SQUARE(sm[ml + jr*m_numElemConstraints]); } /* **************************************************** */ @@ -3237,14 +3211,14 @@ L_END_LOOP: * this program with ne being the species and nc being the * elements!! */ - for (j = 0; j < ncTrial; ++j) { - for (i = 0; i < ncTrial; ++i) { + for (size_t j = 0; j < ncTrial; ++j) { + for (size_t i = 0; i < ncTrial; ++i) { sm[i + j*m_numElemConstraints] = m_formulaMatrix[i][j]; } } - for (i = 0; i < m_numRxnTot; ++i) { + for (size_t i = 0; i < m_numRxnTot; ++i) { k = m_indexRxnToSpecies[i]; - for (j = 0; j < ncTrial; ++j) { + for (size_t j = 0; j < ncTrial; ++j) { m_stoichCoeffRxnMatrix[i][j] = - m_formulaMatrix[j][k]; } } @@ -3265,14 +3239,14 @@ L_END_LOOP: */ juse = npos; jlose = npos; - for (j = 0; j < m_numElemConstraints; j++) { + for (size_t j = 0; j < m_numElemConstraints; j++) { if (!(m_elementActive[j])) { if (!strcmp((m_elementName[j]).c_str(), "E")) { juse = j; } } } - for (j = 0; j < m_numElemConstraints; j++) { + for (size_t j = 0; j < m_numElemConstraints; j++) { if (m_elementActive[j]) { if (!strncmp((m_elementName[j]).c_str(), "cn_", 3)) { jlose = j; @@ -3282,8 +3256,8 @@ L_END_LOOP: for (k = 0; k < m_numSpeciesTot; k++) { if (m_speciesUnknownType[k] == VCS_SPECIES_TYPE_INTERFACIALVOLTAGE) { - for (j = 0; j < ncTrial; ++j) { - for (i = 0; i < ncTrial; ++i) { + for (size_t j = 0; j < ncTrial; ++j) { + for (size_t i = 0; i < ncTrial; ++i) { if (i == jlose) { sm[i + j*m_numElemConstraints] = m_formulaMatrix[juse][j]; } else { @@ -3291,9 +3265,9 @@ L_END_LOOP: } } } - for (i = 0; i < m_numRxnTot; ++i) { + for (size_t i = 0; i < m_numRxnTot; ++i) { k = m_indexRxnToSpecies[i]; - for (j = 0; j < ncTrial; ++j) { + for (size_t j = 0; j < ncTrial; ++j) { if (j == jlose) { aw[j] = - m_formulaMatrix[juse][k]; } else { @@ -3309,8 +3283,8 @@ L_END_LOOP: } ct_dgetrs(ctlapack::NoTranspose, ncTrial, 1, sm, m_numElemConstraints, &ipiv[0], aw, m_numElemConstraints, info); - i = k - ncTrial; - for (j = 0; j < ncTrial; j++) { + size_t i = k - ncTrial; + for (size_t j = 0; j < ncTrial; j++) { m_stoichCoeffRxnMatrix[i][j] = aw[j]; } } @@ -3320,9 +3294,9 @@ L_END_LOOP: /* * Calculate the szTmp array for each formation reaction */ - for (i = 0; i < m_numRxnTot; i++) { + for (size_t i = 0; i < m_numRxnTot; i++) { double szTmp = 0.0; - for (j = 0; j < ncTrial; j++) { + for (size_t j = 0; j < ncTrial; j++) { szTmp += fabs(m_stoichCoeffRxnMatrix[i][j]); } m_scSize[i] = szTmp; @@ -3332,11 +3306,11 @@ L_END_LOOP: #ifdef DEBUG_MODE if (m_debug_print_lvl >= 2) { plogf(" --- Components:"); - for (j = 0; j < ncTrial; j++) { + for (size_t j = 0; j < ncTrial; j++) { plogf(" %3d", j); } plogf("\n --- Components Moles:"); - for (j = 0; j < ncTrial; j++) { + for (size_t j = 0; j < ncTrial; j++) { if (m_speciesUnknownType[j] == VCS_SPECIES_TYPE_INTERFACIALVOLTAGE) { plogf(" % -10.3E", 0.0); } else { @@ -3344,12 +3318,12 @@ L_END_LOOP: } } plogf("\n --- NonComponent| Moles |"); - for (j = 0; j < ncTrial; j++) { + for (size_t j = 0; j < ncTrial; j++) { plogf(" %10.10s", m_speciesName[j].c_str()); } //plogf("| m_scSize"); plogf("\n"); - for (i = 0; i < m_numRxnTot; i++) { + for (size_t i = 0; i < m_numRxnTot; i++) { plogf(" --- %3d ", m_indexRxnToSpecies[i]); plogf("%-10.10s", m_speciesName[m_indexRxnToSpecies[i]].c_str()); if (m_speciesUnknownType[m_indexRxnToSpecies[i]] == VCS_SPECIES_TYPE_INTERFACIALVOLTAGE) { @@ -3357,7 +3331,7 @@ L_END_LOOP: } else { plogf("|% -10.3E|", m_molNumSpecies_old[m_indexRxnToSpecies[i]]); } - for (j = 0; j < ncTrial; j++) { + for (size_t j = 0; j < ncTrial; j++) { plogf(" %+7.3f", m_stoichCoeffRxnMatrix[i][j]); } //plogf(" | %6.2f", m_scSize[i]); @@ -3369,24 +3343,23 @@ L_END_LOOP: * Manual check on the satisfaction of the reaction matrix's ability * to conserve elements */ - double sum; double sumMax = -1.0; int iMax = -1; int jMax = -1; - size_t n; - for (i = 0; i < m_numRxnTot; ++i) { + for (size_t i = 0; i < m_numRxnTot; ++i) { k = m_indexRxnToSpecies[i]; - for (j = 0; j < ncTrial; ++j) { + double sum; + for (size_t j = 0; j < ncTrial; ++j) { if (j == jlose) { sum = m_formulaMatrix[juse][k]; - for (n = 0; n < ncTrial; n++) { + for (size_t n = 0; n < ncTrial; n++) { double numElements = m_formulaMatrix[juse][n]; double coeff = m_stoichCoeffRxnMatrix[i][n]; sum += coeff * numElements; } } else { sum = m_formulaMatrix[j][k]; - for (n = 0; n < ncTrial; n++) { + for (size_t n = 0; n < ncTrial; n++) { double numElements = m_formulaMatrix[j][n]; double coeff = m_stoichCoeffRxnMatrix[i][n]; sum += coeff * numElements; @@ -3412,7 +3385,7 @@ L_END_LOOP: plogf("%-5.5s", m_elementName[jMax].c_str()); plogf("\n"); plogf(" "); - for (i=0; i<77; i++) { + for (size_t i=0; i<77; i++) { plogf("-"); } plogf("\n"); @@ -3436,15 +3409,15 @@ L_END_LOOP: * array, m_deltaMolNumPhase[irxn][iphase], * and the phase participation array, PhaseParticipation[irxn][iphase] */ - for (irxn = 0; irxn < m_numRxnTot; ++irxn) { + for (size_t irxn = 0; irxn < m_numRxnTot; ++irxn) { scrxn_ptr = m_stoichCoeffRxnMatrix[irxn]; dptr = m_deltaMolNumPhase[irxn]; - kspec = m_indexRxnToSpecies[irxn]; + size_t kspec = m_indexRxnToSpecies[irxn]; size_t iph = m_phaseID[kspec]; int* pp_ptr = m_phaseParticipation[irxn]; dptr[iph] = 1.0; pp_ptr[iph]++; - for (j = 0; j < ncTrial; ++j) { + for (size_t j = 0; j < ncTrial; ++j) { iph = m_phaseID[j]; if (fabs(scrxn_ptr[j]) <= 1.0e-6) { scrxn_ptr[j] = 0.0; @@ -3722,20 +3695,11 @@ void VCS_SOLVE::vcs_chemPotPhase(const int stateCalc, double* const ac, double* const mu_i, const bool do_deleted) { - vcs_VolPhase* Vphase = m_VolPhaseList[iph]; size_t nkk = Vphase->nSpecies(); - size_t k, kspec; - -#ifdef DEBUG_MODE - //if (m_debug_print_lvl >= 2) { - // plogf(" --- Subroutine vcs_chemPotPhase called for phase %d\n", - // iph); - //} -#endif double tMoles = TPhInertMoles[iph]; - for (k = 0; k < nkk; k++) { - kspec = Vphase->spGlobalIndexVCS(k); + for (size_t k = 0; k < nkk; k++) { + size_t kspec = Vphase->spGlobalIndexVCS(k); tMoles += molNum[kspec]; } double tlogMoles = 0.0; @@ -3749,8 +3713,8 @@ void VCS_SOLVE::vcs_chemPotPhase(const int stateCalc, double phi = Vphase->electricPotential(); double Faraday_phi = m_Faraday_dim * phi; - for (k = 0; k < nkk; k++) { - kspec = Vphase->spGlobalIndexVCS(k); + for (size_t k = 0; k < nkk; k++) { + size_t kspec = Vphase->spGlobalIndexVCS(k); if (kspec >= m_numComponents) { if (!do_deleted && (m_speciesStatus[kspec] == VCS_SPECIES_DELETED)) { @@ -3786,13 +3750,8 @@ void VCS_SOLVE::vcs_chemPotPhase(const int stateCalc, void VCS_SOLVE::vcs_dfe(const int stateCalc, const int ll, const size_t lbot, const size_t ltop) { - size_t l1, l2, iph, kspec, irxn; - size_t iphase; double* tPhMoles_ptr=0; double* actCoeff_ptr=0; - double* tlogMoles=0; - vcs_VolPhase* Vphase; - double* feSpecies=0; double* molNum=0; if (stateCalc == VCS_STATECALC_OLD) { @@ -3843,24 +3802,24 @@ void VCS_SOLVE::vcs_dfe(const int stateCalc, } #endif - tlogMoles = VCS_DATA_PTR(m_TmpPhase); + double* tlogMoles = VCS_DATA_PTR(m_TmpPhase); /* * Might as well recalculate the phase mole vector * and compare to the stored one. They should be correct. */ double* tPhInertMoles = VCS_DATA_PTR(TPhInertMoles); - for (iph = 0; iph < m_numPhases; iph++) { + for (size_t iph = 0; iph < m_numPhases; iph++) { tlogMoles[iph] = tPhInertMoles[iph]; } - for (kspec = 0; kspec < m_numSpeciesTot; kspec++) { + for (size_t kspec = 0; kspec < m_numSpeciesTot; kspec++) { if (m_speciesUnknownType[kspec] != VCS_SPECIES_TYPE_INTERFACIALVOLTAGE) { - iph = m_phaseID[kspec]; + size_t iph = m_phaseID[kspec]; tlogMoles[iph] += molNum[kspec]; } } #ifdef DEBUG_MODE - for (iph = 0; iph < m_numPhases; iph++) { + for (size_t iph = 0; iph < m_numPhases; iph++) { if (! vcs_doubleEqual(tlogMoles[iph], tPhMoles_ptr[iph])) { plogf("phase Moles may be off, iph = %d, %20.14g %20.14g \n", iph, tlogMoles[iph], tPhMoles_ptr[iph]); @@ -3869,13 +3828,13 @@ void VCS_SOLVE::vcs_dfe(const int stateCalc, } #endif vcs_dzero(tlogMoles, m_numPhases); - for (iph = 0; iph < m_numPhases; iph++) { + for (size_t iph = 0; iph < m_numPhases; iph++) { if (tPhMoles_ptr[iph] > 0.0) { tlogMoles[iph] = log(tPhMoles_ptr[iph]); } } - + size_t l1, l2; if (ll != 0) { l1 = lbot; l2 = m_numComponents; @@ -3889,8 +3848,8 @@ void VCS_SOLVE::vcs_dfe(const int stateCalc, * not current. Here we also trigger an update check for each * VolPhase to see if its mole numbers are current with vcs */ - for (iphase = 0; iphase < m_numPhases; iphase++) { - Vphase = m_VolPhaseList[iphase]; + for (size_t iphase = 0; iphase < m_numPhases; iphase++) { + vcs_VolPhase* Vphase = m_VolPhaseList[iphase]; Vphase->updateFromVCS_MoleNumbers(stateCalc); if (!Vphase->m_singleSpecies) { Vphase->sendToVCS_ActCoeff(stateCalc, VCS_DATA_PTR(actCoeff_ptr)); @@ -3906,8 +3865,8 @@ void VCS_SOLVE::vcs_dfe(const int stateCalc, * finish up below with loops over either the major noncomponent * species or the minor noncomponent species. */ - for (kspec = l1; kspec < l2; ++kspec) { - iphase = m_phaseID[kspec]; + for (size_t kspec = l1; kspec < l2; ++kspec) { + size_t iphase = m_phaseID[kspec]; if (m_speciesUnknownType[kspec] == VCS_SPECIES_TYPE_INTERFACIALVOLTAGE) { #ifdef DEBUG_MODE if (molNum[kspec] != m_phasePhi[iphase]) { @@ -3931,7 +3890,7 @@ void VCS_SOLVE::vcs_dfe(const int stateCalc, + m_chargeSpecies[kspec] * m_Faraday_dim * m_phasePhi[iphase]; } else { if (molNum[kspec] <= VCS_DELETE_MINORSPECIES_CUTOFF) { - iph = m_phaseID[kspec]; + size_t iph = m_phaseID[kspec]; if (tPhMoles_ptr[iph] > 0.0) { feSpecies[kspec] = m_SSfeSpecies[kspec] + log(actCoeff_ptr[kspec] * VCS_DELETE_MINORSPECIES_CUTOFF) @@ -3954,10 +3913,10 @@ void VCS_SOLVE::vcs_dfe(const int stateCalc, /* **** MAJORS ONLY ******************************* */ /* ************************************************ */ if (ll < 0) { - for (irxn = 0; irxn < m_numRxnRdc; ++irxn) { - kspec = m_indexRxnToSpecies[irxn]; + for (size_t irxn = 0; irxn < m_numRxnRdc; ++irxn) { + size_t kspec = m_indexRxnToSpecies[irxn]; if (m_speciesStatus[kspec] != VCS_SPECIES_MINOR) { - iphase = m_phaseID[kspec]; + size_t iphase = m_phaseID[kspec]; if (m_speciesUnknownType[kspec] == VCS_SPECIES_TYPE_INTERFACIALVOLTAGE) { #ifdef DEBUG_MODE if (molNum[kspec] != m_phasePhi[iphase]) { @@ -3981,7 +3940,7 @@ void VCS_SOLVE::vcs_dfe(const int stateCalc, + m_chargeSpecies[kspec] * m_Faraday_dim * m_phasePhi[iphase]; } else { if (molNum[kspec] <= VCS_DELETE_MINORSPECIES_CUTOFF) { - iph = m_phaseID[kspec]; + size_t iph = m_phaseID[kspec]; if (tPhMoles_ptr[iph] > 0.0) { feSpecies[kspec] = m_SSfeSpecies[kspec] + log(actCoeff_ptr[kspec] * VCS_DELETE_MINORSPECIES_CUTOFF) @@ -4005,10 +3964,10 @@ void VCS_SOLVE::vcs_dfe(const int stateCalc, /* **** MINORS ONLY ******************************* */ /* ************************************************ */ } else if (ll > 0) { - for (irxn = 0; irxn < m_numRxnRdc; ++irxn) { - kspec = m_indexRxnToSpecies[irxn]; + for (size_t irxn = 0; irxn < m_numRxnRdc; ++irxn) { + size_t kspec = m_indexRxnToSpecies[irxn]; if (m_speciesStatus[kspec] == VCS_SPECIES_MINOR) { - iphase = m_phaseID[kspec]; + size_t iphase = m_phaseID[kspec]; if (m_speciesUnknownType[kspec] == VCS_SPECIES_TYPE_INTERFACIALVOLTAGE) { #ifdef DEBUG_MODE if (molNum[kspec] != m_phasePhi[iphase]) { @@ -4032,7 +3991,7 @@ void VCS_SOLVE::vcs_dfe(const int stateCalc, + m_chargeSpecies[kspec] * m_Faraday_dim * m_phasePhi[iphase]; } else { if (molNum[kspec] <= VCS_DELETE_MINORSPECIES_CUTOFF) { - iph = m_phaseID[kspec]; + size_t iph = m_phaseID[kspec]; if (tPhMoles_ptr[iph] > 0.0) { feSpecies[kspec] = m_SSfeSpecies[kspec] + log(actCoeff_ptr[kspec] * VCS_DELETE_MINORSPECIES_CUTOFF) @@ -4061,7 +4020,6 @@ void VCS_SOLVE::vcs_printSpeciesChemPot(const int stateCalc) const { double mfValue = 1.0; bool zeroedPhase = false; - size_t kspec; const double* molNum = VCS_DATA_PTR(m_molNumSpecies_old); const double* actCoeff_ptr = VCS_DATA_PTR(m_actCoeffSpecies_old); @@ -4075,7 +4033,7 @@ void VCS_SOLVE::vcs_printSpeciesChemPot(const int stateCalc) const for (size_t iph = 0; iph < m_numPhases; iph++) { tMoles[iph] = tPhInertMoles[iph]; } - for (kspec = 0; kspec < m_numSpeciesTot; kspec++) { + for (size_t kspec = 0; kspec < m_numSpeciesTot; kspec++) { if (m_speciesUnknownType[kspec] != VCS_SPECIES_TYPE_INTERFACIALVOLTAGE) { size_t iph = m_phaseID[kspec]; tMoles[iph] += molNum[kspec]; @@ -4088,7 +4046,7 @@ void VCS_SOLVE::vcs_printSpeciesChemPot(const int stateCalc) const printf(" "); vcs_print_line("-", 132); - for (kspec = 0; kspec < m_numSpeciesTot; ++kspec) { + for (size_t kspec = 0; kspec < m_numSpeciesTot; ++kspec) { mfValue = 1.0; size_t iphase = m_phaseID[kspec]; const vcs_VolPhase* Vphase = m_VolPhaseList[iphase]; @@ -4136,22 +4094,20 @@ void VCS_SOLVE::vcs_printSpeciesChemPot(const int stateCalc) const #ifdef DEBUG_MODE void VCS_SOLVE::prneav() const { - size_t j; - bool kerr; std::vector eav(m_numElemConstraints, 0.0); - for (j = 0; j < m_numElemConstraints; ++j) { + for (size_t j = 0; j < m_numElemConstraints; ++j) { for (size_t i = 0; i < m_numSpeciesTot; ++i) { if (m_speciesUnknownType[i] == VCS_SPECIES_TYPE_INTERFACIALVOLTAGE) { eav[j] += m_formulaMatrix[j][i] * m_molNumSpecies_old[i]; } } } - kerr = false; + bool kerr = false; plogf("--------------------------------------------------"); plogf("ELEMENT ABUNDANCE VECTOR:\n"); plogf(" Element Now Orignal Deviation Type\n"); - for (j = 0; j < m_numElemConstraints; ++j) { + for (size_t j = 0; j < m_numElemConstraints; ++j) { plogf(" "); plogf("%-2.2s", (m_elementName[j]).c_str()); plogf(" = %15.6E %15.6E %15.6E %3d\n", @@ -4176,12 +4132,11 @@ void VCS_SOLVE::prneav() const double VCS_SOLVE::l2normdg(double dgLocal[]) const { - double tmp; - size_t irxn; if (m_numRxnRdc <= 0) { return 0.0; } - for (irxn = 0, tmp = 0.0; irxn < m_numRxnRdc; ++irxn) { + double tmp = 0; + for (size_t irxn = 0; irxn < m_numRxnRdc; ++irxn) { size_t kspec = irxn + m_numComponents; if (m_speciesStatus[kspec] == VCS_SPECIES_MAJOR || m_speciesStatus[kspec] == VCS_SPECIES_MINOR || dgLocal[irxn] < 0.0) { @@ -4195,8 +4150,6 @@ double VCS_SOLVE::l2normdg(double dgLocal[]) const double VCS_SOLVE::vcs_tmoles() { - double sum; - vcs_VolPhase* Vphase; for (size_t i = 0; i < m_numPhases; i++) { m_tPhaseMoles_old[i] = TPhInertMoles[i]; } @@ -4205,10 +4158,10 @@ double VCS_SOLVE::vcs_tmoles() m_tPhaseMoles_old[m_phaseID[i]] += m_molNumSpecies_old[i]; } } - sum = 0.0; + double sum = 0.0; for (size_t i = 0; i < m_numPhases; i++) { sum += m_tPhaseMoles_old[i]; - Vphase = m_VolPhaseList[i]; + vcs_VolPhase* Vphase = m_VolPhaseList[i]; // Took out because we aren't updating mole fractions in Vphase // Vphase->TMoles = m_tPhaseMoles_old[i]; if (m_tPhaseMoles_old[i] == 0.0) { @@ -4224,9 +4177,8 @@ double VCS_SOLVE::vcs_tmoles() #ifdef DEBUG_MODE void VCS_SOLVE::check_tmoles() const { - size_t i; double sum = 0.0; - for (i = 0; i < m_numPhases; i++) { + for (size_t i = 0; i < m_numPhases; i++) { double m_tPhaseMoles_old_a = TPhInertMoles[i]; for (size_t k = 0; k < m_numSpeciesTot; k++) { @@ -4249,9 +4201,8 @@ void VCS_SOLVE::check_tmoles() const void VCS_SOLVE::vcs_updateVP(const int vcsState) { - vcs_VolPhase* Vphase; for (size_t i = 0; i < m_numPhases; i++) { - Vphase = m_VolPhaseList[i]; + vcs_VolPhase* Vphase = m_VolPhaseList[i]; if (vcsState == VCS_STATECALC_OLD) { Vphase->setMolesFromVCSCheck(VCS_STATECALC_OLD, VCS_DATA_PTR(m_molNumSpecies_old), @@ -4355,22 +4306,20 @@ bool VCS_SOLVE::vcs_evaluate_speciesType() void VCS_SOLVE::vcs_switch2D(double* const* const Jac, const size_t k1, const size_t k2) const { - size_t i; if (k1 == k2) { return; } - for (i = 0; i < m_numSpeciesTot; i++) { + for (size_t i = 0; i < m_numSpeciesTot; i++) { std::swap(Jac[k1][i], Jac[k2][i]); } - for (i = 0; i < m_numSpeciesTot; i++) { + for (size_t i = 0; i < m_numSpeciesTot; i++) { std::swap(Jac[i][k1], Jac[i][k2]); } } static void print_space(size_t num) { - size_t j; - for (j = 0; j < num; j++) { + for (size_t j = 0; j < num; j++) { plogf(" "); } } @@ -4378,10 +4327,6 @@ static void print_space(size_t num) void VCS_SOLVE::vcs_deltag(const int l, const bool doDeleted, const int vcsState, const bool alterZeroedPhases) { - size_t iph; - size_t irxn, kspec; - bool lneed; - double* dtmp_ptr; int icase = 0; size_t irxnl = m_numRxnRdc; if (doDeleted) { @@ -4423,12 +4368,12 @@ void VCS_SOLVE::vcs_deltag(const int l, const bool doDeleted, /* **** MAJORS and ZEROED SPECIES ONLY ************* */ /* ************************************************* */ if (l < 0) { - for (irxn = 0; irxn < m_numRxnRdc; ++irxn) { - kspec = irxn + m_numComponents; + for (size_t irxn = 0; irxn < m_numRxnRdc; ++irxn) { + size_t kspec = irxn + m_numComponents; if (m_speciesStatus[kspec] != VCS_SPECIES_MINOR) { icase = 0; deltaGRxn[irxn] = feSpecies[m_indexRxnToSpecies[irxn]]; - dtmp_ptr = m_stoichCoeffRxnMatrix[irxn]; + double* dtmp_ptr = m_stoichCoeffRxnMatrix[irxn]; for (kspec = 0; kspec < m_numComponents; ++kspec) { deltaGRxn[irxn] += dtmp_ptr[kspec] * feSpecies[kspec]; if (molNumSpecies[kspec] < VCS_DELETE_MINORSPECIES_CUTOFF && dtmp_ptr[kspec] < 0.0) { @@ -4444,11 +4389,11 @@ void VCS_SOLVE::vcs_deltag(const int l, const bool doDeleted, /* ************************************************* */ /* **** ALL REACTIONS ****************************** */ /* ************************************************* */ - for (irxn = 0; irxn < irxnl; ++irxn) { + for (size_t irxn = 0; irxn < irxnl; ++irxn) { icase = 0; deltaGRxn[irxn] = feSpecies[m_indexRxnToSpecies[irxn]]; - dtmp_ptr = m_stoichCoeffRxnMatrix[irxn]; - for (kspec = 0; kspec < m_numComponents; ++kspec) { + double* dtmp_ptr = m_stoichCoeffRxnMatrix[irxn]; + for (size_t kspec = 0; kspec < m_numComponents; ++kspec) { deltaGRxn[irxn] += dtmp_ptr[kspec] * feSpecies[kspec]; if (molNumSpecies[kspec] < VCS_DELETE_MINORSPECIES_CUTOFF && dtmp_ptr[kspec] < 0.0) { @@ -4463,12 +4408,12 @@ void VCS_SOLVE::vcs_deltag(const int l, const bool doDeleted, /* ************************************************* */ /* **** MINORS AND ZEROED SPECIES ****************** */ /* ************************************************* */ - for (irxn = 0; irxn < m_numRxnRdc; ++irxn) { - kspec = irxn + m_numComponents; + for (size_t irxn = 0; irxn < m_numRxnRdc; ++irxn) { + size_t kspec = irxn + m_numComponents; if (m_speciesStatus[kspec] <= VCS_SPECIES_MINOR) { icase = 0; deltaGRxn[irxn] = feSpecies[m_indexRxnToSpecies[irxn]]; - dtmp_ptr = m_stoichCoeffRxnMatrix[irxn]; + double* dtmp_ptr = m_stoichCoeffRxnMatrix[irxn]; for (kspec = 0; kspec < m_numComponents; ++kspec) { deltaGRxn[irxn] += dtmp_ptr[kspec] * feSpecies[kspec]; if (m_molNumSpecies_old[kspec] < VCS_DELETE_MINORSPECIES_CUTOFF && @@ -4526,13 +4471,13 @@ void VCS_SOLVE::vcs_deltag(const int l, const bool doDeleted, */ //alterZeroedPhases = false; if (alterZeroedPhases && false) { - for (iph = 0; iph < m_numPhases; iph++) { - lneed = false; + for (size_t iph = 0; iph < m_numPhases; iph++) { + bool lneed = false; vcs_VolPhase* Vphase = m_VolPhaseList[iph]; if (! Vphase->m_singleSpecies) { double sum = 0.0; for (size_t k = 0; k < Vphase->nSpecies(); k++) { - kspec = Vphase->spGlobalIndexVCS(k); + size_t kspec = Vphase->spGlobalIndexVCS(k); if (m_speciesUnknownType[kspec] != VCS_SPECIES_TYPE_INTERFACIALVOLTAGE) { sum += molNumSpecies[kspec]; } @@ -4548,10 +4493,10 @@ void VCS_SOLVE::vcs_deltag(const int l, const bool doDeleted, if (lneed) { double poly = 0.0; for (size_t k = 0; k < Vphase->nSpecies(); k++) { - kspec = Vphase->spGlobalIndexVCS(k); + size_t kspec = Vphase->spGlobalIndexVCS(k); // We may need to look at deltaGRxn for components! if (kspec >= m_numComponents) { - irxn = kspec - m_numComponents; + size_t irxn = kspec - m_numComponents; if (deltaGRxn[irxn] > 50.0) { deltaGRxn[irxn] = 50.0; } @@ -4567,9 +4512,9 @@ void VCS_SOLVE::vcs_deltag(const int l, const bool doDeleted, * negative, then the phase will come back into existence. */ for (size_t k = 0; k < Vphase->nSpecies(); k++) { - kspec = Vphase->spGlobalIndexVCS(k); + size_t kspec = Vphase->spGlobalIndexVCS(k); if (kspec >= m_numComponents) { - irxn = kspec - m_numComponents; + size_t irxn = kspec - m_numComponents; deltaGRxn[irxn] = 1.0 - poly; } } @@ -4588,7 +4533,6 @@ void VCS_SOLVE::vcs_deltag(const int l, const bool doDeleted, void VCS_SOLVE::vcs_printDeltaG(const int stateCalc) { - size_t j; double* deltaGRxn = VCS_DATA_PTR(m_deltaGRxn_old); double* feSpecies = VCS_DATA_PTR(m_feSpecies_old); double* molNumSpecies = VCS_DATA_PTR(m_molNumSpecies_old); @@ -4605,15 +4549,15 @@ void VCS_SOLVE::vcs_printDeltaG(const int stateCalc) bool zeroedPhase = false; if (m_debug_print_lvl >= 2) { plogf(" --- DELTA_G TABLE Components:"); - for (j = 0; j < m_numComponents; j++) { + for (size_t j = 0; j < m_numComponents; j++) { plogf(" %3d ", j); } plogf("\n --- Components Moles:"); - for (j = 0; j < m_numComponents; j++) { + for (size_t j = 0; j < m_numComponents; j++) { plogf("%10.3g", m_molNumSpecies_old[j]); } plogf("\n --- NonComponent| Moles | "); - for (j = 0; j < m_numComponents; j++) { + for (size_t j = 0; j < m_numComponents; j++) { plogf("%-10.10s", m_speciesName[j].c_str()); } //plogf("| m_scSize"); @@ -4626,7 +4570,7 @@ void VCS_SOLVE::vcs_printDeltaG(const int stateCalc) } else { plogf("|%10.3g|", m_molNumSpecies_old[m_indexRxnToSpecies[i]]); } - for (j = 0; j < m_numComponents; j++) { + for (size_t j = 0; j < m_numComponents; j++) { plogf(" %6.2f", m_stoichCoeffRxnMatrix[i][j]); } //plogf(" | %6.2f", m_scSize[i]); @@ -4723,10 +4667,6 @@ void VCS_SOLVE::vcs_printDeltaG(const int stateCalc) void VCS_SOLVE::vcs_deltag_Phase(const size_t iphase, const bool doDeleted, const int stateCalc, const bool alterZeroedPhases) { - size_t iph; - size_t irxn, kspec, kcomp; - double* dtmp_ptr; - double* feSpecies=0; double* deltaGRxn=0; double* actCoeffSpecies=0; @@ -4764,7 +4704,7 @@ void VCS_SOLVE::vcs_deltag_Phase(const size_t iphase, const bool doDeleted, * Single species Phase */ if (vPhase->m_singleSpecies) { - kspec = vPhase->spGlobalIndexVCS(0); + size_t kspec = vPhase->spGlobalIndexVCS(0); #ifdef DEBUG_MODE if (iphase != m_phaseID[kspec]) { plogf("vcs_deltag_Phase index error\n"); @@ -4772,10 +4712,10 @@ void VCS_SOLVE::vcs_deltag_Phase(const size_t iphase, const bool doDeleted, } #endif if (kspec >= m_numComponents) { - irxn = kspec - m_numComponents; + size_t irxn = kspec - m_numComponents; deltaGRxn[irxn] = feSpecies[kspec]; - dtmp_ptr = m_stoichCoeffRxnMatrix[irxn]; - for (kcomp = 0; kcomp < m_numComponents; ++kcomp) { + double* dtmp_ptr = m_stoichCoeffRxnMatrix[irxn]; + for (size_t kcomp = 0; kcomp < m_numComponents; ++kcomp) { deltaGRxn[irxn] += dtmp_ptr[kcomp] * feSpecies[kcomp]; } } @@ -4786,17 +4726,16 @@ void VCS_SOLVE::vcs_deltag_Phase(const size_t iphase, const bool doDeleted, else { bool zeroedPhase = true; - for (irxn = 0; irxn < irxnl; ++irxn) { - kspec = m_indexRxnToSpecies[irxn]; + for (size_t irxn = 0; irxn < irxnl; ++irxn) { + size_t kspec = m_indexRxnToSpecies[irxn]; if (m_speciesUnknownType[kspec] != VCS_SPECIES_TYPE_INTERFACIALVOLTAGE) { - iph = m_phaseID[kspec]; - if (iph == iphase) { + if (m_phaseID[kspec] == iphase) { if (m_molNumSpecies_old[kspec] > 0.0) { zeroedPhase = false; } deltaGRxn[irxn] = feSpecies[kspec]; - dtmp_ptr = m_stoichCoeffRxnMatrix[irxn]; - for (kcomp = 0; kcomp < m_numComponents; ++kcomp) { + double* dtmp_ptr = m_stoichCoeffRxnMatrix[irxn]; + for (size_t kcomp = 0; kcomp < m_numComponents; ++kcomp) { deltaGRxn[irxn] += dtmp_ptr[kcomp] * feSpecies[kcomp]; } } @@ -4846,10 +4785,9 @@ void VCS_SOLVE::vcs_deltag_Phase(const size_t iphase, const bool doDeleted, if (alterZeroedPhases) { if (zeroedPhase) { double phaseDG = 1.0; - for (irxn = 0; irxn < irxnl; ++irxn) { - kspec = m_indexRxnToSpecies[irxn]; - iph = m_phaseID[kspec]; - if (iph == iphase) { + for (size_t irxn = 0; irxn < irxnl; ++irxn) { + size_t kspec = m_indexRxnToSpecies[irxn]; + if (m_phaseID[kspec] == iphase) { if (deltaGRxn[irxn] > 50.0) { deltaGRxn[irxn] = 50.0; } @@ -4862,10 +4800,9 @@ void VCS_SOLVE::vcs_deltag_Phase(const size_t iphase, const bool doDeleted, /* * Overwrite the individual dg's with the phase DG. */ - for (irxn = 0; irxn < irxnl; ++irxn) { - kspec = m_indexRxnToSpecies[irxn]; - iph = m_phaseID[kspec]; - if (iph == iphase) { + for (size_t irxn = 0; irxn < irxnl; ++irxn) { + size_t kspec = m_indexRxnToSpecies[irxn]; + if (m_phaseID[kspec] == iphase) { deltaGRxn[irxn] = 1.0 - phaseDG; } } @@ -4876,8 +4813,6 @@ void VCS_SOLVE::vcs_deltag_Phase(const size_t iphase, const bool doDeleted, void VCS_SOLVE::vcs_switch_pos(const bool ifunc, const size_t k1, const size_t k2) { - size_t i1, i2, iph, kp1, kp2; - vcs_VolPhase* pv1, *pv2; if (k1 == k2) { return; } @@ -4891,11 +4826,11 @@ void VCS_SOLVE::vcs_switch_pos(const bool ifunc, const size_t k1, const size_t k /* * Handle the index pointer in the phase structures first */ - pv1 = m_VolPhaseList[m_phaseID[k1]]; - pv2 = m_VolPhaseList[m_phaseID[k2]]; + vcs_VolPhase* pv1 = m_VolPhaseList[m_phaseID[k1]]; + vcs_VolPhase* pv2 = m_VolPhaseList[m_phaseID[k2]]; - kp1 = m_speciesLocalPhaseIndex[k1]; - kp2 = m_speciesLocalPhaseIndex[k2]; + size_t kp1 = m_speciesLocalPhaseIndex[k1]; + size_t kp2 = m_speciesLocalPhaseIndex[k2]; #ifdef DEBUG_MODE if (pv1->spGlobalIndexVCS(kp1) != k1) { plogf("Indexing error in program\n"); @@ -4948,8 +4883,8 @@ void VCS_SOLVE::vcs_switch_pos(const bool ifunc, const size_t k1, const size_t k /* * Find the Rxn indices corresponding to the two species */ - i1 = k1 - m_numComponents; - i2 = k2 - m_numComponents; + size_t i1 = k1 - m_numComponents; + size_t i2 = k2 - m_numComponents; #ifdef DEBUG_MODE if (i1 > (m_numRxnTot - 1) || i2 > (m_numRxnTot - 1)) { @@ -4961,7 +4896,7 @@ void VCS_SOLVE::vcs_switch_pos(const bool ifunc, const size_t k1, const size_t k std::swap(m_stoichCoeffRxnMatrix[i1][j], m_stoichCoeffRxnMatrix[i2][j]); } std::swap(m_scSize[i1], m_scSize[i2]); - for (iph = 0; iph < m_numPhases; iph++) { + for (size_t iph = 0; iph < m_numPhases; iph++) { std::swap(m_deltaMolNumPhase[i1][iph], m_deltaMolNumPhase[i2][iph]); std::swap(m_phaseParticipation[i1][iph], m_phaseParticipation[i2][iph]); @@ -5053,16 +4988,13 @@ double VCS_SOLVE::vcs_birthGuess(const int kspec) void VCS_SOLVE::vcs_setFlagsVolPhases(const bool upToDate, const int stateCalc) { - vcs_VolPhase* Vphase; if (!upToDate) { for (size_t iph = 0; iph < m_numPhases; iph++) { - Vphase = m_VolPhaseList[iph]; - Vphase->setMolesOutOfDate(stateCalc); + m_VolPhaseList[iph]->setMolesOutOfDate(stateCalc); } } else { for (size_t iph = 0; iph < m_numPhases; iph++) { - Vphase = m_VolPhaseList[iph]; - Vphase->setMolesCurrent(stateCalc); + m_VolPhaseList[iph]->setMolesCurrent(stateCalc); } } } @@ -5070,20 +5002,17 @@ void VCS_SOLVE::vcs_setFlagsVolPhases(const bool upToDate, const int stateCalc) void VCS_SOLVE::vcs_setFlagsVolPhase(const size_t iph, const bool upToDate, const int stateCalc) { - vcs_VolPhase* Vphase = m_VolPhaseList[iph]; if (!upToDate) { - Vphase->setMolesOutOfDate(stateCalc); + m_VolPhaseList[iph]->setMolesOutOfDate(stateCalc); } else { - Vphase->setMolesCurrent(stateCalc); + m_VolPhaseList[iph]->setMolesCurrent(stateCalc); } } void VCS_SOLVE::vcs_updateMolNumVolPhases(const int stateCalc) { - vcs_VolPhase* Vphase; for (size_t iph = 0; iph < m_numPhases; iph++) { - Vphase = m_VolPhaseList[iph]; - Vphase->updateFromVCS_MoleNumbers(stateCalc); + m_VolPhaseList[iph]->updateFromVCS_MoleNumbers(stateCalc); } } diff --git a/src/equil/vcs_species_thermo.cpp b/src/equil/vcs_species_thermo.cpp index ac0da403e..eb8b2858b 100644 --- a/src/equil/vcs_species_thermo.cpp +++ b/src/equil/vcs_species_thermo.cpp @@ -124,9 +124,8 @@ double VCS_SPECIES_THERMO::GStar_R_calc(size_t kglob, double TKelvin, double pres) { char yo[] = "VCS_SPECIES_THERMO::GStar_R_calc "; - double fe, T; - fe = G0_R_calc(kglob, TKelvin); - T = TKelvin; + double fe = G0_R_calc(kglob, TKelvin); + double T = TKelvin; if (UseCanteraCalls) { if (m_VCS_UnitsFormat != VCS_UNITS_MKS) { throw Cantera::CanteraError("VCS_SPECIES_THERMO::GStar_R_calc", @@ -157,9 +156,9 @@ double VCS_SPECIES_THERMO::VolStar_calc(size_t kglob, double TKelvin, double presPA) { char yo[] = "VCS_SPECIES_THERMO::VStar_calc "; - double vol, T; + double vol; - T = TKelvin; + double T = TKelvin; if (UseCanteraCalls) { if (m_VCS_UnitsFormat != VCS_UNITS_MKS) { throw Cantera::CanteraError("VCS_SPECIES_THERMO::VolStar_calc", diff --git a/src/equil/vcs_util.cpp b/src/equil/vcs_util.cpp index 5c58ac840..e8de089ed 100644 --- a/src/equil/vcs_util.cpp +++ b/src/equil/vcs_util.cpp @@ -18,8 +18,7 @@ namespace VCSnonideal #ifndef USE_MEMSET void vcs_dzero(double* vector, int length) { - int i; - for (i = 0; i < length; i++) { + for (int i = 0; i < length; i++) { vector[i] = 0.0; } } @@ -28,8 +27,7 @@ void vcs_dzero(double* vector, int length) #ifndef USE_MEMSET void vcs_izero(int* vector, int length) { - int i; - for (i = 0; i < length; i++) { + for (int i = 0; i < length; i++) { vector[i] = 0; } } @@ -38,8 +36,7 @@ void vcs_izero(int* vector, int length) #ifndef USE_MEMSET void vcs_dcopy(double* const vec_to, const double* const vec_from, int length) { - int i; - for (i = 0; i < length; i++) { + for (int i = 0; i < length; i++) { vec_to[i] = vec_from[i]; } } @@ -48,8 +45,7 @@ void vcs_dcopy(double* const vec_to, const double* const vec_from, int length) #ifndef USE_MEMSET void vcs_icopy(int* vec_to, int* vec_from, int length) { - int i; - for (i = 0; i < length; i++) { + for (int i = 0; i < length; i++) { vec_to[i] = vec_from[i]; } } @@ -109,13 +105,12 @@ void vcs_vicopy(std::vector &vec_to, size_t vcs_optMax(const double* x, const double* xSize, size_t j, size_t n) { - size_t i; size_t largest = j; double big = x[j]; if (xSize) { assert(xSize[j] > 0.0); big *= xSize[j]; - for (i = j + 1; i < n; ++i) { + for (size_t i = j + 1; i < n; ++i) { assert(xSize[i] > 0.0); if ((x[i] * xSize[i]) > big) { largest = i; @@ -123,7 +118,7 @@ size_t vcs_optMax(const double* x, const double* xSize, size_t j, size_t n) } } } else { - for (i = j + 1; i < n; ++i) { + for (size_t i = j + 1; i < n; ++i) { if (x[i] > big) { largest = i; big = x[i]; @@ -135,12 +130,12 @@ size_t vcs_optMax(const double* x, const double* xSize, size_t j, size_t n) int vcs_max_int(const int* vector, int length) { - int i, retn; + int retn; if (vector == NULL || length <= 0) { return 0; } retn = vector[0]; - for (i = 1; i < length; i++) { + for (int i = 1; i < length; i++) { retn = std::max(retn, vector[i]); } return retn; @@ -148,30 +143,23 @@ int vcs_max_int(const int* vector, int length) double vcsUtil_gasConstant(int mu_units) { - double r; switch (mu_units) { case VCS_UNITS_KCALMOL: - r = Cantera::GasConst_cal_mol_K * 1e-3; - break; + return Cantera::GasConst_cal_mol_K * 1e-3; case VCS_UNITS_UNITLESS: - r = 1.0; - break; + return 1.0; case VCS_UNITS_KJMOL: - r = Cantera::GasConstant * 1e-6; - break; + return Cantera::GasConstant * 1e-6; case VCS_UNITS_KELVIN: - r = 1.0; - break; + return 1.0; case VCS_UNITS_MKS: /* joules / kg-mol K = kg m2 / s2 kg-mol K */ - r = Cantera::GasConstant; - break; + return Cantera::GasConstant; default: plogf("vcs_gasConstant error: uknown units: %d\n", mu_units); exit(EXIT_FAILURE); } - return r; } void vcs_print_line(const char* string, int num) @@ -189,77 +177,66 @@ const char* vcs_speciesType_string(int speciesStatus, int length) const char* sss; switch (speciesStatus) { case VCS_SPECIES_COMPONENT: - sss = "Component Species"; - break; + return "Component Species"; case VCS_SPECIES_MAJOR: - sss = "Major Species"; - break; + return "Major Species"; case VCS_SPECIES_MINOR: - sss = "Minor Species"; - break; + return "Minor Species"; case VCS_SPECIES_ZEROEDPHASE: if (length < 48) { - sss = "Set Zeroed-Phase"; + return "Set Zeroed-Phase"; } else { - sss = "Purposely Zeroed-Phase Species (not in problem)"; + return "Purposely Zeroed-Phase Species (not in problem)"; } - break; case VCS_SPECIES_ZEROEDMS: if (length < 23) { - sss = "Zeroed-MS Phase"; + return "Zeroed-MS Phase"; } else { - sss = "Zeroed-MS Phase Species"; + return "Zeroed-MS Phase Species"; } - break; case VCS_SPECIES_ZEROEDSS: if (length < 23) { - sss = "Zeroed-SS Phase"; + return "Zeroed-SS Phase"; } else { - sss = "Zeroed-SS Phase Species"; + return "Zeroed-SS Phase Species"; } - break; case VCS_SPECIES_DELETED: if (length < 22) { - sss = "Deleted Species"; + return "Deleted Species"; } else if (length < 40) { - sss = "Deleted-Small Species"; + return "Deleted-Small Species"; } else { - sss = "Deleted-Small Species in a MS phase"; + return "Deleted-Small Species in a MS phase"; } - break; case VCS_SPECIES_ACTIVEBUTZERO: if (length < 47) { - sss = "Tmp Zeroed in MS"; + return "Tmp Zeroed in MS"; } else { - sss = "Zeroed Species in an active MS phase (tmp)"; + return "Zeroed Species in an active MS phase (tmp)"; } - break; case VCS_SPECIES_STOICHZERO: if (length < 56) { - sss = "Stoich Zeroed in MS"; + return "Stoich Zeroed in MS"; } else { - sss = "Zeroed Species in an active MS phase (Stoich Constraint)"; + return "Zeroed Species in an active MS phase (Stoich Constraint)"; } - break; case VCS_SPECIES_INTERFACIALVOLTAGE: if (length < 29) { - sss = "InterfaceVoltage"; + return "InterfaceVoltage"; } else { - sss = "InterfaceVoltage Species"; + return "InterfaceVoltage Species"; } - break; default: - sss = "unknown species type"; + return "unknown species type"; } - return sss; } void vcs_print_stringTrunc(const char* str, size_t space, int alignment) { - size_t i, ls = 0, rs = 0; + size_t ls = 0, rs = 0; size_t len = strlen(str); if ((len) >= space) { - for (i = 0; i < space; i++) { + for (size_t i = 0; i < space; i++) { plogf("%c", str[i]); } } else { @@ -272,13 +249,13 @@ void vcs_print_stringTrunc(const char* str, size_t space, int alignment) rs = space - len - ls; } if (ls != 0) { - for (i = 0; i < ls; i++) { + for (size_t i = 0; i < ls; i++) { plogf(" "); } } plogf("%s", str); if (rs != 0) { - for (i = 0; i < rs; i++) { + for (size_t i = 0; i < rs; i++) { plogf(" "); } }