diff --git a/include/cantera/kinetics/ReactionStoichMgr.h b/include/cantera/kinetics/ReactionStoichMgr.h index a2aeab9a7..14280a4a3 100644 --- a/include/cantera/kinetics/ReactionStoichMgr.h +++ b/include/cantera/kinetics/ReactionStoichMgr.h @@ -239,10 +239,6 @@ protected: StoichManagerN* m_revproducts; StoichManagerN* m_irrevproducts; vector_fp m_dummy; - -#ifdef INCL_STOICH_WRITER - StoichWriter* m_rwriter; -#endif }; } diff --git a/include/cantera/numerics/DenseMatrix.h b/include/cantera/numerics/DenseMatrix.h index b503f95c4..f7946b573 100644 --- a/include/cantera/numerics/DenseMatrix.h +++ b/include/cantera/numerics/DenseMatrix.h @@ -195,10 +195,6 @@ public: friend int solve(DenseMatrix& A, double* b); friend int solve(DenseMatrix& A, DenseMatrix& b); friend int invert(DenseMatrix& A, int nn); -#ifdef INCL_LEAST_SQUARES - friend int leastSquares(DenseMatrix& A, double* b); -#endif - }; //================================================================================================================== @@ -228,16 +224,6 @@ int solve(DenseMatrix& A, double* b); */ int solve(DenseMatrix& A, DenseMatrix& b); -#ifdef INCL_LEAST_SQUARES -//! Solve Ax = b in the least squares sense -/*! - * @param A Matrix to be inverted in the least squares sense - * @param b Vector b to be solved for - * @todo fix lwork - */ -int leastSquares(DenseMatrix& A, double* b); -#endif - //! Multiply \c A*b and return the result in \c prod. Uses BLAS routine DGEMV. /*! * \f[ diff --git a/include/cantera/numerics/Func1.h b/include/cantera/numerics/Func1.h index 50e229d2e..e2b2e68c3 100644 --- a/include/cantera/numerics/Func1.h +++ b/include/cantera/numerics/Func1.h @@ -7,8 +7,6 @@ #ifndef CT_FUNC1_H #define CT_FUNC1_H -#undef DEBUG_FUNC - #include "cantera/base/ct_defs.h" #include @@ -794,16 +792,6 @@ public: Func1* d3 = &newCompositeFunction(*d1, m_f2->duplicate()); Func1* d2 = &m_f2->derivative(); Func1* p = &newProdFunction(*d3, *d2); -#ifdef DEBUG_FUNC - cout << "Composite1::derivative: \n"; - cout << "f1 = " << m_f1->write("x") << endl; - cout << "f2 = " << m_f2->write("x") << endl; - cout << "d1 = " << d1 << " " << d1->write("x") << endl; - cout << "d3 = " << d3->write("x") << endl; - cout << "d2 = " << d2->write("x") << endl; - cout << "function = \'" + write("x") + "\'\n"; - cout << "derivative = \'" + p->write("x") + "\'\n"; -#endif return *p; } diff --git a/include/cantera/numerics/ctlapack.h b/include/cantera/numerics/ctlapack.h index 327dece0f..d5c3367e8 100644 --- a/include/cantera/numerics/ctlapack.h +++ b/include/cantera/numerics/ctlapack.h @@ -26,7 +26,6 @@ typedef int ftnlen; #define _DGETRF_ dgetrf #define _DGETRS_ dgetrs #define _DGETRI_ dgetri -#define _DGELSS_ dgelss #define _DGBCON_ dgbcon #define _DGBSV_ dgbsv #define _DGBTRF_ dgbtrf @@ -49,7 +48,6 @@ typedef int ftnlen; #define _DGETRF_ dgetrf_ #define _DGETRS_ dgetrs_ #define _DGETRI_ dgetri_ -#define _DGELSS_ dgelss_ #define _DGBCON_ dgbcon_ #define _DGBSV_ dgbsv_ #define _DGBTRF_ dgbtrf_ @@ -128,12 +126,6 @@ extern "C" { int _DGETRI_(const integer* n, doublereal* a, const integer* lda, integer* ipiv, doublereal* work, integer* lwork, integer* info); - int _DGELSS_(integer* m, integer* n, integer* nrhs, - doublereal* a, integer* lda, doublereal* b, integer* ldb, doublereal * - s, doublereal* rcond, integer* rank, doublereal* work, integer* lwork, - integer* info); - - int _DGBSV_(integer* n, integer* kl, integer* ku, integer* nrhs, doublereal* a, integer* lda, integer* ipiv, doublereal* b, integer* ldb, integer* info); @@ -372,21 +364,6 @@ inline void ct_dgetri(int n, doublereal* a, int lda, integer* ipiv, integer f_n = n, f_lda = lda, f_lwork = lwork, f_info = info; _DGETRI_(&f_n, a, &f_lda, ipiv, work, &f_lwork, &f_info); } -//==================================================================================================================== -inline void ct_dgelss(int m, int n, int nrhs, doublereal* a, - int lda, doublereal* b, int ldb, doublereal* s, - doublereal rcond, int& rank, doublereal* work, int lwork, - int& info) -{ - doublereal f_rcond = rcond; - integer f_m = m, f_n = n, f_nrhs = nrhs, f_lda = lda, f_ldb = ldb, - f_rank = rank, f_info = info, f_lwork = lwork; - //f_lwork = 2*(3*min(m,n) + max(2*min(m,n), max(m,n))); - _DGELSS_(&f_m, &f_n, &f_nrhs, a, &f_lda, b, &f_ldb, s, &f_rcond, - &f_rank, work, &f_lwork, &f_info); - info = f_info; - rank = f_rank; -} inline void ct_dscal(int n, doublereal da, doublereal* dx, int incx) { diff --git a/include/cantera/thermo/NasaPoly1.h b/include/cantera/thermo/NasaPoly1.h index ddd5eba50..c5643d0c7 100644 --- a/include/cantera/thermo/NasaPoly1.h +++ b/include/cantera/thermo/NasaPoly1.h @@ -252,13 +252,6 @@ public: for (int i = 2; i < 7; i++) { coeffs[i-2] = m_coeff[i]; } -#ifdef WARN_ABOUT_CHANGES_FROM_VERSION_1_6 - cout << "************************************************\n" - cout << "Warning: NasaPoly1::reportParameters now returns \n" - << "the coefficient array in the same order as in\n" - << "the input file. See file NasaPoly1.h" << endl; - cout << "************************************************\n" -#endif } //! Modify parameters for the standard state diff --git a/src/apps/csvdiff.cpp b/src/apps/csvdiff.cpp index 4a1af1099..8b2f408ad 100644 --- a/src/apps/csvdiff.cpp +++ b/src/apps/csvdiff.cpp @@ -915,11 +915,7 @@ int main(int argc, char* argv[]) /* * Compare the solutions in each file */ - int method = 0; -#define DGG_MODS -#ifdef DGG_MODS - method = 1; -#endif + int method = 1; double slope1, slope2, xatol; int notOK; for (k = 0; k < nColcomparisons; k++) { diff --git a/src/base/application.cpp b/src/base/application.cpp index 35bd35e16..060c4d5ff 100644 --- a/src/base/application.cpp +++ b/src/base/application.cpp @@ -495,20 +495,21 @@ XML_Node* Application::get_XML_File(std::string file, int debug) } else { ff = string("./")+path.substr(0,idot) + ".xml"; } -#ifdef DEBUG_PATHS - cout << "get_XML_File(): Expected location of xml file = " - << ff << endl; -#endif + if (debug > 0) { + writelog("get_XML_File(): Expected location of xml file = " + + ff + "\n"); + } /* * Do a search of the existing XML trees to determine if we have * already processed this file. If we have, return a pointer to * the processed xml tree. */ if (xmlfiles.find(ff) != xmlfiles.end()) { -#ifdef DEBUG_PATHS - cout << "get_XML_File(): File, " << ff << ", was previously read." - << " Retrieving the stored xml tree." << endl; -#endif + if (debug > 0) { + writelog("get_XML_File(): File, " + ff + + ", was previously read." + + " Retrieving the stored xml tree.\n"); + } return xmlfiles[ff]; } /* diff --git a/src/base/ct2ctml.cpp b/src/base/ct2ctml.cpp index c0a800e88..970b75f2b 100644 --- a/src/base/ct2ctml.cpp +++ b/src/base/ct2ctml.cpp @@ -106,9 +106,6 @@ void ct2ctml(const char* file, const int debug) #else string cmd = "sleep " + sleep() + "; " + "\"" + pypath() + "\"" + " " + "\"" + path + "\"" + " &> " + logfile; -#endif -#ifdef DEBUG_PATHS - writelog("ct2ctml: executing the command " + cmd + "\n"); #endif if (debug > 0) { writelog("ct2ctml: executing the command " + cmd + "\n"); @@ -201,14 +198,10 @@ void ct2ctml(const char* file, const int debug) */ void get_CTML_Tree(Cantera::XML_Node* rootPtr, const std::string file, const int debug) { - std::string ff, ext = ""; // find the input file on the Cantera search path std::string inname = findInputFile(file); -#ifdef DEBUG_PATHS - writelog("Found file: "+inname+"\n"); -#endif if (debug > 0) { writelog("Found file: "+inname+"\n"); } @@ -234,20 +227,16 @@ void get_CTML_Tree(Cantera::XML_Node* rootPtr, const std::string file, const int } string ffull = inname.substr(0,idot) + ".xml"; ff = "./" + getBaseName(ffull) + ".xml"; -#ifdef DEBUG_PATHS - writelogf("ffull name = %s\n", ffull.c_str()); - writelogf("ff name = %s\n", ff.c_str()); -#endif + if (debug > 0) { + writelogf("ffull name = %s\n", ffull.c_str()); + writelogf("ff name = %s\n", ff.c_str()); + } } else { ff = inname; } -#ifdef DEBUG_PATHS - writelog("Attempting to parse xml file " + ff + "\n"); -#else if (debug > 0) { writelog("Attempting to parse xml file " + ff + "\n"); } -#endif ifstream fin(ff.c_str()); if (!fin) { throw diff --git a/src/converters/CKParser.cpp b/src/converters/CKParser.cpp index 6aa2276ca..177e4aea2 100644 --- a/src/converters/CKParser.cpp +++ b/src/converters/CKParser.cpp @@ -359,12 +359,6 @@ void CKParser::getCKLine(std::string& s, std::string& comment) // if an end-of-line character is seen, then break. // Check for all common end-of-line characters. if (ch == char13 || (ch == char10 && (m_last_eol != char13))) { -#undef DEBUG_EOL -#ifdef DEBUG_EOL - *m_log << "EOL: found character " << int(ch) << " ending line:" << endl; - *m_log << line << endl; - *m_log << int(m_last_eol) << " " << int('\n') << " " << int(ch) << endl; -#endif m_last_eol = ch; break; } diff --git a/src/converters/ck2ct.cpp b/src/converters/ck2ct.cpp index 104928da0..6e79be5f3 100644 --- a/src/converters/ck2ct.cpp +++ b/src/converters/ck2ct.cpp @@ -202,8 +202,6 @@ static void addTransportParams(FILE* f, string name) throw CanteraError("addTransportParams", "Unrecognized geometry flag for species " + name); } -#define FULL_TRANSPORT_PARAMETER_PRECISION -#ifdef FULL_TRANSPORT_PARAMETER_PRECISION fprintf(f," diam = %g,\n",td.diam); fprintf(f," well_depth = %g",td.welldepth); if (td.polar != 0.0) { @@ -215,19 +213,6 @@ static void addTransportParams(FILE* f, string name) if (td.rot != 0.0) { fprintf(f,",\n rot_relax = %g",td.rot); } -#else - fprintf(f," diam = %8.2f,\n",td.diam); - fprintf(f," well_depth = %8.2f",td.welldepth); - if (td.polar != 0.0) { - fprintf(f,",\n polar = %8.2f",td.polar); - } - if (td.dipole != 0.0) { - fprintf(f,",\n dipole = %8.2f",td.dipole); - } - if (td.rot != 0.0) { - fprintf(f,",\n rot_relax = %8.2f",td.rot); - } -#endif fprintf(f,")"); } diff --git a/src/equil/vcs_VolPhase.cpp b/src/equil/vcs_VolPhase.cpp index 55611aaa1..270fb2acb 100644 --- a/src/equil/vcs_VolPhase.cpp +++ b/src/equil/vcs_VolPhase.cpp @@ -1004,7 +1004,6 @@ void vcs_VolPhase::_updateLnActCoeffJac() if (!m_UpToDate_AC) { _updateActCoeff(); } -#ifndef NOOLD if (!TP_ptr) { return; } @@ -1019,8 +1018,6 @@ void vcs_VolPhase::_updateLnActCoeffJac() lnActCoeffCol[k] /= moles_j_base; } } -#endif - double deltaMoles_j = 0.0; // Make copies of ActCoeff and Xmol_ for use in taking differences diff --git a/src/equil/vcs_rxnadj.cpp b/src/equil/vcs_rxnadj.cpp index 41fb27412..692a9acd7 100644 --- a/src/equil/vcs_rxnadj.cpp +++ b/src/equil/vcs_rxnadj.cpp @@ -339,27 +339,6 @@ size_t VCS_SOLVE::vcs_RxnStepSizes(int& forceComponentCalc, size_t& kSpecial) /* * Delete the single species phase */ -#ifdef OLDSTUFF - m_molNumSpecies_old[kspec] += dss; - m_tPhaseMoles_old[m_phaseID[kspec]] += dss; - for (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]; - } - m_molNumSpecies_old[k] = 0.0; - iph = m_phaseID[k]; - m_tPhaseMoles_old[iph] = 0.0; - Vphase = m_VolPhaseList[iph]; - Vphase->setTotalMoles(0.0); - if (k == kspec) { - m_speciesStatus[kspec] = VCS_SPECIES_ZEROEDSS; - if (m_SSPhase[kspec] != 1) { - printf("vcs_RxnStepSizes:: we shouldn't be here!\n"); - exit(EXIT_FAILURE); - } - } -#else - for (size_t j = 0; j < m_numSpeciesTot; j++) { m_deltaMolNumSpecies[j] = 0.0; } @@ -404,7 +383,6 @@ size_t VCS_SOLVE::vcs_RxnStepSizes(int& forceComponentCalc, size_t& kSpecial) } #endif return iphDel; -#endif } } } /* End of regular processing */ diff --git a/src/kinetics/ReactionStoichMgr.cpp b/src/kinetics/ReactionStoichMgr.cpp index 83877b782..97ff9d78e 100644 --- a/src/kinetics/ReactionStoichMgr.cpp +++ b/src/kinetics/ReactionStoichMgr.cpp @@ -22,16 +22,10 @@ ReactionStoichMgr::ReactionStoichMgr() : m_reactants(0), m_revproducts(0), m_irrevproducts(0) -#ifdef INCL_STOICH_WRITER - , m_rwriter(0) -#endif { m_reactants = new StoichManagerN; m_revproducts = new StoichManagerN; m_irrevproducts = new StoichManagerN; -#ifdef INCL_STOICH_WRITER - m_rwriter = new StoichWriter; -#endif m_dummy.resize(10,1.0); } //==================================================================================================================== @@ -41,26 +35,17 @@ ReactionStoichMgr::~ReactionStoichMgr() delete m_reactants; delete m_revproducts; delete m_irrevproducts; -#ifdef INCL_STOICH_WRITER - delete m_rwriter; -#endif } //==================================================================================================================== ReactionStoichMgr::ReactionStoichMgr(const ReactionStoichMgr& right) : m_reactants(0), m_revproducts(0), m_irrevproducts(0) -#ifdef INCL_STOICH_WRITER - , m_rwriter(0) -#endif { m_reactants = new StoichManagerN(*right.m_reactants); m_revproducts = new StoichManagerN(*right.m_revproducts); m_irrevproducts = new StoichManagerN(*right.m_irrevproducts); m_dummy = right.m_dummy; -#ifdef INCL_STOICH_WRITER - m_rwriter = new StoichManagerN(right.m_writer); -#endif } //==================================================================================================================== ReactionStoichMgr& ReactionStoichMgr::operator=(const ReactionStoichMgr& right) @@ -80,12 +65,6 @@ ReactionStoichMgr& ReactionStoichMgr::operator=(const ReactionStoichMgr& right) m_revproducts = new StoichManagerN(*right.m_revproducts); m_irrevproducts = new StoichManagerN(*right.m_irrevproducts); m_dummy = right.m_dummy; -#ifdef INCL_STOICH_WRITER - if (m_writer) { - delete(m_writer); - } - m_rwriter = new StoichManagerN(right.m_writer); -#endif } return *this; } @@ -128,18 +107,8 @@ add(size_t rxn, const ReactionData& r) // or specified reaction orders, then add it in a general reaction if (isfrac || r.global || rk.size() > 3) { m_reactants->add(rxn, r.reactants, r.rorder, r.rstoich); -#ifdef INCL_STOICH_WRITER - if (m_rwriter) { - m_rwriter->add(rxn, r.reactants, r.order, r.rstoich); - } -#endif } else { m_reactants->add(rxn, rk); -#ifdef INCL_STOICH_WRITER - if (m_rwriter) { - m_rwriter->add(rxn, rk); - } -#endif } std::vector pk; diff --git a/src/kinetics/StoichManager.h b/src/kinetics/StoichManager.h index 314ebd6cb..c54753959 100644 --- a/src/kinetics/StoichManager.h +++ b/src/kinetics/StoichManager.h @@ -925,9 +925,7 @@ public: _writeMultiply(m_cn_list.begin(), m_cn_list.end(), r, out); } - private: - std::vector m_c1_list; std::vector m_c2_list; std::vector m_c3_list; @@ -944,101 +942,6 @@ private: std::map m_loc; }; -#undef INCL_STOICH_WRITER -#ifdef INCL_STOICH_WRITER - -class StoichWriter -{ -public: - - StoichWriter() { - } - - StoichWriter(const StoichWriter& right) : - m_mult(right.m_mult), - m_ir(right.m_ir), - m_dr(right.m_dr), - m_is(right.m_is), - m_ds(right.m_ds) { - } - - StoichWriter& operator=(const StoichWriter& right) { - if (this != &right) { - m_mult = right.m_mult; - m_ir = right.m_ir; - m_dr = right.m_dr; - m_is = right.m_is; - m_ds = right.m_ds; - } - return *this; - } - - - void add(int rxn, const vector_int& k) { - int n, nn = k.size(); - for (n = 0; n < nn; n++) { - if (m_mult[rxn] != "") { - m_mult[rxn] += " * "; - } - m_mult[rxn] += "c[" + int2str(k[n]) + "]"; - m_is[k[n]] += " + rop[" + int2str(rxn) + "]"; - m_ds[k[n]] += " - rop[" + int2str(rxn) + "]"; - m_ir[rxn] += " + grt[" + int2str(k[n]) + "]"; - m_dr[rxn] += " - grt[" + int2str(k[n]) + "]"; - } - } - - void add(int rxn, const vector_int& k, const vector_fp& order, - const vector_fp& stoich) { - int n, nn = k.size(); - std::string s; - for (n = 0; n < nn; n++) { - if (order[n] == 1.0) { - m_mult[rxn] += "*c[" + int2str(k[n]) + "]"; - } else { - m_mult[rxn] += "*pow(c[" _ int2str(k[n]) + "],"+fp2str(order[n])+")"; - } - if (stoich[n] == 1.0) { - m_is[k[n]] += " + r[" + int2str(rxn) + "]"; - m_ds[k[n]] += " - r[" + int2str(rxn) + "]"; - m_ir[rxn] += " + g[" + int2str(k[n]) + "]"; - m_dr[rxn] += " - g[" + int2str(k[n]) + "]"; - } else { - s = fp2str(stoich[n]); - m_is[k[n]] += " + "+s+"*r[" + int2str(rxn) + "]"; - m_ds[k[n]] += " - "+s+"*r[" + int2str(rxn) + "]"; - m_ir[rxn] += " + "+s+"*g[" + int2str(k[n]) + "]"; - m_dr[rxn] += " - "+s+"*g[" + int2str(k[n]) + "]"; - } - } - } - - std::string mult(int rxn) { - return m_mult[rxn]; - } - std::string incrSpec(int k, std::string) { - return m_is[k]; - } - std::string decrSpec(int k) { - return m_ds[k]; - } - std::string incrRxn(int rxn) { - return m_ir[rxn]; - } - std::string decrRxn(int rxn) { - return m_dr[rxn]; - } - -private: - std::map m_mult; - std::map m_ir; - std::map m_dr; - std::map m_is; - std::map m_ds; -}; -#endif - - } #endif diff --git a/src/kinetics/solveSP.cpp b/src/kinetics/solveSP.cpp index d080ab552..e741d7fad 100644 --- a/src/kinetics/solveSP.cpp +++ b/src/kinetics/solveSP.cpp @@ -40,13 +40,6 @@ static doublereal calcWeightedNorm(const doublereal [], const doublereal dx[], s // extern FSUB_TYPE dgetrs_(char *, int *, int *, doublereal *, int *, int [], // doublereal [], int *, int *, unsigned int); // } -/***************************************************************************** - * PROTOTYPES and PREPROC DIRECTIVES FOR MISC. ROUTINES - *****************************************************************************/ - -#ifndef DAMPING -# define DAMPING true -#endif /*************************************************************************** * solveSP Class Definitinos @@ -222,31 +215,15 @@ int solveSP::solveSurfProb(int ifunc, doublereal time_scale, doublereal TKelvin, doublereal damp=1.0, tmp; // Weighted L2 norm of the residual. Currently, this is only // used for IO purposes. It doesn't control convergence. - // Therefore, it is turned off when DEBUG_SOLVESP isn't defined. doublereal resid_norm; doublereal inv_t = 0.0; doublereal t_real = 0.0, update_norm = 1.0E6; bool do_time = false, not_converged = true; -#ifdef DEBUG_SOLVESP -#ifdef DEBUG_SOLVESP_TIME - doublereal t1; -#endif -#else if (m_ioflag > 1) { m_ioflag = 1; } -#endif - -#ifdef DEBUG_SOLVESP -#ifdef DEBUG_SOLVESP_TIME - Cantera::clockWC wc; - if (m_ioflag) { - t1 = wc.secondsWC(); - } -#endif -#endif /* * Set the initial value of the do_time parameter @@ -291,7 +268,7 @@ int solveSP::solveSurfProb(int ifunc, doublereal time_scale, doublereal TKelvin, // m_kin->getNetProductionRates(DATA_PTR(m_netProductionRatesSave)); if (m_ioflag) { - print_header(m_ioflag, ifunc, time_scale, DAMPING, reltol, abstol, + print_header(m_ioflag, ifunc, time_scale, true, reltol, abstol, TKelvin, PGas, DATA_PTR(m_netProductionRatesSave), DATA_PTR(m_XMolKinSpecies)); } @@ -382,17 +359,6 @@ int solveSP::solveSurfProb(int ifunc, doublereal time_scale, doublereal TKelvin, resid_norm = calcWeightedNorm(DATA_PTR(m_wtResid), DATA_PTR(m_resid), m_neq); -#ifdef DEBUG_SOLVESP - if (m_ioflag > 1) { - printIterationHeader(m_ioflag, damp, inv_t, t_real, iter, do_time); - /* - * Print out the residual and jacobian - */ - printResJac(m_ioflag, m_neq, m_Jac, DATA_PTR(m_resid), - DATA_PTR(m_wtResid), resid_norm); - } -#endif - /* * Solve Linear system (with LAPACK). The solution is in resid[] */ @@ -434,11 +400,6 @@ int solveSP::solveSurfProb(int ifunc, doublereal time_scale, doublereal TKelvin, if (do_time) { t_real += time_scale; } -#ifdef DEBUG_SOLVESP - if (m_ioflag) { - printf("\nResidual is small, forcing convergence!\n"); - } -#endif } /* @@ -447,9 +408,7 @@ int solveSP::solveSurfProb(int ifunc, doublereal time_scale, doublereal TKelvin, * in any unknown. */ -#ifdef DAMPING damp = calc_damping(DATA_PTR(m_CSolnSP), DATA_PTR(m_resid), m_neq, &label_d); -#endif /* * Calculate the weighted norm of the update vector @@ -491,11 +450,6 @@ int solveSP::solveSurfProb(int ifunc, doublereal time_scale, doublereal TKelvin, (resid_norm < 1.0e-7 && update_norm*time_scale/t_real < EXTRA_ACCURACY)) { do_time = false; -#ifdef DEBUG_SOLVESP - if (m_ioflag > 1) { - printf("\t\tSwitching to steady solve.\n"); - } -#endif } } else { not_converged = ((update_norm > EXTRA_ACCURACY) || @@ -516,13 +470,6 @@ int solveSP::solveSurfProb(int ifunc, doublereal time_scale, doublereal TKelvin, printf("Continuing anyway\n"); } } -#ifdef DEBUG_SOLVESP -#ifdef DEBUG_SOLVESP_TIME - if (m_ioflag) { - printf("\nEnd of solve, time used: %e\n", wc.secondsWC()-t1); - } -#endif -#endif /* * Decide on what to return in the solution vector @@ -551,9 +498,6 @@ int solveSP::solveSurfProb(int ifunc, doublereal time_scale, doublereal TKelvin, return 1; } -#undef DAMPING - - /* * Update the surface states of the surface phases. */ @@ -1013,89 +957,10 @@ calc_t(doublereal netProdRateSolnSP[], doublereal XMolSolnSP[], *label_factor = 1.0; } inv_timeScale = inv_timeScale / *label_factor; -#ifdef DEBUG_SOLVESP - if (ioflag > 1) { - if (*label_factor > 1.0) { - printf("Delta_t increase due to repeated controlling species = %e\n", - *label_factor); - } - size_t kkin = m_kinSpecIndex[*label]; - InterfaceKinetics* m_kin = m_objects[ispSpecial]; - string sn = m_kin->kineticsSpeciesName(kkin); - printf("calc_t: spec=%d(%s) sf=%e pr=%e dt=%e\n", - *label, sn.c_str(), XMolSolnSP[*label], - netProdRateSolnSP[*label], 1.0/inv_timeScale); - } -#endif - return (inv_timeScale); } /* calc_t */ - -/** - * printResJac(): prints out the residual and Jacobian. - * - */ -#ifdef DEBUG_SOLVESP -void solveSP::printResJac(int ioflag, int neq, const Array2D& Jac, - doublereal resid[], doublereal wtRes[], - doublereal norm) -{ - int i, j, isp, nsp, irowKSI, irowISP; - int kstartKSI; - int kindexSP = 0; - string sname, pname, cname; - if (ioflag > 1) { - printf(" Printout of residual and jacobian\n"); - printf("\t Residual: weighted norm = %10.4e\n", norm); - printf("\t Index Species_Name Residual " - "Resid/wtRes wtRes\n"); - for (isp = 0; isp < m_numSurfPhases; isp++) { - nsp = m_nSpeciesSurfPhase[isp]; - InterfaceKinetics* m_kin = m_objects[isp]; - int surfPhaseIndex = m_kinObjPhaseIDSurfPhase[isp]; - m_kin->getNetProductionRates(DATA_PTR(m_numEqn1)); - kstartKSI = m_kin->kineticsSpeciesIndex(0, surfPhaseIndex); - SurfPhase* sp_ptr = m_ptrsSurfPhase[isp]; - pname = sp_ptr->id(); - - for (int k = 0; k < nsp; k++, kindexSP++) { - sname = sp_ptr->speciesName(k); - cname = pname + ":" + sname; - printf("\t %d: %-24s: %11.3e %11.3e %11.3e\n", kindexSP, - cname.c_str(), resid[kindexSP], - resid[kindexSP]/wtRes[kindexSP], wtRes[kindexSP]); - } - } - if (m_bulkFunc == BULK_DEPOSITION) { - for (isp = 0; isp < m_numBulkPhasesSS; isp++) { - // fill in - } - } - if (ioflag > 2) { - printf("\t Jacobian:\n"); - for (i = 0; i < m_neq; i++) { - irowISP = m_kinObjIndex[i]; - InterfaceKinetics* m_kin = m_objects[irowISP]; - irowKSI = m_kinSpecIndex[i]; - ThermoPhase& THref = m_kin->speciesPhase(irowKSI); - int phaseIndex = m_kin->speciesPhaseIndex(irowKSI); - kstartKSI = m_kin->kineticsSpeciesIndex(0, phaseIndex); - int klocal = i - m_eqnIndexStartSolnPhase[irowISP]; - sname = THref.speciesName(klocal); - printf("\t Row %d:%-16s:\n", i, sname.c_str()); - printf("\t "); - for (j = 0; j < m_neq; j++) { - printf("%10.4e ", Jac(i,j)); - } - printf("\n"); - } - } - } -} /* printResJac */ -#endif - /* * Optional printing at the start of the solveSP problem */ @@ -1143,61 +1008,6 @@ void solveSP::print_header(int ioflag, int ifunc, doublereal time_scale, printf(" Reltol = %9.3e, Abstol = %9.3e\n", reltol, abstol); } - /* - * Print out the initial guess - */ -#ifdef DEBUG_SOLVESP - if (ioflag > 1) { - printf("\n================================ INITIAL GUESS " - "========================================\n"); - int kindexSP = 0; - for (size_t isp = 0; isp < m_numSurfPhases; isp++) { - InterfaceKinetics* m_kin = m_objects[isp]; - int surfIndex = m_kin->surfacePhaseIndex(); - int nPhases = m_kin->nPhases(); - m_kin->getNetProductionRates(netProdRate); - updateMFKinSpecies(XMolKinSpecies, isp); - - printf("\n IntefaceKinetics Object # %d\n\n", isp); - - printf("\t Number of Phases = %d\n", nPhases); - printf("\t Temperature = %10.3e Kelvin\n", TKelvin); - printf("\t Pressure = %10.3g Pa\n\n", PGas); - printf("\t Phase:SpecName Prod_Rate MoleFraction kindexSP\n"); - printf("\t -------------------------------------------------------" - "----------\n"); - - int kspindex = 0; - bool inSurfacePhase = false; - for (int ip = 0; ip < nPhases; ip++) { - if (ip == surfIndex) { - inSurfacePhase = true; - } else { - inSurfacePhase = false; - } - ThermoPhase& THref = m_kin->thermo(ip); - int nsp = THref.nSpecies(); - string pname = THref.id(); - for (int k = 0; k < nsp; k++) { - string sname = THref.speciesName(k); - string cname = pname + ":" + sname; - if (inSurfacePhase) { - printf("\t %-24s %10.3e %10.3e %d\n", cname.c_str(), - netProdRate[kspindex], XMolKinSpecies[kspindex], - kindexSP); - kindexSP++; - } else { - printf("\t %-24s %10.3e %10.3e\n", cname.c_str(), - netProdRate[kspindex], XMolKinSpecies[kspindex]); - } - kspindex++; - } - } - printf("==========================================================" - "=================================\n"); - } - } -#endif if (ioflag == 1) { printf("\n\n\t Iter Time Del_t Damp DelX " " Resid Name-Time Name-Damp\n"); @@ -1252,53 +1062,6 @@ void solveSP::printIteration(int ioflag, doublereal damp, int label_d, } printf("\n"); } -#ifdef DEBUG_SOLVESP - else if (ioflag > 1) { - - updateMFSolnSP(XMolSolnSP); - printf("\n\t Weighted norm of update = %10.4e\n", update_norm); - - printf("\t Name Prod_Rate XMol Conc " - " Conc_Old wtConc"); - if (damp < 1.0) { - printf(" UnDamped_Conc"); - } - printf("\n"); - printf("\t---------------------------------------------------------" - "-----------------------------\n"); - int kindexSP = 0; - for (size_t isp = 0; isp < m_numSurfPhases; isp++) { - int nsp = m_nSpeciesSurfPhase[isp]; - InterfaceKinetics* m_kin = m_objects[isp]; - //int surfPhaseIndex = m_kinObjPhaseIDSurfPhase[isp]; - m_kin->getNetProductionRates(DATA_PTR(m_numEqn1)); - for (int k = 0; k < nsp; k++, kindexSP++) { - int kspIndex = m_kinSpecIndex[kindexSP]; - nm = m_kin->kineticsSpeciesName(kspIndex); - printf("\t%-16s %10.3e %10.3e %10.3e %10.3e %10.3e ", - nm.c_str(), - m_numEqn1[kspIndex], - XMolSolnSP[kindexSP], - CSolnSP[kindexSP], CSolnSP[kindexSP]+damp*resid[kindexSP], - wtSpecies[kindexSP]); - if (damp < 1.0) { - printf("%10.4e ", CSolnSP[kindexSP]+(damp-1.0)*resid[kindexSP]); - if (label_d == kindexSP) { - printf(" Damp "); - } - } - if (label_t == kindexSP) { - printf(" Tctrl"); - } - printf("\n"); - } - - } - - printf("\t--------------------------------------------------------" - "------------------------------\n"); - } -#endif } /* printIteration */ @@ -1349,130 +1112,6 @@ void solveSP::printFinal(int ioflag, doublereal damp, int label_d, int label_t, } printf(" -- success\n"); } -#ifdef DEBUG_SOLVESP - else if (ioflag > 1) { - - - printf("\n================================== FINAL RESULT =========" - "==================================================\n"); - updateMFSolnSP(XMolSolnSP); - printf("\n Weighted norm of solution update = %10.4e\n", update_norm); - printf(" Weighted norm of residual update = %10.4e\n\n", resid_norm); - - printf(" Name Prod_Rate XMol Conc " - " wtConc Resid Resid/wtResid wtResid"); - if (damp < 1.0) { - printf(" UnDamped_Conc"); - } - printf("\n"); - printf("---------------------------------------------------------------" - "---------------------------------------------\n"); - int kindexSP = 0; - for (size_t isp = 0; isp < m_numSurfPhases; isp++) { - int nsp = m_nSpeciesSurfPhase[isp]; - InterfaceKinetics* m_kin = m_objects[isp]; - //int surfPhaseIndex = m_kinObjPhaseIDSurfPhase[isp]; - m_kin->getNetProductionRates(DATA_PTR(m_numEqn1)); - for (int k = 0; k < nsp; k++, kindexSP++) { - int kspIndex = m_kinSpecIndex[kindexSP]; - nm = m_kin->kineticsSpeciesName(kspIndex); - printf("%-16s %10.3e %10.3e %10.3e %10.3e %10.3e %10.3e %10.3e", - nm.c_str(), - m_numEqn1[kspIndex], - XMolSolnSP[kindexSP], - CSolnSP[kindexSP], - wtSpecies[kindexSP], - resid[kindexSP], - resid[kindexSP]/wtRes[kindexSP], wtRes[kindexSP]); - if (damp < 1.0) { - printf("%10.4e ", CSolnSP[kindexSP]+(damp-1.0)*resid[kindexSP]); - if (label_d == kindexSP) { - printf(" Damp "); - } - } - if (label_t == kindexSP) { - printf(" Tctrl"); - } - printf("\n"); - } - - } - printf("---------------------------------------------------------------" - "---------------------------------------------\n"); - doublereal* XMolKinSpecies = DATA_PTR(m_numEqn2); - kindexSP = 0; - for (size_t isp = 0; isp < m_numSurfPhases; isp++) { - InterfaceKinetics* m_kin = m_objects[isp]; - int surfIndex = m_kin->surfacePhaseIndex(); - int nPhases = m_kin->nPhases(); - m_kin->getNetProductionRates(netProdRateKinSpecies); - - updateMFKinSpecies(XMolKinSpecies, isp); - - printf("\n IntefaceKinetics Object # %d\n\n", isp); - - printf("\t Number of Phases = %d\n", nPhases); - printf("\t Temperature = %10.3e Kelvin\n", TKelvin); - printf("\t Pressure = %10.3g Pa\n\n", PGas); - printf("\t Phase:SpecName Prod_Rate MoleFraction kindexSP\n"); - printf("\t--------------------------------------------------------------" - "---\n"); - - int kspindex = 0; - bool inSurfacePhase = false; - for (int ip = 0; ip < nPhases; ip++) { - if (ip == surfIndex) { - inSurfacePhase = true; - } else { - inSurfacePhase = false; - } - ThermoPhase& THref = m_kin->thermo(ip); - int nsp = THref.nSpecies(); - string pname = THref.id(); - for (k = 0; k < nsp; k++) { - string sname = THref.speciesName(k); - string cname = pname + ":" + sname; - if (inSurfacePhase) { - printf("\t%-24s %10.3e %10.3e %d\n", cname.c_str(), - netProdRateKinSpecies[kspindex], XMolKinSpecies[kspindex], kindexSP); - kindexSP++; - } else { - printf("\t%-24s %10.3e %10.3e\n", cname.c_str(), - netProdRateKinSpecies[kspindex], XMolKinSpecies[kspindex]); - } - kspindex++; - } - } - } - printf("\n"); - printf("===============================================================" - "============================================\n\n"); - } -#endif } -#ifdef DEBUG_SOLVESP -void solveSP:: -printIterationHeader(int ioflag, doublereal damp,doublereal inv_t, doublereal t_real, - int iter, bool do_time) -{ - if (ioflag > 1) { - printf("\n===============================Iteration %5d " - "=================================\n", iter); - if (do_time) { - printf(" Transient step with: Real Time_n-1 = %10.4e sec,", t_real); - printf(" Time_n = %10.4e sec\n", t_real + 1.0/inv_t); - printf(" Delta t = %10.4e sec", 1.0/inv_t); - } else { - printf(" Steady Solve "); - } - if (damp < 1.0) { - printf(", Damping value = %10.4e\n", damp); - } else { - printf("\n"); - } - } -} -#endif - } diff --git a/src/kinetics/solveSP.h b/src/kinetics/solveSP.h index 86a466b87..fa46bf578 100644 --- a/src/kinetics/solveSP.h +++ b/src/kinetics/solveSP.h @@ -260,12 +260,6 @@ private: doublereal TKelvin, doublereal PGas, doublereal netProdRate[], doublereal XMolKinSpecies[]); -#ifdef DEBUG_SOLVESP - - void printResJac(int ioflag, int neq, const Array2D& Jac, - doublereal resid[], doublereal wtResid[], doublereal norm); -#endif - //! Printing routine that gets called after every iteration void printIteration(int ioflag, doublereal damp, int label_d, int label_t, doublereal inv_t, doublereal t_real, size_t iter, @@ -338,22 +332,6 @@ private: const Array2D& Jac, const doublereal CSolnSP[], const doublereal abstol, const doublereal reltol); -#ifdef DEBUG_SOLVESP - //! Utility routine to print a header for high lvls of debugging - /*! - * @param ioflag Lvl of debugging - * @param damp lvl of damping - * @param inv_t Inverse of the value of delta T - * @param t_real Value of the time - * @param iter Interation number - * @param do_time boolean indicating whether time stepping is taking - * place - */ - void printIterationHeader(int ioflag, doublereal damp, - doublereal inv_t, doublereal t_real, int iter, - bool do_time); -#endif - /** * Update the surface states of the surface phases. */ diff --git a/src/numerics/DenseMatrix.cpp b/src/numerics/DenseMatrix.cpp index ca28cdd8c..a129e0bcb 100644 --- a/src/numerics/DenseMatrix.cpp +++ b/src/numerics/DenseMatrix.cpp @@ -236,36 +236,6 @@ int solve(DenseMatrix& A, DenseMatrix& b) return info; } //==================================================================================================================== - -#ifdef INCL_LEAST_SQUARES -/** @todo fix lwork */ -int leastSquares(DenseMatrix& A, double* b) -{ - int info = 0; - int rank = 0; - double rcond = -1.0; - // fix this! - int lwork = 6000; // 2*(3*min(m,n) + max(2*min(m,n), max(m,n))); - vector_fp work(lwork); - vector_fp s(min(static_cast(A.nRows()), - static_cast(A.nColumns()))); - ct_dgelss(static_cast(A.nRows()), - static_cast(A.nColumns()), 1, A.ptrColumn(0), - static_cast(A.nRows()), b, - static_cast(A.nColumns()), &s[0], //.begin(), - rcond, rank, &work[0], work.size(), info); - if (info != 0) { - if (A.m_printLevel) { - writelogf("leastSquares(): DGELSS returned INFO = %d\n", info); - } - if (! A.m_useReturnErrorCode) { - throw CELapackError("leastSquares()", "DGELSS returned INFO = " + int2str(info)); - } - } - return info; -} -#endif -//==================================================================================================================== void multiply(const DenseMatrix& A, const double* const b, double* const prod) { ct_dgemv(ctlapack::ColMajor, ctlapack::NoTranspose, diff --git a/src/numerics/Func1.cpp b/src/numerics/Func1.cpp index 9ba3ecd99..fe2cfe923 100644 --- a/src/numerics/Func1.cpp +++ b/src/numerics/Func1.cpp @@ -161,11 +161,6 @@ Func1& Sin1::derivative() const { Func1* c = new Cos1(m_c); Func1* r = &newTimesConstFunction(*c, m_c); -#ifdef DEBUG_FUNC - cout << "Sin1::derivative: \n"; - cout << "function = \'" + write("x") + "\'\n"; - cout << "derivative = \'" + r->write("x") + "\'\n"; -#endif return *r; } /*****************************************************************************/ @@ -174,11 +169,6 @@ Func1& Cos1::derivative() const { Func1* s = new Sin1(m_c); Func1* r = &newTimesConstFunction(*s, -m_c); -#ifdef DEBUG_FUNC - cout << "Cos1::derivative: \n"; - cout << "function = \'" + write("x") + "\'\n"; - cout << "derivative = \'" + r->write("x") + "\'\n"; -#endif return *r; } @@ -225,11 +215,6 @@ Func1& Pow1::derivative() const Func1* f = new Pow1(m_c - 1.0); r = &newTimesConstFunction(*f, m_c); } -#ifdef DEBUG_FUNC - cout << "Pow1::derivative: \n"; - cout << "function = \'" + write("x") + "\'\n"; - cout << "derivative = \'" + r->write("x") + "\'\n"; -#endif return *r; } @@ -567,11 +552,6 @@ Func1& newRatioFunction(Func1& f1, Func1& f2) Func1& newCompositeFunction(Func1& f1, Func1& f2) { - //#ifdef DEBUG_FUNC - //cout << "creating new composite function." << endl; - //cout << "f1 = " << f1.write("x") << " " << f1.ID() << endl; - //cout << "f2 = " << f2.write("x") << " " << f2.ID() << endl; - //#endif if (isZero(f1)) { delete &f1; delete &f2; diff --git a/src/numerics/lapack.h b/src/numerics/lapack.h index 746c9d20b..ba0726359 100644 --- a/src/numerics/lapack.h +++ b/src/numerics/lapack.h @@ -10,7 +10,6 @@ #if defined(NEEDS_F77_TRANSLATION) #if defined(F77EXTERNS_UPPERCASE_NOTRAILINGBAR) -#define dgelss_ DGELSS #define dgetrs_ DGETRS #define dgetrf_ DGETRF #define dgetri_ DGETRI @@ -19,7 +18,6 @@ #define simplx_ SIMPLX #define splin2_ SPLIN2 #define splie2_ SPLIE2 -#define dgelss_ DGELSS #endif #endif diff --git a/src/python/methods.h b/src/python/methods.h index 631565d8a..9ed8f190e 100644 --- a/src/python/methods.h +++ b/src/python/methods.h @@ -297,10 +297,6 @@ static PyMethodDef ct_methods[] = { {"ct_appdelete", pyct_appdelete, METH_VARARGS}, -#ifdef INCL_USER_PYTHON -#include "usermethods.h" -#endif - {NULL, NULL} /* sentinel */ }; diff --git a/src/python/pycantera.cpp b/src/python/pycantera.cpp index 7707b7867..1213bf15e 100644 --- a/src/python/pycantera.cpp +++ b/src/python/pycantera.cpp @@ -54,11 +54,6 @@ static PyObject* ErrorObject; #include "ctonedim_methods.cpp" #include "ctmultiphase_methods.cpp" -#ifdef INCL_USER_PYTHON -#include "ctuser.h" -#include "ctuser_methods.cpp" -#endif - static PyObject* pyct_appdelete(PyObject* self, PyObject* args) { diff --git a/src/python/pyutils.h b/src/python/pyutils.h index 269375558..f0f050edd 100644 --- a/src/python/pyutils.h +++ b/src/python/pyutils.h @@ -1,5 +1,5 @@ -#ifndef CTPY_UTILS -#define CTPY_UTILS +#ifndef CT_PYUTILS_H +#define CT_PYUTILS_H #include "Python.h" diff --git a/src/spectra/spectralUtilities.h b/src/spectra/spectralUtilities.h index 65f984368..3b4bb4130 100644 --- a/src/spectra/spectralUtilities.h +++ b/src/spectra/spectralUtilities.h @@ -1,5 +1,5 @@ -#ifndef CT_SPEC_UTILS -#define CT_SPEC_UTILS +#ifndef CT_SPEC_UTILS_H +#define CT_SPEC_UTILS_H #include "Nuclei.h" diff --git a/src/thermo/MineralEQ3.cpp b/src/thermo/MineralEQ3.cpp index c2f3e3417..5746d7e9f 100644 --- a/src/thermo/MineralEQ3.cpp +++ b/src/thermo/MineralEQ3.cpp @@ -554,18 +554,6 @@ void MineralEQ3::setParametersFromXML(const XML_Node& eosdata) doublereal MineralEQ3::LookupGe(const std::string& elemName) { -#ifdef OLDWAY - int num = sizeof(geDataTable) / sizeof(struct GeData); - string s3 = elemName.substr(0,3); - for (int i = 0; i < num; i++) { - //if (!std::strncmp(elemName.c_str(), aWTable[i].name, 3)) { - if (s3 == geDataTable[i].name) { - return (geDataTable[i].GeValue); - } - } - throw CanteraError("LookupGe", "element " + s + " not found"); - return -1.0; -#else size_t iE = elementIndex(elemName); if (iE == npos) { throw CanteraError("PDSS_HKFT::LookupGe", "element " + elemName + " not found"); @@ -577,7 +565,6 @@ doublereal MineralEQ3::LookupGe(const std::string& elemName) } geValue *= (-298.15); return geValue; -#endif } void MineralEQ3::convertDGFormation() diff --git a/src/thermo/MixtureFugacityTP.cpp b/src/thermo/MixtureFugacityTP.cpp index 16036f975..3ec7f13bd 100644 --- a/src/thermo/MixtureFugacityTP.cpp +++ b/src/thermo/MixtureFugacityTP.cpp @@ -1161,21 +1161,6 @@ doublereal MixtureFugacityTP::calculatePsat(doublereal TKelvin, doublereal& mola int stab; doublereal presLast = pres; - -#ifdef DDDD - double pVec[100]; - int n = 0; - for (int i = 0; i < 50; i++) { - pVec[n++] = 3.40E6 + 0.01E5 * i; - } - - for (int i = 0; i < 50; i++) { - stab = corr0(TKelvin, pVec[i], RhoLiquid, RhoGas, liqGRT, gasGRT); - printf("p = %g, T = %g, stab = %d, Rl = %g Rg = %g, Gl = %g, Gg = %g\n", - pVec[i], TKelvin, stab, RhoLiquid, RhoGas,liqGRT, gasGRT); - } -#endif - /* * First part of the calculation involves finding a pressure at which the * gas and the liquid state coexists. diff --git a/src/thermo/PDSS_HKFT.cpp b/src/thermo/PDSS_HKFT.cpp index 5ac555f6b..01646c5c8 100644 --- a/src/thermo/PDSS_HKFT.cpp +++ b/src/thermo/PDSS_HKFT.cpp @@ -1194,49 +1194,6 @@ doublereal PDSS_HKFT::gstar(const doublereal temp, const doublereal pres, const return res; } - -#ifdef OLDWAY - -/* awData structure */ -/*! - * Database for atomic molecular weights - * - * Values are taken from the 1989 Standard Atomic Weights, CRC - * - * awTable[] is a static function with scope limited to this file. - * It can only be referenced via the static Elements class function, - * LookupWtElements(). - * - * units = kg / kg-mol (or equivalently gm / gm-mol) - * - * (note: this structure was picked because it's simple, compact, - * and extensible). - * - */ -struct GeData { - char name[4]; ///< Null Terminated name, First letter capitalized - doublereal GeValue; /// < Gibbs free energies of elements J kmol-1 -}; - -//! Values of G_elements(T=298.15,1atm) -/*! - * all units are Joules kmol-1 - */ - -static struct GeData geDataTable[] = { - {"H", -19.48112E6}, // NIST Webbook - Cox, Wagman 1984 - {"Na", -15.29509E6}, // NIST Webbook - Cox, Wagman 1984 - {"O", -30.58303E6}, // NIST Webbook - Cox, Wagman 1984 - {"Cl", -33.25580E6}, // NIST Webbook - Cox, Wagman 1984 - {"Si", -5.61118E6}, // Janaf - {"C", -1.71138E6}, // barin, Knack, NBS Bulletin 1971 - {"S", -9.55690E6}, // Yellow - webbook - {"Al", -8.42870E6}, // Webbook polynomial - {"K", -19.26943E6}, // Webbook - {"Fe", -8.142476E6}, // Nist Webbook - Cox, Wagman 1984 - {"E", 0.0} // Don't overcount -}; -#endif //! Static function to look up Element Free Energies /*! * @@ -1253,18 +1210,6 @@ static struct GeData geDataTable[] = { */ doublereal PDSS_HKFT::LookupGe(const std::string& elemName) { -#ifdef OLDWAY - int num = sizeof(geDataTable) / sizeof(struct GeData); - string s3 = elemName.substr(0,3); - for (int i = 0; i < num; i++) { - //if (!std::strncmp(elemName.c_str(), aWTable[i].name, 3)) { - if (s3 == geDataTable[i].name) { - return (geDataTable[i].GeValue); - } - } - throw CanteraError("LookupGe", "element " + s + " not found"); - return -1.0; -#else size_t iE = m_tp->elementIndex(elemName); if (iE == npos) { throw CanteraError("PDSS_HKFT::LookupGe", "element " + elemName + " not found"); @@ -1276,7 +1221,6 @@ doublereal PDSS_HKFT::LookupGe(const std::string& elemName) } geValue *= (-298.15); return geValue; -#endif } void PDSS_HKFT::convertDGFormation() diff --git a/src/thermo/SpeciesThermoFactory.cpp b/src/thermo/SpeciesThermoFactory.cpp index 04fa24b44..2569d8d59 100644 --- a/src/thermo/SpeciesThermoFactory.cpp +++ b/src/thermo/SpeciesThermoFactory.cpp @@ -383,81 +383,6 @@ static void installNasaThermoFromXML(std::string speciesName, sp.install(speciesName, k, NASA, &c[0], tmin, tmax, p0); } -#ifdef INCL_NASA96 - -//! Install a NASA96 polynomial thermodynamic property parameterization for species k into a SpeciesThermo instance. -/*! - * This is called by method installThermoForSpecies if a NASA block is found in the XML input. - * - * @param speciesName String name of the species - * @param sp SpeciesThermo object that will receive the nasa polynomial object - * @param k Species index within the phase - * @param f0ptr Ptr to the first XML_Node for the first NASA polynomial - * @param f1ptr Ptr to the first XML_Node for the first NASA polynomial - */ -static void installNasa96ThermoFromXML(std::string speciesName, SpeciesThermo& sp, int k, - const XML_Node* f0ptr, const XML_Node* f1ptr) -{ - doublereal tmin0, tmax0, tmin1, tmax1, tmin, tmid, tmax; - - const XML_Node& f0 = *f0ptr; - bool dualRange = false; - if (f1ptr) { - dualRange = true; - } - tmin0 = fpValue(f0["Tmin"]); - tmax0 = fpValue(f0["Tmax"]); - tmin1 = tmax0; - tmax1 = tmin1 + 0.0001; - if (dualRange) { - tmin1 = fpValue((*f1ptr)["Tmin"]); - tmax1 = fpValue((*f1ptr)["Tmax"]); - } - - - doublereal p0 = OneAtm; - if (f0.hasAttrib("P0")) { - p0 = fpValue(f0["P0"]); - } - if (f0.hasAttrib("Pref")) { - p0 = fpValue(f0["Pref"]); - } - - vector_fp c0, c1; - if (fabs(tmax0 - tmin1) < 0.01) { - tmin = tmin0; - tmid = tmax0; - tmax = tmax1; - getFloatArray(f0.child("floatArray"), c0, false); - if (dualRange) { - getFloatArray(f1ptr->child("floatArray"), c1, false); - } else { - c1.resize(7,0.0); - copy(c0.begin(), c0.end(), c1.begin()); - } - } else if (fabs(tmax1 - tmin0) < 0.01) { - tmin = tmin1; - tmid = tmax1; - tmax = tmax0; - getFloatArray(f1ptr->child("floatArray"), c0, false); - getFloatArray(f0.child("floatArray"), c1, false); - } else { - throw CanteraError("installNasaThermo", - "non-continuous temperature ranges."); - } - vector_fp c(15); - c[0] = tmid; - c[1] = c0[5]; - c[2] = c0[6]; - copy(c0.begin(), c0.begin()+5, c.begin() + 3); - c[8] = c1[5]; - c[9] = c1[6]; - copy(c1.begin(), c1.begin()+5, c.begin() + 10); - sp.install(speciesName, k, NASA, &c[0], tmin, tmax, p0); -} - -#endif - //! Look up the elemental reference state entropies /*! * @param elemName String name of the element @@ -465,18 +390,6 @@ static void installNasa96ThermoFromXML(std::string speciesName, SpeciesThermo& s */ static doublereal LookupGe(const std::string& elemName, ThermoPhase* th_ptr) { -#ifdef OLDWAY - int num = sizeof(geDataTable) / sizeof(struct GeData); - string s3 = elemName.substr(0,3); - for (int i = 0; i < num; i++) { - //if (!std::strncmp(elemName.c_str(), aWTable[i].name, 3)) { - if (s3 == geDataTable[i].name) { - return (geDataTable[i].GeValue); - } - } - throw CanteraError("LookupGe", "element " + s + " not found"); - return -1.0; -#else size_t iE = th_ptr->elementIndex(elemName); if (iE == npos) { throw CanteraError("PDSS_HKFT::LookupGe", "element " + elemName + " not found"); @@ -488,7 +401,6 @@ static doublereal LookupGe(const std::string& elemName, ThermoPhase* th_ptr) } geValue *= (-298.15); return geValue; -#endif } //! Convert delta G formulation