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