diff --git a/Cantera/src/equil/BasisOptimize.cpp b/Cantera/src/equil/BasisOptimize.cpp index 92b01ed2a..4b3b097fd 100644 --- a/Cantera/src/equil/BasisOptimize.cpp +++ b/Cantera/src/equil/BasisOptimize.cpp @@ -1,4 +1,3 @@ - /** * @file BasisOptimize.cpp * Functions which calculation optimized basis of the @@ -6,9 +5,7 @@ */ /* - * $Author$ - * $Date$ - * $Revision$ + * $Id$ */ #include "ct_defs.h" @@ -18,7 +15,7 @@ using namespace Cantera; using namespace std; -#ifdef DEBUG_BASISOPTIMIZE +#ifdef DEBUG_MODE namespace Cantera { int BasisOptimize_print_lvl = 0; } @@ -120,7 +117,7 @@ int Cantera::BasisOptimize(int *usedZeroedSpecies, bool doFormRxn, } } -#ifdef DEBUG_BASISOPTIMIZE +#ifdef DEBUG_MODE double molSave = 0.0; if (BasisOptimize_print_lvl >= 1) { writelog(" "); for(i=0; i<77; i++) writelog("-"); writelog("\n"); @@ -186,7 +183,7 @@ int Cantera::BasisOptimize(int *usedZeroedSpecies, bool doFormRxn, formRxnMatrix.resize(nspecies*ne, 0.0); } -#ifdef DEBUG_BASISOPTIMIZE +#ifdef DEBUG_MODE /* * For debugging purposes keep an unmodified copy of the array. */ @@ -233,7 +230,7 @@ int Cantera::BasisOptimize(int *usedZeroedSpecies, bool doFormRxn, * Assign a small negative number to the component that we have * just found, in order to take it out of further consideration. */ -#ifdef DEBUG_BASISOPTIMIZE +#ifdef DEBUG_MODE molSave = molNum[kk]; #endif molNum[kk] = USEDBEFORE; @@ -293,7 +290,7 @@ int Cantera::BasisOptimize(int *usedZeroedSpecies, bool doFormRxn, /* **** REARRANGE THE DATA ****************** */ /* ****************************************** */ if (jr != k) { -#ifdef DEBUG_BASISOPTIMIZE +#ifdef DEBUG_MODE if (BasisOptimize_print_lvl >= 1) { kk = orderVectorSpecies[k]; sname = mphase->speciesName(kk); @@ -378,7 +375,7 @@ int Cantera::BasisOptimize(int *usedZeroedSpecies, bool doFormRxn, throw CanteraError("basopt", "mlequ returned an error condition"); } -#ifdef DEBUG_BASISOPTIMIZE +#ifdef DEBUG_MODE if (Cantera::BasisOptimize_print_lvl >= 1) { writelog(" ---\n"); writelogf(" --- Number of Components = %d\n", nComponents); @@ -425,7 +422,7 @@ int Cantera::BasisOptimize(int *usedZeroedSpecies, bool doFormRxn, -#ifdef DEBUG_BASISOPTIMIZE +#ifdef DEBUG_MODE static void print_stringTrunc(const char *str, int space, int alignment) /*********************************************************************** @@ -543,7 +540,7 @@ static int amax(double *x, int j, int n) { for (k = i + 1; k < n; ++k) { if (c[k + i * idem] != 0.0) goto FOUND_PIVOT; } -#ifdef DEBUG_BASISOPTIMIZE +#ifdef DEBUG_MODE writelogf("vcs_mlequ ERROR: Encountered a zero column: %d\n", i); #endif return 1; @@ -624,7 +621,7 @@ int Cantera::ElemRearrange(int nComponents, const vector_fp & elementAbundances, int nspecies = mphase->nSpecies(); double test = -1.0E10; -#ifdef DEBUG_BASISOPTIMIZE +#ifdef DEBUG_MODE if (BasisOptimize_print_lvl > 0) { writelog(" "); for(i=0; i<77; i++) writelog("-"); writelog("\n"); writelog(" --- Subroutine ElemRearrange() called to "); @@ -716,7 +713,7 @@ int Cantera::ElemRearrange(int nComponents, const vector_fp & elementAbundances, // When we are here, there is an error usually. // We haven't found the number of elements necessary. // This is signalled by returning jr != nComponents. -#ifdef DEBUG_BASISOPTIMIZE +#ifdef DEBUG_MODE if (BasisOptimize_print_lvl > 0) { writelogf("Error exit: returning with nComponents = %d\n", jr); } @@ -794,17 +791,17 @@ int Cantera::ElemRearrange(int nComponents, const vector_fp & elementAbundances, /* **** REARRANGE THE DATA ****************** */ /* ****************************************** */ if (jr != k) { -#ifdef DEBUG_BASISOPTIMIZE +#ifdef DEBUG_MODE if (BasisOptimize_print_lvl > 0) { - kk = orderVectorElements[k]; - ename = mphase->elementName(kk); - writelog(" --- "); + kk = orderVectorElements[k]; + ename = mphase->elementName(kk); + writelog(" --- "); writelogf("%-2.2s", ename.c_str()); - writelog("replaces "); - kk = orderVectorElements[jr]; - ename = mphase->elementName(kk); - writelogf("%-2.2s", ename.c_str()); - writelogf(" as element %3d\n", jr); + writelog("replaces "); + kk = orderVectorElements[jr]; + ename = mphase->elementName(kk); + writelogf("%-2.2s", ename.c_str()); + writelogf(" as element %3d\n", jr); } #endif switch_pos(orderVectorElements, jr, k); diff --git a/Cantera/src/equil/ChemEquil.cpp b/Cantera/src/equil/ChemEquil.cpp index 61b156072..8ec384626 100755 --- a/Cantera/src/equil/ChemEquil.cpp +++ b/Cantera/src/equil/ChemEquil.cpp @@ -1,9 +1,11 @@ /** * * @file ChemEquil.cpp - * * Chemical equilibrium. Implementation file for class * ChemEquil. + */ + + /* * * * $Id$ @@ -30,7 +32,7 @@ using namespace std; #include "stringUtils.h" #include "MultiPhase.h" -#ifdef DEBUG_CHEMEQUIL +#ifdef DEBUG_MODE #include "stdio.h" int Cantera::ChemEquil_print_lvl = 0; //static char sbuf[1024]; @@ -276,7 +278,7 @@ namespace Cantera { */ update(s); -#ifdef DEBUG_CHEMEQUIL +#ifdef DEBUG_MODE if (ChemEquil_print_lvl > 0) { writelog("setInitialMoles: Estimated Mole Fractions\n"); writelogf(" Temperature = %g\n", s.temperature()); @@ -373,7 +375,7 @@ namespace Cantera { doublereal rrt = 1.0/(GasConstant* s.temperature()); scale(mu_RT.begin(), mu_RT.end(), mu_RT.begin(), rrt); -#ifdef DEBUG_CHEMEQUIL +#ifdef DEBUG_MODE if (ChemEquil_print_lvl > 0) { for (m = 0; m < m_nComponents; m++) { int isp = m_component[m]; @@ -423,7 +425,7 @@ namespace Cantera { } } -#ifdef DEBUG_CHEMEQUIL +#ifdef DEBUG_MODE if (ChemEquil_print_lvl > 0) { writelog(" id CompSpecies ChemPot EstChemPot Diff\n"); for (m = 0; m < m_nComponents; m++) { @@ -510,7 +512,7 @@ namespace Cantera { "Input ThermoPhase is incompatible with initialization"); } -#ifdef DEBUG_CHEMEQUIL +#ifdef DEBUG_MODE int n; const vector& eNames = s.elementNames(); #endif @@ -832,7 +834,7 @@ namespace Cantera { // Compute the Jacobian matrix equilJacobian(s, x, elMolesGoal, jac, xval, yval); -#ifdef DEBUG_CHEMEQUIL +#ifdef DEBUG_MODE if (ChemEquil_print_lvl > 0) { writelogf("Jacobian matrix %d:\n", iter); for (m = 0; m <= m_mm; m++) { @@ -914,7 +916,7 @@ namespace Cantera { } if (fctr != 1.0) { addLogEntry("WARNING: factor to keep solution in bounds", fctr); -#ifdef DEBUG_CHEMEQUIL +#ifdef DEBUG_MODE if (ChemEquil_print_lvl > 0) { writelogf("WARNING Soln Damping because of bounds: %g\n", fctr); } @@ -1092,7 +1094,7 @@ namespace Cantera { for (m = 0; m < nvar; m++) { x[m] = oldx[m] + damp * step[m]; } -#ifdef DEBUG_CHEMEQUIL +#ifdef DEBUG_MODE if (ChemEquil_print_lvl > 0) { writelogf("Solution Unknowns: damp = %g\n", damp); writelog(" X_new X_old Step\n"); @@ -1146,7 +1148,7 @@ namespace Cantera { +fp2str(elmFracGoal[m])+")"); } -#ifdef DEBUG_CHEMEQUIL +#ifdef DEBUG_MODE if (ChemEquil_print_lvl > 0 && !m_doResPerturb) { writelog("Residual: ElFracGoal ElFracCurrent Resid\n"); for (n = 0; n < m_mm; n++) { @@ -1168,7 +1170,7 @@ namespace Cantera { endLogGroup("ChemEquil::equilResidual"); } -#ifdef DEBUG_CHEMEQUIL +#ifdef DEBUG_MODE if (ChemEquil_print_lvl > 0 && !m_doResPerturb) { writelog(" Goal Xvalue Resid\n"); writelogf(" XX : % -14.7E % -14.7E % -10.5E\n", xval, xx, resid[m_mm]); @@ -1395,7 +1397,7 @@ namespace Cantera { } -#ifdef DEBUG_CHEMEQUIL +#ifdef DEBUG_MODE const vector& eNames = s.elementNames(); if (ChemEquil_print_lvl > 0) { writelog("estimateEP_Brinkley::\n\n"); @@ -1442,7 +1444,7 @@ namespace Cantera { /* * Calculate the mole numbers of species */ -#ifdef DEBUG_CHEMEQUIL +#ifdef DEBUG_MODE if (ChemEquil_print_lvl > 0) { writelogf("START ITERATION %d:\n", iter); } @@ -1458,7 +1460,7 @@ namespace Cantera { } -#ifdef DEBUG_CHEMEQUIL +#ifdef DEBUG_MODE if (ChemEquil_print_lvl > 0) { writelog(" Species: Calculated_Moles Calculated_Mole_Fraction\n"); for (k = 0; k < m_kk; k++) { @@ -1493,7 +1495,7 @@ namespace Cantera { } } } -#ifdef DEBUG_CHEMEQUIL +#ifdef DEBUG_MODE if (ChemEquil_print_lvl > 0) { if (!normalStep) { writelogf(" NOTE: iter(%d) Doing an abnormal step due to row %d\n", iter, iM); @@ -1561,7 +1563,7 @@ namespace Cantera { } nCutoff = 1.0E-9 * n_t_calc; -#ifdef DEBUG_CHEMEQUIL +#ifdef DEBUG_MODE if (ChemEquil_print_lvl > 0) { writelog(" Lump Sum Elements Calculation: \n"); } @@ -1589,7 +1591,7 @@ namespace Cantera { } } } -#ifdef DEBUG_CHEMEQUIL +#ifdef DEBUG_MODE if (ChemEquil_print_lvl > 0) { string nnn = eNames[m]; writelogf(" %5s %3d : %5d %5d\n",nnn.c_str(), lumpSum[m], kMSp, kMSp2); @@ -1647,7 +1649,7 @@ namespace Cantera { for (m = 0; m < m_mm; m++) { if (a1(m,m) < 1.0E-50) { -#ifdef DEBUG_CHEMEQUIL +#ifdef DEBUG_MODE if (ChemEquil_print_lvl > 0) { writelogf(" NOTE: Diagonalizing the analytical Jac row %d\n", m); } @@ -1669,7 +1671,7 @@ namespace Cantera { resid[m_mm] = n_t - n_t_calc; -#ifdef DEBUG_CHEMEQUIL +#ifdef DEBUG_MODE if (ChemEquil_print_lvl > 0) { writelog("Matrix:\n"); for (m = 0; m <= m_mm; m++) { @@ -1684,7 +1686,7 @@ namespace Cantera { tmp = resid[m_mm] /(n_t + 1.0E-15); sum += tmp * tmp; -#ifdef DEBUG_CHEMEQUIL +#ifdef DEBUG_MODE if (ChemEquil_print_lvl > 0) { writelogf("(it %d) Convergence = %g\n", iter, sum); } @@ -1710,7 +1712,7 @@ namespace Cantera { tmp += fabs(a1(m,n)); } if (m < m_mm && tmp < 1.0E-30) { -#ifdef DEBUG_CHEMEQUIL +#ifdef DEBUG_MODE if (ChemEquil_print_lvl > 0) { writelogf(" NOTE: Diagonalizing row %d\n", m); } @@ -1729,7 +1731,7 @@ namespace Cantera { resid[m] *= tmp; } -#ifdef DEBUG_CHEMEQUIL +#ifdef DEBUG_MODE if (ChemEquil_print_lvl > 0) { writelog("Row Summed Matrix:\n"); for (m = 0; m <= m_mm; m++) { @@ -1778,7 +1780,7 @@ namespace Cantera { } } if (sameAsRow >= 0 || lumpSum[m]) { -#ifdef DEBUG_CHEMEQUIL +#ifdef DEBUG_MODE if (ChemEquil_print_lvl > 0) { if (lumpSum[m]) { writelogf("Lump summing row %d, due to rank deficiency analysis\n", m); @@ -1797,7 +1799,7 @@ namespace Cantera { } } -#ifdef DEBUG_CHEMEQUIL +#ifdef DEBUG_MODE if (ChemEquil_print_lvl > 0 && modifiedMatrix) { writelog("Row Summed, MODIFIED Matrix:\n"); for (m = 0; m <= m_mm; m++) { @@ -1815,7 +1817,7 @@ namespace Cantera { } catch (CanteraError) { addLogEntry("estimateEP_Brinkley:Jacobian is singular."); -#ifdef DEBUG_CHEMEQUIL +#ifdef DEBUG_MODE if (ChemEquil_print_lvl > 0) { writelog("Matrix is SINGULAR.ERROR\n"); } @@ -1848,7 +1850,7 @@ namespace Cantera { } } } -#ifdef DEBUG_CHEMEQUIL +#ifdef DEBUG_MODE if (ChemEquil_print_lvl > 0) { if (beta != 1.0) { writelogf("(it %d) Beta = %g\n", iter, beta); @@ -1866,7 +1868,7 @@ namespace Cantera { n_t *= exp(beta * resid[m_mm]); -#ifdef DEBUG_CHEMEQUIL +#ifdef DEBUG_MODE if (ChemEquil_print_lvl > 0) { writelogf("(it %d) OLD_SOLUTION NEW SOLUTION (undamped updated)\n", iter); for (m = 0; m < m_mm; m++) { @@ -1878,7 +1880,7 @@ namespace Cantera { #endif } exit: -#ifdef DEBUG_CHEMEQUIL +#ifdef DEBUG_MODE if (ChemEquil_print_lvl > 0) { double temp = s.temperature(); double pres = s.pressure(); @@ -1905,7 +1907,7 @@ namespace Cantera { s.getMoleFractions(DATA_PTR(m_molefractions)); int k; -#ifdef DEBUG_CHEMEQUIL +#ifdef DEBUG_MODE int maxPosEloc = -1; int maxNegEloc = -1; double maxPosVal = -1.0; @@ -1943,7 +1945,7 @@ namespace Cantera { if (sumPos >= sumNeg) { if ( sumPos <= 0.0) return; double factor = (elMolesGoal[m_eloc] + sumNeg) / sumPos; -#ifdef DEBUG_CHEMEQUIL +#ifdef DEBUG_MODE if (ChemEquil_print_lvl > 0) { if (factor < 0.9999999999) { string nnn = s.speciesName(maxPosEloc); @@ -1960,7 +1962,7 @@ namespace Cantera { } } else { double factor = (-elMolesGoal[m_eloc] + sumPos) / sumNeg; -#ifdef DEBUG_CHEMEQUIL +#ifdef DEBUG_MODE if (ChemEquil_print_lvl > 0) { if (factor < 0.9999999999) { string nnn = s.speciesName(maxNegEloc); diff --git a/Cantera/src/equil/ChemEquil.h b/Cantera/src/equil/ChemEquil.h index bef2aa316..7e6ef2c74 100755 --- a/Cantera/src/equil/ChemEquil.h +++ b/Cantera/src/equil/ChemEquil.h @@ -247,7 +247,7 @@ namespace Cantera { }; -#ifdef DEBUG_CHEMEQUIL +#ifdef DEBUG_MODE extern int ChemEquil_print_lvl; #endif diff --git a/Cantera/src/equil/Makefile.in b/Cantera/src/equil/Makefile.in index 17a75d1f1..c13c582b7 100644 --- a/Cantera/src/equil/Makefile.in +++ b/Cantera/src/equil/Makefile.in @@ -31,20 +31,16 @@ endif #LOCAL_DEFS=-DDEBUG_MODE # -# Local Define to turn on if you want to debug ChemEquil: +# Local Define to turn on if you want to debug ChemEquil: or +# BasisOptimize # -#LOCAL_DEFS=-DDEBUG_CHEMEQUIL -# -# Local define to turn on debug statements for BasisOptimize: -# -#LOCAL_DEFS=-DDEBUG_BASISOPTIMIZE +#LOCAL_DEFS=-DDEBUG_MODE # ifeq ($(do_AltLinProg), 1) DALT_STR=-DALTLINPROG endif -LOCAL_DEFS=-DDEBUG_BASISOPTIMIZE -DDEBUG_CHEMEQUIL \ - -DDEBUG $(DALT_STR) +LOCAL_DEFS= $(DALT_STR) # # PIC_FLAG=@PIC@ diff --git a/Cantera/src/equil/MultiPhase.h b/Cantera/src/equil/MultiPhase.h index 8228cdd24..5a3a91fcd 100644 --- a/Cantera/src/equil/MultiPhase.h +++ b/Cantera/src/equil/MultiPhase.h @@ -728,7 +728,7 @@ namespace Cantera { vector_int & orderVectorSpecies, vector_int & orderVectorElements); -#ifdef DEBUG_BASISOPTIMIZE +#ifdef DEBUG_MODE extern int BasisOptimize_print_lvl; #endif } diff --git a/Cantera/src/equil/MultiPhaseEquil.cpp b/Cantera/src/equil/MultiPhaseEquil.cpp index fdddec117..7c78baeff 100644 --- a/Cantera/src/equil/MultiPhaseEquil.cpp +++ b/Cantera/src/equil/MultiPhaseEquil.cpp @@ -1008,7 +1008,7 @@ namespace Cantera { VolPM[k], VolPhaseVolumes ); } } -#ifdef DEBUG +#ifdef DEBUG_MODE /* * Check consistency: These should be equal */ diff --git a/Cantera/src/equil/vcs_MultiPhaseEquil.cpp b/Cantera/src/equil/vcs_MultiPhaseEquil.cpp index 7036c812f..533e89e35 100644 --- a/Cantera/src/equil/vcs_MultiPhaseEquil.cpp +++ b/Cantera/src/equil/vcs_MultiPhaseEquil.cpp @@ -1,11 +1,16 @@ /** * @file vcs_MultiPhaseEquil.cpp - * + * Driver routines for equilibrium solvers */ /* * * $Id$ */ +/* + * Copywrite (2006) Sandia Corporation. Under the terms of + * Contract DE-AC04-94AL85000 with Sandia Corporation, the + * U.S. Government retains certain rights in this software. + */ #include "vcs_MultiPhaseEquil.h" #include "vcs_prob.h" #include "vcs_internal.h" @@ -734,7 +739,7 @@ namespace Cantera { double Temp = m_mix->temperature(); double pres = m_mix->pressure(); double *mf = VCS_DATA_PTR(m_vprob->mf); -#ifdef DEBUG +#ifdef DEBUG_MODE double *fe = VCS_DATA_PTR(m_vprob->m_gibbsSpecies); #endif std::vector VolPM; @@ -857,7 +862,7 @@ namespace Cantera { } } -#ifdef DEBUG +#ifdef DEBUG_MODE /* * Check consistency: These should be equal */ diff --git a/Cantera/src/equil/vcs_elem.cpp b/Cantera/src/equil/vcs_elem.cpp index e3d6d6168..0af328933 100644 --- a/Cantera/src/equil/vcs_elem.cpp +++ b/Cantera/src/equil/vcs_elem.cpp @@ -210,7 +210,7 @@ int VCS_SOLVE::vcs_elcorr(double aa[], double x[]) int i, j, retn = 0, kspec, goodSpec, its; double xx, par, saveDir, dir; -#ifdef DEBUG +#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(ga), m_numElemConstraints); @@ -316,7 +316,7 @@ int VCS_SOLVE::vcs_elcorr(double aa[], double x[]) double maxPermissible = gai[i] / atomComp; if (soln[kspec] > maxPermissible) { -#ifdef DEBUG +#ifdef DEBUG_MODE if (vcs_debug_print_lvl >= 3) { plogf(" --- vcs_elcorr: Reduced species %s from %g to %g due to %s max bounds constraint\n", SpName[kspec].c_str(), soln[kspec], maxPermissible, ElName[i].c_str()); @@ -331,7 +331,7 @@ int VCS_SOLVE::vcs_elcorr(double aa[], double x[]) } else { spStatus[kspec] = VCS_SPECIES_ZEROEDMS; } -#ifdef DEBUG +#ifdef DEBUG_MODE if (vcs_debug_print_lvl >= 2) { plogf(" --- vcs_elcorr: Zeroed species %s and changed status to %d due to max bounds constraint\n", SpName[kspec].c_str(), spStatus[kspec]); @@ -578,7 +578,7 @@ int VCS_SOLVE::vcs_elcorr(double aa[], double x[]) L_CLEANUP: ; vcs_tmoles(); -#ifdef DEBUG +#ifdef DEBUG_MODE l2after = 0.0; for (i = 0; i < m_numElemConstraints; ++i) { l2after += SQUARE(ga[i] - gai[i]); diff --git a/Cantera/src/equil/vcs_elem_rearrange.cpp b/Cantera/src/equil/vcs_elem_rearrange.cpp index c2adbbca2..bb9fcbcdd 100644 --- a/Cantera/src/equil/vcs_elem_rearrange.cpp +++ b/Cantera/src/equil/vcs_elem_rearrange.cpp @@ -60,7 +60,7 @@ int VCS_SOLVE::vcs_elem_rearrange(double *aw, double *sa, double *sm, int j, k, l, i, jl, ml, jr, lindep, ielem; int ncomponents = m_numComponents; double test = -1.0E10; -#ifdef DEBUG +#ifdef DEBUG_MODE if (vcs_debug_print_lvl >= 2) { plogf(" "); for(i=0; i<77; i++) plogf("-"); plogf("\n"); plogf(" --- Subroutine elem_rearrange() called to "); @@ -179,7 +179,7 @@ int VCS_SOLVE::vcs_elem_rearrange(double *aw, double *sa, double *sm, /* **** REARRANGE THE DATA ****************** */ /* ****************************************** */ if (jr != k) { -#ifdef DEBUG +#ifdef DEBUG_MODE if (vcs_debug_print_lvl >= 2) { plogf(" --- "); plogf("%-2.2s", (ElName[k]).c_str()); plogf("(%9.2g) replaces ", gai[k]); @@ -214,7 +214,7 @@ void VCS_SOLVE::vcs_switch_elem_pos(int ipos, int jpos) { int j; double dtmp; vcs_VolPhase *volPhase; -#ifdef DEBUG +#ifdef DEBUG_MODE if (ipos < 0 || ipos > (m_numElemConstraints - 1) || jpos < 0 || jpos > (m_numElemConstraints - 1) ) { plogf("vcs_switch_elem_pos: ifunc = 0: inappropriate args: %d %d\n", diff --git a/Cantera/src/equil/vcs_funcVtot.cpp b/Cantera/src/equil/vcs_funcVtot.cpp index cdebb5756..2931d66b0 100644 --- a/Cantera/src/equil/vcs_funcVtot.cpp +++ b/Cantera/src/equil/vcs_funcVtot.cpp @@ -1,12 +1,16 @@ -/* ======================================================================= */ -/* -------------------------------------------------- */ -/* | RCS Head Information on zuzax.pchem.sandia.gov | */ -/* -------------------------------------------------- */ -/* $RCSfile$ */ -/* $Author$ */ -/* $Date$ */ -/* $Revision$ */ -/* ======================================================================= */ +/** + * @file vcs_funcVtot.cpp + * Routine to calculate tht total volume of an extrinsic system. + */ +/* + * $Id$ + */ +/* + * Copywrite (2006) Sandia Corporation. Under the terms of + * Contract DE-AC04-94AL85000 with Sandia Corporation, the + * U.S. Government retains certain rights in this software. + */ + #include #include #include @@ -36,7 +40,7 @@ double vcs_funcVtot(double xval, double Vtarget, int varID, void *fptrPassthroug } else if (varID == 1) { vptr->Pres = xval; } -#ifdef DEBUG +#ifdef DEBUG_MODE retn = vptr->vcs_TP(1, 1, 10000, vptr->T, vptr->Pres); #else retn = vptr->vcs_TP(0, 0, 10000, vptr->T, vptr->Pres); @@ -48,7 +52,7 @@ double vcs_funcVtot(double xval, double Vtarget, int varID, void *fptrPassthroug } vol = vptr->vcs_VolTotal(vptr->T, vptr->Pres, VCS_DATA_PTR(vptr->soln), VCS_DATA_PTR(vptr->VolPM)); -#ifdef DEBUG +#ifdef DEBUG_MODE vptr->vcs_report(retn); #endif diff --git a/Cantera/src/equil/vcs_inest.cpp b/Cantera/src/equil/vcs_inest.cpp index 31d2c2d7d..a1afb8601 100644 --- a/Cantera/src/equil/vcs_inest.cpp +++ b/Cantera/src/equil/vcs_inest.cpp @@ -41,7 +41,7 @@ namespace VCSnonideal { * Linear programming module is based on using dbolm. ***************************************************************************/ { - int conv, j, k, lt, ikl, kspec, iph, irxn, jj; + int conv, k, lt, ikl, kspec, iph, irxn; double s, s1, xl, par; int finished; int nspecies = m_numSpeciesTot; @@ -66,7 +66,7 @@ namespace VCSnonideal { #ifdef ALTLINPROG vcs_setMolesLinProg(); #else - + int j, jj; std::vector ax(m_numElemConstraints*nspecies, 0.0); std::vector bb(m_numElemConstraints, 0.0); std::vector cc(nspecies, 0.0); @@ -94,7 +94,7 @@ namespace VCSnonideal { VCS_DATA_PTR(bb), neActive, nspecies, neActive); #endif -#ifdef DEBUG +#ifdef DEBUG_MODE if (vcs_debug_print_lvl >= 2) { plogf("%s Mole Numbers returned from linear programming (vcs_inest initial guess):\n", pprefix); @@ -107,8 +107,8 @@ namespace VCSnonideal { "programming (vcs_inest initial guess):\n", pprefix); plogf("%s Element Goal Actual\n", pprefix); - jj = 0; - for (j = 0; j < m_numElemConstraints; j++) { + int jj = 0; + for (int j = 0; j < m_numElemConstraints; j++) { if (ElActive[j]) { double tmp = 0.0; for (kspec = 0; kspec < nspecies; ++kspec) { @@ -202,7 +202,7 @@ namespace VCSnonideal { } } vcs_deltag(0, true); -#ifdef DEBUG +#ifdef DEBUG_MODE if (vcs_debug_print_lvl >= 2) { for (kspec = 0; kspec < nspecies; ++kspec) { plogf("%s", pprefix); plogf("%-12.12s", SpName[kspec].c_str()); @@ -255,7 +255,7 @@ namespace VCSnonideal { } } } -#ifdef DEBUG +#ifdef DEBUG_MODE if (vcs_debug_print_lvl >= 2) { for (kspec = 0; kspec < nspecies; ++kspec) { if (SpeciesUnknownType[kspec] != VCS_SPECIES_TYPE_INTERFACIALVOLTAGE) { @@ -361,7 +361,7 @@ namespace VCSnonideal { } while (!finished); finished: ; -#ifdef DEBUG +#ifdef DEBUG_MODE if (vcs_debug_print_lvl >= 2) { plogf("%s Final Mole Numbers produced by inest:\n", pprefix); @@ -415,7 +415,7 @@ namespace VCSnonideal { /* * Go get the estimate of the solution */ -#ifdef DEBUG +#ifdef DEBUG_MODE if (vcs_debug_print_lvl >= 2) { plogf("%sGo find an initial estimate for the equilibrium problem\n", pprefix); @@ -439,7 +439,7 @@ namespace VCSnonideal { */ int rangeCheck = vcs_elabcheck(1); if (!vcs_elabcheck(0)) { -#ifdef DEBUG +#ifdef DEBUG_MODE if (vcs_debug_print_lvl >= 2) { plogf("%sInitial guess failed element abundances\n", pprefix); plogf("%sCall vcs_elcorr to attempt fix\n", pprefix); @@ -454,7 +454,7 @@ namespace VCSnonideal { "constraints is probable\n", pprefix); retn = -1; } else { -#ifdef DEBUG +#ifdef DEBUG_MODE if (vcs_debug_print_lvl >= 2) { if (rangeCheck) { plogf("%sInitial guess now satisfies element abundances\n", pprefix); @@ -469,7 +469,7 @@ namespace VCSnonideal { } } else { -#ifdef DEBUG +#ifdef DEBUG_MODE if (vcs_debug_print_lvl >= 2) { if (rangeCheck) { plogf("%sInitial guess satisfies element abundances\n", pprefix); @@ -483,7 +483,7 @@ namespace VCSnonideal { #endif } -#ifdef DEBUG +#ifdef DEBUG_MODE if (vcs_debug_print_lvl >= 2) { plogf("%sTotal Dimensionless Gibbs Free Energy = %15.7E\n", pprefix, vcs_Total_Gibbs(VCS_DATA_PTR(soln), VCS_DATA_PTR(m_gibbsSpecies), diff --git a/Cantera/src/equil/vcs_internal.h b/Cantera/src/equil/vcs_internal.h index 6d0d96193..413c9e5bf 100644 --- a/Cantera/src/equil/vcs_internal.h +++ b/Cantera/src/equil/vcs_internal.h @@ -120,7 +120,7 @@ extern void vcsUtil_stsw(std::vector & vecStrings, int, int); /* Externals for vcs_root1d.c */ extern int vcsUtil_root1d(double, double, int, VCS_FUNC_PTR , void *, - double , int, double *); + double , int, double *, int printLvl = 0); /* Externals defined in vcs_timer_generic.c */ diff --git a/Cantera/src/equil/vcs_linmaxc.cpp b/Cantera/src/equil/vcs_linmaxc.cpp index a209eb56e..a9249b801 100644 --- a/Cantera/src/equil/vcs_linmaxc.cpp +++ b/Cantera/src/equil/vcs_linmaxc.cpp @@ -16,7 +16,7 @@ #ifdef hpux #define dbocls_ dbocls #endif -#ifdef DEBUG +#ifdef DEBUG_MODE //extern int vcs_debug_print_lvl; #endif @@ -129,7 +129,7 @@ int linprogmax(double *XMOLES, double *CC, double *AX, double *BB, for (j = 0; j < NCOLS; j++) { XMOLES[j] = X[j]; } -#ifdef DEBUG +#ifdef DEBUG_MODE //sum = 0.0; //for (j = 0; j < NCOLS; j++) { // sum += XMOLES[j] * CC[j]; diff --git a/Cantera/src/equil/vcs_nasa_poly.cpp b/Cantera/src/equil/vcs_nasa_poly.cpp index a7e5f187e..baa98ae01 100644 --- a/Cantera/src/equil/vcs_nasa_poly.cpp +++ b/Cantera/src/equil/vcs_nasa_poly.cpp @@ -167,7 +167,7 @@ double vcs_G0_NASA(double TKelvin, VCS_NASA_POLY *poly_ptr) * Find the temperature region */ if (TKelvin <= *Tlim) { -#ifdef DEBUG +#ifdef DEBUG_MODE plogf("vcs_G0_NASA error: TKelvin below lowest bounds %g\n", *Tlim); #endif iRegion = 0; @@ -182,7 +182,7 @@ double vcs_G0_NASA(double TKelvin, VCS_NASA_POLY *poly_ptr) Tlim++; if (TKelvin <= *Tlim) goto L_FOUNDREGION; } -#ifdef DEBUG +#ifdef DEBUG_MODE plogf("vcs_G0_NASA error: TKelvin above highest bounds %g\n", *(Tlim)); #endif @@ -238,7 +238,7 @@ double vcs_H0_NASA(double TKelvin, VCS_NASA_POLY *poly_ptr) * Find the temperature region */ if (TKelvin <= *Tlim) { -#ifdef DEBUG +#ifdef DEBUG_MODE plogf("vcs_H0_NASA error: TKelvin below lowest bounds\n"); #endif iRegion = 0; @@ -252,7 +252,7 @@ double vcs_H0_NASA(double TKelvin, VCS_NASA_POLY *poly_ptr) Tlim++; if (TKelvin <= *Tlim) goto L_FOUNDREGION; } -#ifdef DEBUG +#ifdef DEBUG_MODE plogf("vcs_H0_NASA error: TKelvin above highest bounds\n"); #endif iRegion--; @@ -309,7 +309,7 @@ double vcs_Cp0_NASA(double TKelvin, VCS_NASA_POLY *poly_ptr) * Find the temperature region */ if (TKelvin <= *Tlim) { -#ifdef DEBUG +#ifdef DEBUG_MODE plogf("vcs_Cp0_NASA error: TKelvin below lowest bounds\n"); #endif iRegion = 0; @@ -323,7 +323,7 @@ double vcs_Cp0_NASA(double TKelvin, VCS_NASA_POLY *poly_ptr) Tlim++; if (TKelvin <= *Tlim) goto L_FOUNDREGION; } -#ifdef DEBUG +#ifdef DEBUG_MODE plogf("vcs_Cp0_NASA error: TKelvin above highest bounds\n"); #endif @@ -376,7 +376,7 @@ double vcs_S0_NASA(double TKelvin, VCS_NASA_POLY *poly_ptr) * Find the temperature region */ if (TKelvin <= *Tlim) { -#ifdef DEBUG +#ifdef DEBUG_MODE plogf("vcs_S0_NASA error: TKelvin below lowest bounds\n"); #endif iRegion = 0; @@ -390,7 +390,7 @@ double vcs_S0_NASA(double TKelvin, VCS_NASA_POLY *poly_ptr) Tlim++; if (TKelvin <= *Tlim) goto L_FOUNDREGION; } -#ifdef DEBUG +#ifdef DEBUG_MODE plogf("vcs_S0_NASA error: TKelvin above highest bounds\n"); #endif iRegion--; diff --git a/Cantera/src/equil/vcs_prep.cpp b/Cantera/src/equil/vcs_prep.cpp index 71b7849c7..3031a46ea 100644 --- a/Cantera/src/equil/vcs_prep.cpp +++ b/Cantera/src/equil/vcs_prep.cpp @@ -123,7 +123,7 @@ int VCS_SOLVE::vcs_prep_oneTime(int printLvl) double *aw, *sa, *sm, *ss; bool modifiedSoln = false; -#ifdef DEBUG +#ifdef DEBUG_MODE vcs_debug_print_lvl = printLvl; #endif diff --git a/Cantera/src/equil/vcs_prob.cpp b/Cantera/src/equil/vcs_prob.cpp index 83efc2709..ccbc97cc3 100644 --- a/Cantera/src/equil/vcs_prob.cpp +++ b/Cantera/src/equil/vcs_prob.cpp @@ -54,7 +54,7 @@ VCS_PROB::VCS_PROB(int nsp, int nel, int nph) : m_Iterations(0), m_NumBasisOptimizations(0), m_printLvl(0) -#ifdef DEBUG +#ifdef DEBUG_MODE , vcs_debug_print_lvl(0) #endif @@ -603,7 +603,7 @@ void VCS_PROB::reportCSV(const std::string &reportFile) { } } -#ifdef DEBUG +#ifdef DEBUG_MODE /* * Check consistency: These should be equal */ @@ -623,7 +623,7 @@ void VCS_PROB::reportCSV(const std::string &reportFile) { } -#ifdef DEBUG +#ifdef DEBUG_MODE void VCS_PROB::setDebugPrintLvl(int lvl) { vcs_debug_print_lvl = lvl; } diff --git a/Cantera/src/equil/vcs_prob.h b/Cantera/src/equil/vcs_prob.h index 3a9bf419b..2fa99eee0 100644 --- a/Cantera/src/equil/vcs_prob.h +++ b/Cantera/src/equil/vcs_prob.h @@ -239,7 +239,7 @@ namespace VCSnonideal { //! Print level for print routines int m_printLvl; -#ifdef DEBUG +#ifdef DEBUG_MODE //! Debug print lvl int vcs_debug_print_lvl; #endif @@ -361,7 +361,7 @@ namespace VCSnonideal { void reportCSV(const std::string &reportFile); -#ifdef DEBUG +#ifdef DEBUG_MODE //! Set the debug level /*! * @param vcs_debug_print_lvl input debug level diff --git a/Cantera/src/equil/vcs_root1d.cpp b/Cantera/src/equil/vcs_root1d.cpp index f3ca72d0a..ea19a613d 100644 --- a/Cantera/src/equil/vcs_root1d.cpp +++ b/Cantera/src/equil/vcs_root1d.cpp @@ -1,12 +1,16 @@ -/* ======================================================================= */ -/* -------------------------------------------------- */ -/* | RCS Head Information on zuzax.pchem.sandia.gov | */ -/* -------------------------------------------------- */ -/* $RCSfile$ */ -/* $Author$ */ -/* $Date$ */ -/* $Revision$ */ -/* ======================================================================= */ +/** + * @file vcs_root1d.cpp + * Code for a one dimensional root finder program. + */ +/* + * $Id$ + */ +/* + * Copywrite (2006) Sandia Corporation. Under the terms of + * Contract DE-AC04-94AL85000 with Sandia Corporation, the + * U.S. Government retains certain rights in this software. + */ + #include #include #include @@ -19,9 +23,8 @@ namespace VCSnonideal { /*****************************************************************************/ /*****************************************************************************/ /*****************************************************************************/ -#ifdef DEBUG_ROOT1D -static void print_funcEval(FILE *fp, double xval, double fval, int its) - +#ifdef DEBUG_MODE +static void print_funcEval(FILE *fp, double xval, double fval, int its) { fprintf(fp,"\n"); fprintf(fp,"...............................................................\n"); @@ -36,24 +39,21 @@ static void print_funcEval(FILE *fp, double xval, double fval, int its) /*****************************************************************************/ /*****************************************************************************/ /*****************************************************************************/ - -int vcsUtil_root1d(double xmin, double xmax, int itmax, - VCS_FUNC_PTR func, void *fptrPassthrough, - double FuncTargVal, int varID, - double *xbest) - - /************************************************************************** + /************************************************************************** * * vcs_root1d: * * Driver for solving a 1D function. ***************************************************************************/ -{ + int vcsUtil_root1d(double xmin, double xmax, int itmax, + VCS_FUNC_PTR func, void *fptrPassthrough, + double FuncTargVal, int varID, + double *xbest, int printLvl) { static int callNum = 0; const char *stre = "vcs_root1d ERROR: "; const char *strw = "vcs_root1d WARNING: "; int converged = FALSE, err = FALSE; -#ifdef DEBUG_ROOT1D +#ifdef DEBUG_MODE char fileName[80]; FILE *fp; #endif @@ -68,12 +68,18 @@ int vcsUtil_root1d(double xmin, double xmax, int itmax, double c[9], f[3], xn1, xn2, x0 = 0.0, f0 = 0.0, root, theta, xquad; callNum++; -#ifdef DEBUG_ROOT1D - sprintf(fileName, "rootfd_%d.log", callNum); - fp = fopen(fileName, "w"); - fprintf(fp, " Iter TP_its xval Func_val | Reasoning\n"); - fprintf(fp, "-----------------------------------------------------" - "-------------------------------\n"); +#ifdef DEBUG_MODE + if (printLvl >= 3) { + sprintf(fileName, "rootfd_%d.log", callNum); + fp = fopen(fileName, "w"); + fprintf(fp, " Iter TP_its xval Func_val | Reasoning\n"); + fprintf(fp, "-----------------------------------------------------" + "-------------------------------\n"); + } +#else + if (printLvl >= 3) { + plogf("WARNING: vcsUtil_root1d: printlvl >= 3, but debug mode not turned on\n"); + } #endif if (xmax <= xmin) { plogf("%sxmin and xmax are bad: %g %g\n", stre, xmin, xmax); @@ -84,10 +90,11 @@ int vcsUtil_root1d(double xmin, double xmax, int itmax, x1 = (xmin + xmax) / 2.0; } f1 = func(x1, FuncTargVal, varID, fptrPassthrough, &err); -#ifdef DEBUG_ROOT1D - print_funcEval(x1, f1, its); - - fprintf(fp, "%-5d %-5d %-15.5E %-15.5E\n", -2, VCount->Its, x1, f1); +#ifdef DEBUG_MODE + if (printLvl >= 3) { + print_funcEval(fp, x1, f1, its); + fprintf(fp, "%-5d %-5d %-15.5E %-15.5E\n", -2, 0, x1, f1); + } #endif if (f1 == 0.0) { *xbest = x1; @@ -103,9 +110,11 @@ int vcsUtil_root1d(double xmin, double xmax, int itmax, x2 = x1 * 1.1; if (x2 > xmax) x2 = x1 - (xmax - xmin) / 100.; f2 = func(x2, FuncTargVal, varID, fptrPassthrough, &err); -#ifdef DEBUG_ROOT1D - print_funcEval(x2, f2, its); - fprintf(fp, "%-5d %-5d %-15.5E %-15.5E", -1, VCount->Its, x2, f2); +#ifdef DEBUG_MODE + if (printLvl >= 3) { + print_funcEval(fp, x2, f2, its); + fprintf(fp, "%-5d %-5d %-15.5E %-15.5E", -1, 0, x2, f2); + } #endif if (FuncTargVal != 0.0) { @@ -146,8 +155,10 @@ int vcsUtil_root1d(double xmin, double xmax, int itmax, } else { xnew = x2 - f2 / slope; } -#ifdef DEBUG_ROOT1D - fprintf(fp, " | xlin = %-9.4g", xnew); +#ifdef DEBUG_MODE + if (printLvl >= 3) { + fprintf(fp, " | xlin = %-9.4g", xnew); + } #endif /* @@ -171,9 +182,11 @@ int vcsUtil_root1d(double xmin, double xmax, int itmax, theta = fabs(xquad - xnew) / fabs(xnew - x2); theta = MIN(1.0, theta); xnew = theta * xnew + (1.0 - theta) * xquad; -#ifdef DEBUG_ROOT1D - if (theta != 1.0) { - fprintf(fp, " | xquad = %-9.4g", xnew); +#ifdef DEBUG_MODE + if (printLvl >= 3) { + if (theta != 1.0) { + fprintf(fp, " | xquad = %-9.4g", xnew); + } } #endif } else { @@ -184,8 +197,10 @@ int vcsUtil_root1d(double xmin, double xmax, int itmax, if ((DSIGN(xnew - x2) == DSIGN(x2 - x1)) && (DSIGN(x2 - x1) == DSIGN(x1 - x0)) ) { xnew += xnew - x2; -#ifdef DEBUG_ROOT1D - fprintf(fp, " | xquada = %-9.4g", xnew); +#ifdef DEBUG_MODE + if (printLvl >= 3) { + fprintf(fp, " | xquada = %-9.4g", xnew); + } #endif } } @@ -207,14 +222,18 @@ int vcsUtil_root1d(double xmin, double xmax, int itmax, slope = fabs(x2 - x1) / 10.; if (fabs(xnew - x1) < slope) { xnew = x1 + DSIGN(xnew-x1) * slope; -#ifdef DEBUG_ROOT1D - fprintf(fp, " | x10%% = %-9.4g", xnew); +#ifdef DEBUG_MODE + if (printLvl >= 3) { + fprintf(fp, " | x10%% = %-9.4g", xnew); + } #endif } if (fabs(xnew - x2) < slope) { xnew = x2 + DSIGN(xnew-x2) * slope; -#ifdef DEBUG_ROOT1D - fprintf(fp, " | x10%% = %-9.4g", xnew); +#ifdef DEBUG_MODE + if (printLvl >= 3) { + fprintf(fp, " | x10%% = %-9.4g", xnew); + } #endif } } else { @@ -225,8 +244,10 @@ int vcsUtil_root1d(double xmin, double xmax, int itmax, slope = 2.0 * fabs(x2 - x1); if (fabs(slope) < fabs(xnew - x2)) { xnew = x2 + DSIGN(xnew-x2) * slope; -#ifdef DEBUG_ROOT1D - fprintf(fp, " | xlimitsize = %-9.4g", xnew); +#ifdef DEBUG_MODE + if (printLvl >= 3) { + fprintf(fp, " | xlimitsize = %-9.4g", xnew); + } #endif } } @@ -234,18 +255,22 @@ int vcsUtil_root1d(double xmin, double xmax, int itmax, if (xnew > xmax) { xnew = x2 + (xmax - x2) / 2.0; -#ifdef DEBUG_ROOT1D - fprintf(fp, " | xlimitmax = %-9.4g", xnew); +#ifdef DEBUG_MODE + if (printLvl >= 3) { + fprintf(fp, " | xlimitmax = %-9.4g", xnew); + } #endif } if (xnew < xmin) { xnew = x2 + (x2 - xmin) / 2.0; -#ifdef DEBUG_ROOT1D - fprintf(fp, " | xlimitmin = %-9.4g", xnew); +#ifdef DEBUG_MODE + if (printLvl >= 3) { + fprintf(fp, " | xlimitmin = %-9.4g", xnew); + } #endif } if (foundStraddle) { -#ifdef DEBUG_ROOT1D +#ifdef DEBUG_MODE slope = xnew; #endif if (posStraddle) { @@ -265,18 +290,22 @@ int vcsUtil_root1d(double xmin, double xmax, int itmax, if (xnew < xPosF) xnew = (xPosF + x2)/2; } } -#ifdef DEBUG_ROOT1D - if (slope != xnew) { - fprintf(fp, " | xstraddle = %-9.4g", xnew); +#ifdef DEBUG_MODE + if (printLvl >= 3) { + if (slope != xnew) { + fprintf(fp, " | xstraddle = %-9.4g", xnew); + } } #endif } fnew = func(xnew, FuncTargVal, varID, fptrPassthrough, &err); -#ifdef DEBUG_ROOT1D - fprintf(fp,"\n"); - print_funcEval(xnew, fnew, its); - fprintf(fp, "%-5d %-5d %-15.5E %-15.5E", its, VCount->Its, xnew, fnew); +#ifdef DEBUG_MODE + if (printLvl >= 3) { + fprintf(fp,"\n"); + print_funcEval(fp, xnew, fnew, its); + fprintf(fp, "%-5d %-5d %-15.5E %-15.5E", its, 0, xnew, fnew); + } #endif if (foundStraddle) { @@ -327,20 +356,30 @@ int vcsUtil_root1d(double xmin, double xmax, int itmax, its++; } while (! converged && its < itmax); if (converged) { -#ifdef DEBUG_ROOT1D - plogf("vcs_root1d success: convergence achieved\n"); - fprintf(fp, " | vcs_root1d success in %d its, fnorm = %g\n", its, fnorm); + if (printLvl >= 1) { + plogf("vcs_root1d success: convergence achieved\n"); + } +#ifdef DEBUG_MODE + if (printLvl >= 3) { + fprintf(fp, " | vcs_root1d success in %d its, fnorm = %g\n", its, fnorm); + } #endif } else { retn = VCS_FAILED_CONVERGENCE; - plogf("vcs_root1d ERROR: maximum iterations exceeded without convergence\n"); -#ifdef DEBUG_ROOT1D - fprintf(fp, "\nvcs_root1d failure in %d its\n", its); + if (printLvl >= 1) { + plogf("vcs_root1d ERROR: maximum iterations exceeded without convergence\n"); + } +#ifdef DEBUG_MODE + if (printLvl >= 3) { + fprintf(fp, "\nvcs_root1d failure in %d its\n", its); + } #endif } *xbest = x2; -#ifdef DEBUG_ROOT1D - fclose(fp); +#ifdef DEBUG_MODE + if (printLvl >= 3) { + fclose(fp); + } #endif return retn; } diff --git a/Cantera/src/equil/vcs_rxnadj.cpp b/Cantera/src/equil/vcs_rxnadj.cpp index c538afb12..4066ee3ce 100644 --- a/Cantera/src/equil/vcs_rxnadj.cpp +++ b/Cantera/src/equil/vcs_rxnadj.cpp @@ -1,12 +1,15 @@ -/* ======================================================================= */ -/* -------------------------------------------------- */ -/* | RCS Head Information on zuzax.pchem.sandia.gov | */ -/* -------------------------------------------------- */ -/* $RCSfile$ */ -/* $Author$ */ -/* $Date$ */ -/* $Revision$ */ -/* ======================================================================= */ +/** + * @file vcs_rxnadj.cpp + * routines for carrying out various line adjustments + */ +/* + * $Id$ + */ +/* + * Copywrite (2006) Sandia Corporation. Under the terms of + * Contract DE-AC04-94AL85000 with Sandia Corporation, the + * U.S. Government retains certain rights in this software. + */ #include #include #include @@ -52,7 +55,7 @@ int VCS_SOLVE::vcs_rxn_adj_cg(void) int irxn, j, k, kspec, soldel = 0; double s, xx, dss; double *dnPhase_irxn; -#ifdef DEBUG +#ifdef DEBUG_MODE char ANOTE[128]; plogf(" "); for (j = 0; j < 77; j++) plogf("-"); plogf("\n --- Subroutine rxn_adj_cg() called\n"); @@ -66,7 +69,7 @@ int VCS_SOLVE::vcs_rxn_adj_cg(void) * this algorithm. If not, we bail out. */ for (irxn = 0; irxn < m_numRxnRdc; ++irxn) { -#ifdef DEBUG +#ifdef DEBUG_MODE sprintf(ANOTE,"Normal Calc"); #endif @@ -82,14 +85,14 @@ int VCS_SOLVE::vcs_rxn_adj_cg(void) * should be replaced with something more relativistic */ if (dg[irxn] < -1.0e-4) { -#ifdef DEBUG +#ifdef DEBUG_MODE (void) sprintf(ANOTE, "MultSpec: come alive DG = %11.3E", dg[irxn]); #endif ds[kspec] = 1.0e-10; spStatus[irxn] = VCS_SPECIES_MAJOR; --(m_numRxnMinorZeroed); } else { -#ifdef DEBUG +#ifdef DEBUG_MODE (void) sprintf(ANOTE, "MultSpec: still dead DG = %11.3E", dg[irxn]); #endif ds[kspec] = 0.0; @@ -106,7 +109,7 @@ int VCS_SOLVE::vcs_rxn_adj_cg(void) * in this mode. */ if (fabs(dg[irxn]) <= tolmaj2) { -#ifdef DEBUG +#ifdef DEBUG_MODE sprintf(ANOTE,"Skipped: converged DG = %11.3E\n", dg[irxn]); plogf(" --- "); plogf("%-12.12s", SpName[kspec].c_str()); plogf(" %12.4E %12.4E | %s\n", soln[kspec], ds[kspec], ANOTE); @@ -118,7 +121,7 @@ int VCS_SOLVE::vcs_rxn_adj_cg(void) * their values are to be decreasing anyway. */ if (spStatus[irxn] <= VCS_SPECIES_MINOR && dg[irxn] >= 0.0) { -#ifdef DEBUG +#ifdef DEBUG_MODE sprintf(ANOTE,"Skipped: IC = %3d and DG >0: %11.3E\n", spStatus[irxn], dg[irxn]); plogf(" --- "); plogf("%-12.12s", SpName[kspec].c_str()); @@ -196,7 +199,7 @@ int VCS_SOLVE::vcs_rxn_adj_cg(void) } soln[k] = 0.0; TPhMoles[PhaseID[k]] = 0.0; -#ifdef DEBUG +#ifdef DEBUG_MODE plogf(" --- vcs_st2 Special section to delete "); plogf("%-12.12s", SpName[k].c_str()); plogf("\n --- Immediate return - Restart iteration\n"); @@ -212,7 +215,7 @@ int VCS_SOLVE::vcs_rxn_adj_cg(void) } } } /* End of regular processing */ -#ifdef DEBUG +#ifdef DEBUG_MODE plogf(" --- "); plogf("%-12.12s", SpName[kspec].c_str()); plogf(" %12.4E %12.4E | %s\n", soln[kspec], ds[kspec], ANOTE); #endif @@ -231,7 +234,7 @@ int VCS_SOLVE::vcs_rxn_adj_cg(void) */ -#ifdef DEBUG +#ifdef DEBUG_MODE plogf(" "); for (j = 0; j < 77; j++) plogf("-"); plogf("\n"); #endif return soldel; @@ -385,7 +388,7 @@ double VCS_SOLVE::deltaG_Recalc_Rxn(int irxn, const double *const molNum, /*****************************************************************************/ /*****************************************************************************/ /*****************************************************************************/ -#ifdef DEBUG +#ifdef DEBUG_MODE double VCS_SOLVE::vcs_line_search(int irxn, double dx_orig, char *ANOTE) #else double VCS_SOLVE::vcs_line_search(int irxn, double dx_orig) @@ -418,7 +421,7 @@ double VCS_SOLVE::vcs_line_search(int irxn, double dx_orig) if (deltaGOrig > 0.0) { if (dx_orig > 0.0) { dx = 0.0; -#ifdef DEBUG +#ifdef DEBUG_MODE if (vcs_debug_print_lvl >= 2) { //plogf(" --- %s :Warning possible error dx>0 dg > 0\n", SpName[kspec]); } @@ -429,7 +432,7 @@ double VCS_SOLVE::vcs_line_search(int irxn, double dx_orig) } else if (deltaGOrig < 0.0) { if (dx_orig < 0.0) { dx = 0.0; -#ifdef DEBUG +#ifdef DEBUG_MODE if (vcs_debug_print_lvl >= 2) { //plogf(" --- %s :Warning possible error dx<0 dg < 0\n", SpName[kspec]); } @@ -511,13 +514,13 @@ double VCS_SOLVE::vcs_line_search(int irxn, double dx_orig) finalize: if (its >= MAXITS) { -#ifdef DEBUG +#ifdef DEBUG_MODE sprintf(ANOTE,"Rxn reduced to zero step size from %g to %g (MAXITS)", dx_orig, dx); return dx; #endif } -#ifdef DEBUG +#ifdef DEBUG_MODE if (dx != dx_orig) { sprintf(ANOTE,"Line Search reduced step size from %g to %g", dx_orig, dx); diff --git a/Cantera/src/equil/vcs_setMolesLinProg.cpp b/Cantera/src/equil/vcs_setMolesLinProg.cpp index 00d6849c4..f7673206b 100644 --- a/Cantera/src/equil/vcs_setMolesLinProg.cpp +++ b/Cantera/src/equil/vcs_setMolesLinProg.cpp @@ -26,7 +26,7 @@ using namespace std; namespace VCSnonideal { -#ifdef DEBUG +#ifdef DEBUG_MODE static void printProgress(const vector &spName, const vector &soln, const vector &ff) { @@ -60,7 +60,7 @@ int VCS_SOLVE::vcs_setMolesLinProg() { int ik, irxn; double test = -1.0E-10; -#ifdef DEBUG +#ifdef DEBUG_MODE std::string pprefix(" --- seMolesLinProg "); if (vcs_debug_print_lvl >= 2) { plogf(" --- call setInitialMoles\n"); @@ -96,7 +96,7 @@ int VCS_SOLVE::vcs_setMolesLinProg() { } } -#ifdef DEBUG +#ifdef DEBUG_MODE if (vcs_debug_print_lvl >= 2) { printProgress(SpName, soln, ff); } @@ -105,7 +105,7 @@ int VCS_SOLVE::vcs_setMolesLinProg() { while (redo) { if (!vcs_elabcheck(0)) { -#ifdef DEBUG +#ifdef DEBUG_MODE if (vcs_debug_print_lvl >= 2) { plogf("%s Mole numbers failing element abundances\n", pprefix.c_str()); plogf("%sCall vcs_elcorr to attempt fix\n", pprefix.c_str()); @@ -130,7 +130,7 @@ int VCS_SOLVE::vcs_setMolesLinProg() { test, &usedZeroedSpecies); if (retn != VCS_SUCCESS) return retn; -#ifdef DEBUG +#ifdef DEBUG_MODE if (vcs_debug_print_lvl >= 2) { plogf("iteration %d\n", iter); } @@ -169,7 +169,7 @@ int VCS_SOLVE::vcs_setMolesLinProg() { // with a new set of components if (!redo) { if (delta_xi < 1.0e-10 && (soln[ik] >= 1.0E-10)) { -#ifdef DEBUG +#ifdef DEBUG_MODE if (vcs_debug_print_lvl >= 2) { plogf(" --- Component too small: %s\n", SpName[jcomp].c_str()); } @@ -206,14 +206,14 @@ int VCS_SOLVE::vcs_setMolesLinProg() { // updateMixMoles(); // Update the phase objects with the contents of the soln vector // vcs_updateVP(0); -#ifdef DEBUG +#ifdef DEBUG_MODE if (vcs_debug_print_lvl >= 2) { printProgress(SpName, soln, ff); } #endif } -#ifdef DEBUG +#ifdef DEBUG_MODE if (vcs_debug_print_lvl == 1) { printProgress(SpName, soln, ff); plogf(" --- setInitialMoles end\n"); @@ -229,4 +229,3 @@ int VCS_SOLVE::vcs_setMolesLinProg() { } #endif } - diff --git a/Cantera/src/equil/vcs_solve.cpp b/Cantera/src/equil/vcs_solve.cpp index 0d7059a38..63c58adae 100644 --- a/Cantera/src/equil/vcs_solve.cpp +++ b/Cantera/src/equil/vcs_solve.cpp @@ -50,7 +50,7 @@ namespace VCSnonideal { Vol(0.0), Faraday_dim(1.602e-19 * 6.022136736e26), m_VCount(0), -#ifdef DEBUG +#ifdef DEBUG_MODE vcs_debug_print_lvl(0), #endif m_VCS_UnitsFormat(VCS_UNITS_UNITLESS) @@ -491,7 +491,7 @@ namespace VCSnonideal { */ NPhase = nph; -#ifdef DEBUG +#ifdef DEBUG_MODE vcs_debug_print_lvl = pub->vcs_debug_print_lvl; #endif diff --git a/Cantera/src/equil/vcs_solve.h b/Cantera/src/equil/vcs_solve.h index c6fd4d963..f139a96ab 100644 --- a/Cantera/src/equil/vcs_solve.h +++ b/Cantera/src/equil/vcs_solve.h @@ -202,7 +202,7 @@ public: double vcs_Hessian_diag_adj(int, double); double vcs_Hessian_actCoeff_diag(int irxn); void vcs_CalcLnActCoeffJac(const double * const moleSpeciesVCS); -#ifdef DEBUG +#ifdef DEBUG_MODE double vcs_line_search(int irxn, double dx_orig, char *ANOTE); #else double vcs_line_search(int irxn, double dx_orig); @@ -294,7 +294,7 @@ private: * @param do_delete */ double minor_alt_calc(int kspec, int irxn, int *do_delete -#ifdef DEBUG +#ifdef DEBUG_MODE , char *ANOTE #endif ); @@ -302,7 +302,7 @@ private: int force(int iti); void vcs_switch2D(double * const * const Jac, int k1, int k2); double l2normdg(double dg[]); -#ifdef DEBUG +#ifdef DEBUG_MODE void prneav(void); void checkDelta1(double * const ds, double * const delTPhMoles, int kspec); #endif @@ -779,7 +779,7 @@ public: VCS_COUNTERS *m_VCount; -#ifdef DEBUG +#ifdef DEBUG_MODE int vcs_debug_print_lvl; #endif diff --git a/Cantera/src/equil/vcs_solve_TP.cpp b/Cantera/src/equil/vcs_solve_TP.cpp index 0ff61dd06..a5aee834f 100644 --- a/Cantera/src/equil/vcs_solve_TP.cpp +++ b/Cantera/src/equil/vcs_solve_TP.cpp @@ -1,4 +1,4 @@ -/*! +/** * @file vcs_solve_TP.cpp Implementation file that contains the * main algorithm for finding an equilibrium */ @@ -34,13 +34,13 @@ namespace VCSnonideal { static void print_space(int num); -#ifdef DEBUG +#ifdef DEBUG_MODE //static double minor_alt_calc(int, int, int *, char *); #else //static double minor_alt_calc(int, int, int *); #endif -#ifdef DEBUG -# ifdef DEBUG_MORE +#ifdef DEBUG_MODE +# ifdef DEBUG_NOT static void prneav(void); static int prnfm(void); # endif @@ -49,7 +49,7 @@ static int prnfm(void); /*****************************************************************************/ /*****************************************************************************/ -#ifdef DEBUG +#ifdef DEBUG_MODE void VCS_SOLVE::checkDelta1(double * const dsLocal, double * const delTPhMoles, int kspec) { std::vector dchange(NPhase, 0.0); @@ -172,7 +172,7 @@ int VCS_SOLVE::vcs_solve_TP(int print_lvl, int printDetails, int maxit) vcs_VolPhase *Vphase; double *sc_irxn = NULL; /* Stoichiometric coefficients for cur rxn */ double *dnPhase_irxn; -#ifdef DEBUG +#ifdef DEBUG_MODE char ANOTE[128]; /* * Set the debug print lvl to the same as the print lvl. @@ -359,7 +359,7 @@ int VCS_SOLVE::vcs_solve_TP(int print_lvl, int printDetails, int maxit) spStatus[irxn] = vcs_species_type(kspec); if (spStatus[irxn] == VCS_SPECIES_MINOR) { spStatus[irxn] = VCS_SPECIES_MAJOR; -#ifdef DEBUG +#ifdef DEBUG_MODE if (vcs_debug_print_lvl >= 2) { plogf(" --- Minor species changed to major: "); plogf("%-12s\n", SpName[kspec].c_str()); @@ -373,7 +373,7 @@ int VCS_SOLVE::vcs_solve_TP(int print_lvl, int printDetails, int maxit) im = (m_numRxnMinorZeroed == m_numRxnRdc); lec = FALSE; if (! vcs_elabcheck(0)) { -#ifdef DEBUG +#ifdef DEBUG_MODE if (vcs_debug_print_lvl >= 2) { plogf(" --- Element Abundance check failed\n"); } @@ -381,7 +381,7 @@ int VCS_SOLVE::vcs_solve_TP(int print_lvl, int printDetails, int maxit) vcs_elcorr(VCS_DATA_PTR(sm), VCS_DATA_PTR(wx)); vcs_dfe(VCS_DATA_PTR(soln), 0, 0, 0, m_numSpeciesRdc); } -#ifdef DEBUG +#ifdef DEBUG_MODE else { if (vcs_debug_print_lvl >= 2) { plogf(" --- Element Abundance check passed\n"); @@ -474,7 +474,7 @@ int VCS_SOLVE::vcs_solve_TP(int print_lvl, int printDetails, int maxit) /* - we don't need to start over. */ } } else { -#ifdef DEBUG +#ifdef DEBUG_MODE if (vcs_debug_print_lvl >= 2) { if (im) { plogf(" --- vcs_RxnStepSizes not called because all" @@ -513,7 +513,7 @@ int VCS_SOLVE::vcs_solve_TP(int print_lvl, int printDetails, int maxit) */ goto L_RETURN_BLOCK; } -#ifdef DEBUG +#ifdef DEBUG_MODE if (vcs_debug_print_lvl >= 2) { plogf(" --- Main Loop Treatment of each non-component species "); if (iti == 0) plogf("- Full Calculation:\n"); @@ -528,14 +528,14 @@ int VCS_SOLVE::vcs_solve_TP(int print_lvl, int printDetails, int maxit) sc_irxn = sc[irxn]; iph = PhaseID[kspec]; Vphase = VPhaseList[iph]; -#ifdef DEBUG +#ifdef DEBUG_MODE ANOTE[0] = '\0'; #endif /********************************************************************/ /********************** VOLTAGE SPECIES **************************/ /********************************************************************/ if (spStatus[irxn] == VCS_SPECIES_INTERFACIALVOLTAGE) { -#ifdef DEBUG +#ifdef DEBUG_MODE dx = minor_alt_calc(kspec, irxn, &soldel, ANOTE); #else dx = minor_alt_calc(kspec, irxn, &soldel); @@ -548,7 +548,7 @@ int VCS_SOLVE::vcs_solve_TP(int print_lvl, int printDetails, int maxit) /********************** ZEROED OUT SPECIES **************************/ /********************************************************************/ bool resurrect = true; -#ifdef DEBUG +#ifdef DEBUG_MODE if (vcs_debug_print_lvl >= 3) { plogf(" --- %s currently zeroed (SpStatus=%-2d):", SpName[kspec].c_str(), spStatus[irxn]); @@ -567,7 +567,7 @@ int VCS_SOLVE::vcs_solve_TP(int print_lvl, int printDetails, int maxit) wt[kspec] = soln[kspec]; ds[kspec] = 0.0; resurrect = false; -#ifdef DEBUG +#ifdef DEBUG_MODE sprintf(ANOTE, "Species stays zeroed: DG = %11.4E", dg[irxn]); if (dg[irxn] < 0.0) { @@ -589,7 +589,7 @@ int VCS_SOLVE::vcs_solve_TP(int print_lvl, int printDetails, int maxit) if (atomComp > 0.0) { double maxPermissible = gai[j] / atomComp; if (maxPermissible < VCS_DELETE_MINORSPECIES_CUTOFF) { -#ifdef DEBUG +#ifdef DEBUG_MODE sprintf(ANOTE, "Species stays zeroed even though dG neg, because of %s elemAbund", ElName[j].c_str()); #endif @@ -606,7 +606,7 @@ int VCS_SOLVE::vcs_solve_TP(int print_lvl, int printDetails, int maxit) if (resurrect) { if (Vphase->Existence == 0) Vphase->Existence = 1; --m_numRxnMinorZeroed; -#ifdef DEBUG +#ifdef DEBUG_MODE if (vcs_debug_print_lvl >= 2) { plogf(" --- Zeroed species changed to major: "); plogf("%-12s\n", SpName[kspec].c_str()); @@ -624,7 +624,7 @@ int VCS_SOLVE::vcs_solve_TP(int print_lvl, int printDetails, int maxit) dx = wt[kspec] - soln[kspec]; } ds[kspec] = dx; -#ifdef DEBUG +#ifdef DEBUG_MODE sprintf(ANOTE, "Born:IC=-1 to IC=1:DG=%11.4E", dg[irxn]); #endif } else { @@ -644,7 +644,7 @@ int VCS_SOLVE::vcs_solve_TP(int print_lvl, int printDetails, int maxit) wt[kspec] = soln[kspec]; ds[kspec] = 0.0; dx = 0.0; -#ifdef DEBUG +#ifdef DEBUG_MODE sprintf(ANOTE,"minor species not considered"); if (vcs_debug_print_lvl >= 2) { plogf(" --- "); plogf("%-12s", SpName[kspec].c_str()); @@ -669,7 +669,7 @@ int VCS_SOLVE::vcs_solve_TP(int print_lvl, int printDetails, int maxit) * If soldel is true on return, then we branch to the section * that deletes a species from the current set of active species. */ -#ifdef DEBUG +#ifdef DEBUG_MODE dx = minor_alt_calc(kspec, irxn, &soldel, ANOTE); #else dx = minor_alt_calc(kspec, irxn, &soldel); @@ -680,7 +680,7 @@ int VCS_SOLVE::vcs_solve_TP(int print_lvl, int printDetails, int maxit) /***** DELETE MINOR SPECIES LESS THAN VCS_DELETE_SPECIES_CUTOFF */ /***** MOLE NUMBER */ /*******************************************************************/ -#ifdef DEBUG +#ifdef DEBUG_MODE if (vcs_debug_print_lvl >= 2) { plogf(" --- Delete minor species in multispec phase: %-12s\n", SpName[kspec].c_str()); @@ -704,7 +704,7 @@ int VCS_SOLVE::vcs_solve_TP(int print_lvl, int printDetails, int maxit) * Therefore, we will decrement the species counter, here. */ --irxn; -#ifdef DEBUG +#ifdef DEBUG_MODE goto L_MAIN_LOOP_END_NO_PRINT; #else goto L_MAIN_LOOP_END; @@ -714,7 +714,7 @@ int VCS_SOLVE::vcs_solve_TP(int print_lvl, int printDetails, int maxit) /********************************************************************/ /*********************** MAJOR SPECIES ******************************/ /********************************************************************/ -#ifdef DEBUG +#ifdef DEBUG_MODE sprintf(ANOTE, "Normal Major Calc"); #endif /* @@ -726,7 +726,7 @@ int VCS_SOLVE::vcs_solve_TP(int print_lvl, int printDetails, int maxit) wt[kspec] = soln[kspec]; ds[kspec] = 0.0; dx = 0.0; -#ifdef DEBUG +#ifdef DEBUG_MODE sprintf(ANOTE, "major species is converged"); if (vcs_debug_print_lvl >= 2) { plogf(" --- "); plogf("%-12s", SpName[kspec].c_str()); @@ -751,7 +751,7 @@ int VCS_SOLVE::vcs_solve_TP(int print_lvl, int printDetails, int maxit) } else { dx = 0.0; ds[kspec] = 0.0; -#ifdef DEBUG +#ifdef DEBUG_MODE sprintf(ANOTE, "dx set to 0, DG flipped sign due to " "changed initial point"); #endif @@ -767,7 +767,7 @@ int VCS_SOLVE::vcs_solve_TP(int print_lvl, int printDetails, int maxit) * or we restart the entire iteration. */ if (wt[kspec] <= 0.0) { -#ifdef DEBUG +#ifdef DEBUG_MODE sprintf(ANOTE, "initial nonpos moles= %11.3E", wt[kspec]); #endif @@ -803,7 +803,7 @@ int VCS_SOLVE::vcs_solve_TP(int print_lvl, int printDetails, int maxit) if (wt[kspec] < 0.005 * TMoles) { iph = PhaseID[kspec]; if (wt[kspec] < (TPhMoles[iph] * 0.01)) { -#ifdef DEBUG +#ifdef DEBUG_MODE if (vcs_debug_print_lvl >= 2) { plogf(" --- Major species changed to minor: "); plogf("%-12s\n", SpName[kspec].c_str()); @@ -840,7 +840,7 @@ int VCS_SOLVE::vcs_solve_TP(int print_lvl, int printDetails, int maxit) wt[kspec] = soln[kspec] + dx; if (wt[kspec] > 0.0) { ds[kspec] = dx; -#ifdef DEBUG +#ifdef DEBUG_MODE sprintf(ANOTE, "zeroing SS phase created a neg component species " "-> reducing step size instead"); @@ -853,7 +853,7 @@ int VCS_SOLVE::vcs_solve_TP(int print_lvl, int printDetails, int maxit) iph = PhaseID[kspec]; Vphase = VPhaseList[iph]; Vphase->Existence = 0; -#ifdef DEBUG +#ifdef DEBUG_MODE sprintf(ANOTE, "zero SS phase: moles went neg"); #endif /* @@ -896,7 +896,7 @@ int VCS_SOLVE::vcs_solve_TP(int print_lvl, int printDetails, int maxit) for (ll = irxn+1; ll < m_numRxnRdc; ++ll) { dgl[ll] = dg[ll]; } -#ifdef DEBUG +#ifdef DEBUG_MODE if (vcs_debug_print_lvl >= 2) { if (spStatus[irxn] >= 0) { plogf(" --- SS species changed to zeroedss: "); @@ -925,7 +925,7 @@ int VCS_SOLVE::vcs_solve_TP(int print_lvl, int printDetails, int maxit) if (dx != 0.0 && (soln[kspec] > 0.0) && (SpeciesUnknownType[kspec] != VCS_SPECIES_TYPE_INTERFACIALVOLTAGE)) { double dx_old = dx; -#ifdef DEBUG +#ifdef DEBUG_MODE dx = vcs_line_search(irxn, dx_old, ANOTE); #else dx = vcs_line_search(irxn, dx_old); @@ -943,7 +943,7 @@ int VCS_SOLVE::vcs_solve_TP(int print_lvl, int printDetails, int maxit) * to the reaction delta that we just computed. * This should keep the amount of material constant. */ -#ifdef DEBUG +#ifdef DEBUG_MODE if (ds[kspec] != dx) { plogf("we have a problem!\n"); exit(-1); @@ -962,16 +962,16 @@ int VCS_SOLVE::vcs_solve_TP(int print_lvl, int printDetails, int maxit) DelTPhMoles[iph] += dx * dnPhase_irxn[iph]; } } -#ifdef DEBUG +#ifdef DEBUG_MODE checkDelta1(VCS_DATA_PTR(ds), VCS_DATA_PTR(DelTPhMoles), kspec+1); #endif /* * Branch point for returning - */ -#ifndef DEBUG +#ifndef DEBUG_MODE L_MAIN_LOOP_END: ; #endif -#ifdef DEBUG +#ifdef DEBUG_MODE if (vcs_debug_print_lvl >= 2) { wt[kspec] = soln[kspec] + ds[kspec]; plogf(" --- "); plogf("%-12.12s", SpName[kspec].c_str()); @@ -983,7 +983,7 @@ int VCS_SOLVE::vcs_solve_TP(int print_lvl, int printDetails, int maxit) #endif /**************** END OF MAIN LOOP OVER FORMATION REACTIONS ************/ } -#ifdef DEBUG +#ifdef DEBUG_MODE if (vcs_debug_print_lvl >= 2) { for (k = 0; k < m_numComponents; k++) { plogf(" --- "); plogf("%-12.12s", SpName[k].c_str()); @@ -1007,7 +1007,7 @@ int VCS_SOLVE::vcs_solve_TP(int print_lvl, int printDetails, int maxit) xx = -ds[k] / soln[k]; if (par < xx) { par = xx; -#ifdef DEBUG +#ifdef DEBUG_MODE ll = k; #endif } @@ -1027,7 +1027,7 @@ int VCS_SOLVE::vcs_solve_TP(int print_lvl, int printDetails, int maxit) if (par <= 1.01 && par > 0.0) { /* Reduce the size of the step by the multiplicative factor, par */ par *= 0.99; -#ifdef DEBUG +#ifdef DEBUG_MODE if (vcs_debug_print_lvl >= 2) { plogf(" --- Reduction in step size due to component "); plogf("%s", SpName[ll].c_str()); @@ -1043,7 +1043,7 @@ int VCS_SOLVE::vcs_solve_TP(int print_lvl, int printDetails, int maxit) } else { par = 1.0; } -#ifdef DEBUG +#ifdef DEBUG_MODE checkDelta1(VCS_DATA_PTR(ds), VCS_DATA_PTR(DelTPhMoles), m_numSpeciesTot); #endif @@ -1091,13 +1091,13 @@ int VCS_SOLVE::vcs_solve_TP(int print_lvl, int printDetails, int maxit) vcs_dfe(VCS_DATA_PTR(wt), 1, 1, 0, m_numSpeciesRdc); if (printDetails) { if (iti != 0) { -#ifdef DEBUG +#ifdef DEBUG_MODE if (vcs_debug_print_lvl >= 2) { plogf(" *** vcs_dfe for printout only:"); } #endif vcs_dfe(VCS_DATA_PTR(wt), 1, 1, 0, m_numSpeciesRdc); -#ifdef DEBUG +#ifdef DEBUG_MODE if (vcs_debug_print_lvl >= 2) { plogf(" *** vcs_deltag for printout only:"); } @@ -1153,8 +1153,8 @@ int VCS_SOLVE::vcs_solve_TP(int print_lvl, int printDetails, int maxit) if (m_VCount->Its > 150) { plogf(" --- Troublesome solve\n"); } -#ifdef DEBUG -#ifdef DEBUG_MORE +#ifdef DEBUG_MODE +#ifdef DEBUG_NOT if (vcs_debug_print_lvl >= 3) { prneav(); } @@ -1181,14 +1181,14 @@ int VCS_SOLVE::vcs_solve_TP(int print_lvl, int printDetails, int maxit) if (printDetails && forced) { if (iti != 0) { -#ifdef DEBUG +#ifdef DEBUG_MODE if (vcs_debug_print_lvl >= 3) { plogf(" *** vcs_dfe for printout only:"); } #endif vcs_updateVP(0); vcs_dfe(VCS_DATA_PTR(soln), 0, 1, 0, m_numSpeciesRdc); -#ifdef DEBUG +#ifdef DEBUG_MODE if (vcs_debug_print_lvl >= 3) { plogf(" *** vcs_deltag call for printouts only;"); } @@ -1249,7 +1249,7 @@ int VCS_SOLVE::vcs_solve_TP(int print_lvl, int printDetails, int maxit) */ ++(m_VCount->Its); ++it1; -#ifdef DEBUG +#ifdef DEBUG_MODE if (vcs_debug_print_lvl >= 2) { plogf(" --- Increment counter increased, step is accepted: %4d\n", m_VCount->Its); @@ -1293,7 +1293,7 @@ int VCS_SOLVE::vcs_solve_TP(int print_lvl, int printDetails, int maxit) } } if (soldel) { -#ifdef DEBUG +#ifdef DEBUG_MODE if (vcs_debug_print_lvl >= 1) { plogf(" --- Setting microscopic phase %d to zero\n", iph); } @@ -1334,14 +1334,14 @@ int VCS_SOLVE::vcs_solve_TP(int print_lvl, int printDetails, int maxit) /*************************************************************************/ /***************** CHECK FOR ELEMENT ABUNDANCE****************************/ /*************************************************************************/ -#ifdef DEBUG +#ifdef DEBUG_MODE if (vcs_debug_print_lvl >= 2) { plogf(" --- Normal element abundance check"); } #endif vcs_elab(); if (! vcs_elabcheck(0)) { -#ifdef DEBUG +#ifdef DEBUG_MODE if (vcs_debug_print_lvl >= 2) { plogf(" - failed -> redoing element abundances.\n"); } @@ -1351,7 +1351,7 @@ int VCS_SOLVE::vcs_solve_TP(int print_lvl, int printDetails, int maxit) vcs_deltag(0, true); uptodate_minors = TRUE; } -#ifdef DEBUG +#ifdef DEBUG_MODE else { if (vcs_debug_print_lvl >= 2) { plogf(" - passed\n"); @@ -1386,7 +1386,7 @@ int VCS_SOLVE::vcs_solve_TP(int print_lvl, int printDetails, int maxit) for (j = m_numComponents - 1; j >= 0; j--) { if (soln[l] > soln[j]) { if (sc[i][j] != 0.0) { -#ifdef DEBUG +#ifdef DEBUG_MODE if (vcs_debug_print_lvl >= 2) { plogf(" --- Get a new basis because %s", SpName[l].c_str()); plogf(" is larger than comp %s", SpName[j].c_str()); @@ -1399,12 +1399,12 @@ int VCS_SOLVE::vcs_solve_TP(int print_lvl, int printDetails, int maxit) } else { break; } -#ifdef DEBUG_HKM +#ifdef DEBUG_NOT if (spStatus[i] == VCS_SPECIES_ZEROEDMS) { if (soln[j] == 0.0) { if (sc[i][j] != 0.0) { if (dg[i] < 0.0) { -#ifdef DEBUG +#ifdef DEBUG_MODE if (vcs_debug_print_lvl >= 2) { plogf(" --- Get a new basis because %s", SpName[l].c_str()); plogf(" has dg < 0.0 and comp %s has zero mole num", SpName[j].c_str()); @@ -1426,7 +1426,7 @@ int VCS_SOLVE::vcs_solve_TP(int print_lvl, int printDetails, int maxit) for (j = 0; j < m_numComponents; ++j) { if (soln[l] > soln[j]) { if (sc[i][j] != 0.0) { -#ifdef DEBUG +#ifdef DEBUG_MODE if (vcs_debug_print_lvl >= 2) { plogf(" --- Get a new basis because "); plogf("%s", SpName[l].c_str()); @@ -1439,12 +1439,12 @@ int VCS_SOLVE::vcs_solve_TP(int print_lvl, int printDetails, int maxit) goto L_COMPONENT_CALC; } } -#ifdef DEBUG_HKM +#ifdef DEBUG_NOT if (spStatus[i] == VCS_SPECIES_ZEROEDMS) { if (soln[j] == 0.0) { if (sc[i][j] != 0.0) { if (dg[i] < 0.0) { -#ifdef DEBUG +#ifdef DEBUG_MODE if (vcs_debug_print_lvl >= 2) { plogf(" --- Get a new basis because %s", SpName[l].c_str()); plogf(" has dg < 0.0 and comp %s has zero mole num", SpName[j].c_str()); @@ -1461,7 +1461,7 @@ int VCS_SOLVE::vcs_solve_TP(int print_lvl, int printDetails, int maxit) } } } -#ifdef DEBUG +#ifdef DEBUG_MODE if (vcs_debug_print_lvl >= 2) { plogf(" --- Check for an optimum basis passed\n"); } @@ -1474,7 +1474,7 @@ int VCS_SOLVE::vcs_solve_TP(int print_lvl, int printDetails, int maxit) * Go right to the check equilibrium section */ if (iti == 0) { -#ifdef DEBUG +#ifdef DEBUG_MODE if (vcs_debug_print_lvl >= 2) { plogf(" --- Reevaluate major-minor status of noncomponents:\n"); } @@ -1485,7 +1485,7 @@ int VCS_SOLVE::vcs_solve_TP(int print_lvl, int printDetails, int maxit) int speciesType = vcs_species_type(kspec); if (speciesType < VCS_SPECIES_MINOR) { -#ifdef DEBUG +#ifdef DEBUG_MODE if (vcs_debug_print_lvl >= 2) { if (spStatus[irxn] >= VCS_SPECIES_MINOR) { plogf(" --- major/minor species is now zeroed out: %s\n", @@ -1495,7 +1495,7 @@ int VCS_SOLVE::vcs_solve_TP(int print_lvl, int printDetails, int maxit) #endif ++m_numRxnMinorZeroed; } else if (speciesType == VCS_SPECIES_MINOR) { -#ifdef DEBUG +#ifdef DEBUG_MODE if (vcs_debug_print_lvl >= 2) { if (spStatus[irxn] != VCS_SPECIES_MINOR) { if (spStatus[irxn] == VCS_SPECIES_MAJOR) { @@ -1513,7 +1513,7 @@ int VCS_SOLVE::vcs_solve_TP(int print_lvl, int printDetails, int maxit) ++m_numRxnMinorZeroed; } else if (speciesType == VCS_SPECIES_MAJOR) { if (spStatus[irxn] != VCS_SPECIES_MAJOR) { -#ifdef DEBUG +#ifdef DEBUG_MODE if (vcs_debug_print_lvl >= 2) { if (spStatus[irxn] == VCS_SPECIES_MINOR) { plogf(" --- Noncomponent turned from minor to major: "); @@ -1548,7 +1548,7 @@ int VCS_SOLVE::vcs_solve_TP(int print_lvl, int printDetails, int maxit) /*************************************************************************/ L_EQUILIB_CHECK: ; if (! im) { -#ifdef DEBUG +#ifdef DEBUG_MODE if (vcs_debug_print_lvl >= 2) { plogf(" --- Equilibrium check for major species: "); } @@ -1563,7 +1563,7 @@ int VCS_SOLVE::vcs_solve_TP(int print_lvl, int printDetails, int maxit) */ goto L_RETURN_BLOCK; } else { -#ifdef DEBUG +#ifdef DEBUG_MODE if (vcs_debug_print_lvl >= 2) { plogf("%s failed\n", SpName[ir[irxn]].c_str()); } @@ -1581,13 +1581,13 @@ int VCS_SOLVE::vcs_solve_TP(int print_lvl, int printDetails, int maxit) } } } -#ifdef DEBUG +#ifdef DEBUG_MODE if (vcs_debug_print_lvl >= 2) { plogf(" MAJOR SPECIES CONVERGENCE achieved\n"); } #endif } -#ifdef DEBUG +#ifdef DEBUG_MODE else { if (vcs_debug_print_lvl >= 2) { plogf(" MAJOR SPECIES CONVERGENCE achieved " @@ -1613,7 +1613,7 @@ int VCS_SOLVE::vcs_solve_TP(int print_lvl, int printDetails, int maxit) vcs_deltag(1, false); uptodate_minors = TRUE; } -#ifdef DEBUG +#ifdef DEBUG_MODE if (vcs_debug_print_lvl >= 2) { plogf(" --- Equilibrium check for minor species: "); } @@ -1628,7 +1628,7 @@ int VCS_SOLVE::vcs_solve_TP(int print_lvl, int printDetails, int maxit) */ goto L_RETURN_BLOCK; } -#ifdef DEBUG +#ifdef DEBUG_MODE if (vcs_debug_print_lvl >= 2) { plogf("%s failed\n", SpName[ir[irxn]].c_str()); } @@ -1641,7 +1641,7 @@ int VCS_SOLVE::vcs_solve_TP(int print_lvl, int printDetails, int maxit) goto L_MAINLOOP_ALL_SPECIES; } } -#ifdef DEBUG +#ifdef DEBUG_MODE if (vcs_debug_print_lvl >= 2) { plogf(" CONVERGENCE achieved\n"); } @@ -1659,7 +1659,7 @@ int VCS_SOLVE::vcs_solve_TP(int print_lvl, int printDetails, int maxit) /* LEC is only true when we are near the end game */ if (lec) { if (!giveUpOnElemAbund) { -#ifdef DEBUG +#ifdef DEBUG_MODE if (vcs_debug_print_lvl >= 2) { plogf(" --- Check the Full Element Abundances: "); } @@ -1670,7 +1670,7 @@ int VCS_SOLVE::vcs_solve_TP(int print_lvl, int printDetails, int maxit) * the element abundances, and then go do a major step */ if (! vcs_elabcheck(1) ) { -#ifdef DEBUG +#ifdef DEBUG_MODE if (vcs_debug_print_lvl >= 2) { if (! vcs_elabcheck(0)) { plogf(" failed\n"); @@ -1682,7 +1682,7 @@ int VCS_SOLVE::vcs_solve_TP(int print_lvl, int printDetails, int maxit) // delete? goto L_ELEM_ABUND_CHECK; } -#ifdef DEBUG +#ifdef DEBUG_MODE if (vcs_debug_print_lvl >= 2) { plogf(" passed\n"); } @@ -1762,7 +1762,7 @@ int VCS_SOLVE::vcs_solve_TP(int print_lvl, int printDetails, int maxit) /* * Probably an unrecoverable range error */ -#ifdef DEBUG +#ifdef DEBUG_MODE if (vcs_debug_print_lvl >= 2) { plogf(" --- vcs_solve_tp: RANGE SPACE ERROR ENCOUNTERED\n"); plogf(" --- vcs_solve_tp: - Giving up on NE Element Abundance satisfaction \n"); @@ -1912,7 +1912,7 @@ int VCS_SOLVE::vcs_solve_TP(int print_lvl, int printDetails, int maxit) /*****************************************************************************/ double VCS_SOLVE::minor_alt_calc(int kspec, int irxn, int *do_delete -#ifdef DEBUG +#ifdef DEBUG_MODE , char *ANOTE #endif ) @@ -1973,7 +1973,7 @@ double VCS_SOLVE::minor_alt_calc(int kspec, int irxn, int *do_delete if (dg_irxn < -20.) { dg_irxn = -20.; } -#ifdef DEBUG +#ifdef DEBUG_MODE sprintf(ANOTE,"minor species alternative calc"); #endif if (dg_irxn >= 82.0) { @@ -2049,7 +2049,7 @@ double VCS_SOLVE::minor_alt_calc(int kspec, int irxn, int *do_delete * HKM -> Need to check the sign */ dx = dg[irxn]/ Faraday_dim; -#ifdef DEBUG +#ifdef DEBUG_MODE sprintf(ANOTE,"voltage species alternative calc"); #endif } @@ -2279,7 +2279,7 @@ void VCS_SOLVE::vcs_reinsert_deleted(int kspec) { int i, k, irxn = kspec - m_numComponents; int *phaseID = VCS_DATA_PTR(PhaseID); double dx; -#ifdef DEBUG +#ifdef DEBUG_MODE if (vcs_debug_print_lvl >= 2) { plogf(" --- Add back a deleted species: %-12s\n", SpName[kspec].c_str()); } @@ -2359,7 +2359,7 @@ void VCS_SOLVE::delete_multiphase(int iph) { * set the phase existence flag to dead */ Vphase->Existence = 0; -#ifdef DEBUG +#ifdef DEBUG_MODE if (vcs_debug_print_lvl >= 2) { plogf(" --- delete_multiphase %d, %s\n", iph, Vphase->PhaseName.c_str()); } @@ -2427,7 +2427,7 @@ void VCS_SOLVE::delete_multiphase(int iph) { ++(m_numRxnRdc); ++(m_numSpeciesRdc); -#ifdef DEBUG +#ifdef DEBUG_MODE if (vcs_debug_print_lvl >= 2) { plogf(" --- Make %s", SpName[kspec].c_str()); plogf(" an active but zeroed species because its phase " @@ -2461,7 +2461,7 @@ int VCS_SOLVE::recheck_deleted(void) { int iph, kspec, irxn, npb; double *xtcutoff = VCS_DATA_PTR(TmpPhase); -#ifdef DEBUG +#ifdef DEBUG_MODE if (vcs_debug_print_lvl >= 2) { plogf(" --- Start rechecking deleted species in multispec phases\n"); } @@ -2639,7 +2639,7 @@ int VCS_SOLVE::force(int iti) return FALSE; } #endif -#ifdef DEBUG +#ifdef DEBUG_MODE if (vcs_debug_print_lvl >= 2) { plogf(" --- subroutine FORCE: End Slope = %g\n", s2); } @@ -2654,7 +2654,7 @@ int VCS_SOLVE::force(int iti) } #ifdef DEBUG_NOT if (s1 >= 0.0) { -#ifdef DEBUG +#ifdef DEBUG_MODE if (vcs_debug_print_lvl >= 2) { plogf(" --- subroutine FORCE produced no adjustments,"); plogf(" failed s1 test -PROBLEM!!\n"); @@ -2663,7 +2663,7 @@ int VCS_SOLVE::force(int iti) return FALSE; } #endif -#ifdef DEBUG +#ifdef DEBUG_MODE if (vcs_debug_print_lvl >= 2) { plogf(" --- subroutine FORCE: Beginning Slope = %g\n", s1); } @@ -2676,14 +2676,14 @@ int VCS_SOLVE::force(int iti) al = s1 / (s1 - s2); } if (al >= 0.95 || al < 0.0) { -#ifdef DEBUG +#ifdef DEBUG_MODE if (vcs_debug_print_lvl >= 2) { plogf(" --- subroutine FORCE produced no adjustments (al = %g)\n", al); } #endif return FALSE; } -#ifdef DEBUG +#ifdef DEBUG_MODE if (vcs_debug_print_lvl >= 2) { plogf(" --- subroutine FORCE produced a damping factor = %g\n", al); } @@ -2700,7 +2700,7 @@ int VCS_SOLVE::force(int iti) } vcs_updateVP(0); -#ifdef DEBUG +#ifdef DEBUG_MODE if (vcs_debug_print_lvl >= 2) { plogf(" --- subroutine FORCE adjusted the mole " "numbers, AL = %10.3f\n", al); @@ -2746,7 +2746,7 @@ int VCS_SOLVE::vcs_RxnStepSizes() { double s, xx, dss; vcs_VolPhase *Vphase = 0; double *dnPhase_irxn; -#ifdef DEBUG +#ifdef DEBUG_MODE char ANOTE[128]; if (vcs_debug_print_lvl >= 2) { plogf(" "); for (j = 0; j < 82; j++) plogf("-"); plogf("\n"); @@ -2768,7 +2768,7 @@ int VCS_SOLVE::vcs_RxnStepSizes() { ************************************************************************/ for (irxn = 0; irxn < m_numRxnRdc; ++irxn) { -#ifdef DEBUG +#ifdef DEBUG_MODE sprintf(ANOTE,"Normal Calc"); #endif @@ -2797,13 +2797,13 @@ int VCS_SOLVE::vcs_RxnStepSizes() { double trphmoles = tphmoles / TMoles; if (trphmoles > VCS_DELETE_PHASE_CUTOFF) { ds[kspec] = TMoles * VCS_SMALL_MULTIPHASE_SPECIES; -#ifdef DEBUG +#ifdef DEBUG_MODE sprintf(ANOTE, "MultSpec: small species born again DG = %11.3E", dg[irxn]); #endif } else { -#ifdef DEBUG +#ifdef DEBUG_MODE sprintf(ANOTE, "MultSpec: phase come alive DG = %11.3E", dg[irxn]); #endif Vphase = VPhaseList[iph]; @@ -2812,7 +2812,7 @@ int VCS_SOLVE::vcs_RxnStepSizes() { } --(m_numRxnMinorZeroed); } else { -#ifdef DEBUG +#ifdef DEBUG_MODE sprintf(ANOTE, "MultSpec: still dead DG = %11.3E", dg[irxn]); #endif ds[kspec] = 0.0; @@ -2829,7 +2829,7 @@ int VCS_SOLVE::vcs_RxnStepSizes() { * in this mode. */ if (fabs(dg[irxn]) <= tolmaj2) { -#ifdef DEBUG +#ifdef DEBUG_MODE sprintf(ANOTE,"Skipped: superconverged DG = %11.3E", dg[irxn]); if (vcs_debug_print_lvl >= 2) { plogf(" --- %-12.12s", SpName[kspec].c_str()); @@ -2844,7 +2844,7 @@ int VCS_SOLVE::vcs_RxnStepSizes() { * their values are to be decreasing anyway. */ if ((spStatus[irxn] != VCS_SPECIES_MAJOR) && (dg[irxn] >= 0.0)) { -#ifdef DEBUG +#ifdef DEBUG_MODE sprintf(ANOTE,"Skipped: IC = %3d and DG >0: %11.3E", spStatus[irxn], dg[irxn]); if (vcs_debug_print_lvl >= 2) { @@ -2886,7 +2886,7 @@ int VCS_SOLVE::vcs_RxnStepSizes() { if (UseActCoeffJac) { double s_old = s; s = vcs_Hessian_diag_adj(irxn, s_old); -#ifdef DEBUG +#ifdef DEBUG_MODE if (s_old != s) { sprintf(ANOTE, "Normal calc: diag adjusted from %g " "to %g due to act coeff", s_old, s); @@ -2905,14 +2905,14 @@ int VCS_SOLVE::vcs_RxnStepSizes() { double negChangeComp = - stoicC * ds[kspec]; if (negChangeComp > soln[j]) { if (soln[j] > 0.0) { -#ifdef DEBUG +#ifdef DEBUG_MODE sprintf(ANOTE, "Delta damped from %g " "to %g due to component %d (%10s) going neg", ds[kspec], -soln[j]/stoicC, j, SpName[j].c_str()); #endif ds[kspec] = - soln[j] / stoicC; } else { -#ifdef DEBUG +#ifdef DEBUG_MODE sprintf(ANOTE, "Delta damped from %g " "to %g due to component %d (%10s) zero", ds[kspec], -soln[j]/stoicC, j, SpName[j].c_str()); @@ -2924,7 +2924,7 @@ int VCS_SOLVE::vcs_RxnStepSizes() { } // Implement a damping term that limits ds to the size of the mole number if (-ds[kspec] > soln[kspec]) { -#ifdef DEBUG +#ifdef DEBUG_MODE sprintf(ANOTE, "Delta damped from %g " "to %g due to %s going negative", ds[kspec], -soln[kspec], SpName[kspec].c_str()); @@ -2990,7 +2990,7 @@ int VCS_SOLVE::vcs_RxnStepSizes() { Vphase = VPhaseList[iph]; Vphase->Existence = 0; TPhMoles[iph] = 0.0; -#ifdef DEBUG +#ifdef DEBUG_MODE if (vcs_debug_print_lvl >= 2) { plogf(" --- vcs_RxnStepSizes Special section to delete %s\n", SpName[k].c_str()); @@ -3008,7 +3008,7 @@ int VCS_SOLVE::vcs_RxnStepSizes() { } } } /* End of regular processing */ -#ifdef DEBUG +#ifdef DEBUG_MODE if (vcs_debug_print_lvl >= 2) { plogf(" --- %-12.12s", SpName[kspec].c_str()); plogf(" %12.4E %12.4E %12.4E | %s\n", @@ -3017,7 +3017,7 @@ int VCS_SOLVE::vcs_RxnStepSizes() { #endif } /* End of loop over SpeciesUnknownType */ } /* End of loop over non-component stoichiometric formation reactions */ -#ifdef DEBUG +#ifdef DEBUG_MODE if (vcs_debug_print_lvl >= 2) { plogf(" "); vcs_print_line("-", 82); } @@ -3063,7 +3063,7 @@ void VCS_SOLVE::vcs_deltag(int l, bool doDeleted) { irxnl = m_numRxnTot; } -#ifdef DEBUG +#ifdef DEBUG_MODE if (vcs_debug_print_lvl >= 2) { plogf(" --- Subroutine vcs_deltag called for "); if (l < 0) { @@ -3277,7 +3277,7 @@ int VCS_SOLVE::vcs_basopt(int ifirst, double aw[], double sa[], double sm[], int jlose = -1; double *dptr, *scrxn_ptr; double tsecond = vcs_second(); -#ifdef DEBUG +#ifdef DEBUG_MODE if (vcs_debug_print_lvl >= 2) { plogf(" "); for(i=0; i<77; i++) plogf("-"); plogf("\n"); plogf(" --- Subroutine BASOPT called to "); @@ -3465,7 +3465,7 @@ int VCS_SOLVE::vcs_basopt(int ifirst, double aw[], double sa[], double sm[], for (i = 0; i < m_numSpeciesTot; ++i) { ir[i] = ncTrial + i; } -#ifdef DEBUG +#ifdef DEBUG_MODE if (vcs_debug_print_lvl >= 2) { plogf(" --- Total number of components found = %3d (ne = %d)\n ", ncTrial, m_numElemConstraints); @@ -3531,7 +3531,7 @@ int VCS_SOLVE::vcs_basopt(int ifirst, double aw[], double sa[], double sm[], /* **** REARRANGE THE DATA ****************** */ /* ****************************************** */ if (jr != k) { -#ifdef DEBUG +#ifdef DEBUG_MODE if (vcs_debug_print_lvl >= 2) { plogf(" --- %-12.12s", (SpName[k]).c_str()); plogf("(%9.2g) replaces %-12.12s", soln[k], SpName[jr].c_str()); @@ -3541,7 +3541,7 @@ int VCS_SOLVE::vcs_basopt(int ifirst, double aw[], double sa[], double sm[], vcs_switch_pos(FALSE, jr, k); vcsUtil_dsw(aw, jr, k); } -#ifdef DEBUG +#ifdef DEBUG_MODE else { if (vcs_debug_print_lvl >= 2) { plogf(" --- %-12.12s", SpName[k].c_str()); @@ -3682,7 +3682,7 @@ int VCS_SOLVE::vcs_basopt(int ifirst, double aw[], double sa[], double sm[], } -#ifdef DEBUG +#ifdef DEBUG_MODE if (vcs_debug_print_lvl >= 2) { plogf(" --- Components:"); for (j = 0; j < ncTrial; j++) { @@ -3805,7 +3805,7 @@ int VCS_SOLVE::vcs_species_type(int kspec) double negChangeComp = - stoicC; if (negChangeComp > 0.0) { if (soln[j] < 1.0E-60) { -#ifdef DEBUG +#ifdef DEBUG_MODE if (vcs_debug_print_lvl >= 2) { plogf(" --- %s would have popped back into existance but" " needed component %s is zero\n", @@ -3829,7 +3829,7 @@ int VCS_SOLVE::vcs_species_type(int kspec) if (atomComp > 0.0) { double maxPermissible = gai[j] / atomComp; if (maxPermissible < VCS_DELETE_MINORSPECIES_CUTOFF) { -#ifdef DEBUG +#ifdef DEBUG_MODE if (vcs_debug_print_lvl >= 2) { plogf(" --- %s would have popped back into existance but" " needed element %s is zero\n", @@ -3948,7 +3948,7 @@ void VCS_SOLVE::vcs_chemPotPhase(int iph, const double *const molNum, int nkk = Vphase->NVolSpecies; int k, kspec; -#ifdef DEBUG +#ifdef DEBUG_MODE //if (vcs_debug_print_lvl >= 2) { // plogf(" --- Subroutine vcs_chemPotPhase called for phase %d\n", // iph); @@ -3980,7 +3980,7 @@ void VCS_SOLVE::vcs_chemPotPhase(int iph, const double *const molNum, } } if (SpeciesUnknownType[kspec] == VCS_SPECIES_TYPE_INTERFACIALVOLTAGE) { -#ifdef DEBUG +#ifdef DEBUG_MODE if (molNum[kspec] != phi) { plogf("We have an inconsistency!\n"); exit(-1); @@ -4107,7 +4107,7 @@ void VCS_SOLVE::vcs_dfe(double *z, int kk, int ll, int lbot, int ltop) vcs_VolPhase *Vphase; VCS_SPECIES_THERMO *st_ptr; -#ifdef DEBUG +#ifdef DEBUG_MODE if (vcs_debug_print_lvl >= 2) { if (ll == 0) { if (lbot != 0) { @@ -4146,7 +4146,7 @@ void VCS_SOLVE::vcs_dfe(double *z, int kk, int ll, int lbot, int ltop) tlogMoles[iph] += z[kspec]; } } -#ifdef DEBUG +#ifdef DEBUG_MODE for (iph = 0; iph < NPhase; iph++) { if (! vcs_doubleEqual(tlogMoles[iph], tPhMoles_ptr[iph])) { plogf("phase Moles may be off, iph = %d, %20.14g %20.14g \n", @@ -4202,7 +4202,7 @@ void VCS_SOLVE::vcs_dfe(double *z, int kk, int ll, int lbot, int ltop) for (kspec = l1; kspec < l2; ++kspec) { iphase = PhaseID[kspec]; if (SpeciesUnknownType[kspec] == VCS_SPECIES_TYPE_INTERFACIALVOLTAGE) { -#ifdef DEBUG +#ifdef DEBUG_MODE if (z[kspec] != phasePhi[iphase]) { plogf("We have an inconsistency!\n"); exit(-1); @@ -4244,7 +4244,7 @@ void VCS_SOLVE::vcs_dfe(double *z, int kk, int ll, int lbot, int ltop) kspec = ir[irxn]; iphase = PhaseID[kspec]; if (SpeciesUnknownType[kspec] == VCS_SPECIES_TYPE_INTERFACIALVOLTAGE) { -#ifdef DEBUG +#ifdef DEBUG_MODE if (z[kspec] != phasePhi[iphase]) { plogf("We have an inconsistency!\n"); exit(-1); @@ -4287,7 +4287,7 @@ void VCS_SOLVE::vcs_dfe(double *z, int kk, int ll, int lbot, int ltop) kspec = ir[irxn]; iphase = PhaseID[kspec]; if (SpeciesUnknownType[kspec] == VCS_SPECIES_TYPE_INTERFACIALVOLTAGE) { -#ifdef DEBUG +#ifdef DEBUG_MODE if (z[kspec] != phasePhi[iphase]) { plogf("We have an inconsistency!\n"); exit(-1); @@ -4478,7 +4478,7 @@ void vcsUtil_stsw(std::vector & vstr, int i1, int i2) /*****************************************************************************/ /*****************************************************************************/ /*****************************************************************************/ -#ifdef DEBUG +#ifdef DEBUG_MODE void VCS_SOLVE::prneav(void) @@ -4672,7 +4672,7 @@ void VCS_SOLVE::vcs_switch_pos(int ifunc, int k1, int k2) vcs_VolPhase *pv1, *pv2; VCS_SPECIES_THERMO *st_tmp; if (k1 == k2) return; -#ifdef DEBUG +#ifdef DEBUG_MODE if (k1 < 0 || k1 > (m_numSpeciesTot - 1) || k2 < 0 || k2 > (m_numSpeciesTot - 1) ) { plogf("vcs_switch_pos: ifunc = 0: inappropriate args: %d %d\n", @@ -4687,7 +4687,7 @@ void VCS_SOLVE::vcs_switch_pos(int ifunc, int k1, int k2) kp1 = indPhSp[k1]; kp2 = indPhSp[k2]; -#ifdef DEBUG +#ifdef DEBUG_MODE if (pv1->IndSpecies[kp1] != k1) { plogf("Indexing error in program\n"); exit(-1); @@ -4740,7 +4740,7 @@ void VCS_SOLVE::vcs_switch_pos(int ifunc, int k1, int k2) */ i1 = k1 - m_numComponents; i2 = k2 - m_numComponents; -#ifdef DEBUG +#ifdef DEBUG_MODE if (i1 < 0 || i1 > (m_numRxnTot - 1) || i2 < 0 || i2 > (m_numRxnTot - 1) ) { plogf("switch_pos: ifunc = 1: inappropriate noncomp values: %d %d\n", @@ -4795,7 +4795,7 @@ void VCS_SOLVE::vcs_deltag_Phase(int iphase, bool doDeleted) { if (doDeleted) irxnl = m_numRxnTot; vcs_VolPhase *vPhase = VPhaseList[iphase]; -#ifdef DEBUG +#ifdef DEBUG_MODE if (vcs_debug_print_lvl >= 2) { plogf(" --- Subroutine vcs_deltag_Phase called for phase %d\n", iphase); @@ -4807,7 +4807,7 @@ void VCS_SOLVE::vcs_deltag_Phase(int iphase, bool doDeleted) { */ if (vPhase->SingleSpecies) { kspec = vPhase->IndSpecies[0]; -#ifdef DEBUG +#ifdef DEBUG_MODE if (iphase != PhaseID[kspec]) { plogf("vcs_deltag_Phase index error\n"); exit(-1); @@ -4938,7 +4938,7 @@ double VCS_SOLVE::vcs_birthGuess(int kspec) { /* * Logic to handle species in multiple species phases */ -#ifdef DEBUG +#ifdef DEBUG_MODE char ANOTE[32]; double dxm = minor_alt_calc(kspec, irxn, &soldel, ANOTE); #else diff --git a/Cantera/src/equil/vcs_species_thermo.cpp b/Cantera/src/equil/vcs_species_thermo.cpp index 63bf82efc..8b9b890a6 100644 --- a/Cantera/src/equil/vcs_species_thermo.cpp +++ b/Cantera/src/equil/vcs_species_thermo.cpp @@ -291,7 +291,7 @@ VolStar_calc(int kglob, double TKelvin, double pres) */ double VCS_SPECIES_THERMO::G0_R_calc(int kglob, double TKelvin) { -#ifdef DEBUG +#ifdef DEBUG_MODE char yo[] = "VS_SPECIES_THERMO::G0_R_calc "; #endif double fe, H, S; @@ -323,7 +323,7 @@ double VCS_SPECIES_THERMO::G0_R_calc(int kglob, double TKelvin) fe = vcs_G0_NASA(TKelvin, (VCS_NASA_POLY *) SS0_Params); break; default: -#ifdef DEBUG +#ifdef DEBUG_MODE plogf("%sERROR: unknown model\n", yo); #endif exit(-1); @@ -355,7 +355,7 @@ double VCS_SPECIES_THERMO::G0_R_calc(int kglob, double TKelvin) */ double VCS_SPECIES_THERMO::eval_ac(int kglob) { -#ifdef DEBUG +#ifdef DEBUG_MODE char yo[] = "VCS_SPECIES_THERMO::eval_ac "; #endif double ac; @@ -388,7 +388,7 @@ double VCS_SPECIES_THERMO::eval_ac(int kglob) plogf("Not implemented Yet\n"); exit(-1); default: -#ifdef DEBUG +#ifdef DEBUG_MODE plogf("%sERROR: unknown model\n", yo); #endif exit(-1);