diff --git a/Cantera/cxx/include/Cantera.mak.in b/Cantera/cxx/include/Cantera.mak.in index 67b6b1885..ffbc46197 100644 --- a/Cantera/cxx/include/Cantera.mak.in +++ b/Cantera/cxx/include/Cantera.mak.in @@ -23,6 +23,8 @@ # in_CanteraBuildTree = 0 +CANTERA_VERSION=@ctversion@ + ############################################################################### # CANTERA CORE ############################################################################### @@ -122,6 +124,10 @@ endif # CANTERA_TOTAL_INCLUDES= $(CANTERA_CORE_INCLUDES) $(CANTERA_BOOST_INCLUDES) $(CANTERA_CVODE_INCLUDE) # +# You can add this into the compilation environment to identify the version number +# +CANTERA_DEFINES = -DCANTERA_VERSION=@ctversion@ +# # LIBS and LIBS should be the same ... # CANTERA_TOTAL_LIBS2 = -L$(CANTERA_LIBSDIR) @LOCAL_LIBS@ diff --git a/Cantera/cxx/include/Cantera_bt.mak.in b/Cantera/cxx/include/Cantera_bt.mak.in index 54c124c6e..88cdc0633 100644 --- a/Cantera/cxx/include/Cantera_bt.mak.in +++ b/Cantera/cxx/include/Cantera_bt.mak.in @@ -23,6 +23,8 @@ # in_CanteraBuildTree = 0 +CANTERA_VERSION=@ctversion@ + ############################################################################### # CANTERA CORE ############################################################################### @@ -117,6 +119,10 @@ endif #################################################################### CANTERA_TOTAL_INCLUDES= $(CANTERA_CORE_INCLUDES) $(CANTERA_BOOST_INCLUDES) $(CANTERA_CVODE_INCLUDE) +# +# You can add this into the compilation environment to identify the version number +# +CANTERA_DEFINES = -DCANTERA_VERSION=@ctversion@ CANTERA_TOTAL_LIBS2 = @LOCAL_LIB_DIRS@ @LOCAL_LIBS@ diff --git a/Cantera/src/kinetics/InterfaceKinetics.cpp b/Cantera/src/kinetics/InterfaceKinetics.cpp index db03c7cdc..2f5cd0cd7 100644 --- a/Cantera/src/kinetics/InterfaceKinetics.cpp +++ b/Cantera/src/kinetics/InterfaceKinetics.cpp @@ -52,6 +52,10 @@ namespace Cantera { m_has_coverage_dependence(false), m_has_electrochem_rxns(false), m_has_exchange_current_density_formulation(false), + m_phaseExistsCheck(false), + m_phaseExists(0), + m_rxnPhaseIsReactant(0), + m_rxnPhaseIsProduct(0), m_ioFlag(0) { if (thermo != 0) addPhase(*thermo); @@ -67,6 +71,10 @@ namespace Cantera { if (m_integrator) { delete m_integrator; } + for (int i = 0; i < m_ii; i++) { + delete m_rxnPhaseIsReactant[i]; + delete m_rxnPhaseIsProduct[i]; + } } //==================================================================================================================== // Copy Constructor for the %InterfaceKinetics object. @@ -92,6 +100,10 @@ namespace Cantera { m_has_coverage_dependence(false), m_has_electrochem_rxns(false), m_has_exchange_current_density_formulation(false), + m_phaseExistsCheck(false), + m_phaseExists(0), + m_rxnPhaseIsReactant(0), + m_rxnPhaseIsProduct(0), m_ioFlag(0) { m_kdata = new InterfaceKineticsData; @@ -111,13 +123,21 @@ namespace Cantera { */ InterfaceKinetics& InterfaceKinetics:: operator=(const InterfaceKinetics &right) { + int i; /* * Check for self assignment. */ if (this == &right) return *this; + for (i = 0; i < m_ii; i++) { + delete (m_rxnPhaseIsReactant[i]); + delete (m_rxnPhaseIsProduct[i]); + } + Kinetics::operator=(right); + + m_kk = right.m_kk; m_revindex = right.m_revindex; m_rates = right.m_rates; @@ -138,7 +158,6 @@ namespace Cantera { m_E = right.m_E; m_surf = right.m_surf; //DANGER - shallow copy m_integrator = right.m_integrator; //DANGER - shallow copy - m_beta = right.m_beta; m_ctrxn = right.m_ctrxn; m_ctrxn_ecdf = right.m_ctrxn_ecdf; @@ -149,6 +168,24 @@ namespace Cantera { m_has_coverage_dependence = right.m_has_coverage_dependence; m_has_electrochem_rxns = right.m_has_electrochem_rxns; m_has_exchange_current_density_formulation = right.m_has_exchange_current_density_formulation; + m_phaseExistsCheck = right.m_phaseExistsCheck; + m_phaseExists = right.m_phaseExists; + + + m_rxnPhaseIsReactant.resize(m_ii, 0); + m_rxnPhaseIsProduct.resize(m_ii, 0); + int np = nPhases(); + for (i = 0; i < m_ii; i++) { + m_rxnPhaseIsReactant[i] = new bool[np]; + m_rxnPhaseIsProduct[i] = new bool[np]; + for (int p = 0; p < np; p++) { + m_rxnPhaseIsReactant[i][p] = right.m_rxnPhaseIsReactant[i][p]; + m_rxnPhaseIsProduct[i][p] = right.m_rxnPhaseIsProduct[i][p]; + } + } + + m_rxnPhaseIsProduct = right.m_rxnPhaseIsProduct; + m_ioFlag = right.m_ioFlag; return *this; @@ -639,7 +676,57 @@ namespace Cantera { for (int j = 0; j != m_ii; ++j) { ropnet[j] = ropf[j] - ropr[j]; - } + } + + + /* + * For reactions involving multiple phases, we must check that the phase + * being consumed actually exists. This is particularly important for + * phases that are stoichiometric phases containing one species with a unity activity + */ + if (m_phaseExistsCheck) { + for (int j = 0; j != m_ii; ++j) { + if ((ropr[j] > ropf[j]) && (ropr[j] > 0.0)) { + for (int p = 0; p < nPhases(); p++) { + if (m_rxnPhaseIsProduct[j][p]) { + if (! m_phaseExists[p]) { + ropnet[j] = 0.0; + ropr[j] = ropf[j]; + if (ropf[j] > 0.0) { + for (int rp = 0; rp < nPhases(); rp++) { + if (m_rxnPhaseIsReactant[j][rp]) { + if (! m_phaseExists[rp]) { + ropnet[j] = 0.0; + ropr[j] = ropf[j] = 0.0;; + } + } + } + } + } + } + } + } else if ((ropf[j] > ropr[j]) && (ropf[j] > 0.0)) { + for (int p = 0; p < nPhases(); p++) { + if (m_rxnPhaseIsReactant[j][p]) { + if (! m_phaseExists[p]) { + ropnet[j] = 0.0; + ropf[j] = ropr[j]; + if (ropf[j] > 0.0) { + for (int rp = 0; rp < nPhases(); rp++) { + if (m_rxnPhaseIsProduct[j][rp]) { + if (! m_phaseExists[rp]) { + ropnet[j] = 0.0; + ropf[j] = ropr[j] = 0.0; + } + } + } + } + } + } + } + } + } + } m_kdata->m_ROP_ok = true; } @@ -851,11 +938,37 @@ namespace Cantera { addElementaryReaction(r); // operations common to all reaction types - installReagents( r ); + installReagents(r); //installGroups(reactionNumber(), r.rgroups, r.pgroups); incrementRxnCount(); m_rxneqn.push_back(r.equation); - } + + m_rxnPhaseIsReactant.resize(m_ii, 0); + m_rxnPhaseIsProduct.resize(m_ii, 0); + + int np = nPhases(); + int i = m_ii -1; + m_rxnPhaseIsReactant[i] = new bool[np]; + m_rxnPhaseIsProduct[i] = new bool[np]; + + for (int p = 0; p < np; p++) { + m_rxnPhaseIsReactant[i][p] = false; + m_rxnPhaseIsProduct[i][p] = false; + } + + const vector_int& vr = reactants(i); + for (int ik = 0; ik < (int) vr.size(); ik++) { + int k = vr[ik]; + int p = speciesPhaseIndex(k); + m_rxnPhaseIsReactant[i][p] = true; + } + const vector_int& vp = products(i); + for (int ik = 0; ik < (int) vp.size(); ik++) { + int k = vp[ik]; + int p = speciesPhaseIndex(k); + m_rxnPhaseIsProduct[i][p] = true; + } + } //==================================================================================================================== void InterfaceKinetics::addElementaryReaction(const ReactionData& r) { int iloc; @@ -1008,7 +1121,7 @@ namespace Cantera { * calculates rates of species production from reaction rates of * progress. */ - m_rxnstoich.add( reactionNumber(), r); + m_rxnstoich.add(reactionNumber(), r); /* * register reaction in lists of reversible and irreversible rxns. */ @@ -1016,11 +1129,16 @@ namespace Cantera { m_revindex.push_back(reactionNumber()); m_nrev++; } else { - m_irrev.push_back( reactionNumber() ); + m_irrev.push_back(reactionNumber()); m_nirrev++; } } - + //=============================================================================================== + void InterfaceKinetics::addPhase(thermo_t &thermo) { + Kinetics::addPhase(thermo); + m_phaseExists.push_back(true); + } + //================================================================================================ /** * Prepare the class for the addition of reactions. This function * must be called after instantiation of the class, but before @@ -1044,7 +1162,7 @@ namespace Cantera { m_pot.resize(m_kk, 0.0); m_phi.resize(np, 0.0); } - + //================================================================================================ /** * Finish adding reactions and prepare for use. This function * must be called after all reactions are entered into the mechanism @@ -1071,10 +1189,13 @@ namespace Cantera { m_deltaG0.resize(m_ii, 0.0); m_ProdStanConcReac.resize(m_ii, 0.0); + if (m_thermo.size() != m_phaseExists.size()) { + throw CanteraError("InterfaceKinetics::finalize", "internal error"); + } + m_finalized = true; } - doublereal InterfaceKinetics::electrochem_beta(int irxn) const{ int n = m_ctrxn.size(); for (int i = 0; i < n; i++) { @@ -1089,7 +1210,7 @@ namespace Cantera { bool InterfaceKinetics::ready() const { return (m_finalized); } - + //================================================================================================ // Advance the surface coverages in time /* * @param tstep Time value to advance the surface coverages @@ -1106,7 +1227,7 @@ namespace Cantera { delete m_integrator; m_integrator = 0; } - + //================================================================================================ // Solve for the pseudo steady-state of the surface problem /* * Solve for the steady state of the surface problem. @@ -1135,7 +1256,25 @@ namespace Cantera { */ m_integrator->solvePseudoSteadyStateProblem(ifuncOverride, timeScaleOverride); } + //================================================================================================ + void InterfaceKinetics::setPhaseExistence(const int iphase, const bool exists) { + if (iphase < 0 || iphase >= (int) m_thermo.size()) { + throw CanteraError("InterfaceKinetics:setPhaseExistence", "out of bounds"); + } + if (exists) { + if (!m_phaseExists[iphase]) { + m_phaseExistsCheck--; + m_phaseExists[iphase] = true; + } + } else { + if (m_phaseExists[iphase]) { + m_phaseExistsCheck++; + m_phaseExists[iphase] = false; + } + } + } + //================================================================================================ void EdgeKinetics::finalize() { m_rwork.resize(nReactions()); int ks = reactionPhaseIndex(); @@ -1147,7 +1286,8 @@ namespace Cantera { "expected interface dimension = 1, but got dimension = " +int2str(m_surf->nDim())); m_finalized = true; - } + } + //================================================================================================ } diff --git a/Cantera/src/kinetics/InterfaceKinetics.h b/Cantera/src/kinetics/InterfaceKinetics.h index b569afc7f..05e4bc73f 100644 --- a/Cantera/src/kinetics/InterfaceKinetics.h +++ b/Cantera/src/kinetics/InterfaceKinetics.h @@ -374,6 +374,19 @@ namespace Cantera { */ //@{ + //! Add a phase to the kinetics manager object. + /*! + * This must be done before the function init() is called or + * before any reactions are input. + * + * This function calls the Kinetics operator addPhase. + * It also sets the following functions + * + * m_phaseExists[] + * + * @param thermo Reference to the ThermoPhase to be added. + */ + virtual void addPhase(thermo_t& thermo); //! Prepare the class for the addition of reactions. /*! @@ -403,7 +416,10 @@ namespace Cantera { virtual bool ready() const; - + //! Internal routine that updates the Rates of Progress of the reactions + /*! + * This is actually the guts of the functionality of the object + */ void updateROP(); @@ -507,6 +523,16 @@ namespace Cantera { */ void applyExchangeCurrentDensityFormulation(doublereal* const kfwd); + //! Set the existence of a phase in the reaction object + /*! + * Tell the kinetics object whether a phase in the object exists. + * This is actually an extrinsic specification that must be carried out on top of the + * intrinsic calculation of the reaction rate + * + * @param iphase Index of the phase. This is the order within the internal thermo vector object + * @param exists Boolean indicating whether the phase exists or not + */ + void setPhaseExistence(const int iphase, const bool exists); protected: @@ -515,7 +541,7 @@ namespace Cantera { //! m_kk is the number of species in all of the phases //! that participate in this kinetics mechanism. - int m_kk; + int m_kk; //! List of reactions numbers which are reversible reactions /*! @@ -626,7 +652,7 @@ namespace Cantera { */ vector_fp m_mu0; - //! Vector of phase potentials + //! Vector of phase electric potentials /*! * Temporary vector containing the potential of each phase * in the kinetics object @@ -720,6 +746,27 @@ namespace Cantera { */ bool m_has_exchange_current_density_formulation; + //! Int flag to indicate that some phases in the kinetics mechanism are + //! non-existent. + /*! + * We change the ROP vectors to make sure that non-existent phases are treated + * correctly in the kinetics operator. The value of this is equal to the number + * of phases which don't exist. + */ + int m_phaseExistsCheck; + + //! Vector of booleans indicating whether phases exist or not + /*! + * Vector of booleans indicating whether a phase exists or not. + * We use this to set the ROP's so that unphysical things don't happen + * + * length = number of phases in the object + * By default all phases exist. + */ + std::vector m_phaseExists; + + std::vector m_rxnPhaseIsReactant; + std::vector m_rxnPhaseIsProduct; int m_ioFlag; private: diff --git a/Cantera/src/kinetics/Kinetics.cpp b/Cantera/src/kinetics/Kinetics.cpp index 1c45823a2..225ba887d 100644 --- a/Cantera/src/kinetics/Kinetics.cpp +++ b/Cantera/src/kinetics/Kinetics.cpp @@ -224,7 +224,8 @@ namespace Cantera { throw CanteraError("speciesPhase", "unknown species "+nm); } - /** + //============================================================================================== + /* * This function takes as an argument the kineticsSpecies index * (i.e., the list index in the list of species in the kinetics * manager) and returns the index of the phase owning the diff --git a/Cantera/src/thermo/IdealGasPhase.h b/Cantera/src/thermo/IdealGasPhase.h index b83cec970..915e7bfa0 100644 --- a/Cantera/src/thermo/IdealGasPhase.h +++ b/Cantera/src/thermo/IdealGasPhase.h @@ -330,6 +330,7 @@ namespace Cantera { */ IdealGasPhase(const IdealGasPhase &right); + //! Asignment operator /*! * Assignment operator for the object. Constructed diff --git a/config.h.in b/config.h.in index 9cc487cbe..60a4d8e8b 100755 --- a/config.h.in +++ b/config.h.in @@ -4,6 +4,19 @@ #ifndef CT_CONFIG_H #define CT_CONFIG_H +//---------------------------- Version Flags ------------------// +// Cantera version -> this will be a double-quoted string value +// refering to branch number within svn +#undef CANTERA_VERSION + +// Integer for major number of Cantera +#define CANTERA_VERSION_MAJORNUMBER 18 +// Flag indicating it's part of major version 18 +#define CANTERA_VERSION_18 1 +// Flag indicating it's a development version +#define CANTERA_VERSION_18_XXX 1 +// Flag indictaing that its part of 1.8_LiquidTransportDevelop branch +#define CANTERA_VERSION_18_LTD 1 //------------------------ Development flags ------------------// // @@ -64,9 +77,6 @@ typedef int ftnlen; // Fortran hidden string length type #undef DARWIN #undef HAS_SSTREAM -// Cantera version -#undef CANTERA_VERSION - // Identify whether the operating system is cygwin's overlay of // windows, with gcc being used as the compiler. #undef CYGWIN diff --git a/preconfig b/preconfig index ee7b62fea..b6ef3e80f 100755 --- a/preconfig +++ b/preconfig @@ -488,8 +488,21 @@ CT_SHARED_LIB=${CT_SHARED_LIB:=clib} # lowercase 'helvetica'. RPFONT=${RPFONT:="Helvetica"} -# Don't change this. -CANTERA_VERSION=${CANTERA_VERSION:="1.8.0"} +# +# Version Control: +# +# We define a new string for every branch created in svn. +# +# Currently expected versions: +# "1.7" Old version of Cantera +# "1.8.0" Released version of 1.8 +# "1.8.x" Development version of 1.8 +# "1.8_liquidTransportDevelop" branched version of 1.8 having to deal with liquid transport +# +# This field gets written into config.h and is also an autoconf variable. +# +# +CANTERA_VERSION=${CANTERA_VERSION:="1.8_liquidTransportDevelop"} #----------------------------------------------------------------------- #------------------- don't change anything below!! ---------------------