Cantera  2.5.1
IDA_Solver.h
Go to the documentation of this file.
1 /**
2  * @file IDA_Solver.h
3  * Header file for class IDA_Solver
4  */
5 
6 // This file is part of Cantera. See License.txt in the top-level directory or
7 // at https://cantera.org/license.txt for license and copyright information.
8 
9 #ifndef CT_IDA_SOLVER_H
10 #define CT_IDA_SOLVER_H
11 
12 #include "DAE_Solver.h"
13 
14 #include "sundials/sundials_nvector.h"
15 
16 // These constants are defined internally in the IDA package, ida.c
17 #define IDA_NN 0
18 #define IDA_SS 1
19 #define IDA_SV 2
20 #define IDA_WF 3
21 
22 #define REAL_WORKSPACE_SIZE 0
23 
24 namespace Cantera
25 {
26 
27 class ResidData;
28 
29 /**
30  * Wrapper for Sundials IDA solver
31  *
32  * @attention This class currently does not have any test cases or examples. Its
33  * implementation may be incomplete, and future changes to Cantera may
34  * unexpectedly cause this class to stop working. If you use this class,
35  * please consider contributing examples or test cases. In the absence of
36  * new tests or examples, this class may be deprecated and removed in a
37  * future version of Cantera. See
38  * https://github.com/Cantera/cantera/issues/267 for additional information.
39  */
40 class IDA_Solver : public DAE_Solver
41 {
42 public:
43  //! Constructor.
44  /*!
45  * Default settings: dense Jacobian, no user-supplied Jacobian function,
46  * Newton iteration.
47  *
48  * @param f Function that will supply the time dependent residual to be solved
49  */
51 
52  virtual ~IDA_Solver();
53 
54  virtual void setTolerances(doublereal reltol,
55  doublereal* abstol);
56 
57  virtual void setTolerances(doublereal reltol, doublereal abstol);
58 
59  virtual void setLinearSolverType(int solverType);
60 
61  //! Set up the problem to use a dense linear direct solver
62  virtual void setDenseLinearSolver();
63 
64  //! Set up the problem to use a band solver
65  /*!
66  * @param m_upper upper band width of the matrix
67  * @param m_lower lower band width of the matrix
68  */
69  virtual void setBandedLinearSolver(int m_upper, int m_lower);
70 
71  virtual void setMaxOrder(int n);
72 
73  //! Set the maximum number of time steps
74  /*!
75  * @param n input of maximum number of time steps
76  */
77  virtual void setMaxNumSteps(int n);
78 
79  //! Set the initial step size
80  /*!
81  * @param h0 initial step size value
82  */
83  virtual void setInitialStepSize(doublereal h0);
84 
85  //! Set the stop time
86  /*!
87  * @param tstop the independent variable value past which the solution is
88  * not to proceed.
89  */
90  virtual void setStopTime(doublereal tstop);
91 
92  //! Get the current step size from IDA via a call
93  /*!
94  * @returns the current step size.
95  */
96  virtual double getCurrentStepFromIDA();
97 
98  //! Set the form of the Jacobian
99  /*!
100  * @param formJac Form of the Jacobian
101  * 0 numerical Jacobian
102  * 1 analytical Jacobian given by the evalJacobianDP() function
103  */
104  virtual void setJacobianType(int formJac);
105 
106 
107  virtual void setMaxErrTestFailures(int n);
108 
109  //! Set the maximum number of nonlinear iterations on a timestep
110  /*!
111  * @param n Set the max iterations. The default is 4, which seems awfully low to me.
112  */
113  virtual void setMaxNonlinIterations(int n);
114 
115  //! Set the maximum number of nonlinear solver convergence failures
116  /*!
117  * @param n Value of nonlin failures. If value is exceeded, the calculation terminates.
118  */
119  virtual void setMaxNonlinConvFailures(int n);
120 
121 
122  virtual void inclAlgebraicInErrorTest(bool yesno);
123 
124  /**
125  * Get the value of a solver-specific output parameter.
126  */
127  virtual doublereal getOutputParameter(int flag) const;
128 
129  virtual void correctInitial_Y_given_Yp(doublereal* y, doublereal* yp,
130  doublereal tout);
131 
132  virtual void correctInitial_YaYp_given_Yd(doublereal* y, doublereal* yp, doublereal tout);
133 
134  //! Step the system to a final value of the time
135  /*!
136  * @param tout Final value of the time
137  * @returns the IDASolve() return flag
138  *
139  * The return values for IDASolve are described below. (The numerical return
140  * values are defined above in this file.) All unsuccessful returns give a
141  * negative return value.
142  *
143  * IDA_SUCCESS
144  * IDASolve succeeded and no roots were found.
145  *
146  * IDA_ROOT_RETURN: IDASolve succeeded, and found one or more roots.
147  * If nrtfn > 1, call IDAGetRootInfo to see which g_i were found
148  * to have a root at (*tret).
149  *
150  * IDA_TSTOP_RETURN:
151  * IDASolve returns computed results for the independent variable
152  * value tstop. That is, tstop was reached.
153  *
154  * IDA_MEM_NULL:
155  * The IDA_mem argument was NULL.
156  *
157  * IDA_ILL_INPUT:
158  * One of the inputs to IDASolve is illegal. This includes the
159  * situation when a component of the error weight vectors
160  * becomes < 0 during internal stepping. It also includes the
161  * situation where a root of one of the root functions was found
162  * both at t0 and very near t0. The ILL_INPUT flag
163  * will also be returned if the linear solver function IDA---
164  * (called by the user after calling IDACreate) failed to set one
165  * of the linear solver-related fields in ida_mem or if the linear
166  * solver's init routine failed. In any case, the user should see
167  * the printed error message for more details.
168  *
169  * IDA_TOO_MUCH_WORK:
170  * The solver took mxstep internal steps but could not reach tout.
171  * The default value for mxstep is MXSTEP_DEFAULT = 500.
172  *
173  * IDA_TOO_MUCH_ACC:
174  * The solver could not satisfy the accuracy demanded by the user
175  * for some internal step.
176  *
177  * IDA_ERR_FAIL:
178  * Error test failures occurred too many times (=MXETF = 10) during
179  * one internal step.
180  *
181  * IDA_CONV_FAIL:
182  * Convergence test failures occurred too many times (= MXNCF = 10)
183  * during one internal step.
184  *
185  * IDA_LSETUP_FAIL:
186  * The linear solver's setup routine failed
187  * in an unrecoverable manner.
188  *
189  * IDA_LSOLVE_FAIL:
190  * The linear solver's solve routine failed
191  * in an unrecoverable manner.
192  *
193  * IDA_CONSTR_FAIL:
194  * The inequality constraints were violated,
195  * and the solver was unable to recover.
196  *
197  * IDA_REP_RES_ERR:
198  * The user's residual function repeatedly returned a recoverable
199  * error flag, but the solver was unable to recover.
200  *
201  * IDA_RES_FAIL:
202  * The user's residual function returned a nonrecoverable error
203  * flag.
204  */
205  virtual int solve(doublereal tout);
206 
207  virtual doublereal step(doublereal tout);
208 
209  virtual void init(doublereal t0);
210 
211  virtual doublereal solution(int k) const;
212 
213  virtual const doublereal* solutionVector() const;
214 
215  virtual doublereal derivative(int k) const;
216 
217  virtual const doublereal* derivativeVector() const;
218 
219  void* IDAMemory() {
220  return m_ida_mem;
221  }
222 
223 protected:
224  //! Pointer to the IDA memory for the problem
225  void* m_ida_mem;
226  void* m_linsol; //!< Sundials linear solver object
227  void* m_linsol_matrix; //!< matrix used by Sundials
228 
229  //! Initial value of the time
230  doublereal m_t0;
231 
232  //! Current value of the solution vector
233  N_Vector m_y;
234 
235  //! Current value of the derivative of the solution vector
236  N_Vector m_ydot;
237  N_Vector m_id;
238  N_Vector m_constraints;
239  N_Vector m_abstol;
240  int m_type;
241 
242  int m_itol;
243  int m_iter;
244  doublereal m_reltol;
245  doublereal m_abstols;
246  int m_nabs;
247 
248  //! Maximum value of the timestep allowed
249  doublereal m_hmax;
250 
251  //! Minimum value of the timestep allowed
252  doublereal m_hmin;
253 
254  //! Value of the initial time step
255  doublereal m_h0;
256 
257  //! Maximum number of time steps allowed
259 
260  //! maximum time step order of the method
261  int m_maxord;
262 
263  //! Form of the Jacobian
264  /*!
265  * 0 numerical Jacobian created by IDA
266  * 1 analytical Jacobian. Must have populated the evalJacobianDP()
267  * function in the ResidJacEval class.
268  * 2 numerical Jacobian formed by the ResidJacEval class (unimplemented)
269  */
271 
272  //! maximum time
273  doublereal m_tstop;
274 
275  //! Value of the previous, previous time
276  doublereal m_told_old;
277 
278  //! Value of the previous time
279  doublereal m_told;
280 
281  //! Value of the current time
282  doublereal m_tcurrent;
283 
284  //! Value of deltaT for the current step
285  doublereal m_deltat;
286 
287  //! maximum number of error test failures
289 
290  //! Maximum number of nonlinear solver iterations at one solution
291  /*!
292  * If zero, this is the default of 4.
293  */
295 
296  //! Maximum number of nonlinear convergence failures
298 
299  //! If true, the algebraic variables don't contribute to error tolerances
301 
302  std::unique_ptr<ResidData> m_fdata;
303  int m_mupper;
304  int m_mlower;
305 };
306 
307 }
308 
309 #endif
Header file for class DAE_Solver.
Wrapper for DAE solvers.
Definition: DAE_Solver.h:76
Wrapper for Sundials IDA solver.
Definition: IDA_Solver.h:41
virtual doublereal step(doublereal tout)
Take one internal step.
Definition: IDA_Solver.cpp:638
virtual void setStopTime(doublereal tstop)
Set the stop time.
Definition: IDA_Solver.cpp:289
virtual void correctInitial_Y_given_Yp(doublereal *y, doublereal *yp, doublereal tout)
Calculate consistent value of the starting solution given the starting solution derivatives.
Definition: IDA_Solver.cpp:547
virtual void setMaxNumSteps(int n)
Set the maximum number of time steps.
Definition: IDA_Solver.cpp:279
virtual void setBandedLinearSolver(int m_upper, int m_lower)
Set up the problem to use a band solver.
Definition: IDA_Solver.cpp:267
virtual doublereal solution(int k) const
the current value of solution component k.
Definition: IDA_Solver.cpp:204
N_Vector m_y
Current value of the solution vector.
Definition: IDA_Solver.h:233
virtual void setJacobianType(int formJac)
Set the form of the Jacobian.
Definition: IDA_Solver.cpp:301
virtual void correctInitial_YaYp_given_Yd(doublereal *y, doublereal *yp, doublereal tout)
Calculate consistent value of the algebraic constraints and derivatives at the start of the problem.
Definition: IDA_Solver.cpp:575
void * m_linsol
Sundials linear solver object.
Definition: IDA_Solver.h:226
int m_maxNonlinConvFails
Maximum number of nonlinear convergence failures.
Definition: IDA_Solver.h:297
virtual void setInitialStepSize(doublereal h0)
Set the initial step size.
Definition: IDA_Solver.cpp:284
virtual int solve(doublereal tout)
Step the system to a final value of the time.
Definition: IDA_Solver.cpp:604
int m_maxNonlinIters
Maximum number of nonlinear solver iterations at one solution.
Definition: IDA_Solver.h:294
IDA_Solver(ResidJacEval &f)
Constructor.
Definition: IDA_Solver.cpp:150
doublereal m_told
Value of the previous time.
Definition: IDA_Solver.h:279
virtual void setDenseLinearSolver()
Set up the problem to use a dense linear direct solver.
Definition: IDA_Solver.cpp:262
doublereal m_deltat
Value of deltaT for the current step.
Definition: IDA_Solver.h:285
int m_maxErrTestFails
maximum number of error test failures
Definition: IDA_Solver.h:288
virtual void setTolerances(doublereal reltol, doublereal *abstol)
Set error tolerances.
Definition: IDA_Solver.cpp:224
virtual void setMaxNonlinIterations(int n)
Set the maximum number of nonlinear iterations on a timestep.
Definition: IDA_Solver.cpp:322
doublereal m_tcurrent
Value of the current time.
Definition: IDA_Solver.h:282
int m_maxord
maximum time step order of the method
Definition: IDA_Solver.h:261
virtual doublereal getOutputParameter(int flag) const
Get the value of a solver-specific output parameter.
Definition: IDA_Solver.cpp:662
virtual void init(doublereal t0)
initialize.
Definition: IDA_Solver.cpp:341
int m_formJac
Form of the Jacobian.
Definition: IDA_Solver.h:270
doublereal m_hmin
Minimum value of the timestep allowed.
Definition: IDA_Solver.h:252
doublereal m_hmax
Maximum value of the timestep allowed.
Definition: IDA_Solver.h:249
doublereal m_told_old
Value of the previous, previous time.
Definition: IDA_Solver.h:276
virtual void setMaxNonlinConvFailures(int n)
Set the maximum number of nonlinear solver convergence failures.
Definition: IDA_Solver.cpp:327
virtual double getCurrentStepFromIDA()
Get the current step size from IDA via a call.
Definition: IDA_Solver.cpp:294
virtual doublereal derivative(int k) const
the current value of the derivative of solution component k.
Definition: IDA_Solver.cpp:214
doublereal m_h0
Value of the initial time step.
Definition: IDA_Solver.h:255
void * m_linsol_matrix
matrix used by Sundials
Definition: IDA_Solver.h:227
N_Vector m_ydot
Current value of the derivative of the solution vector.
Definition: IDA_Solver.h:236
doublereal m_t0
Initial value of the time.
Definition: IDA_Solver.h:230
int m_maxsteps
Maximum number of time steps allowed.
Definition: IDA_Solver.h:258
doublereal m_tstop
maximum time
Definition: IDA_Solver.h:273
int m_setSuppressAlg
If true, the algebraic variables don't contribute to error tolerances.
Definition: IDA_Solver.h:300
void * m_ida_mem
Pointer to the IDA memory for the problem.
Definition: IDA_Solver.h:225
Wrappers for the function evaluators for Nonlinear solvers and Time steppers.
Definition: ResidJacEval.h:56
Namespace for the Cantera kernel.
Definition: AnyMap.cpp:264