diff --git a/Cantera/src/kinetics/Kinetics.h b/Cantera/src/kinetics/Kinetics.h index 124ad846d..ca0a47bbb 100755 --- a/Cantera/src/kinetics/Kinetics.h +++ b/Cantera/src/kinetics/Kinetics.h @@ -149,7 +149,7 @@ namespace Cantera { /// mix_defs.h. virtual int type() { return 0; } - /// Number of reactions in the reaction mechanism. + //! Number of reactions in the reaction mechanism. int nReactions() const {return m_ii;} //@} @@ -297,18 +297,20 @@ namespace Cantera { return m_start[n] + k; } - /** - * Return the std::string name of the kth species in the kinetics - * manager. k is an integer from 0 to ktot - 1, where ktot is - * the number of species in the kinetics manager, which is the - * sum of the number of species in all phases participating in - * the kinetics manager. If k is out of bounds, the std::string - * "" is returned. - * - * @param k species index - */ - std::string kineticsSpeciesName(int k) const; - + + //! Return the std::string name of the kth species in the kinetics + //! manager. + /*! + * k is an integer from 0 to ktot - 1, where ktot is + * the number of species in the kinetics manager, which is the + * sum of the number of species in all phases participating in + * the kinetics manager. If k is out of bounds, the std::string + * "" is returned. + * + * @param k species index + */ + std::string kineticsSpeciesName(int k) const; + /** * This routine will look up a species number based on * the input std::string nm. The lookup of species will @@ -721,9 +723,9 @@ namespace Cantera { */ //@{ - /** - * Add a phase to the kinetics manager object. This must - * be done before the function init() is called or + //! Add a phase to the kinetics manager object. + /*! + * This must be done before the function init() is called or * before any reactions are input. * The following fields are updated: * m_start -> vector of integers, containing the @@ -740,7 +742,7 @@ namespace Cantera { * * @param thermo Reference to the ThermoPhase to be added. */ - void addPhase(thermo_t& thermo); + virtual void addPhase(thermo_t& thermo); /** * Prepare the class for the addition of reactions. This diff --git a/Cantera/src/kinetics/importKinetics.cpp b/Cantera/src/kinetics/importKinetics.cpp index 8e0be3003..0afe47722 100644 --- a/Cantera/src/kinetics/importKinetics.cpp +++ b/Cantera/src/kinetics/importKinetics.cpp @@ -63,29 +63,18 @@ public: std::vector m_typ; //! vector of bools. std::vector m_rev; - ~rxninfo() { - m_eqn.clear(); - m_dup.clear(); - m_nr.clear(); - m_typ.clear(); - m_rdata.clear(); - } - bool installReaction(int i, const XML_Node& r, Kinetics* k, - std::string default_phase, int rule, - bool validate_rxn) ; + ~rxninfo() { + m_eqn.clear(); + m_dup.clear(); + m_nr.clear(); + m_typ.clear(); + m_rdata.clear(); + } + bool installReaction(int i, const XML_Node& r, Kinetics* k, + std::string default_phase, int rule, + bool validate_rxn) ; }; -//! Temporary storage of rxninfo -//rxninfo* _rxns = 0; -//! @name utilitydefines. -//@{ -#define _reactiondata m_rdata -#define _eqn m_eqn -#define _dup m_dup -#define _nr m_nr -#define _typ m_typ -#define _rev m_rev -//@} /* @@ -193,12 +182,11 @@ public: * allowing the calling routine to skip this reaction * and continue. */ - static bool getReagents(const XML_Node& rxn, kinetics_t& kin, int rp, - string default_phase, - vector_int& spnum, vector_fp& stoich, vector_fp& order, - int rule) { + bool getReagents(const XML_Node& rxn, kinetics_t& kin, int rp, + std::string default_phase, + vector_int& spnum, vector_fp& stoich, vector_fp& order, + int rule) { - string rptype; /* @@ -605,11 +593,17 @@ public: * * @ingroup kineticsmgr */ - bool rxninfo::installReaction(int i, const XML_Node& r, Kinetics* k, + bool rxninfo::installReaction(int i, const XML_Node& r, Kinetics* k, string default_phase, int rule, bool validate_rxn) { Kinetics& kin = *k; + + /* Check to see that we are in fact at a reaction node */ + if (r.name() != "reaction") { + throw CanteraError(" rxninfo::installReaction", + " expected xml node reaction, got " + r.name()); + } /* * We use the ReactionData object to store initial values read * in from the xml data. Then, when we have collected everything @@ -618,52 +612,56 @@ public: * the ReactionData object). */ ReactionData rdata; - rdata.reactionType = ELEMENTARY_RXN; // default - vector_int reac, prod; - string eqn, type; - int nn, eqlen; - vector_fp dummy; - - // check to see if the reaction is specified to be a duplicate - // of another reaction, or to allow a negative pre-exponential. + + // Check to see if the reaction is specified to be a duplicate + // of another reaction. It's an error if the reaction is a + // duplicate and this is not set. int dup = 0; if (r.hasAttrib("duplicate")) dup = 1; + + // Check to see if the reaction rate constant can be negative + // It's an error if a negative rate constant is found and + // this is not set. int negA = 0; if (r.hasAttrib("negative_A")) negA = 1; - /* * This seemingly simple expression goes and finds the child element, * "equation". Then it treats all of the contents of the "equation" * as a string, and returns it the variable eqn. We post process * the string to convert [ and ] characters into < and >, which * cannot be stored in an XML file. + * The string eqn is just used for IO purposes. It isn't parsed + * for the identities of reactants or products. */ - if (r.hasChild("equation")) + string eqn = ""; + if (r.hasChild("equation")) { eqn = r("equation"); - else - eqn = ""; - - eqlen = static_cast(eqn.size()); + } + int eqlen = static_cast(eqn.size()); + int nn; for (nn = 0; nn < eqlen; nn++) { if (eqn[nn] == '[') eqn[nn] = '<'; if (eqn[nn] == ']') eqn[nn] = '>'; } - bool ok; + // get the reactants - ok = getReagents(r, kin, 1, default_phase, rdata.reactants, + + bool ok = getReagents(r, kin, 1, default_phase, rdata.reactants, rdata.rstoich, rdata.order, rule); //cout << "Reactants: " << endl; - int npp = rdata.reactants.size(); - int nj; + //int npp = rdata.reactants.size(); + //int nj; //for (nj = 0; nj < npp; nj++) { // cout << rdata.reactants[nj] << " " << rdata.rstoich[nj] << endl; //} + /* * Get the products. We store the id of products in rdata.products - */ + */ + vector_fp dummy; ok = ok && getReagents(r, kin, -1, default_phase, rdata.products, rdata.pstoich, dummy, rule); //cout << "Products: " << endl;npp = rdata.products.size(); @@ -679,11 +677,9 @@ public: // reversible. Default is irreversible. rdata.reversible = false; string isrev = r["reversible"]; - if (isrev == "yes" || isrev == "true") + if (isrev == "yes" || isrev == "true") { rdata.reversible = true; - - - string typ = r["type"]; + } /* * If reaction orders are specified, then this reaction @@ -702,13 +698,13 @@ public: rdata.global = true; } - /* - * Seaarch the reaction element for the attribute "type". + * Search the reaction element for the attribute "type". * If found, then branch on the type, to fill in appropriate * fields in rdata. - */ - + */ + rdata.reactionType = ELEMENTARY_RXN; + string typ = r["type"]; if (typ == "falloff") { rdata.reactionType = FALLOFF_RXN; rdata.falloffType = SIMPLE_FALLOFF; @@ -726,10 +722,10 @@ public: else if (typ == "edge") { rdata.reactionType = EDGE_RXN; } - else if (typ != "") + else if (typ != "") { throw CanteraError("installReaction", "Unknown reaction type: " + typ); - + } /* * Look for undeclared duplicate reactions. */ @@ -746,35 +742,29 @@ public: for (nn = 0; nn < np; nn++) { rxnstoich[rdata.products[nn]+1] += rdata.pstoich[nn]; } - int nrxns = static_cast(_reactiondata.size()); + int nrxns = static_cast(m_rdata.size()); for (nn = 0; nn < nrxns; nn++) { - if ((int(rdata.reactants.size()) == _nr[nn]) - && (rdata.reactionType == _typ[nn])) { - c = isDuplicateReaction(rxnstoich, _reactiondata[nn]); + if ((int(rdata.reactants.size()) == m_nr[nn]) + && (rdata.reactionType == m_typ[nn])) { + c = isDuplicateReaction(rxnstoich, m_rdata[nn]); if (c > 0.0 || (c < 0.0 && rdata.reversible) - || (c < 0.0 && _rev[nn])) { - if ((!dup || !_dup[nn])) { + || (c < 0.0 && m_rev[nn])) { + if ((!dup || !m_dup[nn])) { string msg = string("Undeclared duplicate reactions detected: \n") - +"Reaction "+int2str(nn+1)+": "+_eqn[nn] + +"Reaction "+int2str(nn+1)+": "+m_eqn[nn] +"\nReaction "+int2str(i+1)+": "+eqn+"\n"; - //_reactiondata.clear(); - //_eqn.clear(); - //_rev.clear(); - //_nr.clear(); - //_typ.clear(); - //_dup.clear(); - throw CanteraError("installReaction",msg); + throw CanteraError("installReaction", msg); } } } } - _dup.push_back(dup); - _rev.push_back(rdata.reversible); - _eqn.push_back(eqn); - _nr.push_back(rdata.reactants.size()); - _typ.push_back(rdata.reactionType); - _reactiondata.push_back(rxnstoich); + m_dup.push_back(dup); + m_rev.push_back(rdata.reversible); + m_eqn.push_back(eqn); + m_nr.push_back(rdata.reactants.size()); + m_typ.push_back(rdata.reactionType); + m_rdata.push_back(rxnstoich); } rdata.equation = eqn; @@ -821,10 +811,7 @@ public: bool installReactionArrays(const XML_Node& p, Kinetics& kin, std::string default_phase, bool check_for_duplicates) { - std::auto_ptr< rxninfo > _rxns( new rxninfo ) ; - //if (_rxns == 0) { - //_rxns = new rxninfo; - //} + const std::auto_ptr< rxninfo > _rxns( new rxninfo ) ; //_eqn.clear(); //_dup.clear(); //_nr.clear(); @@ -860,8 +847,6 @@ public: * of the current xml node. */ const XML_Node* rdata = get_XML_Node(rxns["datasrc"], &rxns.root()); - //const XML_Node* rdata = find_XML(rxns["datasrc"],&rxns.root(), - // "","","reactionData"); /* * If the reactionArray element has a child element named * "skip", and if the attribute of skip called "species" has @@ -945,14 +930,7 @@ public: * the true number of reactions in the mechanism, itot. */ kin.finalize(); - //writer = 0; - //_eqn.clear(); - //_dup.clear(); - //_nr.clear(); - //_typ.clear(); - //_reactiondata.clear(); - //delete _rxns; - //_rxns = 0; + return true; } @@ -969,13 +947,15 @@ public: * phaseArray containing a listing of other phases * that participate in the kinetics mechanism. * - * @param th This is a list of ThermoPhase pointers containing + * @param th This is a list of ThermoPhase pointers which must + * include all of * the phases that participate in the kinetics - * reactions. All of the phases must have already + * operator. All of the phases must have already * been initialized and formed within Cantera. * However, their pointers should not have been * added to the Kinetics object; this addition - * is carried out here. + * is carried out here. Additional phases may + * be include; these have no effect. * * @param k This is a pointer to the kinetics manager class * that will be initialized with a kinetics @@ -988,8 +968,10 @@ public: Kinetics& kin = *k; - // This phase will be the default one - string default_phase = phase["id"]; + // This phase will be the owning phase for the kinetics operator + // For interfaces, it is the surface phase between two volumes. + // For homogeneous kinetics, it's the current volumetric phase. + string owning_phase = phase["id"]; bool check_for_duplicates = false; if (phase.parent()->hasChild("validate")) { @@ -1007,7 +989,7 @@ public: const XML_Node& pa = phase.child("phaseArray"); getStringArray(pa, phase_ids); } - phase_ids.push_back(default_phase); + phase_ids.push_back(owning_phase); int np = static_cast(phase_ids.size()); int nt = static_cast(th.size()); @@ -1048,7 +1030,7 @@ public: kin.init(); // Install the reactions. - return installReactionArrays(phase, kin, default_phase, check_for_duplicates); + return installReactionArrays(phase, kin, owning_phase, check_for_duplicates); } /* diff --git a/Cantera/src/kinetics/importKinetics.h b/Cantera/src/kinetics/importKinetics.h index 2bd7ef58a..404937bc7 100644 --- a/Cantera/src/kinetics/importKinetics.h +++ b/Cantera/src/kinetics/importKinetics.h @@ -69,6 +69,40 @@ namespace Cantera { const ReactionData &rdata, doublereal errorTolerance = 1.0e-3); + /** + * Get the reactants or products of a reaction. The information + * is returned in the spnum, stoich, and order vectors. The + * length of the vectors is the number of different types of + * reactants or products found for the reaction. + * + * Input + * -------- + * rxn -> xml node pointing to the reaction element + * in the xml tree. + * kin -> Reference to the kinetics object to install + * the information into. + * rp = 1 -> Go get the reactants for a reaction + * -1 -> Go get the products for a reaction + * default_phase = String name for the default phase + * to loop up species in. + * Output + * ----------- + * spnum = vector of species numbers found. + * Length is number of reactants or products. + * stoich = stoichiometric coefficient of the reactant or product + * Length is number of reactants or products. + * order = Order of the reactant and product in the reaction + * rate expression + * @param rule If we fail to find a species, we will throw an error + * if rule != 1. If rule = 1, we simply return false, + * allowing the calling routine to skip this reaction + * and continue. + */ + bool getReagents(const XML_Node& rxn, kinetics_t& kin, int rp, + std::string default_phase, + vector_int& spnum, vector_fp& stoich, vector_fp& order, + int rule); + //! Read the rate coefficient data from the XML file. /*! * Extract the rate coefficient for a reaction from the xml node, kf. @@ -128,11 +162,12 @@ namespace Cantera { * For example, if phase I is an interface phase between bulk * phases A and B. Then, the XML_Node for phase I should be * the first argument. - * The vector of %ThermoPhase objects should be consist of pointers + * The vector of %ThermoPhase objects should consist of pointers * to phases I, A, and B. * * @param phase This is an xml node containing a description - * of a phase. Within the phase is a XML element + * of the owning phase for the kinetics object. + * Within the phase is a XML element * called reactionArray containing the location * of the description of the reactions that make * up the kinetics object. @@ -140,21 +175,24 @@ namespace Cantera { * phaseArray containing a listing of other phases * that participate in the kinetics mechanism. * - * @param th This is a list of ThermoPhase pointers containing + * @param th This is a list of ThermoPhase pointers which must + * include all of * the phases that participate in the kinetics - * reactions. All of the phases must have already + * operator. All of the phases must have already * been initialized and formed within Cantera. * However, their pointers should not have been * added to the Kinetics object; this addition - * is carried out here. + * is carried out here. Additional phases may + * be include in the list; these have no effect. * - * @param kin This is a pointer to a bare kinetics manager class + * @param kin This is a pointer to a kinetics manager class * that will be initialized with the kinetics - * mechanism. + * mechanism. Inherited Kinetics classes may be + * used here. * * @ingroup kineticsmgr * - */ + */ bool importKinetics(const XML_Node& phase, std::vector th, Kinetics* kin);