diff --git a/Cantera/src/kinetics/ImplicitSurfChem.cpp b/Cantera/src/kinetics/ImplicitSurfChem.cpp index 82991afc8..2c2122b0c 100755 --- a/Cantera/src/kinetics/ImplicitSurfChem.cpp +++ b/Cantera/src/kinetics/ImplicitSurfChem.cpp @@ -1,8 +1,11 @@ /** * @file ImplicitSurfChem.cpp - * - * Implicit integration of surface site density equations - * + * Definitions for the implicit integration of surface site density equations + * (see \ref kineticsmgr and class + * \link Cantera::ImplicitSurfChem ImplicitSurfChem\endlink). + */ + +/* * $Author$ * $Revision$ * $Date$ @@ -23,6 +26,9 @@ using namespace std; namespace Cantera { + + + ImplicitSurfChem::ImplicitSurfChem(vector k) : FuncEval(), m_nv(0), m_integ(0), m_atol(1.e-14), m_rtol(1.e-7), m_maxstep(0.0) @@ -54,6 +60,12 @@ namespace Cantera { m_work.resize(ntmax); } + /** + * Destructor. Deletes the integrator. + */ + ImplicitSurfChem::~ImplicitSurfChem(){ + delete m_integ; + } // overloaded method of FuncEval. Called by the integrator to // get the initial conditions. @@ -76,6 +88,34 @@ namespace Cantera { m_integ->initialize(t0, *this); } + // Integrate from t0 to t1. The integrator is reinitialized first. + /* + * This routine does a time accurate solve from t = t0 to t = t1. + * of the surface problem. + * + * @param t0 Initial Time -> this is an input + * @param t1 Final Time -> This is an input + */ + void ImplicitSurfChem::integrate(doublereal t0, doublereal t1) { + m_integ->initialize(t0, *this); + m_integ->setMaxStepSize(t1 - t0); + m_integ->integrate(t1); + updateState(m_integ->solution()); + } + + // Integrate from t0 to t1 without reinitializing the integrator. + /* + * Use when the coverages have not changed from + * their values on return from the last call to integrate or + * integrate0. + * + * @param t0 Initial Time -> this is an input + * @param t1 Final Time -> This is an input + */ + void ImplicitSurfChem::integrate0(doublereal t0, doublereal t1) { + m_integ->integrate(t1); + updateState(m_integ->solution()); + } void ImplicitSurfChem::updateState(doublereal* c) { int loc = 0; @@ -86,29 +126,29 @@ namespace Cantera { } - /** - * Called by the integrator to evaluate ydot given y at time 'time'. - */ - void ImplicitSurfChem::eval(doublereal time, doublereal* y, - doublereal* ydot, doublereal* p) - { - int n; - updateState(y); // synchronize the surface state(s) with y - doublereal rs0, sum; - int loc, k, kstart; - for (n = 0; n < m_nsurf; n++) { - rs0 = 1.0/m_surf[n]->siteDensity(); - m_kin[n]->getNetProductionRates(DATA_PTR(m_work)); - kstart = m_kin[n]->kineticsSpeciesIndex(0,m_surfindex[n]); - sum = 0.0; - loc = 0; - for (k = 1; k < m_nsp[n]; k++) { - ydot[k + loc] = m_work[kstart + k] * rs0 * m_surf[n]->size(k); - sum -= ydot[k]; - } - ydot[loc] = sum; - loc += m_nsp[n]; - } + /** + * Called by the integrator to evaluate ydot given y at time 'time'. + */ + void ImplicitSurfChem::eval(doublereal time, doublereal* y, + doublereal* ydot, doublereal* p) + { + int n; + updateState(y); // synchronize the surface state(s) with y + doublereal rs0, sum; + int loc, k, kstart; + for (n = 0; n < m_nsurf; n++) { + rs0 = 1.0/m_surf[n]->siteDensity(); + m_kin[n]->getNetProductionRates(DATA_PTR(m_work)); + kstart = m_kin[n]->kineticsSpeciesIndex(0,m_surfindex[n]); + sum = 0.0; + loc = 0; + for (k = 1; k < m_nsp[n]; k++) { + ydot[k + loc] = m_work[kstart + k] * rs0 * m_surf[n]->size(k); + sum -= ydot[k]; + } + ydot[loc] = sum; + loc += m_nsp[n]; } + } } diff --git a/Cantera/src/kinetics/ImplicitSurfChem.h b/Cantera/src/kinetics/ImplicitSurfChem.h index b8f4fa417..bcf834071 100755 --- a/Cantera/src/kinetics/ImplicitSurfChem.h +++ b/Cantera/src/kinetics/ImplicitSurfChem.h @@ -1,8 +1,11 @@ /** * @file ImplicitSurfChem.h - * - * Implicit integration of surface site density equations. - * + * Declarations for the implicit integration of surface site density equations + * (see \ref kineticsmgr and class + * \link Cantera::ImplicitSurfChem ImplicitSurfChem\endlink). + */ + +/* * $Author$ * $Revision$ * $Date$ @@ -26,91 +29,153 @@ namespace Cantera { - /** - * Advances the surface coverages of an associated SurfacePhase - * object in time by implicitly integrating \f[ \dot \theta_k = - * \dot s_k (\sigma_k / s_0)\f] + + //! Advances the surface coverages of the associated set of SurfacePhase + //! objects in time + /*! + * This function advances a set of SurfacePhase objects, each + * associated with one InterfaceKinetics object, in time. + * The following equation is used for each surface phase, i. + * + * \f[ + * \dot \theta_k = \dot s_k (\sigma_k / s_0) + * \f] + * + * In this equation, + * \f$ \theta_k \f$ is the site coverage for the kth species. + * \f$ \dot s_k \f$ is the source term for the kth species + * \f$ \sigma_k \f$ is the number of surface sites covered by + * each species k. + * \f$ s_0 \f$ is the total site density of the interfacial phase. + * + * Additionally, the 0'th equation in the set is discarded. Instead the + * alternate equation is solved for + * + * \f[ + * \sum_{k=0}^{N-1} \dot \theta_k = 0 + * \f] + * + * This last equation serves to ensure that sum of the \f$ \theta_k \f$ + * values stays constant. + * + * The object uses the CVODE software to advance the surface equations. + * + * The solution vector used by this object is as follows. + * For each surface phase with \f$ N_s \f$ surface sites, + * it consists of the surface coverages + * \f$ \theta_k \f$ for \f$ k = 0, N_s - 1 \f$ + * + * @ingroup kineticsmgr + * + */ + class ImplicitSurfChem : public FuncEval { + + public: + + + //! Constructor for multiple surfaces. + /*! + * @param k Vector of pointers to InterfaceKinetics objects + * Each object consists of a surface or an edge containing + * internal degrees of freedom representing the concentration + * of surface adsorbates. */ - class ImplicitSurfChem : public FuncEval { + ImplicitSurfChem(std::vector k); - public: + /** + * Destructor. Deletes the integrator. + */ + virtual ~ImplicitSurfChem(); - /** - * Constructor. - */ - //ImplicitSurfChem(InterfaceKinetics& kin); + /** + * Overloads the virtual function + * declared in FuncEval. + */ + virtual void initialize(doublereal t0 = 0.0); - /** - * Constructor for multiple surfaces. - */ - ImplicitSurfChem(std::vector k); - /** - * Destructor. Deletes the integrator. - */ - virtual ~ImplicitSurfChem(){ delete m_integ; } + //! Integrate from t0 to t1. The integrator is reinitialized first. + /*! + * This routine does a time accurate solve from t = t0 to t = t1. + * of the surface problem. + * + * @param t0 Initial Time -> this is an input + * @param t1 Final Time -> This is an input + */ + void integrate(doublereal t0, doublereal t1); + + + //! Integrate from t0 to t1 without reinitializing the integrator. + /*! + * Use when the coverages have not changed from + * their values on return from the last call to integrate or + * integrate0. + * + * @param t0 Initial Time -> this is an input + * @param t1 Final Time -> This is an input + */ + void integrate0(doublereal t0, doublereal t1); + + + // overloaded methods of class FuncEval + + //! Return the number of equations + virtual int neq() { return m_nv; } + + //! Evaluate the value of ydot[k] at the current conditions + /*! + * @param t Time (seconds) + * @param y Vector containing the current solution vector + * @param ydot Output vector containing the value of the + * derivative of the surface coverages. + * @param p Unused parameter pass-through parameter vector + */ + virtual void eval(doublereal t, doublereal* y, doublereal* ydot, + doublereal* p); + + //! Set the initial conditions for the solution vector + /*! + * @param t0 Initial time + * @param leny Length of the solution vector + * @param y Value of the solution vector to be used. + * On output, this contains the initial value + * of the solution. + */ + virtual void getInitialConditions(doublereal t0, + size_t leny, doublereal* y); + + + protected: + + //! Set the mixture to a state consistent with solution + //! vector y. + /*! + * This function will set the surface site factions + * in the underlying %SurfPhase objects to the current + * value of the solution vector. + * + * @param y Current value of the solution vector. + * The lenth is equal to the sum of the number of surface + * sites in all the surface phases + */ + void updateState(doublereal* y); - /** - * Overloads the virtual function - * declared in FuncEval. - */ - virtual void initialize(doublereal t0 = 0.0); + std::vector m_surf; + std::vector m_kin; + vector_int m_nsp; + vector_int m_surfindex; + int m_nsurf; + int m_nv; + //int m_nsp, m_surfindex; + Integrator* m_integ; // pointer to integrator + doublereal m_atol, m_rtol; // tolerances + doublereal m_maxstep; // max step size + vector_fp m_work; + private: - /** - * Integrate from t0 to t1. The integrator is reinitialized - * first. - */ - void integrate(doublereal t0, doublereal t1) { - m_integ->initialize(t0, *this); - m_integ->setMaxStepSize(t1 - t0); - m_integ->integrate(t1); - updateState(m_integ->solution()); - } - - /** - * Integrate from t0 to t1 without reinitializing the - * integrator. Use when the coverages have not changed from - * their values on return from the last call to integrate or - * integrate0. - */ - void integrate0(doublereal t0, doublereal t1) { - m_integ->integrate(t1); - updateState(m_integ->solution()); - } - - // overloaded methods of class FuncEval - virtual int neq() { return m_nv; } - virtual void eval(doublereal t, doublereal* y, doublereal* ydot, - doublereal* p); - virtual void getInitialConditions(doublereal t0, - size_t leny, doublereal* y); - - - protected: - - /** - * Set the mixture to a state consistent with solution - * vector y. - */ - void updateState(doublereal* y); - - std::vector m_surf; - std::vector m_kin; - vector_int m_nsp; - vector_int m_surfindex; - int m_nsurf; - int m_nv; - //int m_nsp, m_surfindex; - Integrator* m_integ; // pointer to integrator - doublereal m_atol, m_rtol; // tolerances - doublereal m_maxstep; // max step size - vector_fp m_work; - - private: - - }; + }; } #endif