00001 /** 00002 * @file Kinetics.h 00003 * Base class for kinetics managers and also contains the kineticsmgr 00004 * module documentation (see \ref kineticsmgr and class 00005 * \link Cantera::Kinetics Kinetics\endlink). 00006 * 00007 * $Date: 2009-12-10 18:20:44 -0500 (Thu, 10 Dec 2009) $ 00008 * $Revision: 309 $ 00009 */ 00010 00011 // Copyright 2001-2004 California Institute of Technology 00012 00013 #ifndef CT_KINETICS_H 00014 #define CT_KINETICS_H 00015 00016 #include "ctexceptions.h" 00017 #include "ThermoPhase.h" 00018 #include "mix_defs.h" 00019 00020 namespace Cantera { 00021 00022 // forward references 00023 class ReactionData; 00024 00025 /** 00026 * @defgroup chemkinetics Chemical Kinetics 00027 */ 00028 00029 /// @defgroup kineticsmgr Kinetics Managers 00030 /// @section kinmodman Models and Managers 00031 /// 00032 /// A kinetics manager is a C++ class that implements a kinetics 00033 /// model; a kinetics model is a set of mathematical equation 00034 /// describing how various kinetic quanities are to be computed -- 00035 /// reaction rates, species production rates, etc. Many different 00036 /// kinetics models might be defined to handle different types of 00037 /// kinetic processes. For example, one kinetics model might use 00038 /// expressions valid for elementary reactions in ideal gas 00039 /// mixtures. It might, for example, require the reaction orders 00040 /// to be integral and equal to the forward stoichiometric 00041 /// coefficients, require that each reaction be reversible with a 00042 /// reverse rate satisfying detailed balance, include 00043 /// pressure-dependent unimolecular reactions, etc. Another 00044 /// kinetics model might be designed for heterogeneous chemistry 00045 /// at interfaces, and might allow empirical reaction orders, 00046 /// coverage-dependent activation energies, irreversible 00047 /// reactions, and include effects of potential differences across 00048 /// the interface on reaction rates. 00049 /// 00050 /// A kinetics manager implements a kinetics model. Since the 00051 /// model equations may be complex and expensive to evaluate, a 00052 /// kinetics manager may adopt various strategies to 'manage' the 00053 /// computation and evaluate the expressions efficiently. For 00054 /// example, if there are rate coefficients or other quantities 00055 /// that depend only on temperature, a manager class may choose to 00056 /// store these quantities internally, and re-evaluate them only 00057 /// when the temperature has actually changed. Or a manager 00058 /// designed for use with reaction mechanisms with a few repeated 00059 /// activation energies might precompute the terms \f$ exp(-E/RT) 00060 /// \f$, instead of evaluating the exponential repeatedly for each 00061 /// reaction. There are many other possible 'management styles', 00062 /// each of which might be better suited to some reaction 00063 /// mechanisms than others. 00064 /// 00065 /// But however a manager structures the internal computation, the 00066 /// tasks the manager class must perform are, for the most part, 00067 /// the same. It must be able to compute reaction rates, species 00068 /// production rates, equilibrium constants, etc. Therefore, all 00069 /// kinetics manager classes should have a common set of public 00070 /// methods, but differ in how they implement these methods. 00071 /// 00072 /// A kinetics manager computes reaction rates of progress, 00073 /// species production rates, equilibrium constants, and similar 00074 /// quantities for a reaction mechanism. All kinetics manager 00075 /// classes derive from class Kinetics, which defines a common 00076 /// public interface for all kinetics managers. Each derived class 00077 /// overloads the virtual methods of Kinetics to implement a 00078 /// particular kinetics model. 00079 /// 00080 /// For example, class GasKinetics implements reaction rate 00081 /// expressions appropriate for homogeneous reactions in ideal gas 00082 /// mixtures, and class InterfaceKinetics implements expressions 00083 /// appropriate for heterogeneous mechanisms at interfaces, 00084 /// including how to handle reactions involving charged species of 00085 /// phases with different electric potentials --- something that 00086 /// class GasKinetics doesn't deal with at all. 00087 /// 00088 /// Kinetics managers may be also created that hard-wire a 00089 /// particular reaction mechanism in C++ code. This can often 00090 /// result in faster performance. An example of this is the 00091 /// kinetics manager GRI30_Kinetics that hard-wires the rate 00092 /// expressions for the natural gas combustion mechanism GRI-3.0. 00093 /// 00094 /// Many of the methods of class Kinetics write into arrays the 00095 /// values of some quantity for each species, for example the net 00096 /// production rate. These methods always write the results into 00097 /// flat arrays, ordered by phase in the order the phase was 00098 /// added, and within a phase in the order the species were added 00099 /// to the phase (which is the same ordering as in the input 00100 /// file). Example: suppose a heterogeneous mechanism involves 00101 /// three phases -- a bulk phase 'a', another bulk phase 'b', and 00102 /// the surface phase 'a:b' at the a/b interface. Phase 'a' 00103 /// contains 12 species, phase 'b' contains 3, and at the 00104 /// interface there are 5 adsorbed species defined in phase 00105 /// 'a:b'. Then methods like getNetProductionRates(doublereal* net) 00106 /// will write and output array of length 20, beginning at the location 00107 /// pointed to by 'net'. The first 12 values will be the net production 00108 /// rates for all 12 species of phase 'a' (even if some do not participate 00109 /// in the reactions), the next 3 will be for phase 'b', and finally the 00110 /// net production rates for the surface species will occupy the last 00111 /// 5 locations. 00112 /// @ingroup chemkinetics 00113 00114 00115 //! Public interface for kinetics managers. 00116 /*! 00117 * This class serves as a base class to derive 'kinetics 00118 * managers', which are classes that manage homogeneous chemistry 00119 * within one phase, or heterogeneous chemistry at one 00120 * interface. The virtual methods of this class are meant to be 00121 * overloaded in subclasses. The non-virtual methods perform 00122 * generic functions and are implemented in Kinetics. They should 00123 * not be overloaded. Only those methods required by a subclass 00124 * need to be overloaded; the rest will throw exceptions if 00125 * called. 00126 * 00127 * When the nomenclature "kinetics species index" is used below, 00128 * this means that the species index ranges over all species in 00129 * all phases handled by the kinetics manager. 00130 * 00131 * @ingroup kineticsmgr 00132 */ 00133 class Kinetics { 00134 00135 public: 00136 00137 //! typedef for ThermoPhase 00138 typedef ThermoPhase thermo_t; 00139 00140 /** 00141 * @name Constructors and General Information about Mechanism 00142 */ 00143 //@{ 00144 00145 /// Default constructor. 00146 Kinetics(); 00147 00148 /// This constructor initializes with a starting phase. 00149 /// @deprecated 00150 // Kinetics(thermo_t* thermo); 00151 00152 /// Destructor. 00153 virtual ~Kinetics(); 00154 00155 //!Copy Constructor for the %Kinetics object. 00156 /*! 00157 * Currently, this is not fully implemented. If called it will 00158 * throw an exception. 00159 */ 00160 Kinetics(const Kinetics &); 00161 00162 //! Assignment operator 00163 /*! 00164 * This is NOT a virtual function. 00165 * 00166 * @param right Reference to %Kinetics object to be copied into the 00167 * current one. 00168 */ 00169 Kinetics& operator=(const Kinetics &right); 00170 00171 00172 //! Duplication routine for objects which inherit from 00173 //! Kinetics 00174 /*! 00175 * This virtual routine can be used to duplicate %Kinetics objects 00176 * inherited from %Kinetics even if the application only has 00177 * a pointer to %Kinetics to work with. 00178 * 00179 * These routines are basically wrappers around the derived copy 00180 * constructor. 00181 */ 00182 virtual Kinetics *duplMyselfAsKinetics() const; 00183 00184 00185 //! Identifies the subclass of the Kinetics manager type. 00186 /*! 00187 * These are listed in mix_defs.h. 00188 */ 00189 virtual int ID() const; 00190 00191 //! Identifies the kinetics manager type. 00192 /*! 00193 * Each class derived from Kinetics should overload this method to 00194 * return a unique integer. Standard values are defined in file 00195 * mix_defs.h. 00196 */ 00197 virtual int type() const; 00198 00199 //! Number of reactions in the reaction mechanism. 00200 int nReactions() const {return m_ii;} 00201 00202 //@} 00203 00204 00205 /** 00206 * @name Information/Lookup Functions about Phases and Species 00207 */ 00208 //@{ 00209 00210 /** 00211 * The number of phases participating in the reaction 00212 * mechanism. For a homogeneous reaction mechanism, this will 00213 * always return 1, but for a heterogeneous mechanism it will 00214 * return the total number of phases in the mechanism. 00215 */ 00216 int nPhases() const { return static_cast<int>(m_thermo.size()); } 00217 00218 /** 00219 * Return the phase index of a phase in the list of phases 00220 * defined within the object. 00221 * 00222 * @param ph std::string name of the phase 00223 * 00224 * If a -1 is returned, then the phase is not defined in 00225 * the Kinetics object. 00226 */ 00227 int phaseIndex(std::string ph) { 00228 if (m_phaseindex.find(ph) == m_phaseindex.end()) { 00229 return -1; 00230 } 00231 else { 00232 return m_phaseindex[ph] - 1; 00233 } 00234 } 00235 00236 /** 00237 * This returns the integer index of the phase which has 00238 * ThermoPhase type cSurf. For heterogeneous mechanisms, this 00239 * identifies the one surface phase. For homogeneous 00240 * mechanisms, this reurns -1. 00241 */ 00242 int surfacePhaseIndex() { return m_surfphase; } 00243 00244 /** 00245 * Phase where the reactions occur. For heterogeneous 00246 * mechanisms, one of the phases in the list of phases 00247 * represents the 2D interface or 1D edge at which the 00248 * reactions take place. This method returns the index of the 00249 * phase with the smallest spatial dimension (1, 2, or 3) 00250 * among the list of phases. If there is more than one, the 00251 * index of the first one is returned. For homogeneous 00252 * mechanisms, the value 0 is returned. 00253 */ 00254 int reactionPhaseIndex() { return m_rxnphase; } 00255 00256 00257 /** 00258 * This method returns a reference to the nth ThermoPhase 00259 * object defined in this kinetics mechanism. It is typically 00260 * used so that member functions of the ThermoPhase object may 00261 * be called. For homogeneous mechanisms, there is only one 00262 * object, and this method can be called without an argument 00263 * to access it. 00264 * 00265 * @param n Index of the ThermoPhase being sought. 00266 */ 00267 thermo_t& thermo(int n=0) { return *m_thermo[n]; } 00268 const thermo_t& thermo(int n=0) const { return *m_thermo[n]; } 00269 00270 /** 00271 * This method returns a reference to the nth ThermoPhase 00272 * defined in this kinetics mechanism. 00273 * It is typically used so that member functions of the 00274 * ThermoPhase may be called. @deprecated This method is redundant. 00275 * 00276 * @param n Index of the ThermoPhase being sought. 00277 */ 00278 thermo_t& phase(int n=0) { 00279 deprecatedMethod("Kinetics","phase","thermo"); 00280 return *m_thermo[n]; 00281 } 00282 /** 00283 * This method returns a reference to the nth ThermoPhase 00284 * defined in this kinetics mechanism. 00285 * It is typically used so that member functions of the 00286 * ThermoPhase may be called. @deprecated This method is redundant. 00287 * 00288 * @param n Index of the ThermoPhase being sought. 00289 */ 00290 const thermo_t& phase(int n=0) const { 00291 deprecatedMethod("Kinetics","phase","thermo"); 00292 return *m_thermo[n]; 00293 } 00294 00295 /** 00296 * The total number of species in all phases participating in 00297 * the kinetics mechanism. This is useful to dimension arrays 00298 * for use in calls to methods that return the species 00299 * production rates, for example. 00300 */ 00301 int nTotalSpecies() const { 00302 int n=0, np; 00303 np = nPhases(); 00304 for (int p = 0; p < np; p++) n += thermo(p).nSpecies(); 00305 return n; 00306 } 00307 00308 /** 00309 * Returns the starting index of the species in the nth phase 00310 * associated with the reaction mechanism. 00311 * 00312 * @param n Return the index of first species in the nth phase 00313 * associated with the reaction mechanism. 00314 */ 00315 int start(int n) { 00316 deprecatedMethod("Kinetics","start","kineticsSpeciesIndex(0,n)"); 00317 return m_start[n]; 00318 } 00319 00320 00321 /** 00322 * The location of species k of phase n in species arrays. 00323 * Kinetics manager classes return species production rates in 00324 * flat arrays, with the species of each phases following one 00325 * another, in the order the phases were added. This method 00326 * is useful to find the value for a particular species of a 00327 * particular phase in arrrays returned from methods like 00328 * getCreationRates that return an array of species-specific 00329 * quantities. 00330 * 00331 * Example: suppose a heterogeneous mechanism involves three 00332 * phases. The first contains 12 species, the second 26, and 00333 * the third 3. Then species arrays must have size at least 00334 * 41, and positions 0 - 11 are the values for the species in 00335 * the first phase, positions 12 - 37 are the values for the 00336 * species in the second phase, etc. Then 00337 * kineticsSpeciesIndex(7, 0) = 7, kineticsSpeciesIndex(4, 1) 00338 * = 16, and kineticsSpeciesIndex(2, 2) = 40. 00339 * 00340 * @param k species index 00341 * @param n phase index for the species 00342 */ 00343 int kineticsSpeciesIndex(int k, int n) const { 00344 return m_start[n] + k; 00345 } 00346 00347 00348 //! Return the std::string name of the kth species in the kinetics 00349 //! manager. 00350 /*! 00351 * k is an integer from 0 to ktot - 1, where ktot is 00352 * the number of species in the kinetics manager, which is the 00353 * sum of the number of species in all phases participating in 00354 * the kinetics manager. If k is out of bounds, the std::string 00355 * "<unknown>" is returned. 00356 * 00357 * @param k species index 00358 */ 00359 std::string kineticsSpeciesName(int k) const; 00360 00361 /** 00362 * This routine will look up a species number based on 00363 * the input std::string nm. The lookup of species will 00364 * occur for all phases listed in the kinetics object, 00365 * unless the std::string ph refers to a specific phase of 00366 * the object. 00367 * 00368 * return 00369 * - If a match is found, the position in the species list 00370 * is returned. 00371 * - If a specific phase is specified and no match is found, 00372 * the value -1 is returned. 00373 * - If no match is found in any phase, the value -2 is returned. 00374 * 00375 * @param nm Input string name of the species 00376 * @param ph Input string name of the phase. Defaults to "<any>" 00377 */ 00378 int kineticsSpeciesIndex(std::string nm, std::string ph = "<any>") const; 00379 00380 /** 00381 * This function looks up the std::string name of a species and 00382 * returns a reference to the ThermoPhase object of the 00383 * phase where the species resides. 00384 * Will throw an error if the species std::string doesn't match. 00385 * 00386 * @param nm String containing the name of the species. 00387 */ 00388 thermo_t& speciesPhase(std::string nm); 00389 00390 /** 00391 * This function takes as an argument the kineticsSpecies index 00392 * (i.e., the list index in the list of species in the kinetics 00393 * manager) and returns the species' owning ThermoPhase object. 00394 * 00395 * @param k Species index 00396 */ 00397 thermo_t& speciesPhase(int k) { 00398 return thermo(speciesPhaseIndex(k)); 00399 } 00400 00401 /** 00402 * This function takes as an argument the kineticsSpecies index 00403 * (i.e., the list index in the list of species in the kinetics 00404 * manager) and returns the index of the phase owning the 00405 * species. 00406 * 00407 * @param k Species index 00408 */ 00409 int speciesPhaseIndex(int k); 00410 00411 //@} 00412 00413 00414 00415 /** 00416 * @name Reaction Rates Of Progress 00417 */ 00418 //@{ 00419 00420 //! Return the forward rates of progress of the reactions 00421 /*! 00422 * Forward rates of progress. Return the forward rates of 00423 * progress in array fwdROP, which must be dimensioned at 00424 * least as large as the total number of reactions. 00425 * 00426 * @param fwdROP Output vector containing forward rates 00427 * of progress of the reactions. Length: m_ii. 00428 */ 00429 virtual void getFwdRatesOfProgress(doublereal* fwdROP) { 00430 err("getFwdRatesOfProgress"); 00431 } 00432 00433 //! Return the Reverse rates of progress of the reactions 00434 /*! 00435 * Return the reverse rates of 00436 * progress in array revROP, which must be dimensioned at 00437 * least as large as the total number of reactions. 00438 * 00439 * @param revROP Output vector containing reverse rates 00440 * of progress of the reactions. Length: m_ii. 00441 */ 00442 virtual void getRevRatesOfProgress(doublereal* revROP) { 00443 err("getRevRatesOfProgress"); 00444 } 00445 00446 /** 00447 * Net rates of progress. Return the net (forward - reverse) 00448 * rates of progress in array netROP, which must be 00449 * dimensioned at least as large as the total number of 00450 * reactions. 00451 * 00452 * @param netROP Output vector of the net ROP. Length: m_ii. 00453 */ 00454 virtual void getNetRatesOfProgress(doublereal* netROP) { 00455 err("getNetRatesOfProgress"); 00456 } 00457 00458 00459 00460 //! Return a vector of Equilibrium constants. 00461 /*! 00462 * Return the equilibrium constants of 00463 * the reactions in concentration units in array kc, which 00464 * must be dimensioned at least as large as the total number 00465 * of reactions. 00466 * 00467 * @param kc Output vector containing the equilibrium constants. 00468 * Length: m_ii. 00469 */ 00470 virtual void getEquilibriumConstants(doublereal* kc) { 00471 err("getEquilibriumConstants"); 00472 } 00473 00474 /** 00475 * Change in species properties. Given an array of molar species 00476 * property values \f$ z_k, k = 1, \dots, K \f$, return the 00477 * array of reaction values 00478 * \f[ 00479 * \Delta Z_i = \sum_k \nu_{k,i} z_k, i = 1, \dots, I. 00480 * \f] 00481 * For example, if this method is called with the array of 00482 * standard-state molar Gibbs free energies for the species, 00483 * then the values returned in array \c deltaProperty would be 00484 * the standard-state Gibbs free energies of reaction for each 00485 * reaction. 00486 * 00487 * @param property Input vector of property value. Length: m_kk. 00488 * @param deltaProperty Output vector of deltaRxn. Length: m_ii. 00489 */ 00490 virtual void getReactionDelta(const doublereal* property, 00491 doublereal* deltaProperty) { 00492 err("getReactionDelta"); 00493 } 00494 00495 /** 00496 * Return the vector of values for the reaction gibbs free 00497 * energy change. These values depend upon the concentration 00498 * of the solution. 00499 * 00500 * units = J kmol-1 00501 * 00502 * @param deltaG Output vector of deltaG's for reactions 00503 * Length: m_ii. 00504 */ 00505 virtual void getDeltaGibbs( doublereal* deltaG) { 00506 err("getDeltaGibbs"); 00507 } 00508 00509 /** 00510 * Return the vector of values for the reactions change in 00511 * enthalpy. These values depend upon the concentration of 00512 * the solution. 00513 * 00514 * units = J kmol-1 00515 * 00516 * @param deltaH Output vector of deltaH's for reactions 00517 * Length: m_ii. 00518 */ 00519 virtual void getDeltaEnthalpy( doublereal* deltaH) { 00520 err("getDeltaEnthalpy"); 00521 } 00522 00523 /** 00524 * Return the vector of values for the reactions change in 00525 * entropy. These values depend upon the concentration of the 00526 * solution. 00527 * 00528 * units = J kmol-1 Kelvin-1 00529 * 00530 * @param deltaS Output vector of deltaS's for reactions 00531 * Length: m_ii. 00532 */ 00533 virtual void getDeltaEntropy( doublereal* deltaS) { 00534 err("getDeltaEntropy"); 00535 } 00536 00537 /** 00538 * Return the vector of values for the reaction standard state 00539 * gibbs free energy change. These values don't depend upon 00540 * the concentration of the solution. 00541 * 00542 * units = J kmol-1 00543 * 00544 * @param deltaG Output vector of ss deltaG's for reactions 00545 * Length: m_ii. 00546 */ 00547 virtual void getDeltaSSGibbs( doublereal* deltaG) { 00548 err("getDeltaSSGibbs"); 00549 } 00550 00551 /** 00552 * Return the vector of values for the change in the standard 00553 * state enthalpies of reaction. These values don't depend 00554 * upon the concentration of the solution. 00555 * 00556 * units = J kmol-1 00557 * 00558 * @param deltaH Output vector of ss deltaH's for reactions 00559 * Length: m_ii. 00560 */ 00561 virtual void getDeltaSSEnthalpy( doublereal* deltaH) { 00562 err("getDeltaSSEnthalpy"); 00563 } 00564 00565 /** 00566 * Return the vector of values for the change in the standard 00567 * state entropies for each reaction. These values don't 00568 * depend upon the concentration of the solution. 00569 * 00570 * units = J kmol-1 Kelvin-1 00571 * 00572 * @param deltaS Output vector of ss deltaS's for reactions 00573 * Length: m_ii. 00574 */ 00575 virtual void getDeltaSSEntropy( doublereal* deltaS) { 00576 err("getDeltaSSEntropy"); 00577 } 00578 00579 00580 //@} 00581 /** 00582 * @name Species Production Rates 00583 */ 00584 //@{ 00585 00586 /** 00587 * Species creation rates [kmol/m^3/s or kmol/m^2/s]. Return the 00588 * species creation rates in array cdot, which must be 00589 * dimensioned at least as large as the total number of 00590 * species in all phases. @see nTotalSpecies. 00591 * 00592 * @param cdot Output vector of creation rates. 00593 * Length: m_kk. 00594 */ 00595 virtual void getCreationRates(doublereal* cdot) { 00596 err("getCreationRates"); 00597 } 00598 00599 /** 00600 * Species destruction rates [kmol/m^3/s or kmol/m^2/s]. Return 00601 * the species destruction rates in array ddot, which must be 00602 * dimensioned at least as large as the total number of 00603 * species. @see nTotalSpecies. 00604 * 00605 * @param ddot Output vector of destruction rates. 00606 * Length: m_kk. 00607 */ 00608 virtual void getDestructionRates(doublereal* ddot) { 00609 err("getDestructionRates"); 00610 } 00611 00612 /** 00613 * Species net production rates [kmol/m^3/s or kmol/m^2/s]. Return 00614 * the species net production rates (creation - destruction) 00615 * in array wdot, which must be dimensioned at least as large 00616 * as the total number of species. @see nTotalSpecies. 00617 * 00618 * @param wdot Output vector of net production rates. 00619 * Length: m_kk. 00620 */ 00621 virtual void getNetProductionRates(doublereal* wdot) { 00622 err("getNetProductionRates"); 00623 } 00624 00625 //@} 00626 00627 00628 /** 00629 * @name Reaction Mechanism Informational Query Routines 00630 */ 00631 //@{ 00632 00633 /** 00634 * Stoichiometric coefficient of species k as a reactant in 00635 * reaction i. 00636 * 00637 * @param k kinetic species index 00638 * @param i reaction index 00639 */ 00640 virtual doublereal reactantStoichCoeff(int k, int i) const { 00641 err("reactantStoichCoeff"); 00642 return -1.0; 00643 } 00644 00645 /** 00646 * Stoichiometric coefficient of species k as a product in 00647 * reaction i. 00648 * 00649 * @param k kinetic species index 00650 * @param i reaction index 00651 */ 00652 virtual doublereal productStoichCoeff(int k, int i) const { 00653 err("productStoichCoeff"); 00654 return -1.0; 00655 } 00656 00657 /** 00658 * reactant Order of species k in reaction i. 00659 * 00660 * @param k kinetic species index 00661 * @param i reaction index 00662 */ 00663 virtual doublereal reactantOrder(int k, int i) const { 00664 err("reactantOrder"); 00665 return -1.0; 00666 } 00667 00668 /** 00669 * Returns a read-only reference to the vector of reactant 00670 * index numbers for reaction i. 00671 * 00672 * @param i reaction index 00673 */ 00674 virtual const vector_int& reactants(int i) const { 00675 return m_reactants[i]; 00676 } 00677 00678 /** 00679 * Returns a read-only reference to the vector of product 00680 * index numbers for reaction i. 00681 * 00682 * @param i reaction index 00683 */ 00684 virtual const vector_int& products(int i) const { 00685 return m_products[i]; 00686 } 00687 00688 /** 00689 * Flag specifying the type of reaction. The legal values and 00690 * their meaning are specific to the particular kinetics 00691 * manager. 00692 * 00693 * @param i reaction index 00694 */ 00695 virtual int reactionType(int i) const { 00696 err("reactionType"); 00697 return -1; 00698 } 00699 00700 /** 00701 * True if reaction i has been declared to be reversible. If 00702 * isReversible(i) is false, then the reverse rate of progress 00703 * for reaction i is always zero. 00704 * 00705 * @param i reaction index 00706 */ 00707 virtual bool isReversible(int i){ 00708 err("isReversible"); 00709 return false; 00710 } 00711 00712 /** 00713 * Return a std::string representing the reaction. 00714 * 00715 * @param i reaction index 00716 */ 00717 virtual std::string reactionString(int i) const { 00718 err("reactionStd::String"); return "<null>"; 00719 } 00720 00721 /** 00722 * Return the forward rate constants 00723 * 00724 * length is the number of reactions. units depends 00725 * on many issues. @todo DGG: recommend changing name to 00726 * getFwdRateCoefficients. 00727 * 00728 * @param kfwd Output vector containing the foward reaction rate constants. 00729 * Length: m_ii. 00730 */ 00731 virtual void getFwdRateConstants(doublereal *kfwd) { 00732 err("getFwdRateConstants"); 00733 } 00734 00735 /** 00736 * Return the reverse rate constants. 00737 * 00738 * length is the number of reactions. units depends 00739 * on many issues. Note, this routine will return rate constants 00740 * for irreversible reactions if the default for 00741 * doIrreversible is overridden. @todo DGG: recommend changing name to 00742 * getRevRateCoefficients. 00743 * 00744 * @param krev Output vector of reverse rate constants. 00745 * @param doIrreversible boolean indicating whether irreversible reactions 00746 * should be included. 00747 */ 00748 virtual void getRevRateConstants(doublereal *krev, 00749 bool doIrreversible = false) { 00750 err("getFwdRateConstants"); 00751 } 00752 00753 00754 /** 00755 * Return the activation energies in Kelvin. 00756 * 00757 * length is the number of reactions 00758 * 00759 * @param E Ouptut vector of activation energies. 00760 * Length: m_ii. 00761 */ 00762 virtual void getActivationEnergies(doublereal *E) { 00763 err("getActivationEnergies"); 00764 } 00765 00766 00767 //@} 00768 /** 00769 * @name Reaction Mechanism Construction 00770 */ 00771 //@{ 00772 00773 //! Add a phase to the kinetics manager object. 00774 /*! 00775 * This must be done before the function init() is called or 00776 * before any reactions are input. 00777 * The following fields are updated: 00778 * m_start -> vector of integers, containing the 00779 * starting position of the species for 00780 * each phase in the kinetics mechanism. 00781 * m_surfphase -> index of the surface phase. 00782 * m_thermo -> vector of pointers to ThermoPhase phases 00783 * that participate in the kinetics 00784 * mechanism. 00785 * m_phaseindex -> map containing the std::string id of each 00786 * ThermoPhase phase as a key and the 00787 * index of the phase within the kinetics 00788 * manager object as the value. 00789 * 00790 * @param thermo Reference to the ThermoPhase to be added. 00791 */ 00792 virtual void addPhase(thermo_t& thermo); 00793 00794 /** 00795 * Prepare the class for the addition of reactions. This 00796 * method is called by function importKinetics after all 00797 * phases have been added but before any reactions have 00798 * been. The base class method does nothing, but derived 00799 * classes may use this to perform any initialization 00800 * (allocating arrays, etc.) that requires knowing the phases 00801 * and species, but before any reactions are added. 00802 */ 00803 virtual void init() {} 00804 00805 /** 00806 * Finish adding reactions and prepare for use. This method is 00807 * called by function importKinetics after all reactions have 00808 * been entered into the mechanism and before the mechanism is 00809 * used to calculate reaction rates. The base class method 00810 * does nothing, but derived classes may use this to perform 00811 * any initialization (allocating arrays, etc.) that must be 00812 * done after the reactions are entered. 00813 */ 00814 virtual void finalize() {} 00815 00816 /** 00817 * Add a single reaction to the mechanism. This routine 00818 * must be called after init() and before finalize(). 00819 * 00820 * @param r Reference to the ReactionData object for the reaction 00821 * to be added. 00822 */ 00823 virtual void addReaction(const ReactionData& r) { 00824 err("addReaction"); 00825 } 00826 00827 virtual const std::vector<grouplist_t>& reactantGroups(int i) { 00828 //err("reactantGroups"); 00829 return m_dummygroups; 00830 } 00831 00832 virtual const std::vector<grouplist_t>& productGroups(int i) { 00833 //err("productGroups"); 00834 return m_dummygroups; 00835 } 00836 00837 00838 //@} 00839 /** 00840 * @name Altering Reaction Rates 00841 * 00842 * These methods alter reaction rates. They are designed 00843 * primarily for carrying out sensitivity analysis, but may be 00844 * used for any purpose requiring dynamic alteration of rate 00845 * constants. For each reaction, a real-valued multiplier may 00846 * be defined that multiplies the reaction rate 00847 * coefficient. The multiplier may be set to zero to 00848 * completely remove a reaction from the mechanism. 00849 */ 00850 //@{ 00851 00852 /// The current value of the multiplier for reaction i. 00853 /*! 00854 * @param i index of the reaction 00855 */ 00856 doublereal multiplier(int i) const {return m_perturb[i];} 00857 00858 /// Set the multiplier for reaction i to f. 00859 /*! 00860 * @param i index of the reaction 00861 * @param f value of the multiplier. 00862 */ 00863 void setMultiplier(int i, doublereal f) {m_perturb[i] = f;} 00864 00865 //@} 00866 00867 /** 00868 * Increment the number of reactions in the mechanism by one. 00869 * @todo Should be protected? 00870 */ 00871 void incrementRxnCount() { m_ii++; m_perturb.push_back(1.0); } 00872 00873 /** 00874 * Returns true if the kinetics manager has been properly 00875 * initialized and finalized. 00876 */ 00877 virtual bool ready() const { 00878 return false; 00879 } 00880 00881 00882 /** 00883 * Extract from array \c data the portion pertaining to phase \c phase. 00884 * 00885 * @param data data 00886 * @param phase phase 00887 * @param phase_data phase_data 00888 */ 00889 void selectPhase(const doublereal* data, const thermo_t* phase, 00890 doublereal* phase_data); 00891 00892 /// For internal use. May be removed in a future release. 00893 int index(){ return m_index; } 00894 00895 //! Set the index of the Kinetics Manager 00896 /*! 00897 * @param index input index 00898 */ 00899 void setIndex(int index) { m_index = index; } 00900 00901 00902 protected: 00903 00904 00905 //! Number of reactions in the mechanism 00906 int m_ii; 00907 00908 /// Vector of perturbation factors for each reaction's rate of 00909 /// progress vector. It is initialized to one. 00910 /// 00911 vector_fp m_perturb; 00912 00913 /** 00914 * This is a vector of vectors containing the reactants for 00915 * each reaction. The outer vector is over the number of 00916 * reactions, m_ii. The inner vector is a list of species 00917 * indices. If the stoichiometric coefficient for a reactant 00918 * is greater than one, then the reactant is listed 00919 * contiguously in the vector a number of times equal to its 00920 * stoichiometric coefficient. 00921 * NOTE: These vectors will be wrong if there are real 00922 * stoichiometric coefficients in the expression. 00923 */ 00924 std::vector<vector_int> m_reactants; 00925 00926 /** 00927 * This is a vector of vectors containing the products for 00928 * each reaction. The outer vector is over the number of 00929 * reactions, m_ii. The inner vector is a list of species 00930 * indeces. If the stoichiometric coefficient for a product is 00931 * greater than one, then the reactant is listed contiguously 00932 * in the vector a number of times equal to its stoichiometric 00933 * coefficient. 00934 * NOTE: These vectors will be wrong if there are real 00935 * stoichiometric coefficients in the expression. 00936 */ 00937 std::vector<vector_int> m_products; 00938 00939 //! m_thermo is a vector of pointers to ThermoPhase 00940 //! objects. 00941 /*! 00942 * For homogeneous kinetics applications, this vector 00943 * will only have one entry. For interfacial reactions, this 00944 * vector will consist of multiple entries; some of them will 00945 * be surface phases, and the other ones will be bulk phases. 00946 * The order that the objects are listed determines the order 00947 * in which the species comprising each phase are listed in 00948 * the source term vector, originating from the reaction 00949 * mechanism. 00950 * 00951 * Note that this kinetics object doesn't own these ThermoPhase 00952 * objects and is not responsible for creating or deleting 00953 * them. 00954 */ 00955 std::vector<thermo_t*> m_thermo; 00956 00957 /** 00958 * m_start is a vector of integers specifying the beginning position 00959 * for the species vector for the n'th phase in the kinetics 00960 * class. 00961 */ 00962 vector_int m_start; 00963 00964 /** 00965 * Mapping of the phase id, i.e., the id attribute in the xml 00966 * phase element to the position of the phase within the 00967 * kinetics object. Positions start with the value of 1. The 00968 * member function, phaseIndex() decrements by one before 00969 * returning the index value, so that missing phases return 00970 * -1. 00971 */ 00972 std::map<std::string, int> m_phaseindex; 00973 //! Index of the Kinetics Manager 00974 int m_index; 00975 00976 /** 00977 * Index in the list of phases of the one surface phase. 00978 */ 00979 int m_surfphase; 00980 00981 /** 00982 * Index in the list of phases of the one phase where the reactions 00983 * occur. 00984 */ 00985 int m_rxnphase; 00986 00987 /// number of spatial dimensions of lowest-dimensional phase. 00988 int m_mindim; 00989 00990 private: 00991 00992 //! Vector of group lists 00993 std::vector<grouplist_t> m_dummygroups; 00994 00995 //! Function for unhandled situations 00996 /*! 00997 * @param m String error message 00998 */ 00999 void err(std::string m) const; 01000 01001 }; 01002 01003 //! typedef for the kinetics base class 01004 typedef Kinetics kinetics_t; 01005 01006 } 01007 01008 01009 01010 #endif