1 /*
2  * -----------------------------------------------------------------
3  * $Revision: 1.4 $
4  * $Date: 2007/04/06 20:33:24 $
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 implementation file for the optional input and output
14  * functions for the CPODES solver.
15  * -----------------------------------------------------------------
16  */
17 
18 #include <stdio.h>
19 #include <stdlib.h>
20 
21 #include "cpodes_private.h"
22 
23 #define lrw   (cp_mem->cp_lrw)
24 #define liw   (cp_mem->cp_liw)
25 #define lrw1  (cp_mem->cp_lrw1)
26 #define liw1  (cp_mem->cp_liw1)
27 #define lrw1Q (cp_mem->cp_lrw1Q)
28 #define liw1Q (cp_mem->cp_liw1Q)
29 
30 /*
31  * =================================================================
32  * CPODES optional input functions
33  * =================================================================
34  */
35 
36 /*
37  * -----------------------------------------------------------------
38  * Optional inputs to the main solver
39  * -----------------------------------------------------------------
40  */
41 
42 /*
43  * CPodeSetErrHandlerFn
44  *
45  * Specifies the error handler function
46  */
47 
CPodeSetErrHandlerFn(void * cpode_mem,CPErrHandlerFn ehfun,void * eh_data)48 int CPodeSetErrHandlerFn(void *cpode_mem, CPErrHandlerFn ehfun, void *eh_data)
49 {
50   CPodeMem cp_mem;
51 
52   if (cpode_mem==NULL) {
53     cpProcessError(NULL, CP_MEM_NULL, "CPODES", "CPodeSetErrHandlerFn", MSGCP_NO_MEM);
54     return(CP_MEM_NULL);
55   }
56   cp_mem = (CPodeMem) cpode_mem;
57 
58   cp_mem->cp_ehfun = ehfun;
59   cp_mem->cp_eh_data = eh_data;
60 
61   return(CP_SUCCESS);
62 }
63 
64 /*
65  * CPodeSetErrFile
66  *
67  * Specifies the FILE pointer for output (NULL means no messages)
68  */
69 
CPodeSetErrFile(void * cpode_mem,FILE * errfp)70 int CPodeSetErrFile(void *cpode_mem, FILE *errfp)
71 {
72   CPodeMem cp_mem;
73 
74   if (cpode_mem==NULL) {
75     cpProcessError(NULL, CP_MEM_NULL, "CPODES", "CPodeSetErrFile", MSGCP_NO_MEM);
76     return(CP_MEM_NULL);
77   }
78   cp_mem = (CPodeMem) cpode_mem;
79 
80   cp_mem->cp_errfp = errfp;
81 
82   return(CP_SUCCESS);
83 }
84 
85 /*
86  * CPodeSetMaxOrd
87  *
88  * Specifies the maximum method order
89  */
90 
CPodeSetMaxOrd(void * cpode_mem,int maxord)91 int CPodeSetMaxOrd(void *cpode_mem, int maxord)
92 {
93   CPodeMem cp_mem;
94   int qmax_alloc;
95 
96   if (cpode_mem==NULL) {
97     cpProcessError(NULL, CP_MEM_NULL, "CPODES", "CPodeSetMaxOrd", MSGCP_NO_MEM);
98     return(CP_MEM_NULL);
99   }
100   cp_mem = (CPodeMem) cpode_mem;
101 
102   if (maxord <= 0) {
103     cpProcessError(cp_mem, CP_ILL_INPUT, "CPODES", "CPodeSetMaxOrd", MSGCP_NEG_MAXORD);
104     return(CP_ILL_INPUT);
105   }
106 
107   /* Cannot increase maximum order beyond the value that
108      was used when allocating memory */
109   qmax_alloc = cp_mem->cp_qmax_alloc;
110 
111   if (maxord > qmax_alloc) {
112     cpProcessError(cp_mem, CP_ILL_INPUT, "CPODES", "CPodeSetMaxOrd", MSGCP_BAD_MAXORD);
113     return(CP_ILL_INPUT);
114   }
115 
116   cp_mem->cp_qmax = maxord;
117 
118   return(CP_SUCCESS);
119 }
120 
121 /*
122  * CPodeSetMaxNumSteps
123  *
124  * Specifies the maximum number of integration steps
125  */
126 
CPodeSetMaxNumSteps(void * cpode_mem,long int mxsteps)127 int CPodeSetMaxNumSteps(void *cpode_mem, long int mxsteps)
128 {
129   CPodeMem cp_mem;
130 
131   if (cpode_mem==NULL) {
132     cpProcessError(NULL, CP_MEM_NULL, "CPODES", "CPodeSetMaxNumSteps", MSGCP_NO_MEM);
133     return(CP_MEM_NULL);
134   }
135   cp_mem = (CPodeMem) cpode_mem;
136 
137   if (mxsteps < 0) {
138     cpProcessError(cp_mem, CP_ILL_INPUT, "CPODES", "CPodeSetMaxNumSteps", MSGCP_NEG_MXSTEPS);
139     return(CP_ILL_INPUT);
140   }
141 
142   /* Passing 0 sets the default */
143   if (mxsteps == 0)
144     cp_mem->cp_mxstep = MXSTEP_DEFAULT;
145   else
146     cp_mem->cp_mxstep = mxsteps;
147 
148   return(CP_SUCCESS);
149 }
150 
151 /*
152  * CPodeSetMaxHnilWarns
153  *
154  * Specifies the maximum number of warnings for small h
155  */
156 
CPodeSetMaxHnilWarns(void * cpode_mem,int mxhnil)157 int CPodeSetMaxHnilWarns(void *cpode_mem, int mxhnil)
158 {
159   CPodeMem cp_mem;
160 
161   if (cpode_mem==NULL) {
162     cpProcessError(NULL, CP_MEM_NULL, "CPODES", "CPodeSetMaxHnilWarns", MSGCP_NO_MEM);
163     return(CP_MEM_NULL);
164   }
165   cp_mem = (CPodeMem) cpode_mem;
166 
167   cp_mem->cp_mxhnil = mxhnil;
168 
169   return(CP_SUCCESS);
170 }
171 
172 /*
173  *CPodeSetStabLimDet
174  *
175  * Turns on/off the stability limit detection algorithm
176  */
177 
CPodeSetStabLimDet(void * cpode_mem,booleantype sldet)178 int CPodeSetStabLimDet(void *cpode_mem, booleantype sldet)
179 {
180   CPodeMem cp_mem;
181 
182   if (cpode_mem==NULL) {
183     cpProcessError(NULL, CP_MEM_NULL, "CPODES", "CPodeSetStabLimDet", MSGCP_NO_MEM);
184     return(CP_MEM_NULL);
185   }
186   cp_mem = (CPodeMem) cpode_mem;
187 
188   if( sldet && (cp_mem->cp_lmm_type != CP_BDF) ) {
189     cpProcessError(cp_mem, CP_ILL_INPUT, "CPODES", "CPodeSetStabLimDet", MSGCP_SET_SLDET);
190     return(CP_ILL_INPUT);
191   }
192 
193   cp_mem->cp_sldeton = sldet;
194 
195   return(CP_SUCCESS);
196 }
197 
198 /*
199  * CPodeSetInitStep
200  *
201  * Specifies the initial step size
202  */
203 
CPodeSetInitStep(void * cpode_mem,realtype hin)204 int CPodeSetInitStep(void *cpode_mem, realtype hin)
205 {
206   CPodeMem cp_mem;
207 
208   if (cpode_mem==NULL) {
209     cpProcessError(NULL, CP_MEM_NULL, "CPODES", "CPodeSetInitStep", MSGCP_NO_MEM);
210     return(CP_MEM_NULL);
211   }
212   cp_mem = (CPodeMem) cpode_mem;
213 
214   cp_mem->cp_hin = hin;
215 
216   return(CP_SUCCESS);
217 }
218 
219 /*
220  * CPodeSetMinStep
221  *
222  * Specifies the minimum step size
223  */
224 
CPodeSetMinStep(void * cpode_mem,realtype hmin)225 int CPodeSetMinStep(void *cpode_mem, realtype hmin)
226 {
227   CPodeMem cp_mem;
228 
229   if (cpode_mem==NULL) {
230     cpProcessError(NULL, CP_MEM_NULL, "CPODES", "CPodeSetMinStep", MSGCP_NO_MEM);
231     return(CP_MEM_NULL);
232   }
233   cp_mem = (CPodeMem) cpode_mem;
234 
235   if (hmin<0) {
236     cpProcessError(cp_mem, CP_ILL_INPUT, "CPODES", "CPodeSetMinStep", MSGCP_NEG_HMIN);
237     return(CP_ILL_INPUT);
238   }
239 
240   /* Passing 0 sets hmin = zero */
241   if (hmin == ZERO) {
242     cp_mem->cp_hmin = HMIN_DEFAULT;
243     return(CP_SUCCESS);
244   }
245 
246   if (hmin * cp_mem->cp_hmax_inv > ONE) {
247     cpProcessError(cp_mem, CP_ILL_INPUT, "CPODES", "CPodeSetMinStep", MSGCP_BAD_HMIN_HMAX);
248     return(CP_ILL_INPUT);
249   }
250 
251   cp_mem->cp_hmin = hmin;
252 
253   return(CP_SUCCESS);
254 }
255 
256 /*
257  * CPodeSetMaxStep
258  *
259  * Specifies the maximum step size
260  */
261 
CPodeSetMaxStep(void * cpode_mem,realtype hmax)262 int CPodeSetMaxStep(void *cpode_mem, realtype hmax)
263 {
264   realtype hmax_inv;
265   CPodeMem cp_mem;
266 
267   if (cpode_mem==NULL) {
268     cpProcessError(NULL, CP_MEM_NULL, "CPODES", "CPodeSetMaxStep", MSGCP_NO_MEM);
269     return (CP_MEM_NULL);
270   }
271   cp_mem = (CPodeMem) cpode_mem;
272 
273   if (hmax < 0) {
274     cpProcessError(cp_mem, CP_ILL_INPUT, "CPODES", "CPodeSetMaxStep", MSGCP_NEG_HMAX);
275     return(CP_ILL_INPUT);
276   }
277 
278   /* Passing 0 sets hmax = infinity */
279   if (hmax == ZERO) {
280     cp_mem->cp_hmax_inv = HMAX_INV_DEFAULT;
281     return(CP_SUCCESS);
282   }
283 
284   hmax_inv = ONE/hmax;
285   if (hmax_inv * cp_mem->cp_hmin > ONE) {
286     cpProcessError(cp_mem, CP_ILL_INPUT, "CPODES", "CPodeSetMaxStep", MSGCP_BAD_HMIN_HMAX);
287     return(CP_ILL_INPUT);
288   }
289 
290   cp_mem->cp_hmax_inv = hmax_inv;
291 
292   return(CP_SUCCESS);
293 }
294 
295 /*
296  * CPodeSetStopTime
297  *
298  * Specifies the time beyond which the integration is not to
299  * proceed
300  */
301 
CPodeSetStopTime(void * cpode_mem,realtype tstop)302 int CPodeSetStopTime(void *cpode_mem, realtype tstop)
303 {
304   CPodeMem cp_mem;
305 
306   if (cpode_mem==NULL) {
307     cpProcessError(NULL, CP_MEM_NULL, "CPODES", "CPodeSetStopTime", MSGCP_NO_MEM);
308     return (CP_MEM_NULL);
309   }
310   cp_mem = (CPodeMem) cpode_mem;
311 
312   cp_mem->cp_tstop = tstop;
313   cp_mem->cp_tstopset = TRUE;
314 
315   return(CP_SUCCESS);
316 }
317 
318 /*
319  * CPodeSetMaxErrTestFails
320  *
321  * Specifies the maximum number of error test failures during one
322  * step try.
323  */
324 
CPodeSetMaxErrTestFails(void * cpode_mem,int maxnef)325 int CPodeSetMaxErrTestFails(void *cpode_mem, int maxnef)
326 {
327   CPodeMem cp_mem;
328 
329   if (cpode_mem==NULL) {
330     cpProcessError(NULL, CP_MEM_NULL, "CPODES", "CPodeSetMaxErrTestFails", MSGCP_NO_MEM);
331     return (CP_MEM_NULL);
332   }
333   cp_mem = (CPodeMem) cpode_mem;
334 
335   cp_mem->cp_maxnef = maxnef;
336 
337   return(CP_SUCCESS);
338 }
339 
340 /*
341  * CPodeSetMaxConvFails
342  *
343  * Specifies the maximum number of nonlinear convergence failures
344  * during one step try.
345  */
346 
CPodeSetMaxConvFails(void * cpode_mem,int maxncf)347 int CPodeSetMaxConvFails(void *cpode_mem, int maxncf)
348 {
349   CPodeMem cp_mem;
350 
351   if (cpode_mem==NULL) {
352     cpProcessError(NULL, CP_MEM_NULL, "CPODES", "CPodeSetMaxConvFails", MSGCP_NO_MEM);
353     return (CP_MEM_NULL);
354   }
355   cp_mem = (CPodeMem) cpode_mem;
356 
357   cp_mem->cp_maxncf = maxncf;
358 
359   return(CP_SUCCESS);
360 }
361 
362 /*
363  * CPodeSetMaxNonlinIters
364  *
365  * Specifies the maximum number of nonlinear iterations during
366  * one solve.
367  */
368 
CPodeSetMaxNonlinIters(void * cpode_mem,int maxcor)369 int CPodeSetMaxNonlinIters(void *cpode_mem, int maxcor)
370 {
371   CPodeMem cp_mem;
372 
373   if (cpode_mem==NULL) {
374     cpProcessError(NULL, CP_MEM_NULL, "CPODES", "CPodeSetMaxNonlinIters", MSGCP_NO_MEM);
375     return (CP_MEM_NULL);
376   }
377   cp_mem = (CPodeMem) cpode_mem;
378 
379   cp_mem->cp_maxcor = maxcor;
380 
381   return(CP_SUCCESS);
382 }
383 
384 /*
385  * CPodeSetNonlinConvCoef
386  *
387  * Specifies the coeficient in the nonlinear solver convergence
388  * test
389  */
390 
CPodeSetNonlinConvCoef(void * cpode_mem,realtype nlscoef)391 int CPodeSetNonlinConvCoef(void *cpode_mem, realtype nlscoef)
392 {
393   CPodeMem cp_mem;
394 
395   if (cpode_mem==NULL) {
396     cpProcessError(NULL, CP_MEM_NULL, "CPODES", "CPodeSetNonlinConvCoef", MSGCP_NO_MEM);
397     return(CP_MEM_NULL);
398   }
399   cp_mem = (CPodeMem) cpode_mem;
400 
401   cp_mem->cp_nlscoef = nlscoef;
402 
403   return(CP_SUCCESS);
404 }
405 
406 /*
407  * CPodeSetTolerances
408  *
409  * Changes the integration tolerances between calls to CPode().
410  * Here, only CP_SS or CP_SV are allowed.
411  */
412 
CPodeSetTolerances(void * cpode_mem,int tol_type,realtype reltol,void * abstol)413 int CPodeSetTolerances(void *cpode_mem,
414                        int tol_type, realtype reltol, void *abstol)
415 {
416   CPodeMem cp_mem;
417   booleantype neg_abstol;
418 
419   /* Check CPODES memory */
420   if (cpode_mem==NULL) {
421     cpProcessError(NULL, CP_MEM_NULL, "CPODES", "CPodeSetTolerances", MSGCP_NO_MEM);
422     return(CP_MEM_NULL);
423   }
424   cp_mem = (CPodeMem) cpode_mem;
425 
426   /* Check if cpode_mem was allocated */
427   if (cp_mem->cp_MallocDone == FALSE) {
428     cpProcessError(cp_mem, CP_NO_MALLOC, "CPODES", "CPodeSetTolerances", MSGCP_NO_MALLOC);
429     return(CP_NO_MALLOC);
430   }
431 
432   /* Check inputs for legal values */
433   if ( (tol_type != CP_SS) && (tol_type != CP_SV) ) {
434     cpProcessError(cp_mem, CP_ILL_INPUT, "CPODES", "CPodeSetTolerances", MSGCP_BAD_ITOL);
435     return(CP_ILL_INPUT);
436   }
437   if (abstol == NULL) {
438     cpProcessError(cp_mem, CP_ILL_INPUT, "CPODES", "CPodeSetTolerances", MSGCP_NULL_ABSTOL);
439     return(CP_ILL_INPUT);
440   }
441 
442   /* Check positivity of tolerances */
443   if (reltol < ZERO) {
444     cpProcessError(cp_mem, CP_ILL_INPUT, "CPODES", "CPodeSetTolerances", MSGCP_BAD_RELTOL);
445     return(CP_ILL_INPUT);
446   }
447 
448   if (tol_type == CP_SS)
449     neg_abstol = (*((realtype *)abstol) < ZERO);
450   else
451     neg_abstol = (N_VMin((N_Vector)abstol) < ZERO);
452 
453   if (neg_abstol) {
454     cpProcessError(cp_mem, CP_ILL_INPUT, "CPODES", "CPodeSetTolerances", MSGCP_BAD_ABSTOL);
455     return(CP_ILL_INPUT);
456   }
457 
458   /* Copy tolerances into memory */
459   cp_mem->cp_tol_type = tol_type;
460   cp_mem->cp_reltol   = reltol;
461 
462   if (tol_type == CP_SS) {
463     cp_mem->cp_Sabstol = *((realtype *)abstol);
464   } else {
465     if ( !(cp_mem->cp_VabstolMallocDone) ) {
466       cp_mem->cp_Vabstol = N_VClone(cp_mem->cp_ewt);
467       lrw += lrw1;
468       liw += liw1;
469       cp_mem->cp_VabstolMallocDone = TRUE;
470     }
471     N_VScale(ONE, (N_Vector)abstol, cp_mem->cp_Vabstol);
472   }
473 
474   return(CP_SUCCESS);
475 }
476 
477 /*
478  * CPodeSetEwtFn
479  *
480  * Specifies the user-provided function efun of type EwtSet
481  * and a data pointer for efun.
482  */
483 
CPodeSetEwtFn(void * cpode_mem,CPEwtFn efun,void * e_data)484 int CPodeSetEwtFn(void *cpode_mem, CPEwtFn efun, void *e_data)
485 {
486   CPodeMem cp_mem;
487 
488   if (cpode_mem==NULL) {
489     cpProcessError(NULL, CP_MEM_NULL, "CPODES", "CPodeSetEwtFn", MSGCP_NO_MEM);
490     return(CP_MEM_NULL);
491   }
492   cp_mem = (CPodeMem) cpode_mem;
493 
494   cp_mem->cp_efun   = efun;
495   cp_mem->cp_e_data = e_data;
496 
497   return(CP_SUCCESS);
498 }
499 
500 /*
501  * CPodeSetRootDirection
502  *
503  * Specifies the direction of zero-crossings to be monitored.
504  * The default is to monitor both crossings.
505  */
506 
CPodeSetRootDirection(void * cpode_mem,int * rootdir)507 int CPodeSetRootDirection(void *cpode_mem, int *rootdir)
508 {
509   CPodeMem cp_mem;
510   int i, nrt;
511 
512   if (cpode_mem==NULL) {
513     cpProcessError(NULL, CP_MEM_NULL, "CPODES", "CPodeSetRootDirection", MSGCP_NO_MEM);
514     return(CP_MEM_NULL);
515   }
516 
517   cp_mem = (CPodeMem) cpode_mem;
518 
519   nrt = cp_mem->cp_nrtfn;
520   if (nrt==0) {
521     cpProcessError(NULL, CP_ILL_INPUT, "CPODES", "CPodeSetRootDirection", MSGCP_NO_ROOT);
522     return(CP_ILL_INPUT);
523   }
524 
525   for(i=0; i<nrt; i++) cp_mem->cp_rootdir[i] = rootdir[i];
526 
527   return(CP_SUCCESS);
528 }
529 
530 /*
531  * -----------------------------------------------------------------
532  * Optional inputs for projection step
533  * -----------------------------------------------------------------
534  */
535 
536 /*
537  *
538  * CPodeSetProjUpdateErrEst
539  */
540 
CPodeSetProjUpdateErrEst(void * cpode_mem,booleantype proj_err)541 int CPodeSetProjUpdateErrEst(void *cpode_mem, booleantype proj_err)
542 {
543   CPodeMem cp_mem;
544 
545   if (cpode_mem==NULL) {
546     cpProcessError(NULL, CP_MEM_NULL, "CPODES", "CPodeSetProjUpdateErrEst", MSGCP_NO_MEM);
547     return(CP_MEM_NULL);
548   }
549   cp_mem = (CPodeMem) cpode_mem;
550 
551   cp_mem->cp_project_err = proj_err;
552 
553   return(CP_SUCCESS);
554 }
555 
556 /*
557  * CPodeSetProjFrequency
558  */
559 
CPodeSetProjFrequency(void * cpode_mem,long int proj_freq)560 int CPodeSetProjFrequency(void *cpode_mem, long int proj_freq)
561 {
562   CPodeMem cp_mem;
563 
564   if (cpode_mem==NULL) {
565     cpProcessError(NULL, CP_MEM_NULL, "CPODES", "CPodeSetProjFrequency", MSGCP_NO_MEM);
566     return(CP_MEM_NULL);
567   }
568   cp_mem = (CPodeMem) cpode_mem;
569 
570   if (proj_freq < 0) {
571     cpProcessError(cp_mem, CP_ILL_INPUT, "CPODES", "CPodeSetProjFrequency", MSGCP_BAD_FREQ);
572     return(CP_ILL_INPUT);
573   }
574 
575   cp_mem->cp_proj_freq = proj_freq;
576 
577   return(CP_SUCCESS);
578 }
579 
580 /*
581  * CPodeSetProjTestCnstr
582  */
583 
CPodeSetProjTestCnstr(void * cpode_mem,booleantype test_cnstr)584 int CPodeSetProjTestCnstr(void *cpode_mem, booleantype test_cnstr)
585 {
586   CPodeMem cp_mem;
587 
588   if (cpode_mem==NULL) {
589     cpProcessError(NULL, CP_MEM_NULL, "CPODES", "CPodeSetProjTestCnstr", MSGCP_NO_MEM);
590     return(CP_MEM_NULL);
591   }
592   cp_mem = (CPodeMem) cpode_mem;
593 
594   cp_mem->cp_test_cnstr = test_cnstr;
595 
596   return(CP_SUCCESS);
597 }
598 
599 /*
600  * CPodeSetProjLsetupFreq
601  *
602  * Spcifies the frequency of constraint Jacobian evaluations
603  * (actually, the maximum number of steps allowed between calls
604  * to the linear solver's setup function).
605  */
606 
CPodeSetProjLsetupFreq(void * cpode_mem,long int lset_freq)607 int CPodeSetProjLsetupFreq(void *cpode_mem, long int lset_freq)
608 {
609   CPodeMem cp_mem;
610 
611   if (cpode_mem==NULL) {
612     cpProcessError(NULL, CP_MEM_NULL, "CPODES", "CPodeSetProjLsetupFreq", MSGCP_NO_MEM);
613     return(CP_MEM_NULL);
614   }
615   cp_mem = (CPodeMem) cpode_mem;
616 
617   if (lset_freq <= 0) {
618     cpProcessError(cp_mem, CP_ILL_INPUT, "CPODES", "CPodeSetProjLsetupFreq", MSGCP_BAD_LSFREQ);
619     return(CP_ILL_INPUT);
620   }
621 
622   cp_mem->cp_lsetupP_freq = lset_freq;
623 
624   return(CP_SUCCESS);
625 }
626 
627 /*
628  * CPodeSetProjNonlinConvCoef
629  *
630  * Specifies the coeficient in the projection nonlinear solver convergence
631  * test
632  */
633 
CPodeSetProjNonlinConvCoef(void * cpode_mem,realtype prjcoef)634 int CPodeSetProjNonlinConvCoef(void *cpode_mem, realtype prjcoef)
635 {
636   CPodeMem cp_mem;
637 
638   if (cpode_mem==NULL) {
639     cpProcessError(NULL, CP_MEM_NULL, "CPODES", "CPodeSetProjNonlinConvCoef", MSGCP_NO_MEM);
640     return(CP_MEM_NULL);
641   }
642   cp_mem = (CPodeMem) cpode_mem;
643 
644   cp_mem->cp_prjcoef = prjcoef;
645 
646   return(CP_SUCCESS);
647 }
648 
649 /*
650  * -----------------------------------------------------------------
651  * Optional inputs for quadrature integration
652  * -----------------------------------------------------------------
653  */
654 
CPodeSetQuadErrCon(void * cpode_mem,booleantype errconQ,int tol_typeQ,realtype reltolQ,void * abstolQ)655 int CPodeSetQuadErrCon(void *cpode_mem, booleantype errconQ,
656                        int tol_typeQ, realtype reltolQ, void *abstolQ)
657 {
658   CPodeMem cp_mem;
659   booleantype neg_abstol;
660 
661   if (cpode_mem==NULL) {
662     cpProcessError(NULL, CP_MEM_NULL, "CPODES", "CPodeSetQuadErrCon", MSGCP_NO_MEM);
663     return(CP_MEM_NULL);
664   }
665 
666   cp_mem = (CPodeMem) cpode_mem;
667 
668   cp_mem->cp_errconQ = errconQ;
669 
670   /* Ckeck if quadrature was initialized? */
671 
672   if (cp_mem->cp_quadMallocDone == FALSE) {
673     cpProcessError(cp_mem, CP_NO_QUAD, "CPODES", "CPodeSetQuadErrCon", MSGCP_NO_QUAD);
674     return(CP_NO_QUAD);
675   }
676 
677   /* Check inputs */
678 
679   if(errconQ == FALSE) {
680     if (cp_mem->cp_VabstolQMallocDone) {
681       N_VDestroy(cp_mem->cp_VabstolQ);
682       lrw -= lrw1Q;
683       liw -= liw1Q;
684       cp_mem->cp_VabstolQMallocDone = FALSE;
685     }
686     return(CP_SUCCESS);
687   }
688 
689   if ((tol_typeQ != CP_SS) && (tol_typeQ != CP_SV)) {
690     cpProcessError(cp_mem, CP_ILL_INPUT, "CPODES", "CPodeSetQuadErrCon", MSGCP_BAD_ITOLQ);
691     return(CP_ILL_INPUT);
692   }
693 
694   if (abstolQ == NULL) {
695     cpProcessError(cp_mem, CP_ILL_INPUT, "CPODES", "CPodeSetQuadErrCon", MSGCP_NULL_ABSTOLQ);
696     return(CP_ILL_INPUT);
697   }
698 
699   if (reltolQ < ZERO) {
700     cpProcessError(cp_mem, CP_ILL_INPUT, "CPODES", "CPodeSetQuadErrCon", MSGCP_BAD_RELTOLQ);
701     return(CP_ILL_INPUT);
702   }
703 
704   if (tol_typeQ == CP_SS)
705     neg_abstol = (*((realtype *)abstolQ) < ZERO);
706   else
707     neg_abstol = (N_VMin((N_Vector)abstolQ) < ZERO);
708 
709   if (neg_abstol) {
710     cpProcessError(cp_mem, CP_ILL_INPUT, "CPODES", "CPodeSetQuadErrCon", MSGCP_BAD_ABSTOLQ);
711     return(CP_ILL_INPUT);
712   }
713 
714   /* See if we need to free or allocate memory */
715 
716   if ( (tol_typeQ != CP_SV) && (cp_mem->cp_VabstolQMallocDone) ) {
717     N_VDestroy(cp_mem->cp_VabstolQ);
718     lrw -= lrw1Q;
719     liw -= liw1Q;
720     cp_mem->cp_VabstolQMallocDone = FALSE;
721   }
722 
723   if ( (tol_typeQ == CP_SV) && !(cp_mem->cp_VabstolQMallocDone) ) {
724     cp_mem->cp_VabstolQ = N_VClone(cp_mem->cp_tempvQ);
725     lrw += lrw1Q;
726     liw += liw1Q;
727     cp_mem->cp_VabstolQMallocDone = TRUE;
728   }
729 
730   /* Copy tolerances into memory */
731 
732   cp_mem->cp_tol_typeQ = tol_typeQ;
733   cp_mem->cp_reltolQ   = reltolQ;
734 
735   if (tol_typeQ == CP_SS)
736     cp_mem->cp_SabstolQ = *((realtype *)abstolQ);
737   else
738     N_VScale(ONE, (N_Vector)abstolQ, cp_mem->cp_VabstolQ);
739 
740   return(CP_SUCCESS);
741 }
742 
743 /*
744  * =================================================================
745  * CPODE optional output functions
746  * =================================================================
747  */
748 
749 /*
750  * -----------------------------------------------------------------
751  * Optional outputs for IC calculation
752  * -----------------------------------------------------------------
753  */
754 
CPodeGetConsistentIC(void * cpode_mem,N_Vector yy0,N_Vector yp0)755 int CPodeGetConsistentIC(void *cpode_mem, N_Vector yy0, N_Vector yp0)
756 {
757   CPodeMem cp_mem;
758 
759   if (cpode_mem==NULL) {
760     cpProcessError(NULL, CP_MEM_NULL, "CPODES", "CPodeGetConsistentIC", MSGCP_NO_MEM);
761     return(CP_MEM_NULL);
762   }
763   cp_mem = (CPodeMem) cpode_mem;
764 
765   if (cp_mem->cp_next_q != 0) {
766     cpProcessError(cp_mem, CP_ILL_INPUT, "CPODES", "CPodeGetConsistentIC", MSGCP_TOO_LATE);
767     return(CP_ILL_INPUT);
768   }
769 
770   if (yy0 != NULL) N_VScale(ONE, cp_mem->cp_zn[0], yy0);
771   if (cp_mem->cp_ode_type == CP_IMPL)
772     if (yp0 != NULL) N_VScale(ONE, cp_mem->cp_zn[1], yp0);
773 
774   return(CP_SUCCESS);
775 }
776 
777 /*
778  * -----------------------------------------------------------------
779  * Optional outputs for integration
780  * -----------------------------------------------------------------
781  */
782 
783 /*
784  * CPodeGetNumSteps
785  *
786  * Returns the current number of integration steps
787  */
788 
CPodeGetNumSteps(void * cpode_mem,long int * nsteps)789 int CPodeGetNumSteps(void *cpode_mem, long int *nsteps)
790 {
791   CPodeMem cp_mem;
792 
793   if (cpode_mem==NULL) {
794     cpProcessError(NULL, CP_MEM_NULL, "CPODES", "CPodeGetNumSteps", MSGCP_NO_MEM);
795     return(CP_MEM_NULL);
796   }
797   cp_mem = (CPodeMem) cpode_mem;
798 
799   *nsteps = cp_mem->cp_nst;
800 
801   return(CP_SUCCESS);
802 }
803 
804 /*
805  * CPodeGetNumFctEvals
806  *
807  * Returns the current number of calls to f
808  */
809 
CPodeGetNumFctEvals(void * cpode_mem,long int * nfevals)810 int CPodeGetNumFctEvals(void *cpode_mem, long int *nfevals)
811 {
812   CPodeMem cp_mem;
813 
814   if (cpode_mem==NULL) {
815     cpProcessError(NULL, CP_MEM_NULL, "CPODES", "CPodeGetNumFctEvals", MSGCP_NO_MEM);
816     return(CP_MEM_NULL);
817   }
818   cp_mem = (CPodeMem) cpode_mem;
819 
820   *nfevals = cp_mem->cp_nfe;
821 
822   return(CP_SUCCESS);
823 }
824 
825 /*
826  * CPodeGetNumLinSolvSetups
827  *
828  * Returns the current number of calls to the linear solver setup routine
829  */
830 
CPodeGetNumLinSolvSetups(void * cpode_mem,long int * nlinsetups)831 int CPodeGetNumLinSolvSetups(void *cpode_mem, long int *nlinsetups)
832 {
833   CPodeMem cp_mem;
834 
835   if (cpode_mem==NULL) {
836     cpProcessError(NULL, CP_MEM_NULL, "CPODES", "CPodeGetNumLinSolvSetups", MSGCP_NO_MEM);
837     return(CP_MEM_NULL);
838   }
839   cp_mem = (CPodeMem) cpode_mem;
840 
841   *nlinsetups = cp_mem->cp_nsetups;
842 
843   return(CP_SUCCESS);
844 }
845 
846 /*
847  * CPodeGetNumErrTestFails
848  *
849  * Returns the current number of error test failures
850  */
851 
CPodeGetNumErrTestFails(void * cpode_mem,long int * netfails)852 int CPodeGetNumErrTestFails(void *cpode_mem, long int *netfails)
853 {
854   CPodeMem cp_mem;
855 
856   if (cpode_mem==NULL) {
857     cpProcessError(NULL, CP_MEM_NULL, "CPODES", "CPodeGetNumErrTestFails", MSGCP_NO_MEM);
858     return(CP_MEM_NULL);
859   }
860   cp_mem = (CPodeMem) cpode_mem;
861 
862   *netfails = cp_mem->cp_netf;
863 
864   return(CP_SUCCESS);
865 }
866 
867 /*
868  * CPodeGetLastOrder
869  *
870  * Returns the order on the last successful step
871  */
872 
CPodeGetLastOrder(void * cpode_mem,int * qlast)873 int CPodeGetLastOrder(void *cpode_mem, int *qlast)
874 {
875   CPodeMem cp_mem;
876 
877   if (cpode_mem==NULL) {
878     cpProcessError(NULL, CP_MEM_NULL, "CPODES", "CPodeGetLastOrder", MSGCP_NO_MEM);
879     return(CP_MEM_NULL);
880   }
881   cp_mem = (CPodeMem) cpode_mem;
882 
883   *qlast = cp_mem->cp_qu;
884 
885   return(CP_SUCCESS);
886 }
887 
888 /*
889  * CPodeGetCurrentOrder
890  *
891  * Returns the order to be attempted on the next step
892  */
893 
CPodeGetCurrentOrder(void * cpode_mem,int * qcur)894 int CPodeGetCurrentOrder(void *cpode_mem, int *qcur)
895 {
896   CPodeMem cp_mem;
897 
898   if (cpode_mem==NULL) {
899     cpProcessError(NULL, CP_MEM_NULL, "CPODES", "CPodeGetCurrentOrder", MSGCP_NO_MEM);
900     return(CP_MEM_NULL);
901   }
902   cp_mem = (CPodeMem) cpode_mem;
903 
904   *qcur = cp_mem->cp_next_q;
905 
906   return(CP_SUCCESS);
907 }
908 
909 /*
910  * CPodeGetNumStabLimOrderReds
911  *
912  * Returns the number of order reductions triggered by the stability
913  * limit detection algorithm
914  */
915 
CPodeGetNumStabLimOrderReds(void * cpode_mem,long int * nslred)916 int CPodeGetNumStabLimOrderReds(void *cpode_mem, long int *nslred)
917 {
918   CPodeMem cp_mem;
919 
920   if (cpode_mem==NULL) {
921     cpProcessError(NULL, CP_MEM_NULL, "CPODES", "CPodeGetNumStabLimOrderReds", MSGCP_NO_MEM);
922     return(CP_MEM_NULL);
923   }
924   cp_mem = (CPodeMem) cpode_mem;
925 
926   if (cp_mem->cp_sldeton==FALSE)
927     *nslred = 0;
928   else
929     *nslred = cp_mem->cp_nor;
930 
931   return(CP_SUCCESS);
932 }
933 
934 /*
935  * CPodeGetActualInitStep
936  *
937  * Returns the step size used on the first step
938  */
939 
CPodeGetActualInitStep(void * cpode_mem,realtype * hinused)940 int CPodeGetActualInitStep(void *cpode_mem, realtype *hinused)
941 {
942   CPodeMem cp_mem;
943 
944   if (cpode_mem==NULL) {
945     cpProcessError(NULL, CP_MEM_NULL, "CPODES", "CPodeGetActualInitStep", MSGCP_NO_MEM);
946     return(CP_MEM_NULL);
947   }
948   cp_mem = (CPodeMem) cpode_mem;
949 
950   *hinused = cp_mem->cp_h0u;
951 
952   return(CP_SUCCESS);
953 }
954 
955 /*
956  * CPodeGetLastStep
957  *
958  * Returns the step size used on the last successful step
959  */
960 
CPodeGetLastStep(void * cpode_mem,realtype * hlast)961 int CPodeGetLastStep(void *cpode_mem, realtype *hlast)
962 {
963   CPodeMem cp_mem;
964 
965   if (cpode_mem==NULL) {
966     cpProcessError(NULL, CP_MEM_NULL, "CPODES", "CPodeGetLastStep", MSGCP_NO_MEM);
967     return(CP_MEM_NULL);
968   }
969   cp_mem = (CPodeMem) cpode_mem;
970 
971   *hlast = cp_mem->cp_hu;
972 
973   return(CP_SUCCESS);
974 }
975 
976 /*
977  * CPodeGetCurrentStep
978  *
979  * Returns the step size to be attempted on the next step
980  */
981 
CPodeGetCurrentStep(void * cpode_mem,realtype * hcur)982 int CPodeGetCurrentStep(void *cpode_mem, realtype *hcur)
983 {
984   CPodeMem cp_mem;
985 
986   if (cpode_mem==NULL) {
987     cpProcessError(NULL, CP_MEM_NULL, "CPODES", "CPodeGetCurrentStep", MSGCP_NO_MEM);
988     return(CP_MEM_NULL);
989   }
990   cp_mem = (CPodeMem) cpode_mem;
991 
992   *hcur = cp_mem->cp_next_h;
993 
994   return(CP_SUCCESS);
995 }
996 
997 /*
998  * CPodeGetCurrentTime
999  *
1000  * Returns the current value of the independent variable
1001  */
1002 
CPodeGetCurrentTime(void * cpode_mem,realtype * tcur)1003 int CPodeGetCurrentTime(void *cpode_mem, realtype *tcur)
1004 {
1005   CPodeMem cp_mem;
1006 
1007   if (cpode_mem==NULL) {
1008     cpProcessError(NULL, CP_MEM_NULL, "CPODES", "CPodeGetCurrentTime", MSGCP_NO_MEM);
1009     return(CP_MEM_NULL);
1010   }
1011   cp_mem = (CPodeMem) cpode_mem;
1012 
1013   *tcur = cp_mem->cp_tn;
1014 
1015   return(CP_SUCCESS);
1016 }
1017 
1018 /*
1019  * CPodeGetTolScaleFactor
1020  *
1021  * Returns a suggested factor for scaling tolerances
1022  */
1023 
CPodeGetTolScaleFactor(void * cpode_mem,realtype * tolsfact)1024 int CPodeGetTolScaleFactor(void *cpode_mem, realtype *tolsfact)
1025 {
1026   CPodeMem cp_mem;
1027 
1028   if (cpode_mem==NULL) {
1029     cpProcessError(NULL, CP_MEM_NULL, "CPODES", "CPodeGetTolScaleFactor", MSGCP_NO_MEM);
1030     return(CP_MEM_NULL);
1031   }
1032   cp_mem = (CPodeMem) cpode_mem;
1033 
1034   *tolsfact = cp_mem->cp_tolsf;
1035 
1036   return(CP_SUCCESS);
1037 }
1038 
1039 /*
1040  * CPodeGetErrWeights
1041  *
1042  * This routine returns the current weight vector.
1043  */
1044 
CPodeGetErrWeights(void * cpode_mem,N_Vector eweight)1045 int CPodeGetErrWeights(void *cpode_mem, N_Vector eweight)
1046 {
1047   CPodeMem cp_mem;
1048 
1049   if (cpode_mem==NULL) {
1050     cpProcessError(NULL, CP_MEM_NULL, "CPODES", "CPodeGetErrWeights", MSGCP_NO_MEM);
1051     return(CP_MEM_NULL);
1052   }
1053   cp_mem = (CPodeMem) cpode_mem;
1054 
1055   N_VScale(ONE, cp_mem->cp_ewt, eweight);
1056 
1057   return(CP_SUCCESS);
1058 }
1059 
1060 /*
1061  * CPodeGetEstLocalErrors
1062  *
1063  * Returns an estimate of the local error
1064  */
1065 
CPodeGetEstLocalErrors(void * cpode_mem,N_Vector ele)1066 int CPodeGetEstLocalErrors(void *cpode_mem, N_Vector ele)
1067 {
1068   CPodeMem cp_mem;
1069 
1070   if (cpode_mem==NULL) {
1071     cpProcessError(NULL, CP_MEM_NULL, "CPODES", "CPodeGetEstLocalErrors", MSGCP_NO_MEM);
1072     return(CP_MEM_NULL);
1073   }
1074   cp_mem = (CPodeMem) cpode_mem;
1075 
1076   N_VScale(ONE, cp_mem->cp_acor, ele);
1077 
1078   return(CP_SUCCESS);
1079 }
1080 
1081 /*
1082  * CPodeGetWorkSpace
1083  *
1084  * Returns integrator work space requirements
1085  */
1086 
CPodeGetWorkSpace(void * cpode_mem,long int * lenrw,long int * leniw)1087 int CPodeGetWorkSpace(void *cpode_mem, long int *lenrw, long int *leniw)
1088 {
1089   CPodeMem cp_mem;
1090 
1091   if (cpode_mem==NULL) {
1092     cpProcessError(NULL, CP_MEM_NULL, "CPODES", "CPodeGetWorkSpace", MSGCP_NO_MEM);
1093     return(CP_MEM_NULL);
1094   }
1095   cp_mem = (CPodeMem) cpode_mem;
1096 
1097   *leniw = liw;
1098   *lenrw = lrw;
1099 
1100   return(CP_SUCCESS);
1101 }
1102 
1103 /*
1104  * CPodeGetIntegratorStats
1105  *
1106  * Returns integrator statistics
1107  */
1108 
CPodeGetIntegratorStats(void * cpode_mem,long int * nsteps,long int * nfevals,long int * nlinsetups,long int * netfails,int * qlast,int * qcur,realtype * hinused,realtype * hlast,realtype * hcur,realtype * tcur)1109 int CPodeGetIntegratorStats(void *cpode_mem, long int *nsteps, long int *nfevals,
1110                             long int *nlinsetups, long int *netfails, int *qlast,
1111                             int *qcur, realtype *hinused, realtype *hlast,
1112                             realtype *hcur, realtype *tcur)
1113 {
1114   CPodeMem cp_mem;
1115 
1116   if (cpode_mem==NULL) {
1117     cpProcessError(NULL, CP_MEM_NULL, "CPODES", "CPodeGetIntegratorStats", MSGCP_NO_MEM);
1118     return(CP_MEM_NULL);
1119   }
1120   cp_mem = (CPodeMem) cpode_mem;
1121 
1122   *nsteps = cp_mem->cp_nst;
1123   *nfevals = cp_mem->cp_nfe;
1124   *nlinsetups = cp_mem->cp_nsetups;
1125   *netfails = cp_mem->cp_netf;
1126   *qlast = cp_mem->cp_qu;
1127   *qcur = cp_mem->cp_next_q;
1128   *hinused = cp_mem->cp_h0u;
1129   *hlast = cp_mem->cp_hu;
1130   *hcur = cp_mem->cp_next_h;
1131   *tcur = cp_mem->cp_tn;
1132 
1133   return(CP_SUCCESS);
1134 }
1135 
1136 /*
1137  * CPodeGetNumGEvals
1138  *
1139  * Returns the current number of calls to g (for rootfinding)
1140  */
1141 
CPodeGetNumGEvals(void * cpode_mem,long int * ngevals)1142 int CPodeGetNumGEvals(void *cpode_mem, long int *ngevals)
1143 {
1144   CPodeMem cp_mem;
1145 
1146   if (cpode_mem==NULL) {
1147     cpProcessError(NULL, CP_MEM_NULL, "CPODES", "CPodeGetNumGEvals", MSGCP_NO_MEM);
1148     return(CP_MEM_NULL);
1149   }
1150   cp_mem = (CPodeMem) cpode_mem;
1151 
1152   *ngevals = cp_mem->cp_nge;
1153 
1154   return(CP_SUCCESS);
1155 }
1156 
1157 /*
1158  * CPodeGetRootInfo
1159  *
1160  * Returns pointer to array rootsfound showing roots found
1161  */
1162 
CPodeGetRootInfo(void * cpode_mem,int * rootsfound)1163 int CPodeGetRootInfo(void *cpode_mem, int *rootsfound)
1164 {
1165   CPodeMem cp_mem;
1166   int i, nrt;
1167 
1168   if (cpode_mem==NULL) {
1169     cpProcessError(NULL, CP_MEM_NULL, "CPODES", "CPodeGetRootInfo", MSGCP_NO_MEM);
1170     return(CP_MEM_NULL);
1171   }
1172   cp_mem = (CPodeMem) cpode_mem;
1173 
1174   nrt = cp_mem->cp_nrtfn;
1175 
1176   for (i=0; i<nrt; i++) rootsfound[i] = cp_mem->cp_iroots[i];
1177 
1178   return(CP_SUCCESS);
1179 }
1180 
1181 /*
1182  * CPodeGetRootWindow
1183  *
1184  * Returns the window (tLo,tHi] within which root(s) were found. The result
1185  * is meaningless unless this is called right after a root-found return.
1186  */
1187 
CPodeGetRootWindow(void * cpode_mem,realtype * tLo,realtype * tHi)1188 int CPodeGetRootWindow(void *cpode_mem, realtype *tLo, realtype *tHi)
1189 {
1190   CPodeMem cp_mem;
1191 
1192   if (cpode_mem==NULL) {
1193     cpProcessError(NULL, CP_MEM_NULL, "CPODES", "CPodeGetRootWindow", MSGCP_NO_MEM);
1194     return(CP_MEM_NULL);
1195   }
1196   cp_mem = (CPodeMem) cpode_mem;
1197 
1198   *tLo = cp_mem->cp_tlo;
1199   *tHi = cp_mem->cp_thi;
1200 
1201   return(CP_SUCCESS);
1202 }
1203 
1204 /*
1205  * CPodeGetNumNonlinSolvIters
1206  *
1207  * Returns the current number of iterations in the nonlinear solver
1208  */
1209 
CPodeGetNumNonlinSolvIters(void * cpode_mem,long int * nniters)1210 int CPodeGetNumNonlinSolvIters(void *cpode_mem, long int *nniters)
1211 {
1212   CPodeMem cp_mem;
1213 
1214   if (cpode_mem==NULL) {
1215     cpProcessError(NULL, CP_MEM_NULL, "CPODES", "CPodeGetNumNonlinSolvIters", MSGCP_NO_MEM);
1216     return(CP_MEM_NULL);
1217   }
1218   cp_mem = (CPodeMem) cpode_mem;
1219 
1220   *nniters = cp_mem->cp_nni;
1221 
1222   return(CP_SUCCESS);
1223 }
1224 
1225 /*
1226  * CPodeGetNumNonlinSolvConvFails
1227  *
1228  * Returns the current number of convergence failures in the
1229  * nonlinear solver
1230  */
1231 
CPodeGetNumNonlinSolvConvFails(void * cpode_mem,long int * nncfails)1232 int CPodeGetNumNonlinSolvConvFails(void *cpode_mem, long int *nncfails)
1233 {
1234   CPodeMem cp_mem;
1235 
1236   if (cpode_mem==NULL) {
1237     cpProcessError(NULL, CP_MEM_NULL, "CPODES", "CPodeGetNumNonlinSolvConvFails", MSGCP_NO_MEM);
1238     return(CP_MEM_NULL);
1239   }
1240   cp_mem = (CPodeMem) cpode_mem;
1241 
1242   *nncfails = cp_mem->cp_ncfn;
1243 
1244   return(CP_SUCCESS);
1245 }
1246 
1247 /*
1248  * CPodeGetNonlinSolvStats
1249  *
1250  * Returns nonlinear solver statistics
1251  */
1252 
CPodeGetNonlinSolvStats(void * cpode_mem,long int * nniters,long int * nncfails)1253 int CPodeGetNonlinSolvStats(void *cpode_mem, long int *nniters,
1254                             long int *nncfails)
1255 {
1256   CPodeMem cp_mem;
1257 
1258   if (cpode_mem==NULL) {
1259     cpProcessError(NULL, CP_MEM_NULL, "CPODES", "CPodeGetNonlinSolvStats", MSGCP_NO_MEM);
1260     return(CP_MEM_NULL);
1261   }
1262   cp_mem = (CPodeMem) cpode_mem;
1263 
1264   *nniters = cp_mem->cp_nni;
1265   *nncfails = cp_mem->cp_ncfn;
1266 
1267   return(CP_SUCCESS);
1268 }
1269 
1270 /*
1271  * -----------------------------------------------------------------
1272  * Optional outputs for projection step
1273  * -----------------------------------------------------------------
1274  */
1275 
CPodeGetProjNumProj(void * cpode_mem,long int * nproj)1276 int CPodeGetProjNumProj(void *cpode_mem, long int *nproj)
1277 {
1278   CPodeMem cp_mem;
1279 
1280   if (cpode_mem==NULL) {
1281     cpProcessError(NULL, CP_MEM_NULL, "CPODES", "CPodeGetNonlinSolvStats", MSGCP_NO_MEM);
1282     return(CP_MEM_NULL);
1283   }
1284   cp_mem = (CPodeMem) cpode_mem;
1285 
1286   *nproj    = cp_mem->cp_nproj;
1287 
1288   return(CP_SUCCESS);
1289 }
1290 
CPodeGetProjNumCnstrEvals(void * cpode_mem,long int * nce)1291 int CPodeGetProjNumCnstrEvals(void *cpode_mem, long int *nce)
1292 {
1293   CPodeMem cp_mem;
1294 
1295   if (cpode_mem==NULL) {
1296     cpProcessError(NULL, CP_MEM_NULL, "CPODES", "CPodeGetNonlinSolvStats", MSGCP_NO_MEM);
1297     return(CP_MEM_NULL);
1298   }
1299   cp_mem = (CPodeMem) cpode_mem;
1300 
1301   *nce      = cp_mem->cp_nce;
1302 
1303   return(CP_SUCCESS);
1304 }
1305 
CPodeGetProjNumLinSolvSetups(void * cpode_mem,long int * nsetupsP)1306 int CPodeGetProjNumLinSolvSetups(void *cpode_mem, long int *nsetupsP)
1307 {
1308   CPodeMem cp_mem;
1309 
1310   if (cpode_mem==NULL) {
1311     cpProcessError(NULL, CP_MEM_NULL, "CPODES", "CPodeGetNonlinSolvStats", MSGCP_NO_MEM);
1312     return(CP_MEM_NULL);
1313   }
1314   cp_mem = (CPodeMem) cpode_mem;
1315 
1316   *nsetupsP = cp_mem->cp_nsetupsP;
1317 
1318   return(CP_SUCCESS);
1319 }
1320 
CPodeGetProjNumFailures(void * cpode_mem,long int * nprf)1321 int CPodeGetProjNumFailures(void *cpode_mem, long int *nprf)
1322 {
1323   CPodeMem cp_mem;
1324 
1325   if (cpode_mem==NULL) {
1326     cpProcessError(NULL, CP_MEM_NULL, "CPODES", "CPodeGetNonlinSolvStats", MSGCP_NO_MEM);
1327     return(CP_MEM_NULL);
1328   }
1329   cp_mem = (CPodeMem) cpode_mem;
1330 
1331   *nprf = cp_mem->cp_nprf;
1332 
1333   return(CP_SUCCESS);
1334 }
1335 
CPodeGetProjStats(void * cpode_mem,long int * nproj,long int * nce,long int * nsetupsP,long int * nprf)1336 int CPodeGetProjStats(void *cpode_mem, long int *nproj, long int *nce,
1337                       long int *nsetupsP, long int *nprf)
1338 {
1339   CPodeMem cp_mem;
1340 
1341   if (cpode_mem==NULL) {
1342     cpProcessError(NULL, CP_MEM_NULL, "CPODES", "CPodeGetNonlinSolvStats", MSGCP_NO_MEM);
1343     return(CP_MEM_NULL);
1344   }
1345   cp_mem = (CPodeMem) cpode_mem;
1346 
1347   *nproj    = cp_mem->cp_nproj;
1348   *nce      = cp_mem->cp_nce;
1349   *nsetupsP = cp_mem->cp_nsetupsP;
1350   *nprf     = cp_mem->cp_nprf;
1351 
1352   return(CP_SUCCESS);
1353 }
1354 
1355 /*
1356  * -----------------------------------------------------------------
1357  * Optional outputs for quadrature integration
1358  * -----------------------------------------------------------------
1359  */
1360 
CPodeGetQuadNumFunEvals(void * cpode_mem,long int * nqevals)1361 int CPodeGetQuadNumFunEvals(void *cpode_mem, long int *nqevals)
1362 {
1363   CPodeMem cp_mem;
1364 
1365   if (cpode_mem==NULL) {
1366     cpProcessError(NULL, CP_MEM_NULL, "CPODES", "CPodeGetQuadNumFunEvals", MSGCP_NO_MEM);
1367     return(CP_MEM_NULL);
1368   }
1369 
1370   cp_mem = (CPodeMem) cpode_mem;
1371 
1372   if (cp_mem->cp_quadr==FALSE) {
1373     cpProcessError(cp_mem, CP_NO_QUAD, "CPODES", "CPodeGetQuadNumFunEvals", MSGCP_NO_QUAD);
1374     return(CP_NO_QUAD);
1375   }
1376 
1377   *nqevals = cp_mem->cp_nqe;
1378 
1379   return(CP_SUCCESS);
1380 }
1381 
CPodeGetQuadErrWeights(void * cpode_mem,N_Vector eQweight)1382 int CPodeGetQuadErrWeights(void *cpode_mem, N_Vector eQweight)
1383 {
1384   CPodeMem cp_mem;
1385 
1386   if (cpode_mem==NULL) {
1387     cpProcessError(NULL, CP_MEM_NULL, "CPODES", "CPodeGetQuadErrWeights", MSGCP_NO_MEM);
1388     return(CP_MEM_NULL);
1389   }
1390 
1391   cp_mem = (CPodeMem) cpode_mem;
1392 
1393   if (cp_mem->cp_quadr==FALSE) {
1394     cpProcessError(cp_mem, CP_NO_QUAD, "CPODES", "CPodeGetQuadErrWeights", MSGCP_NO_QUAD);
1395     return(CP_NO_QUAD);
1396   }
1397 
1398   if(cp_mem->cp_errconQ) N_VScale(ONE, cp_mem->cp_ewtQ, eQweight);
1399 
1400   return(CP_SUCCESS);
1401 }
1402 
1403 
1404 /*-----------------------------------------------------------------*/
1405 
CPodeGetReturnFlagName(int flag)1406 char *CPodeGetReturnFlagName(int flag)
1407 {
1408   char *name;
1409 
1410   name = (char *)malloc(24*sizeof(char));
1411 
1412   switch(flag) {
1413   case CP_SUCCESS:
1414     sprintf(name,"CP_SUCCESS");
1415     break;
1416   case CP_TSTOP_RETURN:
1417     sprintf(name,"CP_TSTOP_RETURN");
1418     break;
1419   case CP_ROOT_RETURN:
1420     sprintf(name,"CP_ROOT_RETURN");
1421     break;
1422   case CP_TOO_MUCH_WORK:
1423     sprintf(name,"CP_TOO_MUCH_WORK");
1424     break;
1425   case CP_TOO_MUCH_ACC:
1426     sprintf(name,"CP_TOO_MUCH_ACC");
1427     break;
1428   case CP_ERR_FAILURE:
1429     sprintf(name,"CP_ERR_FAILURE");
1430     break;
1431   case CP_CONV_FAILURE:
1432     sprintf(name,"CP_CONV_FAILURE");
1433     break;
1434   case CP_LINIT_FAIL:
1435     sprintf(name,"CP_LINIT_FAIL");
1436     break;
1437   case CP_LSETUP_FAIL:
1438     sprintf(name,"CP_LSETUP_FAIL");
1439     break;
1440   case CP_LSOLVE_FAIL:
1441     sprintf(name,"CP_LSOLVE_FAIL");
1442     break;
1443   case CP_ODEFUNC_FAIL:
1444     sprintf(name,"CP_ODEFUNC_FAIL");
1445     break;
1446   case CP_FIRST_ODEFUNC_ERR:
1447     sprintf(name,"CP_FIRST_ODEFUNC_ERR");
1448     break;
1449   case CP_REPTD_ODEFUNC_ERR:
1450     sprintf(name,"CP_REPTD_ODEFUNC_ERR");
1451     break;
1452   case CP_UNREC_ODEFUNC_ERR:
1453     sprintf(name,"CP_UNREC_ODEFUNC_ERR");
1454     break;
1455   case CP_RTFUNC_FAIL:
1456     sprintf(name,"CP_RTFUNC_FAIL");
1457     break;
1458   case CP_MEM_FAIL:
1459     sprintf(name,"CP_MEM_FAIL");
1460     break;
1461   case CP_MEM_NULL:
1462     sprintf(name,"CP_MEM_NULL");
1463     break;
1464   case CP_ILL_INPUT:
1465     sprintf(name,"CP_ILL_INPUT");
1466     break;
1467   case CP_NO_MALLOC:
1468     sprintf(name,"CP_NO_MALLOC");
1469     break;
1470   case CP_BAD_K:
1471     sprintf(name,"CP_BAD_K");
1472     break;
1473   case CP_BAD_T:
1474     sprintf(name,"CP_BAD_T");
1475     break;
1476   case CP_BAD_DKY:
1477     sprintf(name,"CP_BAD_DKY");
1478     break;
1479   default:
1480     sprintf(name,"NONE");
1481   }
1482 
1483   return(name);
1484 }
1485 
1486