diff --git a/include/cantera/numerics/CVodesIntegrator.h b/include/cantera/numerics/CVodesIntegrator.h index 68907df7d..324db1653 100644 --- a/include/cantera/numerics/CVodesIntegrator.h +++ b/include/cantera/numerics/CVodesIntegrator.h @@ -33,19 +33,19 @@ public: explicit CVodesErr(const std::string& msg) : CanteraError("CVodesIntegrator", msg) {} }; - /** - * Wrapper class for 'cvodes' integrator from LLNL. + * Wrapper class for 'cvodes' integrator from LLNL. * * @see FuncEval.h. Classes that use CVodeInt: * ImplicitChem, ImplicitSurfChem, Reactor - * */ class CVodesIntegrator : public Integrator { - public: - + /** + * Constructor. Default settings: dense jacobian, no user-supplied + * Jacobian function, Newton iteration. + */ CVodesIntegrator(); virtual ~CVodesIntegrator(); virtual void setTolerances(double reltol, size_t n, double* abstol); @@ -87,13 +87,11 @@ public: virtual std::string getErrorInfo(int N); protected: - //! Applies user-specified options to the underlying CVODES solver. Called //! during integrator initialization or reinitialization. void applyOptions(); private: - void sensInit(double t0, FuncEval& func); size_t m_neq; diff --git a/include/cantera/numerics/DAE_Solver.h b/include/cantera/numerics/DAE_Solver.h index 92f6bc9a4..465298418 100644 --- a/include/cantera/numerics/DAE_Solver.h +++ b/include/cantera/numerics/DAE_Solver.h @@ -146,6 +146,8 @@ public: warn("inclAlgebraicInErrorTest"); } + //! Calculate consistent value of the starting solution given the starting + //! solution derivatives /** * This method may be called if the initial conditions do not * satisfy the residual equation F = 0. Given the derivatives @@ -157,12 +159,20 @@ public: warn("correctInitial_Y_given_Yp"); } + //! Calculate consistent value of the algebraic constraints and + //! derivatives at the start of the problem /** * This method may be called if the initial conditions do not * satisfy the residual equation F = 0. Given the initial * values of all differential variables, it computes the * initial values of all algebraic variables and the initial * derivatives of all differential variables. + * @param y Calculated value of the solution vector after the procedure ends + * @param yp Calculated value of the solution derivative after the procedure + * @param The first value of t at which a soluton will be + * requested (from IDASolve). (This is needed here to + * determine the direction of integration and rough scale + * in the independent variable t. */ virtual void correctInitial_YaYp_given_Yd(doublereal* y, doublereal* yp, doublereal tout) { diff --git a/include/cantera/numerics/FuncEval.h b/include/cantera/numerics/FuncEval.h index 29467acda..5a6a6c6f6 100644 --- a/include/cantera/numerics/FuncEval.h +++ b/include/cantera/numerics/FuncEval.h @@ -1,6 +1,5 @@ /** * @file FuncEval.h - * */ // Copyright 2001 California Institute of Technology @@ -12,8 +11,6 @@ namespace Cantera { - - /** * Virtual base class for ODE right-hand-side function evaluators. * Classes derived from FuncEval evaluate the right-hand-side function @@ -25,9 +22,7 @@ namespace Cantera */ class FuncEval { - public: - FuncEval() {} virtual ~FuncEval() {} @@ -47,20 +42,13 @@ public: */ virtual void getInitialConditions(double t0, size_t leny, double* y)=0; - /** - * Number of equations. - */ + //! Number of equations. virtual size_t neq()=0; //! Number of parameters. virtual size_t nparams() { return 0; } - -protected: - -private: - }; } diff --git a/include/cantera/numerics/IDA_Solver.h b/include/cantera/numerics/IDA_Solver.h index 4b4cf5b3e..28f8fdad7 100644 --- a/include/cantera/numerics/IDA_Solver.h +++ b/include/cantera/numerics/IDA_Solver.h @@ -60,17 +60,9 @@ public: virtual ~IDA_Solver(); - /** - * Set error tolerances. This version specifies a scalar - * relative tolerance, and a vector absolute tolerance. - */ virtual void setTolerances(doublereal reltol, doublereal* abstol); - /** - * Set error tolerances. This version specifies a scalar - * relative tolerance, and a scalar absolute tolerance. - */ virtual void setTolerances(doublereal reltol, doublereal abstol); virtual void setLinearSolverType(int solverType); @@ -145,24 +137,9 @@ public: */ virtual doublereal getOutputParameter(int flag) const; - //! Calculate consistent value of the starting solution given the starting solution derivatives - /*! - * This method may be called if the initial conditions do not - * satisfy the residual equation F = 0. Given the derivatives - * of all variables, this method computes the initial y - * values. - */ virtual void correctInitial_Y_given_Yp(doublereal* y, doublereal* yp, doublereal tout); - //! Calculate consistent value of the algebraic constraints and derivatives at the start of the problem - /*! - * This method may be called if the initial conditions do not - * satisfy the residual equation F = 0. Given the initial - * values of all differential variables, it computes the - * initial values of all algebraic variables and the initial - * derivatives of all differential variables. - */ virtual void correctInitial_YaYp_given_Yd(doublereal* y, doublereal* yp, doublereal tout); //! Step the system to a final value of the time @@ -201,7 +178,6 @@ public: * solver's init routine failed. In any case, the user should see * the printed error message for more details. * - * * IDA_TOO_MUCH_WORK: * The solver took mxstep internal steps but could not reach tout. * The default value for mxstep is MXSTEP_DEFAULT = 500. @@ -237,7 +213,6 @@ public: * IDA_RES_FAIL: * The user's residual function returned a nonrecoverable error * flag. - * */ virtual int solve(doublereal tout); @@ -245,15 +220,10 @@ public: virtual void init(doublereal t0); - //! the current value of solution component k. - /*! - * @param k index of the solution - */ virtual doublereal solution(int k) const; virtual const doublereal* solutionVector() const; - //! the current value of the derivative of solution component k. virtual doublereal derivative(int k) const; virtual const doublereal* derivativeVector() const; @@ -263,7 +233,6 @@ public: } protected: - //! Pointer to the IDA memory for the problem void* m_ida_mem; @@ -280,7 +249,6 @@ protected: void* m_abstol; int m_type; - int m_itol; int m_iter; doublereal m_reltol; @@ -350,4 +318,3 @@ protected: #endif #endif - diff --git a/include/cantera/numerics/Integrator.h b/include/cantera/numerics/Integrator.h index a6c810785..fac74afb6 100644 --- a/include/cantera/numerics/Integrator.h +++ b/include/cantera/numerics/Integrator.h @@ -52,9 +52,7 @@ enum IterType { */ class Integrator { - public: - //! Default Constructor Integrator() { } diff --git a/src/numerics/CVodeInt.cpp b/src/numerics/CVodeInt.cpp index a2e39bd86..0b72565d0 100644 --- a/src/numerics/CVodeInt.cpp +++ b/src/numerics/CVodeInt.cpp @@ -45,7 +45,6 @@ extern "C" { f->eval(t, ydata, ydotdata, NULL); } - /** * Function called by cvode to evaluate the Jacobian matrix. * (temporary) @@ -84,12 +83,6 @@ extern "C" { namespace Cantera { - - -/** - * Constructor. Default settings: dense jacobian, no user-supplied - * Jacobian function, Newton iteration. - */ CVodeInt::CVodeInt() : m_neq(0), m_cvode_mem(0), m_t0(0.0), @@ -111,8 +104,6 @@ CVodeInt::CVodeInt() : m_neq(0), fill(m_iopt, m_iopt+OPT_SIZE,0); } - -/// Destructor. CVodeInt::~CVodeInt() { if (m_cvode_mem) { @@ -261,7 +252,6 @@ void CVodeInt::initialize(double t0, FuncEval& func) } } - void CVodeInt::reinitialize(double t0, FuncEval& func) { m_t0 = t0; @@ -333,6 +323,3 @@ int CVodeInt::nEvals() const return m_iopt[NFE]; } } - - - diff --git a/src/numerics/CVodeInt.h b/src/numerics/CVodeInt.h index dea0a8608..c6f234c60 100644 --- a/src/numerics/CVodeInt.h +++ b/src/numerics/CVodeInt.h @@ -23,20 +23,20 @@ public: explicit CVodeErr(const std::string& msg) : CanteraError("CVodeInt", msg) {} }; - /** * Wrapper class for 'cvode' integrator from LLNL. * The unmodified cvode code is in directory ext/cvode. * * @see FuncEval.h. Classes that use CVodeInt: * ImplicitChem, ImplicitSurfChem, Reactor - * */ class CVodeInt : public Integrator { - public: - + /*! + * Constructor. Default settings: dense jacobian, no user-supplied + * Jacobian function, Newton iteration. + */ CVodeInt(); virtual ~CVodeInt(); virtual void setTolerances(double reltol, size_t n, double* abstol); @@ -63,7 +63,6 @@ public: virtual void setMaxErrTestFails(int nmax) {} private: - int m_neq; void* m_cvode_mem; double m_t0; diff --git a/src/numerics/CVodesIntegrator.cpp b/src/numerics/CVodesIntegrator.cpp index 0e7109088..f172f5013 100644 --- a/src/numerics/CVodesIntegrator.cpp +++ b/src/numerics/CVodesIntegrator.cpp @@ -1,6 +1,5 @@ /** * @file CVodesIntegrator.cpp - * */ // Copyright 2001 California Institute of Technology @@ -104,12 +103,6 @@ extern "C" { namespace Cantera { - - -/** - * Constructor. Default settings: dense jacobian, no user-supplied - * Jacobian function, Newton iteration. - */ CVodesIntegrator::CVodesIntegrator() : m_neq(0), m_cvode_mem(0), @@ -140,8 +133,6 @@ CVodesIntegrator::CVodesIntegrator() : //fill(m_iopt, m_iopt+OPT_SIZE,0); } - -/// Destructor. CVodesIntegrator::~CVodesIntegrator() { if (m_cvode_mem) { diff --git a/src/numerics/IDA_Solver.cpp b/src/numerics/IDA_Solver.cpp index 0278fdc3f..2ea1408eb 100644 --- a/src/numerics/IDA_Solver.cpp +++ b/src/numerics/IDA_Solver.cpp @@ -1,6 +1,5 @@ /** * @file IDA_Solver.cpp - * */ // Copyright 2006 California Institute of Technology @@ -58,7 +57,6 @@ public: }; } -//====================================================================================================================== extern "C" { //! Function called by IDA to evaluate the residual, given y and ydot. /*! @@ -96,8 +94,6 @@ extern "C" { //! Function called by by IDA to evaluate the Jacobian, given y and ydot. /*! - * - * * typedef int (*IDADlsDenseJacFn)(sd_size_t N, realtype t, realtype c_j, * N_Vector y, N_Vector yp, N_Vector r, * DlsMat Jac, void *user_data, @@ -125,18 +121,11 @@ extern "C" { f->evalJacobianDP(t, delta_t, c_j, ydata, ydotdata, colPts, rdata); return 0; } - - } namespace Cantera { -//==================================================================================================================== -/* - * Constructor. Default settings: dense jacobian, no user-supplied - * Jacobian function, Newton iteration. - */ IDA_Solver::IDA_Solver(ResidJacEval& f) : DAE_Solver(f), m_ida_mem(0), @@ -172,7 +161,7 @@ IDA_Solver::IDA_Solver(ResidJacEval& f) : m_mlower(0) { } -//==================================================================================================================== + IDA_Solver::~IDA_Solver() { if (m_ida_mem) { @@ -192,27 +181,26 @@ IDA_Solver::~IDA_Solver() } delete m_fdata; } -//==================================================================================================================== + doublereal IDA_Solver::solution(int k) const { return NV_Ith_S(nv(m_y),k); } -//==================================================================================================================== + const doublereal* IDA_Solver::solutionVector() const { return NV_DATA_S(nv(m_y)); } -//==================================================================================================================== + doublereal IDA_Solver::derivative(int k) const { return NV_Ith_S(nv(m_ydot),k); } -//==================================================================================================================== + const doublereal* IDA_Solver::derivativeVector() const { return NV_DATA_S(nv(m_ydot)); } -//==================================================================================================================== void IDA_Solver::setTolerances(double reltol, double* abstol) { @@ -231,7 +219,7 @@ void IDA_Solver::setTolerances(double reltol, double* abstol) } } } -//==================================================================================================================== + void IDA_Solver::setTolerances(doublereal reltol, doublereal abstol) { m_itol = IDA_SS; @@ -244,51 +232,51 @@ void IDA_Solver::setTolerances(doublereal reltol, doublereal abstol) } } } -//==================================================================================================================== + void IDA_Solver::setLinearSolverType(int solverType) { m_type = solverType; } -//==================================================================================================================== + void IDA_Solver::setDenseLinearSolver() { setLinearSolverType(0); } -//==================================================================================================================== + void IDA_Solver::setBandedLinearSolver(int m_upper, int m_lower) { m_type = 2; m_upper = m_mupper; m_mlower = m_lower; } -//==================================================================================================================== + void IDA_Solver::setMaxOrder(int n) { m_maxord = n; } -//==================================================================================================================== + void IDA_Solver::setMaxNumSteps(int n) { m_maxsteps = n; } -//==================================================================================================================== + void IDA_Solver::setInitialStepSize(doublereal h0) { m_h0 = h0; } -//==================================================================================================================== + void IDA_Solver::setStopTime(doublereal tstop) { m_tstop = tstop; } -//==================================================================================================================== + doublereal IDA_Solver::getCurrentStepFromIDA() { doublereal hcur; IDAGetCurrentStep(m_ida_mem, &hcur); return hcur; } -//==================================================================================================================== + void IDA_Solver::setJacobianType(int formJac) { m_formJac = formJac; @@ -301,22 +289,22 @@ void IDA_Solver::setJacobianType(int formJac) } } } -//==================================================================================================================== + void IDA_Solver::setMaxErrTestFailures(int maxErrTestFails) { m_maxErrTestFails = maxErrTestFails; } -//==================================================================================================================== + void IDA_Solver::setMaxNonlinIterations(int n) { m_maxNonlinIters = n; } -//==================================================================================================================== + void IDA_Solver::setMaxNonlinConvFailures(int n) { m_maxNonlinConvFails = n; } -//==================================================================================================================== + void IDA_Solver::inclAlgebraicInErrorTest(bool yesno) { if (yesno) { @@ -326,10 +314,8 @@ void IDA_Solver::inclAlgebraicInErrorTest(bool yesno) } } -//==================================================================================================================== void IDA_Solver::init(doublereal t0) { - m_t0 = t0; m_told = t0; m_told_old = t0; @@ -524,18 +510,8 @@ void IDA_Solver::init(doublereal t0) throw IDA_Err("IDASetSuppressAlg failed."); } } - - - } -//==================================================================================================================== -// Calculate consistent value of the starting solution given the starting solution derivatives -/* - * This method may be called if the initial conditions do not - * satisfy the residual equation F = 0. Given the derivatives - * of all variables, this method computes the initial y - * values. - */ + void IDA_Solver::correctInitial_Y_given_Yp(doublereal* y, doublereal* yp, doublereal tout) { int icopt = IDA_Y_INIT; @@ -566,24 +542,9 @@ void IDA_Solver::correctInitial_Y_given_Yp(doublereal* y, doublereal* yp, doubl yp[i] = yyp[i]; } } -//==================================================================================================================== -/* - * This method may be called if the initial conditions do not - * satisfy the residual equation F = 0. Given the initial - * values of all differential variables, it computes the - * initial values of all algebraic variables and the initial - * derivatives of all differential variables. - * - * @param y Calculated value of the solution vector after the procedure ends - * @param yp Calculated value of the solution derivative after the procedure - * @param The first value of t at which a soluton will be - * requested (from IDASolve). (This is needed here to - * determine the direction of integration and rough scale - * in the independent variable t. - */ + void IDA_Solver::correctInitial_YaYp_given_Yd(doublereal* y, doublereal* yp, doublereal tout) { - int icopt = IDA_YA_YDP_INIT; doublereal tout1 = tout; if (tout == 0.0) { @@ -612,7 +573,7 @@ void IDA_Solver::correctInitial_YaYp_given_Yd(doublereal* y, doublereal* yp, dou yp[i] = yyp[i]; } } -//==================================================================================================================== + int IDA_Solver::solve(double tout) { double tretn; @@ -646,7 +607,7 @@ int IDA_Solver::solve(double tout) } return flag; } -//==================================================================================================================== + double IDA_Solver::step(double tout) { double t; @@ -670,7 +631,7 @@ double IDA_Solver::step(double tout) m_deltat = m_tcurrent - m_told; return t; } -//==================================================================================================================== + doublereal IDA_Solver::getOutputParameter(int flag) const { long int lenrw, leniw; @@ -682,7 +643,6 @@ doublereal IDA_Solver::getOutputParameter(int flag) const } return 0.0; } -//==================================================================================================================== } #endif