Cantera  2.3.0
IDA_Solver.cpp
Go to the documentation of this file.
1 //! @file IDA_Solver.cpp
2 
3 // This file is part of Cantera. See License.txt in the top-level directory or
4 // at http://www.cantera.org/license.txt for license and copyright information.
5 
8 
9 #include "sundials/sundials_types.h"
10 #include "sundials/sundials_math.h"
11 #include "ida/ida.h"
12 #include "ida/ida_dense.h"
13 #include "ida/ida_spgmr.h"
14 #include "ida/ida_band.h"
15 #include "nvector/nvector_serial.h"
16 
17 using namespace std;
18 
19 #if SUNDIALS_VERSION < 25
20 typedef int sd_size_t;
21 #else
22 typedef long int sd_size_t;
23 #endif
24 
25 namespace Cantera
26 {
27 
28 //! A simple class to hold an array of parameter values and a pointer to an
29 //! instance of a subclass of ResidEval.
30 class ResidData
31 {
32 public:
33  ResidData(ResidJacEval* f, IDA_Solver* s, int npar = 0) {
34  m_func = f;
35  m_solver = s;
36  }
37 
38  virtual ~ResidData() {
39  }
40 
41  ResidJacEval* m_func;
42  IDA_Solver* m_solver;
43 };
44 }
45 
46 extern "C" {
47  //! Function called by IDA to evaluate the residual, given y and ydot.
48  /*!
49  * IDA allows passing in a void* pointer to access external data. Instead of
50  * requiring the user to provide a residual function directly to IDA (which
51  * would require using the sundials data types N_Vector, etc.), we define
52  * this function as the single function that IDA always calls. The real
53  * evaluation of the residual is done by an instance of a subclass of
54  * ResidEval, passed in to this function as a pointer in the parameters.
55  *
56  * FROM IDA WRITEUP -> What the IDA solver expects as a return flag from its
57  * residual routines:
58  *
59  * A IDAResFn res should return a value of 0 if successful, a positive value
60  * if a recoverable error occured (e.g. yy has an illegal value), or a
61  * negative value if a nonrecoverable error occured. In the latter case, the
62  * program halts. If a recoverable error occured, the integrator will
63  * attempt to correct and retry.
64  */
65  static int ida_resid(realtype t, N_Vector y, N_Vector ydot, N_Vector r, void* f_data)
66  {
68  Cantera::ResidJacEval* f = d->m_func;
69  Cantera::IDA_Solver* s = d->m_solver;
70  double delta_t = s->getCurrentStepFromIDA();
71  // TODO evaluate evalType. Assumed to be Base_ResidEval
72  int flag = f->evalResidNJ(t, delta_t, NV_DATA_S(y), NV_DATA_S(ydot),
73  NV_DATA_S(r));
74  if (flag < 0) {
75  // This signals to IDA that a nonrecoverable error has occurred.
76  return flag;
77  } else {
78  return 0;
79  }
80  }
81 
82  //! Function called by by IDA to evaluate the Jacobian, given y and ydot.
83  /*!
84  * typedef int (*IDADlsDenseJacFn)(sd_size_t N, realtype t, realtype c_j,
85  * N_Vector y, N_Vector yp, N_Vector r,
86  * DlsMat Jac, void *user_data,
87  * N_Vector tmp1, N_Vector tmp2, N_Vector tmp3);
88  *
89  * A IDADlsDenseJacFn should return
90  * - 0 if successful,
91  * - a positive int if a recoverable error occurred, or
92  * - a negative int if a nonrecoverable error occurred.
93  *
94  * In the case of a recoverable error return, the integrator will attempt to
95  * recover by reducing the stepsize (which changes cj).
96  */
97  static int ida_jacobian(sd_size_t nrows, realtype t, realtype c_j, N_Vector y, N_Vector ydot, N_Vector r,
98  DlsMat Jac, void* f_data, N_Vector tmp1, N_Vector tmp2, N_Vector tmp3)
99  {
100  Cantera::ResidData* d = (Cantera::ResidData*) f_data;
101  Cantera::ResidJacEval* f = d->m_func;
102  Cantera::IDA_Solver* s = d->m_solver;
103  double delta_t = s->getCurrentStepFromIDA();
104  f->evalJacobianDP(t, delta_t, c_j, NV_DATA_S(y), NV_DATA_S(ydot),
105  Jac->cols, NV_DATA_S(r));
106  return 0;
107  }
108 }
109 
110 namespace Cantera
111 {
112 
113 IDA_Solver::IDA_Solver(ResidJacEval& f) :
114  DAE_Solver(f),
115  m_ida_mem(0),
116  m_t0(0.0),
117  m_y(0),
118  m_ydot(0),
119  m_id(0),
120  m_constraints(0),
121  m_abstol(0),
122  m_type(0),
123  m_itol(IDA_SS),
124  m_iter(0),
125  m_reltol(1.e-9),
126  m_abstols(1.e-15),
127  m_nabs(0),
128  m_hmax(0.0),
129  m_hmin(0.0),
130  m_h0(0.0),
131  m_maxsteps(20000),
132  m_maxord(0),
133  m_formJac(0),
134  m_tstop(0.0),
135  m_told_old(0.0),
136  m_told(0.0),
137  m_tcurrent(0.0),
138  m_deltat(0.0),
139  m_maxErrTestFails(-1),
140  m_maxNonlinIters(0),
141  m_maxNonlinConvFails(-1),
142  m_setSuppressAlg(0),
143  m_mupper(0),
144  m_mlower(0)
145 {
146 }
147 
148 IDA_Solver::~IDA_Solver()
149 {
150  if (m_ida_mem) {
151  IDAFree(&m_ida_mem);
152  }
153  if (m_y) {
154  N_VDestroy_Serial(m_y);
155  }
156  if (m_ydot) {
157  N_VDestroy_Serial(m_ydot);
158  }
159  if (m_abstol) {
160  N_VDestroy_Serial(m_abstol);
161  }
162  if (m_constraints) {
163  N_VDestroy_Serial(m_constraints);
164  }
165 }
166 
167 doublereal IDA_Solver::solution(int k) const
168 {
169  return NV_Ith_S(m_y,k);
170 }
171 
172 const doublereal* IDA_Solver::solutionVector() const
173 {
174  return NV_DATA_S(m_y);
175 }
176 
177 doublereal IDA_Solver::derivative(int k) const
178 {
179  return NV_Ith_S(m_ydot,k);
180 }
181 
182 const doublereal* IDA_Solver::derivativeVector() const
183 {
184  return NV_DATA_S(m_ydot);
185 }
186 
187 void IDA_Solver::setTolerances(double reltol, double* abstol)
188 {
189  m_itol = IDA_SV;
190  if (!m_abstol) {
191  m_abstol = N_VNew_Serial(m_neq);
192  }
193  for (int i = 0; i < m_neq; i++) {
194  NV_Ith_S(m_abstol, i) = abstol[i];
195  }
196  m_reltol = reltol;
197  if (m_ida_mem) {
198  int flag = IDASVtolerances(m_ida_mem, m_reltol, m_abstol);
199  if (flag != IDA_SUCCESS) {
200  throw CanteraError("IDA_Solver::setTolerances",
201  "Memory allocation failed.");
202  }
203  }
204 }
205 
206 void IDA_Solver::setTolerances(doublereal reltol, doublereal abstol)
207 {
208  m_itol = IDA_SS;
209  m_reltol = reltol;
210  m_abstols = abstol;
211  if (m_ida_mem) {
212  int flag = IDASStolerances(m_ida_mem, m_reltol, m_abstols);
213  if (flag != IDA_SUCCESS) {
214  throw CanteraError("IDA_Solver::setTolerances",
215  "Memory allocation failed.");
216  }
217  }
218 }
219 
220 void IDA_Solver::setLinearSolverType(int solverType)
221 {
222  m_type = solverType;
223 }
224 
226 {
227  setLinearSolverType(0);
228 }
229 
230 void IDA_Solver::setBandedLinearSolver(int m_upper, int m_lower)
231 {
232  m_type = 2;
233  m_upper = m_mupper;
234  m_mlower = m_lower;
235 }
236 
237 void IDA_Solver::setMaxOrder(int n)
238 {
239  m_maxord = n;
240 }
241 
243 {
244  m_maxsteps = n;
245 }
246 
248 {
249  m_h0 = h0;
250 }
251 
252 void IDA_Solver::setStopTime(doublereal tstop)
253 {
254  m_tstop = tstop;
255 }
256 
258 {
259  doublereal hcur;
260  IDAGetCurrentStep(m_ida_mem, &hcur);
261  return hcur;
262 }
263 
265 {
266  m_formJac = formJac;
267  if (m_ida_mem && m_formJac == 1) {
268  int flag = IDADlsSetDenseJacFn(m_ida_mem, ida_jacobian);
269  if (flag != IDA_SUCCESS) {
270  throw CanteraError("IDA_Solver::setJacobianType",
271  "IDADlsSetDenseJacFn failed.");
272  }
273  }
274 }
275 
276 void IDA_Solver::setMaxErrTestFailures(int maxErrTestFails)
277 {
278  m_maxErrTestFails = maxErrTestFails;
279 }
280 
282 {
283  m_maxNonlinIters = n;
284 }
285 
287 {
289 }
290 
291 void IDA_Solver::inclAlgebraicInErrorTest(bool yesno)
292 {
293  if (yesno) {
294  m_setSuppressAlg = 0;
295  } else {
296  m_setSuppressAlg = 1;
297  }
298 }
299 
300 void IDA_Solver::init(doublereal t0)
301 {
302  m_t0 = t0;
303  m_told = t0;
304  m_told_old = t0;
305  m_tcurrent = t0;
306  if (m_y) {
307  N_VDestroy_Serial(m_y);
308  }
309  if (m_ydot) {
310  N_VDestroy_Serial(m_ydot);
311  }
312  if (m_id) {
313  N_VDestroy_Serial(m_id);
314  }
315  if (m_constraints) {
316  N_VDestroy_Serial(m_constraints);
317  }
318 
319  m_y = N_VNew_Serial(m_neq);
320  m_ydot = N_VNew_Serial(m_neq);
321  m_constraints = N_VNew_Serial(m_neq);
322 
323  for (int i=0; i<m_neq; i++) {
324  NV_Ith_S(m_y, i) = 0.0;
325  NV_Ith_S(m_ydot, i) = 0.0;
326  NV_Ith_S(m_constraints, i) = 0.0;
327  }
328 
329  // get the initial conditions
330  m_resid.getInitialConditions(m_t0, NV_DATA_S(m_y), NV_DATA_S(m_ydot));
331 
332  if (m_ida_mem) {
333  IDAFree(&m_ida_mem);
334  }
335 
336  /* Call IDACreate */
337  m_ida_mem = IDACreate();
338 
339  int flag = 0;
340  if (m_itol == IDA_SV) {
341  flag = IDAInit(m_ida_mem, ida_resid, m_t0, m_y, m_ydot);
342  if (flag != IDA_SUCCESS) {
343  if (flag == IDA_MEM_FAIL) {
344  throw CanteraError("IDA_Solver::init",
345  "Memory allocation failed.");
346  } else if (flag == IDA_ILL_INPUT) {
347  throw CanteraError("IDA_Solver::init",
348  "Illegal value for IDAInit input argument.");
349  } else {
350  throw CanteraError("IDA_Solver::init", "IDAInit failed.");
351  }
352  }
353  flag = IDASVtolerances(m_ida_mem, m_reltol, m_abstol);
354  if (flag != IDA_SUCCESS) {
355  throw CanteraError("IDA_Solver::init", "Memory allocation failed.");
356  }
357  } else {
358  flag = IDAInit(m_ida_mem, ida_resid, m_t0, m_y, m_ydot);
359  if (flag != IDA_SUCCESS) {
360  if (flag == IDA_MEM_FAIL) {
361  throw CanteraError("IDA_Solver::init",
362  "Memory allocation failed.");
363  } else if (flag == IDA_ILL_INPUT) {
364  throw CanteraError("IDA_Solver::init",
365  "Illegal value for IDAInit input argument.");
366  } else {
367  throw CanteraError("IDA_Solver::init", "IDAInit failed.");
368  }
369  }
370  flag = IDASStolerances(m_ida_mem, m_reltol, m_abstols);
371  if (flag != IDA_SUCCESS) {
372  throw CanteraError("IDA_Solver::init", "Memory allocation failed.");
373  }
374  }
375 
376  // set the linear solver type
377  if (m_type == 1 || m_type == 0) {
378  long int N = m_neq;
379  flag = IDADense(m_ida_mem, N);
380  if (flag) {
381  throw CanteraError("IDA_Solver::init", "IDADense failed");
382  }
383  } else if (m_type == 2) {
384  long int N = m_neq;
385  long int nu = m_mupper;
386  long int nl = m_mlower;
387  IDABand(m_ida_mem, N, nu, nl);
388  } else {
389  throw CanteraError("IDA_Solver::init",
390  "unsupported linear solver type");
391  }
392 
393  if (m_formJac == 1) {
394  flag = IDADlsSetDenseJacFn(m_ida_mem, ida_jacobian);
395  if (flag != IDA_SUCCESS) {
396  throw CanteraError("IDA_Solver::init",
397  "IDADlsSetDenseJacFn failed.");
398  }
399  }
400 
401  // pass a pointer to func in m_data
402  m_fdata.reset(new ResidData(&m_resid, this, m_resid.nparams()));
403  flag = IDASetUserData(m_ida_mem, m_fdata.get());
404  if (flag != IDA_SUCCESS) {
405  throw CanteraError("IDA_Solver::init", "IDASetUserData failed.");
406  }
407 
408  // set options
409  if (m_maxord > 0) {
410  flag = IDASetMaxOrd(m_ida_mem, m_maxord);
411  if (flag != IDA_SUCCESS) {
412  throw CanteraError("IDA_Solver::init", "IDASetMaxOrd failed.");
413  }
414  }
415  if (m_maxsteps > 0) {
416  flag = IDASetMaxNumSteps(m_ida_mem, m_maxsteps);
417  if (flag != IDA_SUCCESS) {
418  throw CanteraError("IDA_Solver::init", "IDASetMaxNumSteps failed.");
419  }
420  }
421  if (m_h0 > 0.0) {
422  flag = IDASetInitStep(m_ida_mem, m_h0);
423  if (flag != IDA_SUCCESS) {
424  throw CanteraError("IDA_Solver::init", "IDASetInitStep failed.");
425  }
426  }
427  if (m_tstop > 0.0) {
428  flag = IDASetStopTime(m_ida_mem, m_tstop);
429  if (flag != IDA_SUCCESS) {
430  throw CanteraError("IDA_Solver::init", "IDASetStopTime failed.");
431  }
432  }
433  if (m_maxErrTestFails >= 0) {
434  flag = IDASetMaxErrTestFails(m_ida_mem, m_maxErrTestFails);
435  if (flag != IDA_SUCCESS) {
436  throw CanteraError("IDA_Solver::init",
437  "IDASetMaxErrTestFails failed.");
438  }
439  }
440  if (m_maxNonlinIters > 0) {
441  flag = IDASetMaxNonlinIters(m_ida_mem, m_maxNonlinIters);
442  if (flag != IDA_SUCCESS) {
443  throw CanteraError("IDA_Solver::init",
444  "IDASetmaxNonlinIters failed.");
445  }
446  }
447  if (m_maxNonlinConvFails >= 0) {
448  flag = IDASetMaxConvFails(m_ida_mem, m_maxNonlinConvFails);
449  if (flag != IDA_SUCCESS) {
450  throw CanteraError("IDA_Solver::init",
451  "IDASetMaxConvFails failed.");
452  }
453  }
454  if (m_setSuppressAlg != 0) {
455  flag = IDASetSuppressAlg(m_ida_mem, m_setSuppressAlg);
456  if (flag != IDA_SUCCESS) {
457  throw CanteraError("IDA_Solver::init", "IDASetSuppressAlg failed.");
458  }
459  }
460 }
461 
462 void IDA_Solver::correctInitial_Y_given_Yp(doublereal* y, doublereal* yp, doublereal tout)
463 {
464  doublereal tout1 = tout;
465  if (tout == 0.0) {
466  double h0 = 1.0E-5;
467  if (m_h0 > 0.0) {
468  h0 = m_h0;
469  }
470  tout1 = m_t0 + h0;
471  }
472 
473  int flag = IDACalcIC(m_ida_mem, IDA_Y_INIT, tout1);
474  if (flag != IDA_SUCCESS) {
475  throw CanteraError("IDA_Solver::correctInitial_Y_given_Yp",
476  "IDACalcIC failed: error = {}", flag);
477  }
478 
479  flag = IDAGetConsistentIC(m_ida_mem, m_y, m_ydot);
480  if (flag != IDA_SUCCESS) {
481  throw CanteraError("IDA_Solver::correctInitial_Y_given_Yp",
482  "IDAGetSolution failed: error = {}", flag);
483  }
484  for (int i = 0; i < m_neq; i++) {
485  y[i] = NV_Ith_S(m_y, i);
486  yp[i] = NV_Ith_S(m_ydot, i);
487  }
488 }
489 
490 void IDA_Solver::correctInitial_YaYp_given_Yd(doublereal* y, doublereal* yp, doublereal tout)
491 {
492  int icopt = IDA_YA_YDP_INIT;
493  doublereal tout1 = tout;
494  if (tout == 0.0) {
495  double h0 = 1.0E-5;
496  if (m_h0 > 0.0) {
497  h0 = m_h0;
498  }
499  tout1 = m_t0 + h0;
500  }
501 
502  int flag = IDACalcIC(m_ida_mem, icopt, tout1);
503  if (flag != IDA_SUCCESS) {
504  throw CanteraError("IDA_Solver::correctInitial_YaYp_given_Yd",
505  "IDACalcIC failed: error = {}", flag);
506  }
507 
508  flag = IDAGetConsistentIC(m_ida_mem, m_y, m_ydot);
509  if (flag != IDA_SUCCESS) {
510  throw CanteraError("IDA_Solver::correctInitial_YaYp_given_Yd",
511  "IDAGetSolution failed: error = {}", flag);
512  }
513  for (int i = 0; i < m_neq; i++) {
514  y[i] = NV_Ith_S(m_y, i);
515  yp[i] = NV_Ith_S(m_ydot, i);
516  }
517 }
518 
519 int IDA_Solver::solve(double tout)
520 {
521  double tretn = tout - 1000;
522  int flag;
523  flag = IDASetStopTime(m_ida_mem, tout);
524  if (flag != IDA_SUCCESS) {
525  throw CanteraError("IDA_Solver::solve", "IDA error encountered.");
526  }
527  while (tretn < tout) {
528  if (tout <= m_tcurrent) {
529  throw CanteraError("IDA_Solver::solve", "tout <= tcurrent");
530  }
531  m_told_old = m_told;
532  m_told = m_tcurrent;
533  flag = IDASolve(m_ida_mem, tout, &tretn, m_y, m_ydot, IDA_ONE_STEP);
534  if (flag < 0) {
535  throw CanteraError("IDA_Solver::solve", "IDA error encountered.");
536  } else if (flag == IDA_TSTOP_RETURN) {
537  // we've reached our goal, and have actually integrated past it
538  } else if (flag == IDA_ROOT_RETURN) {
539  // not sure what to do with this yet
540  } else if (flag == IDA_WARNING) {
541  throw CanteraError("IDA_Solver::solve", "IDA Warning encountered.");
542  }
543  m_tcurrent = tretn;
545  };
546 
547  if (flag != IDA_SUCCESS && flag != IDA_TSTOP_RETURN) {
548  throw CanteraError("IDA_Solver::solve", "IDA error encountered.");
549  }
550  return flag;
551 }
552 
553 double IDA_Solver::step(double tout)
554 {
555  double t;
556  int flag;
557  if (tout <= m_tcurrent) {
558  throw CanteraError("IDA_Solver::step", "tout <= tcurrent");
559  }
560  m_told_old = m_told;
561  m_told = m_tcurrent;
562  flag = IDASolve(m_ida_mem, tout, &t, m_y, m_ydot, IDA_ONE_STEP);
563  if (flag < 0) {
564  throw CanteraError("IDA_Solver::step", "IDA error encountered.");
565  } else if (flag == IDA_TSTOP_RETURN) {
566  // we've reached our goal, and have actually integrated past it
567  } else if (flag == IDA_ROOT_RETURN) {
568  // not sure what to do with this yet
569  } else if (flag == IDA_WARNING) {
570  throw CanteraError("IDA_Solver::step", "IDA Warning encountered.");
571  }
572  m_tcurrent = t;
574  return t;
575 }
576 
577 doublereal IDA_Solver::getOutputParameter(int flag) const
578 {
579  long int lenrw, leniw;
580  switch (flag) {
581  case REAL_WORKSPACE_SIZE:
582  flag = IDAGetWorkSpace(m_ida_mem, &lenrw, &leniw);
583  return doublereal(lenrw);
584  break;
585  }
586  return 0.0;
587 }
588 
589 }
virtual void setStopTime(doublereal tstop)
Set the stop time.
Definition: IDA_Solver.cpp:252
virtual int evalJacobianDP(const doublereal t, const doublereal delta_t, doublereal cj, const doublereal *const y, const doublereal *const ydot, doublereal *const *jacobianColPts, doublereal *const resid)
Calculate an analytical Jacobian and the residual at the current time and values. ...
A simple class to hold an array of parameter values and a pointer to an instance of a subclass of Res...
Definition: IDA_Solver.cpp:30
integer m_neq
Number of total equations in the system.
Definition: DAE_Solver.h:249
Wrapper for Sundials IDA solver.
Definition: IDA_Solver.h:53
int m_maxNonlinConvFails
Maximum number of nonlinear convergence failures.
Definition: IDA_Solver.h:308
int nparams() const
Return the number of parameters in the calculation.
Definition: ResidEval.h:157
virtual void setBandedLinearSolver(int m_upper, int m_lower)
Set up the problem to use a band solver.
Definition: IDA_Solver.cpp:230
doublereal m_told
Value of the previous time.
Definition: IDA_Solver.h:290
int m_maxsteps
Maximum number of time steps allowed.
Definition: IDA_Solver.h:269
virtual double getCurrentStepFromIDA()
Get the current step size from IDA via a call.
Definition: IDA_Solver.cpp:257
virtual void setInitialStepSize(doublereal h0)
Set the initial step size.
Definition: IDA_Solver.cpp:247
virtual void setDenseLinearSolver()
Set up the problem to use a dense linear direct solver.
Definition: IDA_Solver.cpp:225
virtual int evalResidNJ(const doublereal t, const doublereal delta_t, const doublereal *const y, const doublereal *const ydot, doublereal *const resid, const ResidEval_Type_Enum evalType=Base_ResidEval, const int id_x=-1, const doublereal delta_x=0.0)
Evaluate the residual function.
int m_maxErrTestFails
maximum number of error test failures
Definition: IDA_Solver.h:299
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:490
doublereal m_tstop
maximum time
Definition: IDA_Solver.h:284
STL namespace.
doublereal m_tcurrent
Value of the current time.
Definition: IDA_Solver.h:293
virtual int getInitialConditions(const doublereal t0, doublereal *const y, doublereal *const ydot)
Fill in the initial conditions.
virtual void setMaxNonlinConvFailures(int n)
Set the maximum number of nonlinear solver convergence failures.
Definition: IDA_Solver.cpp:286
doublereal m_told_old
Value of the previous, previous time.
Definition: IDA_Solver.h:287
static int ida_resid(realtype t, N_Vector y, N_Vector ydot, N_Vector r, void *f_data)
Function called by IDA to evaluate the residual, given y and ydot.
Definition: IDA_Solver.cpp:65
Wrappers for the function evaluators for Nonlinear solvers and Time steppers.
Definition: ResidJacEval.h:55
void * m_ida_mem
Pointer to the IDA memory for the problem.
Definition: IDA_Solver.h:238
virtual void setMaxNumSteps(int n)
Set the maximum number of time steps.
Definition: IDA_Solver.cpp:242
N_Vector m_y
Current value of the solution vector.
Definition: IDA_Solver.h:244
doublereal m_h0
Value of the initial time step.
Definition: IDA_Solver.h:266
virtual int solve(doublereal tout)
Step the system to a final value of the time.
Definition: IDA_Solver.cpp:519
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:462
doublereal m_deltat
Value of deltaT for the current step.
Definition: IDA_Solver.h:296
virtual doublereal solution(int k) const
the current value of solution component k.
Definition: IDA_Solver.cpp:167
Base class for exceptions thrown by Cantera classes.
Definition: ctexceptions.h:65
virtual doublereal derivative(int k) const
the current value of the derivative of solution component k.
Definition: IDA_Solver.cpp:177
virtual void setMaxNonlinIterations(int n)
Set the maximum number of nonlinear iterations on a timestep.
Definition: IDA_Solver.cpp:281
virtual void setJacobianType(int formJac)
Set the form of the Jacobian.
Definition: IDA_Solver.cpp:264
virtual doublereal getOutputParameter(int flag) const
Get the value of a solver-specific output parameter.
Definition: IDA_Solver.cpp:577
virtual void init(doublereal t0)
initialize.
Definition: IDA_Solver.cpp:300
Contains declarations for string manipulation functions within Cantera.
int m_maxNonlinIters
Maximum number of nonlinear solver iterations at one solution.
Definition: IDA_Solver.h:305
doublereal m_t0
Initial value of the time.
Definition: IDA_Solver.h:241
N_Vector m_ydot
Current value of the derivative of the solution vector.
Definition: IDA_Solver.h:247
Header file for class IDA_Solver.
virtual doublereal step(doublereal tout)
Take one internal step.
Definition: IDA_Solver.cpp:553
Namespace for the Cantera kernel.
Definition: application.cpp:29
Wrapper for DAE solvers.
Definition: DAE_Solver.h:75
int m_maxord
maximum time step order of the method
Definition: IDA_Solver.h:272
virtual void setTolerances(doublereal reltol, doublereal *abstol)
Set error tolerances.
Definition: IDA_Solver.cpp:187
static int ida_jacobian(sd_size_t nrows, realtype t, realtype c_j, N_Vector y, N_Vector ydot, N_Vector r, DlsMat Jac, void *f_data, N_Vector tmp1, N_Vector tmp2, N_Vector tmp3)
Function called by by IDA to evaluate the Jacobian, given y and ydot.
Definition: IDA_Solver.cpp:97
int m_formJac
Form of the Jacobian.
Definition: IDA_Solver.h:281
int m_setSuppressAlg
If true, the algebraic variables don&#39;t contribute to error tolerances.
Definition: IDA_Solver.h:311