1 /*
2  * -----------------------------------------------------------------
3  * $Revision: 1.6 $
4  * $Date: 2007/03/20 14:33:17 $
5  * -----------------------------------------------------------------
6  * Programmer: Radu Serban @ LLNL
7  * -----------------------------------------------------------------
8  * Copyright (c) 2006, The Regents of the University of California.
9  * Produced at the Lawrence Livermore National Laboratory.
10  * All rights reserved.
11  * For details, see the LICENSE file.
12  * -----------------------------------------------------------------
13  * This is the interface file for the main CPODES integrator.
14  * -----------------------------------------------------------------
15  *
16  * CPODES is used to solve numerically the ordinary initial value
17  * problem with invariants:
18  *
19  *    y' = f(t,y)   or     F(t,y,y') = 0
20  *    c(t,y) = 0           c(t,y) = 0
21  *    y(t0) = y0           y(t0) = y0; y'(t0) = yp0
22  *
23  * where t0, y0, yp0 in R^N, f: R x R^N -> R^N, F: R x R^N x R^N -> R^N
24  * and c: R x R^N -> R^M.
25  *
26  * -----------------------------------------------------------------
27  */
28 
29 #ifndef _CPODES_H
30 #define _CPODES_H
31 
32 #ifdef __cplusplus  /* wrapper to enable C++ usage */
33 extern "C" {
34 #endif
35 
36 #include <stdio.h>
37 #include <sundials/sundials_nvector.h>
38 
39 /*
40  * =================================================================
41  *              C P O D E     C O N S T A N T S
42  * =================================================================
43  */
44 
45 /*
46  * -----------------------------------------------------------------
47  * Inputs to CPodeCreate, CPodeInit, CPodeReInit, and CPode.
48  * -----------------------------------------------------------------
49  * Symbolic constants for the lmm_type and nls_type inputs to
50  * CPodeCreate, the ode_type and tol_type inputs to CPodeInit and
51  * CPodeReInit, the cnstr_type input to CPodeProjInit, as well as
52  * the input mode to CPode, are given below.
53  *
54  * lmm_type: The user of the CPODES package specifies whether to use
55  *    the CP_ADAMS (Adams-Moulton) or CP_BDF (Backward Differentiation
56  *    Formula) linear multistep method. The BDF method is recommended
57  *    for stiff problems, and the CP_ADAMS method is recommended for
58  *    nonstiff problems.
59  *
60  * nls_type: At each internal time step, a nonlinear equation must
61  *    be solved. The user can specify either CP_FUNCTIONAL iteration,
62  *    which does not require linear algebra, or a CP_NEWTON iteration,
63  *    which requires the solution of linear systems. In the CP_NEWTON
64  *    case, the user also specifies a CPODES linear solver.
65  *    CP_NEWTON is recommended in case of stiff problems.
66  *
67  * ode_type: The ODE system can be given either in explicit form,
68  *    y' = f(t,y) (in which case ode = CP_EXPL) or in implicit form
69  *    (in which case ode = CP_IMPL).
70  *
71  * tol_type: This parameter specifies the relative and absolute
72  *    tolerance types to be used. The CP_SS tolerance type means a
73  *    scalar relative and absolute tolerance. The CP_SV tolerance type
74  *    means a scalar relative tolerance and a vector absolute tolerance
75  *    (a potentially different absolute tolerance for each vector
76  *    component). The CP_WF tolerance type means that the user provides
77  *    a function (of type CPEwtFn) to set the error weight vector.
78  *
79  * proj_type: If performing projection on the invariant manifold,
80  *    this parameter specifies whether to use the internal algorithm
81  *    or a user-provided projection function. The valid values are
82  *    CP_PROJ_USER and CP_PROJ_INERNAL.
83  *    A value proj_type = CP_PROJ_USER indicates that a function to
84  *    perform the projection is supplied by the user.
85  *    A value proj_type = CP_PROJ_INTERNAL indicates that the internal
86  *    CPODES projection algorithm is to be used. In this case, the
87  *    user must specify the constraint function and a linear solver
88  *    to be used.
89  *
90  * proj_norm: Type of the norm in which projection is to be performed.
91  *    The valid values are CP_PROJ_L2NORM (in which case projection
92  *    is done in Euclidian norm) and CP_PROJ_ERRNORM (in which case
93  *    the projection is done in WRMS norm).
94  *
95  * cnstr_type: If the internal projection algorithm is used, then
96  *    cnstr_type specifies the type of constraints.
97  *    If the constraints are linear (cnstr_type = CP_CNSTR_LIN),
98  *    then an improved error estimate for the projected method
99  *    (obtained by projecting the error estimate for the original
100  *    method) is available at no additional cost and will always be used.
101  *    If the constraints are nonlinear (cnstr_type = CP_CNSTR_NONLIN),
102  *    obtaining the global projection operator requires an additional
103  *    evaluation of the constraint Jacobian and constructing its
104  *    More-Penrose pseudo-inverse. In this case, the default is
105  *    to use the error estimate of the unprojected method (unless
106  *    indicated otherwise by the user).
107  *
108  * mode: The mode input parameter to CPode indicates the job of the
109  *    solver for the next user step. The CP_NORMAL mode is to have
110  *    the solver take internal steps until it has reached or just
111  *    passed the user specified tout parameter. The solver then
112  *    interpolates in order to return an approximate value of y(tout).
113  *    The CP_ONE_STEP option tells the solver to just take one internal
114  *    step and return the solution at the point reached by that step.
115  *    The CP_NORMAL_TSTOP and CP_ONE_STEP_TSTOP modes are similar to
116  *    CP_NORMAL and CP_ONE_STEP, respectively, except that the
117  *    integration never proceeds past the value tstop (specified
118  *    through the routine CPodeSetStopTime).
119  * -----------------------------------------------------------------
120  */
121 
122 /* lmm_type */
123 #define CP_ADAMS          1
124 #define CP_BDF            2
125 
126 /* nls_type */
127 #define CP_FUNCTIONAL     1
128 #define CP_NEWTON         2
129 
130 /* ode_type */
131 #define CP_EXPL           1
132 #define CP_IMPL           2
133 
134 /* tol_type */
135 #define CP_SS             1
136 #define CP_SV             2
137 #define CP_WF             3
138 
139 /* proj_type */
140 #define CP_PROJ_USER      1
141 #define CP_PROJ_INTERNAL  2
142 
143 /* proj_norm */
144 #define CP_PROJ_L2NORM    1
145 #define CP_PROJ_ERRNORM   2
146 
147 /* cnstr_type */
148 #define CP_CNSTR_LIN      1
149 #define CP_CNSTR_NONLIN   2
150 
151 /* mode */
152 #define CP_NORMAL         1
153 #define CP_ONE_STEP       2
154 #define CP_NORMAL_TSTOP   3
155 #define CP_ONE_STEP_TSTOP 4
156 
157 /*
158  * ----------------------------------------
159  * CPODES return flags
160  * ----------------------------------------
161  */
162 
163 #define CP_SUCCESS               0
164 #define CP_TSTOP_RETURN          1
165 #define CP_ROOT_RETURN           2
166 
167 #define CP_WARNING              99
168 
169 #define CP_TOO_MUCH_WORK        -1
170 #define CP_TOO_MUCH_ACC         -2
171 #define CP_ERR_FAILURE          -3
172 #define CP_CONV_FAILURE         -4
173 
174 #define CP_LINIT_FAIL           -5
175 #define CP_LSETUP_FAIL          -6
176 #define CP_LSOLVE_FAIL          -7
177 #define CP_ODEFUNC_FAIL         -8
178 #define CP_FIRST_ODEFUNC_ERR    -9
179 #define CP_REPTD_ODEFUNC_ERR    -10
180 #define CP_UNREC_ODEFUNC_ERR    -11
181 #define CP_RTFUNC_FAIL          -12
182 
183 #define CP_MEM_FAIL             -20
184 #define CP_MEM_NULL             -21
185 #define CP_ILL_INPUT            -22
186 #define CP_NO_MALLOC            -23
187 #define CP_BAD_K                -24
188 #define CP_BAD_T                -25
189 #define CP_BAD_DKY              -26
190 #define CP_TOO_CLOSE            -27
191 
192 #define CP_NO_QUAD              -30
193 #define CP_QUADFUNC_FAIL        -31
194 #define CP_FIRST_QUADFUNC_ERR   -32
195 #define CP_REPTD_QUADFUNC_ERR   -33
196 #define CP_UNREC_QUADFUNC_ERR   -34
197 
198 #define CP_BAD_IS               -40
199 #define CP_NO_SENS              -41
200 #define CP_SENSFUNC_FAIL        -42
201 #define CP_FIRST_SENSFUNC_ERR   -43
202 #define CP_REPTD_SENSFUNC_ERR   -44
203 #define CP_UNREC_SENSFUNC_ERR   -45
204 
205 #define CP_PLINIT_FAIL          -50
206 #define CP_PLSETUP_FAIL         -51
207 #define CP_PLSOLVE_FAIL         -52
208 #define CP_CNSTRFUNC_FAIL       -53
209 #define CP_PROJFUNC_FAIL        -54
210 #define CP_PROJ_FAILURE         -55
211 #define CP_REPTD_CNSTRFUNC_ERR  -56
212 #define CP_REPTD_PROJFUNC_ERR   -57
213 
214 #define CP_FIRST_CNSTRFUNC_ERR  -60
215 #define CP_NO_RECOVERY          -61
216 #define CP_LINESEARCH_FAIL      -62
217 
218 /*
219  * =================================================================
220  *              F U N C T I O N   T Y P E S
221  * =================================================================
222  */
223 
224 /*
225  * -----------------------------------------------------------------
226  * Type : CPRhsFn
227  * -----------------------------------------------------------------
228  * If the ODE is given in explicit form, the function which defines
229  * the right-hand side of the ODE system y'=f(t,y) must have type
230  * CPRhsFn.
231  * Such a function takes as input the independent variable value t,
232  * and the dependent variable vector y.  It stores the result of
233  * f(t,y) in the vector fout.  The y and fout arguments are of type
234  * N_Vector. (Allocation of memory for ydot is handled within CODES)
235  * The f_data parameter is the same as the f_data parameter set by
236  * the user through the CPodeInit or CPodeReInit functions. This
237  * user-supplied pointer is passed to the user's fun function every
238  * time it is called.
239  *
240  * A CPRhsFn should return 0 if successful, a negative value if
241  * an unrecoverable error occurred, and a positive value if a
242  * recoverable error (e.g. invalid y values) occurred.
243  * If an unrecoverable occurred, the integration is halted.
244  * If a recoverable error occurred, then (in most cases) CPODES
245  * will try to correct and retry.
246  * -----------------------------------------------------------------
247  */
248 
249 typedef int (*CPRhsFn)(realtype t, N_Vector y,
250                N_Vector fout, void *f_data);
251 
252 /*
253  * ----------------------------------------------------------------
254  * Type : CPResFn
255  * ----------------------------------------------------------------
256  * If the ODE is given in implicit form, the function which
257  * defines the ODE system F(t,y,y')=0 must have type CPResFn.
258  *
259  * A CPResFn takes as input the independent variable value t,
260  * the dependent variable vector y, and the derivative (with
261  * respect to t) of the y vector, yp.  It stores the result of
262  * F(t,y,y') in the vector fout. The y, yp, and fout arguments are
263  * of type N_Vector. The f_data parameter is the pointer f_data
264  * passed by the user to the CPodeInit or CPodeReInit functions.
265  * This user-supplied pointer is passed to the user's fun function
266  * every time it is called.
267  *
268  * A CPResFn function should return a value of 0 if successful, a
269  * positive value if a recoverable error occurred (e.g. y has an
270  * illegal value), or a negative value if a nonrecoverable error
271  * occurred. In the latter case, the program halts. If a recoverable
272  * error occurred, then (in most cases) the integrator will attempt
273  * to correct and retry.
274  * ----------------------------------------------------------------
275  */
276 
277 typedef int (*CPResFn)(realtype t, N_Vector y, N_Vector yp,
278                N_Vector fout, void *f_data);
279 
280 /*
281  * -----------------------------------------------------------------
282  * Type : CPCnstrFn
283  * -----------------------------------------------------------------
284  * The function cfun defining the invariant constraints c(t,y) = 0
285  * must have type CPCnstrFn.
286  *
287  * A CPCnstrFn takes as input the independent variable value t and
288  * the dependent variable vector y.  It stores the result of c(t,y)
289  * in the vector cout.  The y and cout arguments are of type
290  * N_Vector. (Allocation of memory for cout is handled within CPODES)
291  * The c_data parameter is the same as the c_data parameter set by
292  * the user through the CPodeProjDefineConstraints routine.
293  * This user-supplied pointer is passed to the user's cfun function
294  * every time it is called.
295  *
296  * A CPCnstrFn should return 0 if successful, a negative value if
297  * an unrecoverable error occurred, and a positive value if a
298  * recoverable error (e.g. invalid y values) occurred.
299  * If an unrecoverable occurred, the integration is halted.
300  * If a recoverable error occurred, then (in most cases) CPODES
301  * will try to correct and retry.
302  * -----------------------------------------------------------------
303  */
304 
305 typedef int (*CPCnstrFn)(realtype t, N_Vector y,
306              N_Vector cout, void *c_data);
307 
308 /*
309  * -----------------------------------------------------------------
310  * Type : CPProjFn
311  * -----------------------------------------------------------------
312  * A user-supplied function to performs the projection onto the
313  * invariant manifold c(t,y)=0, must have type CPProjFn.
314  *
315  * A CPProjFn takes as input the independent variable t and the
316  * current (corrected) variable vector ycur.
317  * It must compute a correction vector corr such that
318  * y = ycurr + corr lies on the manifold (i.e. c(t,y)=0).
319  * The value epsProj is provided to be used in the stopping test
320  * of a nonlinear solver iteration (the iterations should be
321  * terminated when the WRMS norm of the current iterate update
322  * is below epsProj).
323  *
324  * Note that, if the projection is orthogonal then it can be written as
325  *    y = P * ycur + alpha(t), with P^2 = P,
326  * then the projected error estimate is
327  *    errP = P * err
328  * and ||errP|| <= ||err||.
329  * The vector err contains a (scaled) version of the current error
330  * estimate. CPProjFn may also compute the projected error estimate
331  * and overwrite the vector err with errP. Otherwise, it should leave
332  * err unchanged.
333  *
334  * If err is NULL, a CPProjFn function should attempt NO projection
335  * of the error estimate (in this case, it was called from within
336  * the computation of consistent initial conditions).
337  *
338  * A CPProjFn should return 0 if successful, a negative value if
339  * an unrecoverable error occurred, and a positive value if a
340  * recoverable error (e.g. invalid y values) occurred.
341  * If an unrecoverable occurred, the integration is halted.
342  * If a recoverable error occurred, then (in most cases) CPODES
343  * will try to correct and retry.
344  * -----------------------------------------------------------------
345  * NOTE: If the user's projection routine needs other quantities,
346  *       they are accessible as follows: the error weight vector
347  *       can be obtained by calling CPodeGetErrWeights. The unit
348  *       roundoff is available as UNIT_ROUNDOFF defined in
349  *       sundials_types.h.
350  * -----------------------------------------------------------------
351  */
352 
353 typedef int (*CPProjFn)(realtype t, N_Vector ycur, N_Vector corr,
354             realtype epsProj, N_Vector err, void *pdata);
355 
356 /*
357  * -----------------------------------------------------------------
358  * Type : CPQuadFn
359  * -----------------------------------------------------------------
360  * The qfun function which defines the right hand side of the
361  * quadrature equations q' = fQ(t,y) must have type CPQuadFn.
362  * It takes as input the value of the independent variable t and
363  * the vector of states y and must store the result of fQ in qout.
364  * (Allocation of memory for qout is handled by CPODES).
365  * The q_data parameter is the same as the q_data parameter
366  * set by the user through the CPodeQuadInit function and is
367  * passed to the qfun function every time it is called.
368  *
369  * A CPQuadFn should return 0 if successful, a negative value if
370  * an unrecoverable error occurred, and a positive value if a
371  * recoverable error (e.g. invalid y values) occurred.
372  * If an unrecoverable occurred, the integration is halted.
373  * If a recoverable error occurred, then (in most cases) CPODES
374  * will try to correct and retry.
375  * -----------------------------------------------------------------
376  */
377 
378 typedef int (*CPQuadFn)(realtype t, N_Vector y,
379             N_Vector qout, void *q_data);
380 
381 /*
382  * -----------------------------------------------------------------
383  * Type : CPRootFn
384  * -----------------------------------------------------------------
385  * A function g, which defines a set of functions g_i(t,y,y') whose
386  * roots are sought during the integration, must have type CPRootFn.
387  * The function g takes as input the independent variable value
388  * t, the dependent variable vector y, and its derivative yp=y'.
389  * It stores the nrtfn values g_i(t,y,y') in the realtype array gout.
390  * (Allocation of memory for gout is handled within CPODES.)
391  * The g_data parameter is the same as that passed by the user
392  * to the CPodeRootInit routine.  This user-supplied pointer is
393  * passed to the user's g function every time it is called.
394  *
395  * A CPRootFn should return 0 if successful or a non-zero value
396  * if an error occurred (in which case the integration will be halted).
397  * -----------------------------------------------------------------
398  */
399 
400 typedef int (*CPRootFn)(realtype t, N_Vector y, N_Vector yp,
401             realtype *gout, void *g_data);
402 
403 /*
404  * -----------------------------------------------------------------
405  * Type : CPEwtFn
406  * -----------------------------------------------------------------
407  * A function e, which sets the error weight vector ewt, must have
408  * type CPEwtFn.
409  * The function e takes as input the current dependent variable y.
410  * It must set the vector of error weights used in the WRMS norm:
411  *
412  *   ||y||_WRMS = sqrt [ 1/N * sum ( ewt_i * y_i)^2 ]
413  *
414  * Typically, the vector ewt has components:
415  *
416  *   ewt_i = 1 / (reltol * |y_i| + abstol_i)
417  *
418  * The e_data parameter is the same as that passed by the user
419  * to the CPodeSetEwtFn routine.  This user-supplied pointer is
420  * passed to the user's e function every time it is called.
421  * A CPEwtFn e must return 0 if the error weight vector has been
422  * successfuly set and a non-zero value otherwise.
423  * -----------------------------------------------------------------
424  */
425 
426 typedef int (*CPEwtFn)(N_Vector y, N_Vector ewt, void *e_data);
427 
428 /*
429  * -----------------------------------------------------------------
430  * Type : CPErrHandlerFn
431  * -----------------------------------------------------------------
432  * A function eh, which handles error messages, must have type
433  * CPErrHandlerFn.
434  * The function eh takes as input the error code, the name of the
435  * module reporting the error, the error message, and a pointer to
436  * user data, the same as that passed to CPodeSetErrHandlerFn.
437  *
438  * All error codes are negative, except CP_WARNING which indicates
439  * a warning (the solver continues).
440  *
441  * A CPErrHandlerFn has no return value.
442  * -----------------------------------------------------------------
443  */
444 
445 typedef void (*CPErrHandlerFn)(int error_code,
446                    const char *module, const char *function,
447                    char *msg, void *eh_data);
448 
449 /*
450  * =================================================================
451  *          U S E R - C A L L A B L E   F U N C T I O N S
452  * =================================================================
453  */
454 
455 /*
456  * -----------------------------------------------------------------
457  * Function : CPodeCreate
458  * -----------------------------------------------------------------
459  * CPodeCreate creates an internal memory block for a problem to
460  * be solved by CPODES.
461  *
462  * ode_type  - form in which the ODE system is provided.
463  *             The legal values are CP_EXPL or CP_IMPL (see above).
464  *
465  * lmm_type  - type of linear multistep method to be used.
466  *             The legal values are CP_ADAMS and CP_BDF (see above).
467  *
468  * nls_type  - type of iteration used to solve the nonlinear
469  *             system that arises during each internal time step.
470  *             The legal values are CP_FUNCTIONAL and CP_NEWTON
471  *             for ode_type = CP_EXPL and only CP_NEWTON for
472  *             ode_type = CP_IMPL.
473  *
474  * If successful, CPodeCreate returns a pointer to initialized
475  * problem memory. This pointer should be passed to CPodeInit.
476  * If an initialization error occurs, CPodeCreate prints an error
477  * message to standard err and returns NULL.
478  * -----------------------------------------------------------------
479  */
480 
481 SUNDIALS_EXPORT void *CPodeCreate(int ode_type, int lmm_type, int nls_type);
482 
483 /*
484  * -----------------------------------------------------------------
485  * Function : CPodeInit
486  * -----------------------------------------------------------------
487  * CPodeInit allocates and initializes memory for a problem to be
488  * solved by CPODES.
489  *
490  * cpode_mem - pointer to CPODES memory returned by CPodeCreate.
491  *
492  * fun       - name of the C function defining the ODE system.
493  *             Depending on the ODE type (see CPodeCreate), fun
494  *             should be of type CPRhsFn (ode_type=CP_EXPL) or of
495  *             type CPResFn (ode_type=CP_IMPL).
496  *
497  * f_data    - pointer to user data that will be passed to the
498  *             fun function every time it is called.
499  *
500  * t0        - initial value of t.
501  *
502  * y0        - initial condition vector y(t0).
503  *
504  * yp0       - initial condition vector y'(t0).
505  *
506  * tol_type  - type of tolerances to be used. The legal values are:
507  *             CP_SS (scalar relative and absolute tolerances),
508  *             CP_SV (scalar relative tolerance and vector
509  *                    absolute tolerance).
510  *             CP_WF (indicates that the user will provide a
511  *                    function to evaluate the error weights.
512  *                    In this case, reltol and abstol are ignored.)
513  *
514  * reltol    - scalar relative tolerance scalar.
515  *
516  * abstol    - pointer to the absolute tolerance scalar or
517  *             an N_Vector of absolute tolerances.
518  *
519  * The parameters tol_type, reltol, and abstol define a vector of
520  * error weights, ewt, with components
521  *   ewt[i] = 1/(reltol*abs(y[i]) + abstol)     if tol_type = CP_SS
522  *   ewt[i] = 1/(reltol*abs(y[i]) + abstol[i])  if tol_type = CP_SV.
523  * This vector is used in all error and convergence tests, which
524  * use a weighted RMS norm on all error-like vectors v:
525  *    WRMSnorm(v) = sqrt( (1/N) sum(i=1..N) (v[i]*ewt[i])^2 ),
526  * where N is the problem dimension.
527  *
528  * Return flag:
529  *  CP_SUCCESS if successful
530  *  CP_MEM_NULL if the CPODES memory was NULL
531  *  CP_MEM_FAIL if a memory allocation failed
532  *  CP_ILL_INPUT f an argument has an illegal value.
533  * -----------------------------------------------------------------
534  */
535 
536 SUNDIALS_EXPORT int CPodeInit(void *cpode_mem,
537                   void *fun, void *f_data,
538                   realtype t0, N_Vector y0, N_Vector yp0,
539                   int tol_type, realtype reltol, void *abstol);
540 
541 /*
542  * -----------------------------------------------------------------
543  * Function : CPodeReInit
544  * -----------------------------------------------------------------
545  * CPodeReInit re-initializes CPODES for the solution of a problem,
546  * where a prior call to CPodeInit has been made with the same
547  * problem size N. CPodeReInit performs the same input checking
548  * and initializations that CPodeInit does.
549  * But it does no memory allocation, assuming that the existing
550  * internal memory is sufficient for the new problem.
551  *
552  * The use of CPodeReInit requires that the maximum method order,
553  * maxord, is no larger for the new problem than for the problem
554  * specified in the last call to CPodeInit.  This condition is
555  * automatically fulfilled if the multistep method parameter lmm_type
556  * is unchanged (or changed from CP_ADAMS to CP_BDF) and the default
557  * value for maxord is specified.
558  *
559  * All of the arguments to CPodeReInit have names and meanings
560  * identical to those of CPodeInit.
561  *
562  * The return value of CPodeReInit is equal to CP_SUCCESS = 0 if
563  * there were no errors; otherwise it is a negative int equal to:
564  *   CP_MEM_NULL      indicating cpode_mem was NULL (i.e.,
565  *                    CPodeCreate has not been called).
566  *   CP_NO_MALLOC     indicating that cpode_mem has not been
567  *                    allocated (i.e., CPodeInit has not been
568  *                    called).
569  *   CP_ILL_INPUT     indicating an input argument was illegal
570  *                    (including an attempt to increase maxord).
571  * In case of an error return, an error message is also printed.
572  * -----------------------------------------------------------------
573  */
574 
575 SUNDIALS_EXPORT int CPodeReInit(void *cpode_mem,
576                 void *fun, void *f_data,
577                 realtype t0, N_Vector y0, N_Vector yp0,
578                 int tol_type, realtype reltol, void *abstol);
579 
580 /*
581  * -----------------------------------------------------------------
582  * Function : CPodeProjInit
583  * -----------------------------------------------------------------
584  * CPodeProjInit initializes the internal CPODES coordinate
585  * projection algorithm. It must be called after CPodeCreate and
586  * CPodeInit.
587  *
588  * The arguments are as follows:
589  *
590  * cpode_mem  - pointer to CPODES memory returned by CPodeCreate.
591  *
592  * proj_norm  - the norm in which projection is to be done. Legal
593  *              values are CP_PROJ_L2NORM and CP_PROJ_ERRNORM.
594  *
595  * cnstr_type - the type of constraints.
596  *                CP_CNSTR_LIN :   linear constraints.
597  *                CP_CNSTR_NONLIN: nonlinear constraints.
598  *
599  * cfun       - name of the user-supplied function (type CPCnstrFn)
600  *              defining the invariant manifold.
601  *
602  * c_data     - pointer to user data that will be passed to the
603  *              cfun function every time it is called.
604  *
605  * ctol       - a vector of "absolute tolerances" for the constraints.
606  *              In default operation, this vector is only used as
607  *              a template for cloning other vectors. However, if
608  *              enabled through CPodeSetProjTestCnstr, these values,
609  *              together with reltol, are used to compute the
610  *              constraint WL2 norm and a projection will be
611  *              perfomed only if ||c(t)|_WL2 >= 1.0
612  *
613  * The return value of CPodeProjInit is equal to CP_SUCCESS = 0 if
614  * there were no errors; otherwise it is a negative int equal to:
615  *   CP_MEM_NULL  - cpode_mem was NULL
616  *                  (i.e., CPodeCreate has not been called).
617  *   CP_NO_MALLOC - cpode_mem has not been allocated
618  *                  (i.e., CPodeInit has not been called).
619  *   CP_MEM_FAIL  - a memory allocation failed.
620  *   CP_ILL_INPUT - an input argument was illegal.
621  * In case of an error return, an error message is also printed.
622  * -----------------------------------------------------------------
623  */
624 
625 SUNDIALS_EXPORT int CPodeProjInit(void *cpode_mem, int proj_norm,
626                   int cnstr_type, CPCnstrFn cfun, void *c_data,
627                   N_Vector ctol);
628 
629 /*
630  * -----------------------------------------------------------------
631  * Function : CPodeProjDefine
632  * -----------------------------------------------------------------
633  * CPodeProjDefine initializes coordinate projection using a funciton
634  * provided by the user. It must be called after CPodeInit.
635  *
636  * The arguments are as follows:
637  *
638  * cpode_mem  - pointer to CPODES memory returned by CPodeCreate.
639  *
640  * pfun       - name of the user-supplied function (type CPProjFn)
641  *              which will perform the projection.
642  *
643  * p_data     - pointer to user data that will be passed to the
644  *              pfun function every time it is called.
645  *
646  * The return value of CPodeProjDefine is CP_SUCCESS if there were
647  * no errors, or CP_MEM_NULL if the cpode_mem argument was NULL
648  * (i.e., CPodeCreate has not been called).
649  * In case of an error return, an error message is also printed.
650  * -----------------------------------------------------------------
651  */
652 
653 SUNDIALS_EXPORT int CPodeProjDefine(void *cpode_mem, CPProjFn pfun, void *p_data);
654 
655 /*
656  * -----------------------------------------------------------------
657  * Function : CPodeQuadInit and CPodeQuadReInit
658  * -----------------------------------------------------------------
659  * CVodeQuadInit allocates and initializes memory related to
660  * quadrature integration.
661  *
662  * cpode_mem - pointer to CPODES memory returned by CPodeCreate
663  *
664  * qfun      - the user-provided integrand routine.
665  *
666  * q_data    - pointer to user data that will be passed to the
667  *             qfun function every time it is called.
668  *
669  * q0        - N_Vector with initial values for quadratures
670  *             (typically q0 has all zero components).
671  *
672  * CPodeQuadReInit re-initializes CPODES's quadrature related memory
673  * for a problem, assuming it has already been allocated in prior calls
674  * to CPodeInit and CPodeQuadInit.  All problem specification inputs
675  * are checked for errors. The number of quadratures Nq is assumed to
676  * be unchanged since the previous call to CPodeQuadInit.
677  *
678  * Return values:
679  *  CP_SUCCESS if successful
680  *  CP_MEM_NULL if the CPODES memory was NULL
681  *  CP_MEM_FAIL if a memory allocation failed
682  * -----------------------------------------------------------------
683  */
684 
685 SUNDIALS_EXPORT int CPodeQuadInit(void *cpode_mem, CPQuadFn qfun, void *q_data, N_Vector q0);
686 SUNDIALS_EXPORT int CPodeQuadReInit(void *cpode_mem, CPQuadFn qfun, void *q_data, N_Vector q0);
687 
688 /*
689  * -----------------------------------------------------------------
690  * Function : CPodeRootInit
691  * -----------------------------------------------------------------
692  * CPodeRootInit initializes a rootfinding problem to be solved
693  * during the integration of the ODE system.  It must be called
694  * after CPodeCreate, and before CPode.  The arguments are:
695  *
696  * cpode_mem - pointer to CPODES memory returned by CPodeCreate.
697  *
698  * nrtfn     - number of functions g_i, an int >= 0.
699  *
700  * gfun      - name of user-supplied function, of type CPRootFn,
701  *             defining the functions g_i whose roots are sought.
702  *
703  * g_data    - pointer to user data that will be passed to the
704  *             gfun function every time it is called.
705  *
706  * If a new problem is to be solved with a call to CPodeReInit,
707  * where the new problem has no root functions but the prior one
708  * did, then call CPodeRootInit with nrtfn = 0.
709  *
710  * The return value of CPodeRootInit is CP_SUCCESS = 0 if there were
711  * no errors; otherwise it is a negative int equal to:
712  *   CP_MEM_NULL  - indicating cpode_mem was NULL.
713  *   CP_MEM_FAIL  - indicating a memory allocation failed.
714  *   CP_ILL_INPUT - indicating nrtfn > 0 but gfun = NULL.
715  * In case of an error return, an error message is also printed.
716  * -----------------------------------------------------------------
717  */
718 
719 SUNDIALS_EXPORT int CPodeRootInit(void *cpode_mem, int nrtfn, CPRootFn gfun, void *g_data);
720 
721 /*
722  * -----------------------------------------------------------------
723  * Function : CPodeCalcIC
724  * -----------------------------------------------------------------
725  * CPodeCalcIC calculates corrected initial conditions that are
726  * consistent with the invariant constraints and (for implicit-form
727  * ODEs) with the ODE system itself. It first projects the initial
728  * guess for the state vector (given by the user through CPodeInit
729  * or CPodeReInit) and then, if necessary, computes a state derivative
730  * vector as solution of F(t0, y0, y0') = 0.
731  *
732  * Note: If the initial conditions satisfy both the constraints and
733  * (for implicit-form ODEs) the ODE itself, a call to CPodeCalcIC
734  * is NOT necessary.
735  *
736  * A call to CpodeCalcIC must be preceded by a successful call to
737  * CPodeInit or CPodeReInit for the given ODE problem, and by a
738  * successful call to the linear system solver specification
739  * routine.
740  *
741  * A call to CPodeCalcIC should precede the call(s) to CPode for
742  * the given problem.
743  *
744  * If successful, CPodeCalcIC stores internally the corrected
745  * initial conditions which will be used to start the integration
746  * at the first call to CPode.
747  *
748  * The only argument to CPodeCalcIC is the pointer to the CPODE
749  * memory block returned by CPodeCreate.
750  *
751  * The return value of CPodeCalcIC is one of the following:
752  *
753  * CP_SUCCESS
754  * CP_MEM_NULL
755  * CP_NO_MALLOC
756  * CP_ILL_INPUT
757  * CP_LINIT_FAIL
758  * CP_PLINIT_FAIL
759  * CP_FIRST_CNSTRFUNC_ERR
760  * CP_PROJ_FAILURE
761  * CP_CNSTRFUNC_FAIL
762  * CP_PROJFUNC_FAIL
763  * CP_PLSETUP_FAIL
764  * CP_PLSOLVE_FAIL
765  * CP_NO_RECOVERY
766  *
767  * -----------------------------------------------------------------
768  */
769 
770 SUNDIALS_EXPORT int CPodeCalcIC(void *cpode_mem);
771 
772 /*
773  * -----------------------------------------------------------------
774  * Integrator optional input specification functions
775  * -----------------------------------------------------------------
776  * The following functions can be called to set optional inputs
777  * to values other than the defaults given below:
778  *
779  * Function                |  Optional input / [ default value ]
780  * -----------------------------------------------------------------
781  *                         |
782  * CPodeSetErrHandlerFn    | user-provided ErrHandler function.
783  *                         | [internal]
784  *                         |
785  * CPodeSetErrFile         | the file pointer for an error file
786  *                         | where all CPODES warning and error
787  *                         | messages will be written if the default
788  *                         | internal error handling function is used.
789  *                         | This parameter can be stdout (standard
790  *                         | output), stderr (standard error), or a
791  *                         | file pointer (corresponding to a user
792  *                         | error file opened for writing) returned
793  *                         | by fopen.
794  *                         | If not called, then all messages will
795  *                         | be written to the standard error stream.
796  *                         | [stderr]
797  *                         |
798  * -----------------------------------------------------------------
799  *                         |
800  * CPodeSetEwtFn           | user-provided EwtSet function efun and
801  *                         | a pointer to user data that will be
802  *                         | passed to the user's efun function
803  *                         | every time efun is called.
804  *                         | [internal]
805  *                         | [NULL]
806  *                         |
807  * CPodeSetMaxOrd          | maximum LMM order to be used by the
808  *                         | solver.
809  *                         | [12 for Adams , 5 for BDF]
810  *                         |
811  * CPodeSetMaxNumSteps     | maximum number of internal steps to be
812  *                         | taken by the solver in its attempt to
813  *                         | reach tout.
814  *                         | [500]
815  *                         |
816  * CPodeSetMaxHnilWarns    | maximum number of warning messages
817  *                         | issued by the solver that t+h==t on the
818  *                         | next internal step. A value of -1 means
819  *                         | no such messages are issued.
820  *                         | [10]
821  *                         |
822  * CPodeSetStabLimDet      | flag to turn on/off stability limit
823  *                         | detection (TRUE = on, FALSE = off).
824  *                         | When BDF is used and order is 3 or
825  *                         | greater, CPsldet is called to detect
826  *                         | stability limit.  If limit is detected,
827  *                         | the order is reduced.
828  *                         | [FALSE]
829  *                         |
830  * CPodeSetInitStep        | initial step size.
831  *                         | [estimated by CPODES]
832  *                         |
833  * CPodeSetMinStep         | minimum absolute value of step size
834  *                         | allowed.
835  *                         | [0.0]
836  *                         |
837  * CPodeSetMaxStep         | maximum absolute value of step size
838  *                         | allowed.
839  *                         | [infinity]
840  *                         |
841  * CPodeSetStopTime        | the independent variable value past
842  *                         | which the solution is not to proceed.
843  *                         | [infinity]
844  *                         |
845  * CPodeSetMaxErrTestFails | Maximum number of error test failures
846  *                         | in attempting one step.
847  *                         | [7]
848  *                         |
849  * -----------------------------------------------------------------
850  *                         |
851  * CPodeSetMaxNonlinIters  | Maximum number of nonlinear solver
852  *                         | iterations at one solution.
853  *                         | [3]
854  *                         |
855  * CPodeSetMaxConvFails    | Maximum number of convergence failures
856  *                         | allowed in attempting one step.
857  *                         | [10]
858  *                         |
859  * CPodeSetNonlinConvCoef  | Coefficient in the nonlinear
860  *                         | convergence test.
861  *                         | [0.1]
862  *                         |
863  * -----------------------------------------------------------------
864  *                         |
865  * CPodeSetProjUpdateErrEst| toggles ON/OFF projection of the
866  *                         | error estimate.
867  *                         | [TRUE]
868  *                         |
869  * CPodeSetProjFrequency   | frequency with which the projection
870  *                         | step is performed. A value of 1
871  *                         | indicates that the projection step
872  *                         | will be performed at every step.
873  *                         | A value of 0 will disable projection.
874  *                         | [1]
875  *                         |
876  * CPodeSetProjTestCnstr   | if TRUE, the internal projection
877  *                         | function will be performed only if
878  *                         | the constraint violation is larger
879  *                         | than the prescribed tolerances.
880  *                         | Otherwise, the tolerances are ignored
881  *                         | and projection is always performed.
882  *                         | [FALSE]
883  *                         |
884  * CPodeSetProjLsetupFreq  | frequency with which the linear
885  *                         | solver setup function is called
886  *                         | (i.e. frequency of constraint
887  *                         | Jacobian evaluations). The default
888  *                         | value of 1 forces a Jacobian
889  *                         | evaluation before every single
890  *                         | projection step.
891  *                         | [1]
892  *                         |
893  * CPodeSetProjNonlinConvCoef | Coefficient in the nonlinear
894  *                         | convergence test (for projection).
895  *                         | [0.1]
896  *                         |
897  * -----------------------------------------------------------------
898  *                         |
899  * CPodeSetQuadErrCon      | are quadrature variables considered in
900  *                         | the error control?
901  *                         | If yes, set tolerances for quadrature
902  *                         | integration.
903  *                         | [errconQ = FALSE]
904  *                         | [no tolerances]
905  *                         |
906  * -----------------------------------------------------------------
907  *                         |
908  * CPodeSetTolerances      | Changes the integration tolerances
909  *                         | between calls to CPode.
910  *                         | [set by CPodeInit/CPodeReInit]
911  *                         |
912  * ----------------------------------------------------------------
913  *                         |
914  * CPodeSetRootDirection   | Specifies the direction of zero
915  *                         | crossings to be monitored
916  *                         | [both directions]
917  *                         |
918  * -----------------------------------------------------------------
919  * Return flag:
920  *   CP_SUCCESS   if successful
921  *   CP_MEM_NULL  if the CPODES memory is NULL
922  *   CP_ILL_INPUT if an argument has an illegal value
923  * -----------------------------------------------------------------
924  */
925 
926 SUNDIALS_EXPORT int CPodeSetErrHandlerFn(void *cpode_mem, CPErrHandlerFn ehfun, void *eh_data);
927 SUNDIALS_EXPORT int CPodeSetErrFile(void *cpode_mem, FILE *errfp);
928 
929 SUNDIALS_EXPORT int CPodeSetEwtFn(void *cpode_mem, CPEwtFn efun, void *e_data);
930 SUNDIALS_EXPORT int CPodeSetMaxOrd(void *cpode_mem, int maxord);
931 SUNDIALS_EXPORT int CPodeSetMaxNumSteps(void *cpode_mem, long int mxsteps);
932 SUNDIALS_EXPORT int CPodeSetMaxHnilWarns(void *cpode_mem, int mxhnil);
933 SUNDIALS_EXPORT int CPodeSetStabLimDet(void *cpode_mem, booleantype stldet);
934 SUNDIALS_EXPORT int CPodeSetInitStep(void *cpode_mem, realtype hin);
935 SUNDIALS_EXPORT int CPodeSetMinStep(void *cpode_mem, realtype hmin);
936 SUNDIALS_EXPORT int CPodeSetMaxStep(void *cpode_mem, realtype hmax);
937 SUNDIALS_EXPORT int CPodeSetStopTime(void *cpode_mem, realtype tstop);
938 SUNDIALS_EXPORT int CPodeSetMaxErrTestFails(void *cpode_mem, int maxnef);
939 
940 SUNDIALS_EXPORT int CPodeSetMaxNonlinIters(void *cpode_mem, int maxcor);
941 SUNDIALS_EXPORT int CPodeSetMaxConvFails(void *cpode_mem, int maxncf);
942 SUNDIALS_EXPORT int CPodeSetNonlinConvCoef(void *cpode_mem, realtype nlscoef);
943 
944 SUNDIALS_EXPORT int CPodeSetProjUpdateErrEst(void *cpode_mem, booleantype proj_err);
945 SUNDIALS_EXPORT int CPodeSetProjFrequency(void *cpode_mem, long int proj_freq);
946 SUNDIALS_EXPORT int CPodeSetProjTestCnstr(void *cpode_mem, booleantype test_cnstr);
947 SUNDIALS_EXPORT int CPodeSetProjLsetupFreq(void *cpode_mem, long int proj_lset_freq);
948 SUNDIALS_EXPORT int CPodeSetProjNonlinConvCoef(void *cpode_mem, realtype prjcoef);
949 
950 SUNDIALS_EXPORT int CPodeSetQuadErrCon(void *cpode_mem, booleantype errconQ,
951                        int tol_typeQ, realtype reltolQ, void *abstolQ);
952 
953 SUNDIALS_EXPORT int CPodeSetTolerances(void *cpode_mem,
954                        int tol_type, realtype reltol, void *abstol);
955 
956 SUNDIALS_EXPORT int CPodeSetRootDirection(void *cpode_mem, int *rootdir);
957 
958 /*
959  * -----------------------------------------------------------------
960  * Function : CPode
961  * -----------------------------------------------------------------
962  * CPode integrates the ODE over an interval in t.
963  * If mode is CP_NORMAL, then the solver integrates from its
964  * current internal t value to a point at or beyond tout, then
965  * interpolates to t = tout and returns y(tout) in the user-
966  * allocated vector yout. If mode is CP_ONE_STEP, then the solver
967  * takes one internal time step and returns in yout the value of
968  * y at the new internal time. In this case, tout is used only
969  * during the first call to CPode to determine the direction of
970  * integration and the rough scale of the t variable.  If mode is
971  * CP_NORMAL_TSTOP or CP_ONE_STEP_TSTOP, then CPode returns the
972  * solution at tstop if that comes sooner than tout or the end of
973  * the next internal step, respectively.  In any case,
974  * the time reached by the solver is placed in (*tret). The
975  * user is responsible for allocating the memory for this value.
976  *
977  * cpode_mem is the pointer to CPODE memory returned by
978  *           CPodeCreate.
979  *
980  * tout  is the next time at which a computed solution is desired.
981  *
982  * tret  is a pointer to a real location. CPode sets (*tret) to
983  *       the time reached by the solver and returns yout=y(*tret)
984  *       and ypout=y'(*tret).
985  *
986  * yout  is the computed solution vector. In CP_NORMAL mode with no
987  *       errors and no roots found, yout=y(tout).
988  *
989  * ypout is the computed derivative vector.
990  *
991  * mode  is CP_NORMAL, CP_ONE_STEP, CP_NORMAL_TSTOP, or
992  *       CP_ONE_STEP_TSTOP. These four modes are described above.
993  *
994  * Here is a brief description of each return value:
995  *
996  * CP_SUCCESS:      CPode succeeded and no roots were found.
997  *
998  * CP_ROOT_RETURN:  CPode succeeded, and found one or more roots.
999  *                  If nrtfn > 1, call CPodeGetRootInfo to see
1000  *                  which g_i were found to have a root at (*tret).
1001  *
1002  * CP_TSTOP_RETURN: CPode succeeded and returned at tstop.
1003  *
1004  * CP_MEM_NULL:     The cpode_mem argument was NULL.
1005  *
1006  * CP_NO_MALLOC:    cpode_mem was not allocated.
1007  *
1008  * CP_ILL_INPUT:    One of the inputs to CPode is illegal. This
1009  *                  includes the situation when a component of the
1010  *                  error weight vectors becomes < 0 during
1011  *                  internal time-stepping.  It also includes the
1012  *                  situation where a root of one of the root
1013  *                  functions was found both at t0 and very near t0.
1014  *                  The ILL_INPUT flag will also be returned if the
1015  *                  linear solver routine CP--- (called by the user
1016  *                  after calling CPodeCreate) failed to set one of
1017  *                  the linear solver-related fields in cpode_mem or
1018  *                  if the linear solver's init routine failed. In
1019  *                  any case, the user should see the printed
1020  *                  error message for more details.
1021  *
1022  * CP_TOO_MUCH_WORK: The solver took mxstep internal steps but
1023  *                  could not reach tout. The default value for
1024  *                  mxstep is MXSTEP_DEFAULT = 500.
1025  *
1026  * CP_TOO_MUCH_ACC: The solver could not satisfy the accuracy
1027  *                  demanded by the user for some internal step.
1028  *
1029  * CP_ERR_FAILURE:  Error test failures occurred too many times
1030  *                  (= MXNEF = 7) during one internal time step or
1031  *                  occurred with |h| = hmin.
1032  *
1033  * CP_CONV_FAILURE: Convergence test failures occurred too many
1034  *                  times (= MXNCF = 10) during one internal time
1035  *                  step or occurred with |h| = hmin.
1036  *
1037  * CP_LINIT_FAIL:   The linear solver's initialization function
1038  *                  failed.
1039  *
1040  * CP_LSETUP_FAIL:  The linear solver's setup routine failed in an
1041  *                  unrecoverable manner.
1042  *
1043  * CP_LSOLVE_FAIL:  The linear solver's solve routine failed in an
1044  *                  unrecoverable manner.
1045  * -----------------------------------------------------------------
1046  */
1047 
1048 SUNDIALS_EXPORT int CPode(void *cpode_mem, realtype tout, realtype *tret,
1049               N_Vector yout, N_Vector ypout, int mode);
1050 
1051 /*
1052  * -----------------------------------------------------------------
1053  * Function : CPodeGetDky
1054  * -----------------------------------------------------------------
1055  * CPodeGetDky computes the kth derivative of the y function at
1056  * time t, where tn-hu <= t <= tn, tn denotes the current
1057  * internal time reached, and hu is the last internal step size
1058  * successfully used by the solver. The user may request
1059  * k=0, 1, ..., qu, where qu is the order last used. The
1060  * derivative vector is returned in dky. This vector must be
1061  * allocated by the caller. It is only legal to call this
1062  * function after a successful return from CPode.
1063  *
1064  * cpode_mem is the pointer to CPODES memory returned by
1065  *           CPodeCreate.
1066  *
1067  * t   is the time at which the kth derivative of y is evaluated.
1068  *     The legal range for t is [tn-hu,tn] as described above.
1069  *
1070  * k   is the order of the derivative of y to be computed. The
1071  *     legal range for k is [0,qu] as described above.
1072  *
1073  * dky is the output derivative vector [((d/dy)^k)y](t).
1074  *
1075  * The return value for CPodeGetDky is one of:
1076  *
1077  *   CP_SUCCESS:  CPodeGetDky succeeded.
1078  *
1079  *   CP_BAD_K:    k is not in the range 0, 1, ..., qu.
1080  *
1081  *   CP_BAD_T:    t is not in the interval [tn-hu,tn].
1082  *
1083  *   CP_BAD_DKY:  The dky argument was NULL.
1084  *
1085  *   CP_MEM_NULL: The cpode_mem argument was NULL.
1086  * -----------------------------------------------------------------
1087  */
1088 
1089 SUNDIALS_EXPORT int CPodeGetDky(void *cpode_mem, realtype t, int k, N_Vector dky);
1090 
1091 /*
1092  * -----------------------------------------------------------------
1093  * Quadrature integration solution extraction routines
1094  * -----------------------------------------------------------------
1095  * The following functions can be called to obtain the quadrature
1096  * variables (or derivatives of them) after a successful integration
1097  * step.  If quadratures were not computed, they return CP_NO_QUAD.
1098  * -----------------------------------------------------------------
1099  */
1100 
1101 SUNDIALS_EXPORT int CPodeGetQuad(void *cpode_mem, realtype t, N_Vector yQout);
1102 SUNDIALS_EXPORT int CPodeGetQuadDky(void *cpode_mem, realtype t, int k, N_Vector dky);
1103 
1104 
1105 /*
1106  * -----------------------------------------------------------------
1107  * IC calculation optional output extraction functions
1108  * -----------------------------------------------------------------
1109  * CPodeGetConsistentIC returns the consistent initial conditions
1110  *       computed by CPodeCalcIC
1111  * -----------------------------------------------------------------
1112  */
1113 
1114 SUNDIALS_EXPORT int CPodeGetConsistentIC(void *cpode_mem, N_Vector yy0, N_Vector yp0);
1115 
1116 /*
1117  * -----------------------------------------------------------------
1118  * Integrator optional output extraction functions
1119  * -----------------------------------------------------------------
1120  * The following functions can be called to get optional outputs
1121  * and statistics related to the main integrator.
1122  * -----------------------------------------------------------------
1123  * CPodeGetWorkSpace returns the CPODES real and integer workspaces
1124  * CPodeGetNumSteps returns the cumulative number of internal
1125  *    steps taken by the solver
1126  * CPodeGetNumFctEvals returns the number of calls to the user's
1127  *    fun function
1128  * CPodeGetNumLinSolvSetups returns the number of calls made to
1129  *    the linear solver's setup routine
1130  * CPodeGetNumErrTestFails returns the number of local error test
1131  *    failures that have occurred
1132  * CPodeGetLastOrder returns the order used during the last
1133  *    internal step
1134  * CPodeGetCurrentOrder returns the order to be used on the next
1135  *    internal step
1136  * CPodeGetNumStabLimOrderReds returns the number of order
1137  *    reductions due to stability limit detection
1138  * CPodeGetActualInitStep returns the actual initial step size
1139  *    used by CPODES
1140  * CPodeGetLastStep returns the step size for the last internal step
1141  * CPodeGetCurrentStep returns the step size to be attempted on
1142  *    the next internal step
1143  * CPodeGetCurrentTime returns the current internal time reached
1144  *    by the solver
1145  * CPodeGetTolScaleFactor returns a suggested factor by which the
1146  *    user's tolerances should be scaled when too much accuracy has
1147  *    been requested for some internal step
1148  * CPodeGetErrWeights returns the current error weight vector.
1149  *    The user must allocate space for eweight.
1150  * CPodeGetEstLocalErrors returns the vector of estimated local
1151  *    errors. The user must allocate space for ele.
1152  * CPodeGetNumGEvals returns the number of calls to the user's
1153  *    gfun function (for rootfinding)
1154  * CPodeGetRootInfo returns the indices for which g_i was found to
1155  *    have a root. The user must allocate space for rootsfound.
1156  *    For i = 0 ... nrtfn-1, rootsfound[i] = 1 if g_i has a root,
1157  *    and = 0 if not.
1158  * CPodeGetRootWindow returns the most recent (tLo,tHi] window in
1159  *    which a g_i was found to have a root.
1160  * CPodeGetIntegratorStats retruns most of the optional outputs as
1161  *    a group.
1162  * CPodeGet* return values:
1163  *   CP_SUCCESS   if successful
1164  *   CP_MEM_NULL  if the CPODES memory was NULL
1165  *   CP_NO_SLDET  if stability limit was not turned on
1166  * -----------------------------------------------------------------
1167  */
1168 
1169 SUNDIALS_EXPORT int CPodeGetWorkSpace(void *cpode_mem, long int *lenrw, long int *leniw);
1170 SUNDIALS_EXPORT int CPodeGetNumSteps(void *cpode_mem, long int *nsteps);
1171 SUNDIALS_EXPORT int CPodeGetNumFctEvals(void *cpode_mem, long int *nfevals);
1172 SUNDIALS_EXPORT int CPodeGetNumLinSolvSetups(void *cpode_mem, long int *nlinsetups);
1173 SUNDIALS_EXPORT int CPodeGetNumErrTestFails(void *cpode_mem, long int *netfails);
1174 SUNDIALS_EXPORT int CPodeGetLastOrder(void *cpode_mem, int *qlast);
1175 SUNDIALS_EXPORT int CPodeGetCurrentOrder(void *cpode_mem, int *qcur);
1176 SUNDIALS_EXPORT int CPodeGetNumStabLimOrderReds(void *cpode_mem, long int *nslred);
1177 SUNDIALS_EXPORT int CPodeGetActualInitStep(void *cpode_mem, realtype *hinused);
1178 SUNDIALS_EXPORT int CPodeGetLastStep(void *cpode_mem, realtype *hlast);
1179 SUNDIALS_EXPORT int CPodeGetCurrentStep(void *cpode_mem, realtype *hcur);
1180 SUNDIALS_EXPORT int CPodeGetCurrentTime(void *cpode_mem, realtype *tcur);
1181 SUNDIALS_EXPORT int CPodeGetTolScaleFactor(void *cpode_mem, realtype *tolsfac);
1182 SUNDIALS_EXPORT int CPodeGetErrWeights(void *cpode_mem, N_Vector eweight);
1183 SUNDIALS_EXPORT int CPodeGetEstLocalErrors(void *cpode_mem, N_Vector ele);
1184 SUNDIALS_EXPORT int CPodeGetNumGEvals(void *cpode_mem, long int *ngevals);
1185 SUNDIALS_EXPORT int CPodeGetRootInfo(void *cpode_mem, int *rootsfound);
1186 SUNDIALS_EXPORT int CPodeGetRootWindow(void *cpode_mem, realtype *tlo, realtype *thi);
1187 SUNDIALS_EXPORT int CPodeGetIntegratorStats(void *cpode_mem, long int *nsteps,
1188                         long int *nfevals, long int *nlinsetups,
1189                         long int *netfails, int *qlast,
1190                         int *qcur, realtype *hinused, realtype *hlast,
1191                         realtype *hcur, realtype *tcur);
1192 
1193 /*
1194  * -----------------------------------------------------------------
1195  * Nonlinear solver optional output extraction functions
1196  * -----------------------------------------------------------------
1197  * The following functions can be called to get optional outputs
1198  * and statistics related to the nonlinear solver.
1199  * -----------------------------------------------------------------
1200  * CPodeGetNumNonlinSolvIters returns the number of nonlinear
1201  *    solver iterations performed.
1202  * CPodeGetNumNonlinSolvConvFails returns the number of nonlinear
1203  *    convergence failures.
1204  * CPodeGetNonlinSolvStats returns the nonlinear solver optional
1205  *    outputs in a group.
1206  * -----------------------------------------------------------------
1207  */
1208 
1209 SUNDIALS_EXPORT int CPodeGetNumNonlinSolvIters(void *cpode_mem, long int *nniters);
1210 SUNDIALS_EXPORT int CPodeGetNumNonlinSolvConvFails(void *cpode_mem, long int *nncfails);
1211 SUNDIALS_EXPORT int CPodeGetNonlinSolvStats(void *cpode_mem, long int *nniters,
1212                         long int *nncfails);
1213 
1214 
1215 /*
1216  * -----------------------------------------------------------------
1217  * Projection optional output extraction functions
1218  * -----------------------------------------------------------------
1219  * The following functions can be called to get optional outputs
1220  * and statistics related to the projection step.
1221  * -----------------------------------------------------------------
1222  * -----------------------------------------------------------------
1223  */
1224 
1225 SUNDIALS_EXPORT int CPodeGetProjNumProj(void *cpode_mem, long int *nproj);
1226 SUNDIALS_EXPORT int CPodeGetProjNumCnstrEvals(void *cpode_mem, long int *nce);
1227 SUNDIALS_EXPORT int CPodeGetProjNumLinSolvSetups(void *cpode_mem, long int *nsetupsP);
1228 SUNDIALS_EXPORT int CPodeGetProjNumFailures(void *cpode_mem, long int *nprf);
1229 SUNDIALS_EXPORT int CPodeGetProjStats(void *cpode_mem, long int *nproj,
1230                       long int *nce, long int *nsetupsP,
1231                       long int *nprf);
1232 
1233 /*
1234  * -----------------------------------------------------------------
1235  * Quadrature integration optional output extraction functions
1236  * -----------------------------------------------------------------
1237  * The following functions can be called to get optional outputs
1238  * and statistics related to the integration of quadratures.
1239  * -----------------------------------------------------------------
1240  * CPodeGetQuadNumFunEvals returns the number of calls to the
1241  *    user function qfun defining the integrand
1242  * CPodeGetQuadErrWeights returns the vector of error weights for
1243  *    the quadrature variables. The user must allocate space for ewtQ.
1244  * -----------------------------------------------------------------
1245  */
1246 
1247 SUNDIALS_EXPORT int CPodeGetQuadNumFunEvals(void *cpode_mem, long int *nqevals);
1248 SUNDIALS_EXPORT int CPodeGetQuadErrWeights(void *cpode_mem, N_Vector eQweight);
1249 
1250 /*
1251  * -----------------------------------------------------------------
1252  * The following function returns the name of the constant
1253  * associated with a CPODES return flag
1254  * -----------------------------------------------------------------
1255  */
1256 
1257 SUNDIALS_EXPORT char *CPodeGetReturnFlagName(int flag);
1258 
1259 /*
1260  * -----------------------------------------------------------------
1261  * Function : CPodeFree
1262  * -----------------------------------------------------------------
1263  * CPodeFree frees the problem memory cpode_mem allocated by
1264  * CPodeCreate and CPodeInit.  Its only argument is the pointer
1265  * cpode_mem returned by CPodeCreate.
1266  * -----------------------------------------------------------------
1267  */
1268 
1269 SUNDIALS_EXPORT void CPodeFree(void **cpode_mem);
1270 
1271 /*
1272  * -----------------------------------------------------------------
1273  * Function : CPodeQuadFree
1274  * -----------------------------------------------------------------
1275  * CPodeQuadFree frees the problem memory in cpode_mem allocated
1276  * for quadrature integration. Its only argument is the pointer
1277  * cpode_mem returned by CPodeCreate.
1278  * Note that CPodeQuadFree is called by CPodeFree.
1279  * -----------------------------------------------------------------
1280  */
1281 
1282 SUNDIALS_EXPORT void CPodeQuadFree(void *cpode_mem);
1283 
1284 #ifdef __cplusplus
1285 }
1286 #endif
1287 
1288 #endif
1289