1 /* -----------------------------------------------------------------
2  * Programmer(s): Radu Serban and Cosmin Petra @ LLNL
3  * -----------------------------------------------------------------
4  * SUNDIALS Copyright Start
5  * Copyright (c) 2002-2021, Lawrence Livermore National Security
6  * and Southern Methodist University.
7  * All rights reserved.
8  *
9  * See the top-level LICENSE and NOTICE files for details.
10  *
11  * SPDX-License-Identifier: BSD-3-Clause
12  * SUNDIALS Copyright End
13  * -----------------------------------------------------------------
14  * This is the implementation file for the optional inputs and
15  * outputs for the IDAS solver.
16  * -----------------------------------------------------------------*/
17 
18 #include <stdio.h>
19 #include <stdlib.h>
20 
21 #include "idas_impl.h"
22 #include "sundials/sundials_types.h"
23 #include "sundials/sundials_math.h"
24 
25 #define ZERO    RCONST(0.0)
26 #define HALF    RCONST(0.5)
27 #define ONE     RCONST(1.0)
28 #define TWOPT5  RCONST(2.5)
29 
30 /*
31  * =================================================================
32  * IDA optional input functions
33  * =================================================================
34  */
35 
IDASetErrHandlerFn(void * ida_mem,IDAErrHandlerFn ehfun,void * eh_data)36 int IDASetErrHandlerFn(void *ida_mem, IDAErrHandlerFn ehfun, void *eh_data)
37 {
38   IDAMem IDA_mem;
39 
40   if (ida_mem==NULL) {
41     IDAProcessError(NULL, IDA_MEM_NULL, "IDAS", "IDASetErrHandlerFn", MSG_NO_MEM);
42     return(IDA_MEM_NULL);
43   }
44 
45   IDA_mem = (IDAMem) ida_mem;
46 
47   IDA_mem->ida_ehfun = ehfun;
48   IDA_mem->ida_eh_data = eh_data;
49 
50   return(IDA_SUCCESS);
51 }
52 
53 
IDASetErrFile(void * ida_mem,FILE * errfp)54 int IDASetErrFile(void *ida_mem, FILE *errfp)
55 {
56   IDAMem IDA_mem;
57 
58   if (ida_mem==NULL) {
59     IDAProcessError(NULL, IDA_MEM_NULL, "IDAS", "IDASetErrFile", MSG_NO_MEM);
60     return(IDA_MEM_NULL);
61   }
62 
63   IDA_mem = (IDAMem) ida_mem;
64 
65   IDA_mem->ida_errfp = errfp;
66 
67   return(IDA_SUCCESS);
68 }
69 
70 /*-----------------------------------------------------------------*/
71 
IDASetUserData(void * ida_mem,void * user_data)72 int IDASetUserData(void *ida_mem, void *user_data)
73 {
74   IDAMem IDA_mem;
75 
76   if (ida_mem==NULL) {
77     IDAProcessError(NULL, IDA_MEM_NULL, "IDAS", "IDASetUserData", MSG_NO_MEM);
78     return(IDA_MEM_NULL);
79   }
80 
81   IDA_mem = (IDAMem) ida_mem;
82 
83   IDA_mem->ida_user_data = user_data;
84 
85   return(IDA_SUCCESS);
86 }
87 
88 /*-----------------------------------------------------------------*/
89 
IDASetMaxOrd(void * ida_mem,int maxord)90 int IDASetMaxOrd(void *ida_mem, int maxord)
91 {
92   IDAMem IDA_mem;
93   int maxord_alloc;
94 
95   if (ida_mem==NULL) {
96     IDAProcessError(NULL, IDA_MEM_NULL, "IDAS", "IDASetMaxOrd", MSG_NO_MEM);
97     return(IDA_MEM_NULL);
98   }
99 
100   IDA_mem = (IDAMem) ida_mem;
101 
102   if (maxord <= 0) {
103     IDAProcessError(IDA_mem, IDA_ILL_INPUT, "IDAS", "IDASetMaxOrd", MSG_NEG_MAXORD);
104     return(IDA_ILL_INPUT);
105   }
106 
107   /* Cannot increase maximum order beyond the value that
108      was used when allocating memory */
109   maxord_alloc = IDA_mem->ida_maxord_alloc;
110 
111   if (maxord > maxord_alloc) {
112     IDAProcessError(IDA_mem, IDA_ILL_INPUT, "IDAS", "IDASetMaxOrd", MSG_BAD_MAXORD);
113     return(IDA_ILL_INPUT);
114   }
115 
116   IDA_mem->ida_maxord = SUNMIN(maxord,MAXORD_DEFAULT);
117 
118   return(IDA_SUCCESS);
119 }
120 
121 /*-----------------------------------------------------------------*/
122 
IDASetMaxNumSteps(void * ida_mem,long int mxsteps)123 int IDASetMaxNumSteps(void *ida_mem, long int mxsteps)
124 {
125   IDAMem IDA_mem;
126 
127   if (ida_mem==NULL) {
128     IDAProcessError(NULL, IDA_MEM_NULL, "IDAS", "IDASetMaxNumSteps", MSG_NO_MEM);
129     return(IDA_MEM_NULL);
130   }
131 
132   IDA_mem = (IDAMem) ida_mem;
133 
134   /* Passing mxsteps=0 sets the default. Passing mxsteps<0 disables the test. */
135 
136   if (mxsteps == 0)
137     IDA_mem->ida_mxstep = MXSTEP_DEFAULT;
138   else
139     IDA_mem->ida_mxstep = mxsteps;
140 
141   return(IDA_SUCCESS);
142 }
143 
144 /*-----------------------------------------------------------------*/
145 
IDASetInitStep(void * ida_mem,realtype hin)146 int IDASetInitStep(void *ida_mem, realtype hin)
147 {
148   IDAMem IDA_mem;
149 
150   if (ida_mem==NULL) {
151     IDAProcessError(NULL, IDA_MEM_NULL, "IDAS", "IDASetInitStep", MSG_NO_MEM);
152     return(IDA_MEM_NULL);
153   }
154 
155   IDA_mem = (IDAMem) ida_mem;
156 
157   IDA_mem->ida_hin = hin;
158 
159   return(IDA_SUCCESS);
160 }
161 
162 /*-----------------------------------------------------------------*/
163 
IDASetMaxStep(void * ida_mem,realtype hmax)164 int IDASetMaxStep(void *ida_mem, realtype hmax)
165 {
166   IDAMem IDA_mem;
167 
168   if (ida_mem==NULL) {
169     IDAProcessError(NULL, IDA_MEM_NULL, "IDAS", "IDASetMaxStep", MSG_NO_MEM);
170     return(IDA_MEM_NULL);
171   }
172 
173   IDA_mem = (IDAMem) ida_mem;
174 
175   if (hmax < 0) {
176     IDAProcessError(IDA_mem, IDA_ILL_INPUT, "IDAS", "IDASetMaxStep", MSG_NEG_HMAX);
177     return(IDA_ILL_INPUT);
178   }
179 
180   /* Passing 0 sets hmax = infinity */
181   if (hmax == ZERO) {
182     IDA_mem->ida_hmax_inv = HMAX_INV_DEFAULT;
183     return(IDA_SUCCESS);
184   }
185 
186   IDA_mem->ida_hmax_inv = ONE/hmax;
187 
188   return(IDA_SUCCESS);
189 }
190 
191 /*-----------------------------------------------------------------*/
192 
IDASetStopTime(void * ida_mem,realtype tstop)193 int IDASetStopTime(void *ida_mem, realtype tstop)
194 {
195   IDAMem IDA_mem;
196 
197   if (ida_mem==NULL) {
198     IDAProcessError(NULL, IDA_MEM_NULL, "IDAS", "IDASetStopTime", MSG_NO_MEM);
199     return(IDA_MEM_NULL);
200   }
201 
202   IDA_mem = (IDAMem) ida_mem;
203 
204   /* If IDASolve was called at least once, test if tstop is legal
205    * (i.e. if it was not already passed).
206    * If IDASetStopTime is called before the first call to IDASolve,
207    * tstop will be checked in IDASolve. */
208   if (IDA_mem->ida_nst > 0) {
209 
210     if ( (tstop - IDA_mem->ida_tn) * IDA_mem->ida_hh < ZERO ) {
211       IDAProcessError(IDA_mem, IDA_ILL_INPUT, "IDA", "IDASetStopTime", MSG_BAD_TSTOP, tstop, IDA_mem->ida_tn);
212       return(IDA_ILL_INPUT);
213     }
214 
215   }
216 
217   IDA_mem->ida_tstop = tstop;
218   IDA_mem->ida_tstopset = SUNTRUE;
219 
220   return(IDA_SUCCESS);
221 }
222 
223 /*-----------------------------------------------------------------*/
224 
IDASetNonlinConvCoef(void * ida_mem,realtype epcon)225 int IDASetNonlinConvCoef(void *ida_mem, realtype epcon)
226 {
227   IDAMem IDA_mem;
228 
229   if (ida_mem==NULL) {
230     IDAProcessError(NULL, IDA_MEM_NULL, "IDAS", "IDASetNonlinConvCoef", MSG_NO_MEM);
231     return(IDA_MEM_NULL);
232   }
233 
234   IDA_mem = (IDAMem) ida_mem;
235 
236   if (epcon <= ZERO) {
237     IDAProcessError(IDA_mem, IDA_ILL_INPUT, "IDAS", "IDASetNonlinConvCoef", MSG_NEG_EPCON);
238     return(IDA_ILL_INPUT);
239   }
240 
241   IDA_mem->ida_epcon = epcon;
242 
243   return(IDA_SUCCESS);
244 }
245 
246 /*-----------------------------------------------------------------*/
247 
IDASetMaxErrTestFails(void * ida_mem,int maxnef)248 int IDASetMaxErrTestFails(void *ida_mem, int maxnef)
249 {
250   IDAMem IDA_mem;
251 
252   if (ida_mem==NULL) {
253     IDAProcessError(NULL, IDA_MEM_NULL, "IDAS", "IDASetMaxErrTestFails", MSG_NO_MEM);
254     return (IDA_MEM_NULL);
255   }
256 
257   IDA_mem = (IDAMem) ida_mem;
258 
259   IDA_mem->ida_maxnef = maxnef;
260 
261   return(IDA_SUCCESS);
262 }
263 
264 /*-----------------------------------------------------------------*/
265 
IDASetMaxConvFails(void * ida_mem,int maxncf)266 int IDASetMaxConvFails(void *ida_mem, int maxncf)
267 {
268   IDAMem IDA_mem;
269 
270   if (ida_mem==NULL) {
271     IDAProcessError(NULL, IDA_MEM_NULL, "IDAS", "IDASetMaxConvFails", MSG_NO_MEM);
272     return (IDA_MEM_NULL);
273   }
274 
275   IDA_mem = (IDAMem) ida_mem;
276 
277   IDA_mem->ida_maxncf = maxncf;
278 
279   return(IDA_SUCCESS);
280 }
281 
282 /*-----------------------------------------------------------------*/
283 
IDASetMaxNonlinIters(void * ida_mem,int maxcor)284 int IDASetMaxNonlinIters(void *ida_mem, int maxcor)
285 {
286   IDAMem IDA_mem;
287   booleantype sensi_sim;
288 
289   if (ida_mem==NULL) {
290     IDAProcessError(NULL, IDA_MEM_NULL, "IDAS",
291                     "IDASetMaxNonlinIters", MSG_NO_MEM);
292     return(IDA_MEM_NULL);
293   }
294 
295   IDA_mem = (IDAMem) ida_mem;
296 
297   /* Are we computing sensitivities with the simultaneous approach? */
298   sensi_sim = (IDA_mem->ida_sensi && (IDA_mem->ida_ism==IDA_SIMULTANEOUS));
299 
300   if (sensi_sim) {
301 
302     /* check that the NLS is non-NULL */
303     if (IDA_mem->NLSsim == NULL) {
304       IDAProcessError(NULL, IDA_MEM_FAIL, "IDAS",
305                       "IDASetMaxNonlinIters", MSG_MEM_FAIL);
306       return(IDA_MEM_FAIL);
307     }
308 
309     return(SUNNonlinSolSetMaxIters(IDA_mem->NLSsim, maxcor));
310 
311   } else {
312 
313     /* check that the NLS is non-NULL */
314     if (IDA_mem->NLS == NULL) {
315       IDAProcessError(NULL, IDA_MEM_FAIL, "IDAS",
316                       "IDASetMaxNonlinIters", MSG_MEM_FAIL);
317       return(IDA_MEM_FAIL);
318     }
319 
320     return(SUNNonlinSolSetMaxIters(IDA_mem->NLS, maxcor));
321   }
322 }
323 
324 /*-----------------------------------------------------------------*/
325 
IDASetSuppressAlg(void * ida_mem,booleantype suppressalg)326 int IDASetSuppressAlg(void *ida_mem, booleantype suppressalg)
327 {
328   IDAMem IDA_mem;
329 
330   if (ida_mem==NULL) {
331     IDAProcessError(NULL, IDA_MEM_NULL, "IDAS", "IDASetSuppressAlg", MSG_NO_MEM);
332     return(IDA_MEM_NULL);
333   }
334 
335   IDA_mem = (IDAMem) ida_mem;
336 
337   IDA_mem->ida_suppressalg = suppressalg;
338 
339   return(IDA_SUCCESS);
340 }
341 
342 /*-----------------------------------------------------------------*/
343 
IDASetId(void * ida_mem,N_Vector id)344 int IDASetId(void *ida_mem, N_Vector id)
345 {
346   IDAMem IDA_mem;
347 
348   if (ida_mem==NULL) {
349     IDAProcessError(NULL, IDA_MEM_NULL, "IDAS", "IDASetId", MSG_NO_MEM);
350     return(IDA_MEM_NULL);
351   }
352 
353   IDA_mem = (IDAMem) ida_mem;
354 
355   if (id == NULL) {
356     if (IDA_mem->ida_idMallocDone) {
357       N_VDestroy(IDA_mem->ida_id);
358       IDA_mem->ida_lrw -= IDA_mem->ida_lrw1;
359       IDA_mem->ida_liw -= IDA_mem->ida_liw1;
360     }
361     IDA_mem->ida_idMallocDone = SUNFALSE;
362     return(IDA_SUCCESS);
363   }
364 
365   if ( !(IDA_mem->ida_idMallocDone) ) {
366     IDA_mem->ida_id = N_VClone(id);
367     IDA_mem->ida_lrw += IDA_mem->ida_lrw1;
368     IDA_mem->ida_liw += IDA_mem->ida_liw1;
369     IDA_mem->ida_idMallocDone = SUNTRUE;
370   }
371 
372   /* Load the id vector */
373 
374   N_VScale(ONE, id, IDA_mem->ida_id);
375 
376   return(IDA_SUCCESS);
377 }
378 
379 /*-----------------------------------------------------------------*/
380 
IDASetConstraints(void * ida_mem,N_Vector constraints)381 int IDASetConstraints(void *ida_mem, N_Vector constraints)
382 {
383   IDAMem IDA_mem;
384   realtype temptest;
385 
386   if (ida_mem==NULL) {
387     IDAProcessError(NULL, IDA_MEM_NULL, "IDAS", "IDASetConstraints", MSG_NO_MEM);
388     return(IDA_MEM_NULL);
389   }
390 
391   IDA_mem = (IDAMem) ida_mem;
392 
393   if (constraints == NULL) {
394     if (IDA_mem->ida_constraintsMallocDone) {
395       N_VDestroy(IDA_mem->ida_constraints);
396       IDA_mem->ida_lrw -= IDA_mem->ida_lrw1;
397       IDA_mem->ida_liw -= IDA_mem->ida_liw1;
398     }
399     IDA_mem->ida_constraintsMallocDone = SUNFALSE;
400     IDA_mem->ida_constraintsSet = SUNFALSE;
401     return(IDA_SUCCESS);
402   }
403 
404   /* Test if required vector ops. are defined */
405 
406   if (constraints->ops->nvdiv         == NULL ||
407       constraints->ops->nvmaxnorm     == NULL ||
408       constraints->ops->nvcompare     == NULL ||
409       constraints->ops->nvconstrmask  == NULL ||
410       constraints->ops->nvminquotient == NULL) {
411     IDAProcessError(IDA_mem, IDA_ILL_INPUT, "IDAS", "IDASetConstraints", MSG_BAD_NVECTOR);
412     return(IDA_ILL_INPUT);
413   }
414 
415   /*  Check the constraints vector */
416 
417   temptest = N_VMaxNorm(constraints);
418   if((temptest > TWOPT5) || (temptest < HALF)){
419     IDAProcessError(IDA_mem, IDA_ILL_INPUT, "IDAS", "IDASetConstraints", MSG_BAD_CONSTR);
420     return(IDA_ILL_INPUT);
421   }
422 
423   if ( !(IDA_mem->ida_constraintsMallocDone) ) {
424     IDA_mem->ida_constraints = N_VClone(constraints);
425     IDA_mem->ida_lrw += IDA_mem->ida_lrw1;
426     IDA_mem->ida_liw += IDA_mem->ida_liw1;
427     IDA_mem->ida_constraintsMallocDone = SUNTRUE;
428   }
429 
430   /* Load the constraints vector */
431 
432   N_VScale(ONE, constraints, IDA_mem->ida_constraints);
433 
434   IDA_mem->ida_constraintsSet = SUNTRUE;
435 
436   return(IDA_SUCCESS);
437 }
438 
439 /*
440  * IDASetRootDirection
441  *
442  * Specifies the direction of zero-crossings to be monitored.
443  * The default is to monitor both crossings.
444  */
445 
IDASetRootDirection(void * ida_mem,int * rootdir)446 int IDASetRootDirection(void *ida_mem, int *rootdir)
447 {
448   IDAMem IDA_mem;
449   int i, nrt;
450 
451   if (ida_mem==NULL) {
452     IDAProcessError(NULL, IDA_MEM_NULL, "IDAS", "IDASetRootDirection", MSG_NO_MEM);
453     return(IDA_MEM_NULL);
454   }
455 
456   IDA_mem = (IDAMem) ida_mem;
457 
458   nrt = IDA_mem->ida_nrtfn;
459   if (nrt==0) {
460     IDAProcessError(NULL, IDA_ILL_INPUT, "IDAS", "IDASetRootDirection", MSG_NO_ROOT);
461     return(IDA_ILL_INPUT);
462   }
463 
464   for(i=0; i<nrt; i++) IDA_mem->ida_rootdir[i] = rootdir[i];
465 
466   return(IDA_SUCCESS);
467 }
468 
469 /*
470  * IDASetNoInactiveRootWarn
471  *
472  * Disables issuing a warning if some root function appears
473  * to be identically zero at the beginning of the integration
474  */
475 
IDASetNoInactiveRootWarn(void * ida_mem)476 int IDASetNoInactiveRootWarn(void *ida_mem)
477 {
478   IDAMem IDA_mem;
479 
480   if (ida_mem==NULL) {
481     IDAProcessError(NULL, IDA_MEM_NULL, "IDAS", "IDASetNoInactiveRootWarn", MSG_NO_MEM);
482     return(IDA_MEM_NULL);
483   }
484 
485   IDA_mem = (IDAMem) ida_mem;
486 
487   IDA_mem->ida_mxgnull = 0;
488 
489   return(IDA_SUCCESS);
490 }
491 
492 
493 /*
494  * =================================================================
495  * IDA IC optional input functions
496  * =================================================================
497  */
498 
IDASetNonlinConvCoefIC(void * ida_mem,realtype epiccon)499 int IDASetNonlinConvCoefIC(void *ida_mem, realtype epiccon)
500 {
501   IDAMem IDA_mem;
502 
503   if (ida_mem==NULL) {
504     IDAProcessError(NULL, IDA_MEM_NULL, "IDAS", "IDASetNonlinConvCoefIC", MSG_NO_MEM);
505     return(IDA_MEM_NULL);
506   }
507 
508   IDA_mem = (IDAMem) ida_mem;
509 
510   if (epiccon <= ZERO) {
511     IDAProcessError(IDA_mem, IDA_ILL_INPUT, "IDAS", "IDASetNonlinConvCoefIC", MSG_BAD_EPICCON);
512     return(IDA_ILL_INPUT);
513   }
514 
515   IDA_mem->ida_epiccon = epiccon;
516 
517   return(IDA_SUCCESS);
518 }
519 
520 /*-----------------------------------------------------------------*/
521 
IDASetMaxNumStepsIC(void * ida_mem,int maxnh)522 int IDASetMaxNumStepsIC(void *ida_mem, int maxnh)
523 {
524   IDAMem IDA_mem;
525 
526   if (ida_mem==NULL) {
527     IDAProcessError(NULL, IDA_MEM_NULL, "IDAS", "IDASetMaxNumStepsIC", MSG_NO_MEM);
528     return(IDA_MEM_NULL);
529   }
530 
531   IDA_mem = (IDAMem) ida_mem;
532 
533   if (maxnh <= 0) {
534     IDAProcessError(IDA_mem, IDA_ILL_INPUT, "IDAS", "IDASetMaxNumStepsIC", MSG_BAD_MAXNH);
535     return(IDA_ILL_INPUT);
536   }
537 
538   IDA_mem->ida_maxnh = maxnh;
539 
540   return(IDA_SUCCESS);
541 }
542 
543 /*-----------------------------------------------------------------*/
544 
IDASetMaxNumJacsIC(void * ida_mem,int maxnj)545 int IDASetMaxNumJacsIC(void *ida_mem, int maxnj)
546 {
547   IDAMem IDA_mem;
548 
549   if (ida_mem==NULL) {
550     IDAProcessError(NULL, IDA_MEM_NULL, "IDAS", "IDASetMaxNumJacsIC", MSG_NO_MEM);
551     return(IDA_MEM_NULL);
552   }
553 
554   IDA_mem = (IDAMem) ida_mem;
555 
556    if (maxnj <= 0) {
557     IDAProcessError(IDA_mem, IDA_ILL_INPUT, "IDAS", "IDASetMaxNumJacsIC", MSG_BAD_MAXNJ);
558     return(IDA_ILL_INPUT);
559   }
560 
561   IDA_mem->ida_maxnj = maxnj;
562 
563   return(IDA_SUCCESS);
564 }
565 
566 /*-----------------------------------------------------------------*/
567 
IDASetMaxNumItersIC(void * ida_mem,int maxnit)568 int IDASetMaxNumItersIC(void *ida_mem, int maxnit)
569 {
570   IDAMem IDA_mem;
571 
572   if (ida_mem==NULL) {
573     IDAProcessError(NULL, IDA_MEM_NULL, "IDAS", "IDASetMaxNumItersIC", MSG_NO_MEM);
574     return(IDA_MEM_NULL);
575   }
576 
577   IDA_mem = (IDAMem) ida_mem;
578 
579   if (maxnit <= 0) {
580     IDAProcessError(IDA_mem, IDA_ILL_INPUT, "IDAS", "IDASetMaxNumItersIC", MSG_BAD_MAXNIT);
581     return(IDA_ILL_INPUT);
582   }
583 
584   IDA_mem->ida_maxnit = maxnit;
585 
586   return(IDA_SUCCESS);
587 }
588 
589 /*-----------------------------------------------------------------*/
590 
IDASetMaxBacksIC(void * ida_mem,int maxbacks)591 int IDASetMaxBacksIC(void *ida_mem, int maxbacks)
592 {
593   IDAMem IDA_mem;
594 
595   if (ida_mem==NULL) {
596     IDAProcessError(NULL, IDA_MEM_NULL, "IDA", "IDASetMaxBacksIC", MSG_NO_MEM);
597     return(IDA_MEM_NULL);
598   }
599 
600   IDA_mem = (IDAMem) ida_mem;
601 
602   if (maxbacks <= 0) {
603     IDAProcessError(IDA_mem, IDA_ILL_INPUT, "IDA", "IDASetMaxBacksIC", MSG_IC_BAD_MAXBACKS);
604     return(IDA_ILL_INPUT);
605   }
606 
607   IDA_mem->ida_maxbacks = maxbacks;
608 
609   return(IDA_SUCCESS);
610 }
611 
612 /*-----------------------------------------------------------------*/
613 
IDASetLineSearchOffIC(void * ida_mem,booleantype lsoff)614 int IDASetLineSearchOffIC(void *ida_mem, booleantype lsoff)
615 {
616   IDAMem IDA_mem;
617 
618   if (ida_mem==NULL) {
619     IDAProcessError(NULL, IDA_MEM_NULL, "IDAS", "IDASetLineSearchOffIC", MSG_NO_MEM);
620     return(IDA_MEM_NULL);
621   }
622 
623   IDA_mem = (IDAMem) ida_mem;
624 
625   IDA_mem->ida_lsoff = lsoff;
626 
627   return(IDA_SUCCESS);
628 }
629 
630 /*-----------------------------------------------------------------*/
631 
IDASetStepToleranceIC(void * ida_mem,realtype steptol)632 int IDASetStepToleranceIC(void *ida_mem, realtype steptol)
633 {
634   IDAMem IDA_mem;
635 
636   if (ida_mem==NULL) {
637     IDAProcessError(NULL, IDA_MEM_NULL, "IDAS", "IDASetStepToleranceIC", MSG_NO_MEM);
638     return(IDA_MEM_NULL);
639   }
640 
641   IDA_mem = (IDAMem) ida_mem;
642 
643   if (steptol <= ZERO) {
644     IDAProcessError(IDA_mem, IDA_ILL_INPUT, "IDAS", "IDASetStepToleranceIC", MSG_BAD_STEPTOL);
645     return(IDA_ILL_INPUT);
646   }
647 
648   IDA_mem->ida_steptol = steptol;
649 
650   return(IDA_SUCCESS);
651 }
652 
653 /*
654  * =================================================================
655  * Quadrature optional input functions
656  * =================================================================
657  */
658 
659 /*-----------------------------------------------------------------*/
660 
IDASetQuadErrCon(void * ida_mem,booleantype errconQ)661 int IDASetQuadErrCon(void *ida_mem, booleantype errconQ)
662 {
663   IDAMem IDA_mem;
664 
665   if (ida_mem==NULL) {
666     IDAProcessError(NULL, IDA_MEM_NULL, "IDAS", "IDASetQuadErrCon", MSG_NO_MEM);
667     return(IDA_MEM_NULL);
668   }
669   IDA_mem = (IDAMem) ida_mem;
670 
671   if (IDA_mem->ida_quadMallocDone == SUNFALSE) {
672     IDAProcessError(NULL, IDA_NO_QUAD, "IDAS", "IDASetQuadErrCon", MSG_NO_QUAD);
673     return(IDA_NO_QUAD);
674   }
675 
676   IDA_mem->ida_errconQ = errconQ;
677 
678   return (IDA_SUCCESS);
679 }
680 
681 /*
682  * =================================================================
683  * FSA optional input functions
684  * =================================================================
685  */
686 
IDASetSensDQMethod(void * ida_mem,int DQtype,realtype DQrhomax)687 int IDASetSensDQMethod(void *ida_mem, int DQtype, realtype DQrhomax)
688 {
689   IDAMem IDA_mem;
690 
691   if (ida_mem==NULL) {
692     IDAProcessError(NULL, IDA_MEM_NULL, "IDAS", "IDASetSensDQMethod", MSG_NO_MEM);
693     return(IDA_MEM_NULL);
694   }
695 
696   IDA_mem = (IDAMem) ida_mem;
697 
698   if ( (DQtype != IDA_CENTERED) && (DQtype != IDA_FORWARD) ) {
699     IDAProcessError(IDA_mem, IDA_ILL_INPUT, "IDAS", "IDASetSensDQMethod", MSG_BAD_DQTYPE);
700     return(IDA_ILL_INPUT);
701   }
702 
703   if (DQrhomax < ZERO ) {
704     IDAProcessError(IDA_mem, IDA_ILL_INPUT, "IDAS", "IDASetSensDQMethod", MSG_BAD_DQRHO);
705     return(IDA_ILL_INPUT);
706   }
707 
708   IDA_mem->ida_DQtype = DQtype;
709   IDA_mem->ida_DQrhomax = DQrhomax;
710 
711   return(IDA_SUCCESS);
712 }
713 
714 /*-----------------------------------------------------------------*/
715 
IDASetSensErrCon(void * ida_mem,booleantype errconS)716 int IDASetSensErrCon(void *ida_mem, booleantype errconS)
717 {
718   IDAMem IDA_mem;
719 
720   if (ida_mem==NULL) {
721     IDAProcessError(NULL, IDA_MEM_NULL, "IDAS", "IDASetSensErrCon", MSG_NO_MEM);
722     return(IDA_MEM_NULL);
723   }
724 
725   IDA_mem = (IDAMem) ida_mem;
726   IDA_mem->ida_errconS = errconS;
727 
728   return(IDA_SUCCESS);
729 }
730 
731 /*-----------------------------------------------------------------*/
732 
IDASetSensMaxNonlinIters(void * ida_mem,int maxcorS)733 int IDASetSensMaxNonlinIters(void *ida_mem, int maxcorS)
734 {
735   IDAMem IDA_mem;
736 
737   if (ida_mem==NULL) {
738     IDAProcessError(NULL, IDA_MEM_NULL, "IDAS",
739                     "IDASetSensMaxNonlinIters", MSG_NO_MEM);
740     return(IDA_MEM_NULL);
741   }
742 
743   IDA_mem = (IDAMem) ida_mem;
744 
745   /* check that the NLS is non-NULL */
746   if (IDA_mem->NLSstg == NULL) {
747     IDAProcessError(NULL, IDA_MEM_FAIL, "IDAS",
748                     "IDASetSensMaxNonlinIters", MSG_MEM_FAIL);
749     return(IDA_MEM_FAIL);
750   }
751 
752   return(SUNNonlinSolSetMaxIters(IDA_mem->NLSstg, maxcorS));
753 }
754 
755 /*-----------------------------------------------------------------*/
756 
IDASetSensParams(void * ida_mem,realtype * p,realtype * pbar,int * plist)757 int IDASetSensParams(void *ida_mem, realtype *p, realtype *pbar, int *plist)
758 {
759   IDAMem IDA_mem;
760   int Ns, is;
761 
762   if (ida_mem==NULL) {
763     IDAProcessError(NULL, IDA_MEM_NULL, "IDAS", "IDASetSensParams", MSG_NO_MEM);
764     return(IDA_MEM_NULL);
765   }
766 
767   IDA_mem = (IDAMem) ida_mem;
768 
769   /* Was sensitivity initialized? */
770 
771   if (IDA_mem->ida_sensMallocDone == SUNFALSE) {
772     IDAProcessError(IDA_mem, IDA_NO_SENS, "IDAS", "IDASetSensParams", MSG_NO_SENSI);
773     return(IDA_NO_SENS);
774   }
775 
776   Ns = IDA_mem->ida_Ns;
777 
778   /* Parameters */
779 
780   IDA_mem->ida_p = p;
781 
782   /* pbar */
783 
784   if (pbar != NULL)
785     for (is=0; is<Ns; is++) {
786       if (pbar[is] == ZERO) {
787         IDAProcessError(IDA_mem, IDA_ILL_INPUT, "IDAS", "IDASetSensParams", MSG_BAD_PBAR);
788         return(IDA_ILL_INPUT);
789       }
790       IDA_mem->ida_pbar[is] = SUNRabs(pbar[is]);
791     }
792   else
793     for (is=0; is<Ns; is++)
794       IDA_mem->ida_pbar[is] = ONE;
795 
796   /* plist */
797 
798   if (plist != NULL)
799     for (is=0; is<Ns; is++) {
800       if ( plist[is] < 0 ) {
801         IDAProcessError(IDA_mem, IDA_ILL_INPUT, "IDAS", "IDASetSensParams", MSG_BAD_PLIST);
802         return(IDA_ILL_INPUT);
803       }
804       IDA_mem->ida_plist[is] = plist[is];
805     }
806   else
807     for (is=0; is<Ns; is++)
808       IDA_mem->ida_plist[is] = is;
809 
810   return(IDA_SUCCESS);
811 }
812 
813 /*
814  * -----------------------------------------------------------------
815  * Function: IDASetQuadSensErrCon
816  * -----------------------------------------------------------------
817  * IDASetQuadSensErrCon specifies if quadrature sensitivity variables
818  * are considered or not in the error control.
819  * -----------------------------------------------------------------
820  */
IDASetQuadSensErrCon(void * ida_mem,booleantype errconQS)821 int IDASetQuadSensErrCon(void *ida_mem, booleantype errconQS)
822 {
823   IDAMem IDA_mem;
824 
825   if (ida_mem==NULL) {
826     IDAProcessError(NULL, IDA_MEM_NULL, "IDAS", "IDASetQuadSensErrCon", MSG_NO_MEM);
827     return(IDA_MEM_NULL);
828   }
829   IDA_mem = (IDAMem) ida_mem;
830 
831   /* Was sensitivity initialized? */
832   if (IDA_mem->ida_sensMallocDone == SUNFALSE) {
833     IDAProcessError(IDA_mem, IDA_NO_SENS, "IDAS", "IDASetQuadSensErrCon", MSG_NO_SENSI);
834     return(IDA_NO_SENS);
835   }
836 
837   /* Was quadrature sensitivity initialized? */
838   if (IDA_mem->ida_quadSensMallocDone == SUNFALSE) {
839     IDAProcessError(IDA_mem, IDA_NO_QUADSENS, "IDAS", "IDASetQuadSensErrCon", MSG_NO_SENSI);
840     return(IDA_NO_QUADSENS);
841   }
842 
843   IDA_mem->ida_errconQS = errconQS;
844 
845   return(IDA_SUCCESS);
846 }
847 
848 /*
849  * =================================================================
850  * IDA optional output functions
851  * =================================================================
852  */
853 
IDAGetNumSteps(void * ida_mem,long int * nsteps)854 int IDAGetNumSteps(void *ida_mem, long int *nsteps)
855 {
856   IDAMem IDA_mem;
857 
858   if (ida_mem==NULL) {
859     IDAProcessError(NULL, IDA_MEM_NULL, "IDAS", "IDAGetNumSteps", MSG_NO_MEM);
860     return(IDA_MEM_NULL);
861   }
862 
863   IDA_mem = (IDAMem) ida_mem;
864 
865   *nsteps = IDA_mem->ida_nst;
866 
867   return(IDA_SUCCESS);
868 }
869 
870 /*-----------------------------------------------------------------*/
871 
IDAGetNumResEvals(void * ida_mem,long int * nrevals)872 int IDAGetNumResEvals(void *ida_mem, long int *nrevals)
873 {
874   IDAMem IDA_mem;
875 
876   if (ida_mem==NULL) {
877     IDAProcessError(NULL, IDA_MEM_NULL, "IDAS", "IDAGetNumResEvals", MSG_NO_MEM);
878     return(IDA_MEM_NULL);
879   }
880 
881   IDA_mem = (IDAMem) ida_mem;
882 
883   *nrevals = IDA_mem->ida_nre;
884 
885   return(IDA_SUCCESS);
886 }
887 
888 /*-----------------------------------------------------------------*/
889 
IDAGetNumLinSolvSetups(void * ida_mem,long int * nlinsetups)890 int IDAGetNumLinSolvSetups(void *ida_mem, long int *nlinsetups)
891 {
892   IDAMem IDA_mem;
893 
894   if (ida_mem==NULL) {
895     IDAProcessError(NULL, IDA_MEM_NULL, "IDAS", "IDAGetNumLinSolvSetups", MSG_NO_MEM);
896     return(IDA_MEM_NULL);
897   }
898 
899   IDA_mem = (IDAMem) ida_mem;
900 
901   *nlinsetups = IDA_mem->ida_nsetups;
902 
903   return(IDA_SUCCESS);
904 }
905 
906 /*-----------------------------------------------------------------*/
907 
IDAGetNumErrTestFails(void * ida_mem,long int * netfails)908 int IDAGetNumErrTestFails(void *ida_mem, long int *netfails)
909 {
910   IDAMem IDA_mem;
911 
912   if (ida_mem==NULL) {
913     IDAProcessError(NULL, IDA_MEM_NULL, "IDAS", "IDAGetNumErrTestFails", MSG_NO_MEM);
914     return(IDA_MEM_NULL);
915   }
916 
917   IDA_mem = (IDAMem) ida_mem;
918 
919   *netfails = IDA_mem->ida_netf;
920 
921   return(IDA_SUCCESS);
922 }
923 
924 /*-----------------------------------------------------------------*/
925 
IDAGetNumBacktrackOps(void * ida_mem,long int * nbacktracks)926 int IDAGetNumBacktrackOps(void *ida_mem, long int *nbacktracks)
927 {
928   IDAMem IDA_mem;
929 
930   if (ida_mem==NULL) {
931     IDAProcessError(NULL, IDA_MEM_NULL, "IDAS", "IDAGetNumBacktrackOps", MSG_NO_MEM);
932     return(IDA_MEM_NULL);
933   }
934 
935   IDA_mem = (IDAMem) ida_mem;
936 
937   *nbacktracks = IDA_mem->ida_nbacktr;
938 
939   return(IDA_SUCCESS);
940 }
941 
942 /*-----------------------------------------------------------------*/
943 
IDAGetConsistentIC(void * ida_mem,N_Vector yy0,N_Vector yp0)944 int IDAGetConsistentIC(void *ida_mem, N_Vector yy0, N_Vector yp0)
945 {
946   IDAMem IDA_mem;
947 
948   if (ida_mem == NULL) {
949     IDAProcessError(NULL, IDA_MEM_NULL, "IDAS", "IDAGetConsistentIC", MSG_NO_MEM);
950     return (IDA_MEM_NULL);
951   }
952 
953   IDA_mem = (IDAMem) ida_mem;
954 
955   if (IDA_mem->ida_kused != 0) {
956     IDAProcessError(IDA_mem, IDA_ILL_INPUT, "IDAS", "IDAGetConsistentIC", MSG_TOO_LATE);
957     return(IDA_ILL_INPUT);
958   }
959 
960   if(yy0 != NULL) N_VScale(ONE, IDA_mem->ida_phi[0], yy0);
961   if(yp0 != NULL) N_VScale(ONE, IDA_mem->ida_phi[1], yp0);
962 
963   return(IDA_SUCCESS);
964 }
965 
966 /*-----------------------------------------------------------------*/
967 
IDAGetLastOrder(void * ida_mem,int * klast)968 int IDAGetLastOrder(void *ida_mem, int *klast)
969 {
970   IDAMem IDA_mem;
971 
972   if (ida_mem==NULL) {
973     IDAProcessError(NULL, IDA_MEM_NULL, "IDAS", "IDAGetLastOrder", MSG_NO_MEM);
974     return(IDA_MEM_NULL);
975   }
976 
977   IDA_mem = (IDAMem) ida_mem;
978 
979   *klast = IDA_mem->ida_kused;
980 
981   return(IDA_SUCCESS);
982 }
983 
984 /*-----------------------------------------------------------------*/
985 
IDAGetCurrentOrder(void * ida_mem,int * kcur)986 int IDAGetCurrentOrder(void *ida_mem, int *kcur)
987 {
988   IDAMem IDA_mem;
989 
990   if (ida_mem==NULL) {
991     IDAProcessError(NULL, IDA_MEM_NULL, "IDAS", "IDAGetCurrentOrder", MSG_NO_MEM);
992     return(IDA_MEM_NULL);
993   }
994 
995   IDA_mem = (IDAMem) ida_mem;
996 
997   *kcur = IDA_mem->ida_kk;
998 
999   return(IDA_SUCCESS);
1000 }
1001 
1002 /*-----------------------------------------------------------------*/
1003 
IDAGetCurrentCj(void * ida_mem,realtype * cj)1004 int IDAGetCurrentCj(void *ida_mem, realtype *cj)
1005 {
1006   IDAMem IDA_mem;
1007 
1008   if (ida_mem==NULL) {
1009     IDAProcessError(NULL, IDA_MEM_NULL, "IDAS", "IDAGetCurrentCjRatio", MSG_NO_MEM);
1010     return(IDA_MEM_NULL);
1011   }
1012 
1013   IDA_mem = (IDAMem) ida_mem;
1014 
1015   *cj = IDA_mem->ida_cj;
1016 
1017   return(IDA_SUCCESS);
1018 }
1019 
1020 /*-----------------------------------------------------------------*/
1021 
IDAGetCurrentY(void * ida_mem,N_Vector * ycur)1022 int IDAGetCurrentY(void *ida_mem, N_Vector *ycur)
1023 {
1024   IDAMem IDA_mem;
1025 
1026   if (ida_mem==NULL) {
1027     IDAProcessError(NULL, IDA_MEM_NULL, "IDAS", "IDAGetCurrentY", MSG_NO_MEM);
1028     return(IDA_MEM_NULL);
1029   }
1030 
1031   IDA_mem = (IDAMem) ida_mem;
1032 
1033   *ycur = IDA_mem->ida_yy;
1034 
1035   return(IDA_SUCCESS);
1036 }
1037 
1038 /*-----------------------------------------------------------------*/
1039 
IDAGetCurrentYSens(void * ida_mem,N_Vector ** yS)1040 int IDAGetCurrentYSens(void *ida_mem, N_Vector **yS)
1041 {
1042   IDAMem IDA_mem;
1043 
1044   if (ida_mem==NULL) {
1045     IDAProcessError(NULL, IDA_MEM_NULL, "IDAS", "IDAGetCurrentYSens", MSG_NO_MEM);
1046     return(IDA_MEM_NULL);
1047   }
1048 
1049   IDA_mem = (IDAMem) ida_mem;
1050 
1051   *yS = IDA_mem->ida_yyS;
1052 
1053   return(IDA_SUCCESS);
1054 }
1055 
1056 /*-----------------------------------------------------------------*/
1057 
IDAGetCurrentYp(void * ida_mem,N_Vector * ypcur)1058 int IDAGetCurrentYp(void *ida_mem, N_Vector *ypcur)
1059 {
1060   IDAMem IDA_mem;
1061 
1062   if (ida_mem==NULL) {
1063     IDAProcessError(NULL, IDA_MEM_NULL, "IDAS", "IDAGetCurrentYp", MSG_NO_MEM);
1064     return(IDA_MEM_NULL);
1065   }
1066 
1067   IDA_mem = (IDAMem) ida_mem;
1068 
1069   *ypcur = IDA_mem->ida_yp;
1070 
1071   return(IDA_SUCCESS);
1072 }
1073 
1074 /*-----------------------------------------------------------------*/
1075 
IDAGetCurrentYpSens(void * ida_mem,N_Vector ** ypS)1076 int IDAGetCurrentYpSens(void *ida_mem, N_Vector **ypS)
1077 {
1078   IDAMem IDA_mem;
1079 
1080   if (ida_mem==NULL) {
1081     IDAProcessError(NULL, IDA_MEM_NULL, "IDAS", "IDAGetCurrentYpSens", MSG_NO_MEM);
1082     return(IDA_MEM_NULL);
1083   }
1084 
1085   IDA_mem = (IDAMem) ida_mem;
1086 
1087   *ypS = IDA_mem->ida_ypS;
1088 
1089   return(IDA_SUCCESS);
1090 }
1091 
1092 /*-----------------------------------------------------------------*/
1093 
IDAGetActualInitStep(void * ida_mem,realtype * hinused)1094 int IDAGetActualInitStep(void *ida_mem, realtype *hinused)
1095 {
1096   IDAMem IDA_mem;
1097 
1098   if (ida_mem==NULL) {
1099     IDAProcessError(NULL, IDA_MEM_NULL, "IDAS", "IDAGetActualInitStep", MSG_NO_MEM);
1100     return(IDA_MEM_NULL);
1101   }
1102 
1103   IDA_mem = (IDAMem) ida_mem;
1104 
1105   *hinused = IDA_mem->ida_h0u;
1106 
1107   return(IDA_SUCCESS);
1108 }
1109 
1110 /*-----------------------------------------------------------------*/
1111 
IDAGetLastStep(void * ida_mem,realtype * hlast)1112 int IDAGetLastStep(void *ida_mem, realtype *hlast)
1113 {
1114   IDAMem IDA_mem;
1115 
1116   if (ida_mem==NULL) {
1117     IDAProcessError(NULL, IDA_MEM_NULL, "IDAS", "IDAGetLastStep", MSG_NO_MEM);
1118     return(IDA_MEM_NULL);
1119   }
1120 
1121   IDA_mem = (IDAMem) ida_mem;
1122 
1123   *hlast = IDA_mem->ida_hused;
1124 
1125   return(IDA_SUCCESS);
1126 }
1127 
1128 /*-----------------------------------------------------------------*/
1129 
IDAGetCurrentStep(void * ida_mem,realtype * hcur)1130 int IDAGetCurrentStep(void *ida_mem, realtype *hcur)
1131 {
1132   IDAMem IDA_mem;
1133 
1134   if (ida_mem==NULL) {
1135     IDAProcessError(NULL, IDA_MEM_NULL, "IDAS", "IDAGetCurrentStep", MSG_NO_MEM);
1136     return(IDA_MEM_NULL);
1137   }
1138 
1139   IDA_mem = (IDAMem) ida_mem;
1140 
1141   *hcur = IDA_mem->ida_hh;
1142 
1143   return(IDA_SUCCESS);
1144 }
1145 
1146 /*-----------------------------------------------------------------*/
1147 
IDAGetCurrentTime(void * ida_mem,realtype * tcur)1148 int IDAGetCurrentTime(void *ida_mem, realtype *tcur)
1149 {
1150   IDAMem IDA_mem;
1151 
1152   if (ida_mem==NULL) {
1153     IDAProcessError(NULL, IDA_MEM_NULL, "IDAS", "IDAGetCurrentTime", MSG_NO_MEM);
1154     return(IDA_MEM_NULL);
1155   }
1156 
1157   IDA_mem = (IDAMem) ida_mem;
1158 
1159   *tcur = IDA_mem->ida_tn;
1160 
1161   return(IDA_SUCCESS);
1162 }
1163 
1164 /*-----------------------------------------------------------------*/
1165 
IDAGetTolScaleFactor(void * ida_mem,realtype * tolsfact)1166 int IDAGetTolScaleFactor(void *ida_mem, realtype *tolsfact)
1167 {
1168   IDAMem IDA_mem;
1169 
1170   if (ida_mem==NULL) {
1171     IDAProcessError(NULL, IDA_MEM_NULL, "IDAS", "IDAGetTolScaleFactor", MSG_NO_MEM);
1172     return(IDA_MEM_NULL);
1173   }
1174 
1175   IDA_mem = (IDAMem) ida_mem;
1176 
1177   *tolsfact = IDA_mem->ida_tolsf;
1178 
1179   return(IDA_SUCCESS);
1180 }
1181 
1182 /*-----------------------------------------------------------------*/
1183 
IDAGetErrWeights(void * ida_mem,N_Vector eweight)1184 int IDAGetErrWeights(void *ida_mem, N_Vector eweight)
1185 {
1186   IDAMem IDA_mem;
1187 
1188   if (ida_mem == NULL) {
1189     IDAProcessError(NULL, IDA_MEM_NULL, "IDAS", "IDAGetErrWeights", MSG_NO_MEM);
1190     return (IDA_MEM_NULL);
1191   }
1192 
1193   IDA_mem = (IDAMem) ida_mem;
1194 
1195   N_VScale(ONE, IDA_mem->ida_ewt, eweight);
1196 
1197   return(IDA_SUCCESS);
1198 }
1199 
1200 /*-----------------------------------------------------------------*/
1201 
IDAGetEstLocalErrors(void * ida_mem,N_Vector ele)1202 int IDAGetEstLocalErrors(void *ida_mem, N_Vector ele)
1203 {
1204   IDAMem IDA_mem;
1205 
1206   if (ida_mem == NULL) {
1207     IDAProcessError(NULL, IDA_MEM_NULL, "IDAS", "IDAGetEstLocalErrors", MSG_NO_MEM);
1208     return(IDA_MEM_NULL);
1209   }
1210   IDA_mem = (IDAMem) ida_mem;
1211 
1212   N_VScale(ONE, IDA_mem->ida_ee, ele);
1213 
1214   return(IDA_SUCCESS);
1215 }
1216 
1217 /*-----------------------------------------------------------------*/
1218 
IDAGetWorkSpace(void * ida_mem,long int * lenrw,long int * leniw)1219 int IDAGetWorkSpace(void *ida_mem, long int *lenrw, long int *leniw)
1220 {
1221   IDAMem IDA_mem;
1222 
1223   if (ida_mem==NULL) {
1224     IDAProcessError(NULL, IDA_MEM_NULL, "IDAS", "IDAGetWorkSpace", MSG_NO_MEM);
1225     return(IDA_MEM_NULL);
1226   }
1227 
1228   IDA_mem = (IDAMem) ida_mem;
1229 
1230   *leniw = IDA_mem->ida_liw;
1231   *lenrw = IDA_mem->ida_lrw;
1232 
1233   return(IDA_SUCCESS);
1234 }
1235 
1236 /*-----------------------------------------------------------------*/
1237 
IDAGetIntegratorStats(void * ida_mem,long int * nsteps,long int * nrevals,long int * nlinsetups,long int * netfails,int * klast,int * kcur,realtype * hinused,realtype * hlast,realtype * hcur,realtype * tcur)1238 int IDAGetIntegratorStats(void *ida_mem, long int *nsteps, long int *nrevals,
1239                           long int *nlinsetups, long int *netfails,
1240                           int *klast, int *kcur, realtype *hinused, realtype *hlast,
1241                           realtype *hcur, realtype *tcur)
1242 {
1243   IDAMem IDA_mem;
1244 
1245   if (ida_mem==NULL) {
1246     IDAProcessError(NULL, IDA_MEM_NULL, "IDAS", "IDAGetIntegratorStats", MSG_NO_MEM);
1247     return(IDA_MEM_NULL);
1248   }
1249 
1250   IDA_mem = (IDAMem) ida_mem;
1251 
1252   *nsteps     = IDA_mem->ida_nst;
1253   *nrevals    = IDA_mem->ida_nre;
1254   *nlinsetups = IDA_mem->ida_nsetups;
1255   *netfails   = IDA_mem->ida_netf;
1256   *klast      = IDA_mem->ida_kused;
1257   *kcur       = IDA_mem->ida_kk;
1258   *hinused    = IDA_mem->ida_h0u;
1259   *hlast      = IDA_mem->ida_hused;
1260   *hcur       = IDA_mem->ida_hh;
1261   *tcur       = IDA_mem->ida_tn;
1262 
1263   return(IDA_SUCCESS);
1264 }
1265 
1266 /*-----------------------------------------------------------------*/
1267 
IDAGetNumGEvals(void * ida_mem,long int * ngevals)1268 int IDAGetNumGEvals(void *ida_mem, long int *ngevals)
1269 {
1270   IDAMem IDA_mem;
1271 
1272   if (ida_mem==NULL) {
1273     IDAProcessError(NULL, IDA_MEM_NULL, "IDAS", "IDAGetNumGEvals", MSG_NO_MEM);
1274     return(IDA_MEM_NULL);
1275   }
1276 
1277   IDA_mem = (IDAMem) ida_mem;
1278 
1279   *ngevals = IDA_mem->ida_nge;
1280 
1281   return(IDA_SUCCESS);
1282 }
1283 
1284 /*-----------------------------------------------------------------*/
1285 
IDAGetRootInfo(void * ida_mem,int * rootsfound)1286 int IDAGetRootInfo(void *ida_mem, int *rootsfound)
1287 {
1288   IDAMem IDA_mem;
1289   int i, nrt;
1290 
1291   if (ida_mem==NULL) {
1292     IDAProcessError(NULL, IDA_MEM_NULL, "IDAS", "IDAGetRootInfo", MSG_NO_MEM);
1293     return(IDA_MEM_NULL);
1294   }
1295 
1296   IDA_mem = (IDAMem) ida_mem;
1297 
1298   nrt = IDA_mem->ida_nrtfn;
1299 
1300   for (i=0; i<nrt; i++)
1301     rootsfound[i] = IDA_mem->ida_iroots[i];
1302 
1303   return(IDA_SUCCESS);
1304 }
1305 
1306 /*-----------------------------------------------------------------*/
1307 
IDAGetNumNonlinSolvIters(void * ida_mem,long int * nniters)1308 int IDAGetNumNonlinSolvIters(void *ida_mem, long int *nniters)
1309 {
1310   IDAMem IDA_mem;
1311 
1312   if (ida_mem==NULL) {
1313     IDAProcessError(NULL, IDA_MEM_NULL, "IDAS",
1314                     "IDAGetNumNonlinSolvIters", MSG_NO_MEM);
1315     return(IDA_MEM_NULL);
1316   }
1317 
1318   IDA_mem = (IDAMem) ida_mem;
1319 
1320   *nniters = IDA_mem->ida_nni;
1321 
1322   return(IDA_SUCCESS);
1323 }
1324 
1325 /*-----------------------------------------------------------------*/
1326 
IDAGetNumNonlinSolvConvFails(void * ida_mem,long int * nncfails)1327 int IDAGetNumNonlinSolvConvFails(void *ida_mem, long int *nncfails)
1328 {
1329   IDAMem IDA_mem;
1330 
1331   if (ida_mem==NULL) {
1332     IDAProcessError(NULL, IDA_MEM_NULL, "IDAS", "IDAGetNumNonlinSolvConvFails", MSG_NO_MEM);
1333     return(IDA_MEM_NULL);
1334   }
1335 
1336   IDA_mem = (IDAMem) ida_mem;
1337 
1338   *nncfails = IDA_mem->ida_ncfn;
1339 
1340   return(IDA_SUCCESS);
1341 }
1342 
1343 /*-----------------------------------------------------------------*/
1344 
IDAGetNonlinSolvStats(void * ida_mem,long int * nniters,long int * nncfails)1345 int IDAGetNonlinSolvStats(void *ida_mem, long int *nniters, long int *nncfails)
1346 {
1347   IDAMem IDA_mem;
1348 
1349   if (ida_mem==NULL) {
1350     IDAProcessError(NULL, IDA_MEM_NULL, "IDAS",
1351                     "IDAGetNonlinSolvStats", MSG_NO_MEM);
1352     return(IDA_MEM_NULL);
1353   }
1354 
1355   IDA_mem = (IDAMem) ida_mem;
1356 
1357   *nniters  = IDA_mem->ida_nni;
1358   *nncfails = IDA_mem->ida_ncfn;
1359 
1360   return(IDA_SUCCESS);
1361 }
1362 
1363 /*
1364  * =================================================================
1365  * Quadrature optional output functions
1366  * =================================================================
1367  */
1368 
1369 /*-----------------------------------------------------------------*/
1370 
IDAGetQuadNumRhsEvals(void * ida_mem,long int * nrQevals)1371 int IDAGetQuadNumRhsEvals(void *ida_mem, long int *nrQevals)
1372 {
1373   IDAMem IDA_mem;
1374 
1375   if (ida_mem==NULL) {
1376     IDAProcessError(NULL, IDA_MEM_NULL, "IDAS", "IDAGetQuadNumRhsEvals", MSG_NO_MEM);
1377     return(IDA_MEM_NULL);
1378   }
1379 
1380   IDA_mem = (IDAMem) ida_mem;
1381 
1382   if (IDA_mem->ida_quadr==SUNFALSE) {
1383     IDAProcessError(IDA_mem, IDA_NO_QUAD, "IDAS", "IDAGetQuadNumRhsEvals", MSG_NO_QUAD);
1384     return(IDA_NO_QUAD);
1385   }
1386 
1387   *nrQevals = IDA_mem->ida_nrQe;
1388 
1389   return(IDA_SUCCESS);
1390 }
1391 
1392 /*-----------------------------------------------------------------*/
1393 
IDAGetQuadNumErrTestFails(void * ida_mem,long int * nQetfails)1394 int IDAGetQuadNumErrTestFails(void *ida_mem, long int *nQetfails)
1395 {
1396   IDAMem IDA_mem;
1397 
1398   if (ida_mem==NULL) {
1399     IDAProcessError(NULL, IDA_MEM_NULL, "IDAS", "IDAGetQuadNumErrTestFails", MSG_NO_MEM);
1400     return(IDA_MEM_NULL);
1401   }
1402 
1403   IDA_mem = (IDAMem) ida_mem;
1404 
1405   if (IDA_mem->ida_quadr==SUNFALSE) {
1406     IDAProcessError(IDA_mem, IDA_NO_QUAD, "IDAS", "IDAGetQuadNumErrTestFails", MSG_NO_QUAD);
1407     return(IDA_NO_QUAD);
1408   }
1409 
1410   *nQetfails = IDA_mem->ida_netfQ;
1411 
1412   return(IDA_SUCCESS);
1413 }
1414 
1415 /*-----------------------------------------------------------------*/
1416 
IDAGetQuadErrWeights(void * ida_mem,N_Vector eQweight)1417 int IDAGetQuadErrWeights(void *ida_mem, N_Vector eQweight)
1418 {
1419   IDAMem IDA_mem;
1420 
1421   if (ida_mem==NULL) {
1422     IDAProcessError(NULL, IDA_MEM_NULL, "IDAS", "IDAGetQuadErrWeights", MSG_NO_MEM);
1423     return(IDA_MEM_NULL);
1424   }
1425 
1426   IDA_mem = (IDAMem) ida_mem;
1427 
1428   if (IDA_mem->ida_quadr==SUNFALSE) {
1429     IDAProcessError(IDA_mem, IDA_NO_QUAD, "IDAS", "IDAGetQuadErrWeights", MSG_NO_QUAD);
1430     return(IDA_NO_QUAD);
1431   }
1432 
1433   if(IDA_mem->ida_errconQ)
1434     N_VScale(ONE, IDA_mem->ida_ewtQ, eQweight);
1435 
1436   return(IDA_SUCCESS);
1437 }
1438 
1439 /*-----------------------------------------------------------------*/
1440 
IDAGetQuadStats(void * ida_mem,long int * nrQevals,long int * nQetfails)1441 int IDAGetQuadStats(void *ida_mem, long int *nrQevals, long int *nQetfails)
1442 {
1443   IDAMem IDA_mem;
1444 
1445   if (ida_mem==NULL) {
1446     IDAProcessError(NULL, IDA_MEM_NULL, "IDAS", "IDAGetQuadStats", MSG_NO_MEM);
1447     return(IDA_MEM_NULL);
1448   }
1449 
1450   IDA_mem = (IDAMem) ida_mem;
1451 
1452   if (IDA_mem->ida_quadr==SUNFALSE) {
1453     IDAProcessError(IDA_mem, IDA_NO_QUAD, "IDAS", "IDAGetQuadStats", MSG_NO_QUAD);
1454     return(IDA_NO_QUAD);
1455   }
1456 
1457   *nrQevals = IDA_mem->ida_nrQe;
1458   *nQetfails = IDA_mem->ida_netfQ;
1459 
1460   return(IDA_SUCCESS);
1461 }
1462 
1463 
1464 /*
1465  * =================================================================
1466  * Quadrature FSA optional output functions
1467  * =================================================================
1468  */
1469 
1470 /*-----------------------------------------------------------------*/
1471 
IDAGetQuadSensNumRhsEvals(void * ida_mem,long int * nrhsQSevals)1472 int IDAGetQuadSensNumRhsEvals(void *ida_mem, long int *nrhsQSevals)
1473 {
1474   IDAMem IDA_mem;
1475 
1476   if (ida_mem==NULL) {
1477     IDAProcessError(NULL, IDA_MEM_NULL, "IDAS", "IDAGetQuadSensNumRhsEvals", MSG_NO_MEM);
1478     return(IDA_MEM_NULL);
1479   }
1480 
1481   IDA_mem = (IDAMem) ida_mem;
1482 
1483   if (IDA_mem->ida_quadr_sensi == SUNFALSE) {
1484     IDAProcessError(IDA_mem, IDA_NO_QUADSENS, "IDAS", "IDAGetQuadSensNumRhsEvals", MSG_NO_QUADSENSI);
1485     return(IDA_NO_QUADSENS);
1486   }
1487 
1488   *nrhsQSevals = IDA_mem->ida_nrQSe;
1489 
1490   return(IDA_SUCCESS);
1491 }
1492 
1493 /*-----------------------------------------------------------------*/
1494 
IDAGetQuadSensNumErrTestFails(void * ida_mem,long int * nQSetfails)1495 int IDAGetQuadSensNumErrTestFails(void *ida_mem, long int *nQSetfails)
1496 {
1497   IDAMem IDA_mem;
1498 
1499   if (ida_mem==NULL) {
1500     IDAProcessError(NULL, IDA_MEM_NULL, "IDAS", "IDAGetQuadSensNumErrTestFails", MSG_NO_MEM);
1501     return(IDA_MEM_NULL);
1502   }
1503 
1504   IDA_mem = (IDAMem) ida_mem;
1505 
1506   if (IDA_mem->ida_quadr_sensi == SUNFALSE) {
1507     IDAProcessError(IDA_mem, IDA_NO_QUADSENS, "IDAS", "IDAGetQuadSensNumErrTestFails", MSG_NO_QUADSENSI);
1508     return(IDA_NO_QUADSENS);
1509   }
1510 
1511   *nQSetfails = IDA_mem->ida_netfQS;
1512 
1513   return(IDA_SUCCESS);
1514 }
1515 
1516 /*-----------------------------------------------------------------*/
1517 
IDAGetQuadSensErrWeights(void * ida_mem,N_Vector * eQSweight)1518 int IDAGetQuadSensErrWeights(void *ida_mem, N_Vector *eQSweight)
1519 {
1520   IDAMem IDA_mem;
1521   int is, Ns;
1522 
1523   if (ida_mem==NULL) {
1524     IDAProcessError(NULL, IDA_MEM_NULL, "IDAS", "IDAGetQuadSensErrWeights", MSG_NO_MEM);
1525     return(IDA_MEM_NULL);
1526   }
1527 
1528   IDA_mem = (IDAMem) ida_mem;
1529 
1530   if (IDA_mem->ida_quadr_sensi == SUNFALSE) {
1531     IDAProcessError(IDA_mem, IDA_NO_QUADSENS, "IDAS", "IDAGetQuadSensErrWeights", MSG_NO_QUADSENSI);
1532     return(IDA_NO_QUADSENS);
1533   }
1534   Ns = IDA_mem->ida_Ns;
1535 
1536   if (IDA_mem->ida_errconQS)
1537     for (is=0; is<Ns; is++)
1538       N_VScale(ONE, IDA_mem->ida_ewtQS[is], eQSweight[is]);
1539 
1540   return(IDA_SUCCESS);
1541 }
1542 
1543 /*-----------------------------------------------------------------*/
1544 
IDAGetQuadSensStats(void * ida_mem,long int * nrhsQSevals,long int * nQSetfails)1545 int IDAGetQuadSensStats(void *ida_mem, long int *nrhsQSevals, long int *nQSetfails)
1546 {
1547   IDAMem IDA_mem;
1548 
1549   if (ida_mem==NULL) {
1550     IDAProcessError(NULL, IDA_MEM_NULL, "IDAS", "IDAGetQuadSensStats", MSG_NO_MEM);
1551     return(IDA_MEM_NULL);
1552   }
1553 
1554   IDA_mem = (IDAMem) ida_mem;
1555 
1556   if (IDA_mem->ida_quadr_sensi == SUNFALSE) {
1557     IDAProcessError(IDA_mem, IDA_NO_QUADSENS, "IDAS", "IDAGetQuadSensStats", MSG_NO_QUADSENSI);
1558     return(IDA_NO_QUADSENS);
1559   }
1560 
1561   *nrhsQSevals = IDA_mem->ida_nrQSe;
1562   *nQSetfails = IDA_mem->ida_netfQS;
1563 
1564   return(IDA_SUCCESS);
1565 }
1566 
1567 
1568 
1569 /*
1570  * =================================================================
1571  * FSA optional output functions
1572  * =================================================================
1573  */
1574 
1575 /*-----------------------------------------------------------------*/
1576 
IDAGetSensConsistentIC(void * ida_mem,N_Vector * yyS0,N_Vector * ypS0)1577 int IDAGetSensConsistentIC(void *ida_mem, N_Vector *yyS0, N_Vector *ypS0)
1578 {
1579   IDAMem IDA_mem;
1580   int is;
1581 
1582   if (ida_mem == NULL) {
1583     IDAProcessError(NULL, IDA_MEM_NULL, "IDAS", "IDAGetSensConsistentIC", MSG_NO_MEM);
1584     return (IDA_MEM_NULL);
1585   }
1586 
1587   IDA_mem = (IDAMem) ida_mem;
1588 
1589   if (IDA_mem->ida_sensi==SUNFALSE) {
1590     IDAProcessError(IDA_mem, IDA_NO_SENS, "IDAS", "IDAGetSensConsistentIC", MSG_NO_SENSI);
1591     return(IDA_NO_SENS);
1592   }
1593 
1594   if (IDA_mem->ida_kused != 0) {
1595     IDAProcessError(IDA_mem, IDA_ILL_INPUT, "IDAS", "IDAGetSensConsistentIC", MSG_TOO_LATE);
1596     return(IDA_ILL_INPUT);
1597   }
1598 
1599   if(yyS0 != NULL) {
1600     for (is=0; is<IDA_mem->ida_Ns; is++)
1601       N_VScale(ONE, IDA_mem->ida_phiS[0][is], yyS0[is]);
1602   }
1603 
1604   if(ypS0 != NULL) {
1605     for (is=0; is<IDA_mem->ida_Ns; is++)
1606       N_VScale(ONE, IDA_mem->ida_phiS[1][is], ypS0[is]);
1607   }
1608 
1609   return(IDA_SUCCESS);
1610 }
1611 
1612 /*-----------------------------------------------------------------*/
1613 
IDAGetSensNumResEvals(void * ida_mem,long int * nrSevals)1614 int IDAGetSensNumResEvals(void *ida_mem, long int *nrSevals)
1615 {
1616   IDAMem IDA_mem;
1617 
1618   if (ida_mem==NULL) {
1619     IDAProcessError(NULL, IDA_MEM_NULL, "IDAS", "IDAGeSensNumResEvals", MSG_NO_MEM);
1620     return(IDA_MEM_NULL);
1621   }
1622 
1623   IDA_mem = (IDAMem) ida_mem;
1624 
1625   if (IDA_mem->ida_sensi==SUNFALSE) {
1626     IDAProcessError(IDA_mem, IDA_NO_SENS, "IDAS", "IDAGetSensNumResEvals", MSG_NO_SENSI);
1627     return(IDA_NO_SENS);
1628   }
1629 
1630   *nrSevals = IDA_mem->ida_nrSe;
1631 
1632   return(IDA_SUCCESS);
1633 }
1634 
1635 /*-----------------------------------------------------------------*/
1636 
IDAGetNumResEvalsSens(void * ida_mem,long int * nrevalsS)1637 int IDAGetNumResEvalsSens(void *ida_mem, long int *nrevalsS)
1638 {
1639   IDAMem IDA_mem;
1640 
1641   if (ida_mem==NULL) {
1642     IDAProcessError(NULL, IDA_MEM_NULL, "IDAS", "IDAGetNumResEvalsSens", MSG_NO_MEM);
1643     return(IDA_MEM_NULL);
1644   }
1645 
1646   IDA_mem = (IDAMem) ida_mem;
1647 
1648   if (IDA_mem->ida_sensi==SUNFALSE) {
1649     IDAProcessError(IDA_mem, IDA_NO_SENS, "IDAS", "IDAGetNumResEvalsSens", MSG_NO_SENSI);
1650     return(IDA_NO_SENS);
1651   }
1652 
1653   *nrevalsS = IDA_mem->ida_nreS;
1654 
1655   return(IDA_SUCCESS);
1656 }
1657 
1658 /*-----------------------------------------------------------------*/
1659 
IDAGetSensNumErrTestFails(void * ida_mem,long int * nSetfails)1660 int IDAGetSensNumErrTestFails(void *ida_mem, long int *nSetfails)
1661 {
1662   IDAMem IDA_mem;
1663 
1664   if (ida_mem==NULL) {
1665     IDAProcessError(NULL, IDA_MEM_NULL, "IDAS", "IDAGetSensNumErrTestFails", MSG_NO_MEM);
1666     return(IDA_MEM_NULL);
1667   }
1668 
1669   IDA_mem = (IDAMem) ida_mem;
1670 
1671   if (IDA_mem->ida_sensi==SUNFALSE) {
1672     IDAProcessError(IDA_mem, IDA_NO_SENS, "IDAS", "IDAGetSensNumErrTestFails", MSG_NO_SENSI);
1673     return(IDA_NO_SENS);
1674   }
1675 
1676   *nSetfails = IDA_mem->ida_netfS;
1677 
1678   return(IDA_SUCCESS);
1679 }
1680 
1681 /*-----------------------------------------------------------------*/
1682 
IDAGetSensNumLinSolvSetups(void * ida_mem,long int * nlinsetupsS)1683 int IDAGetSensNumLinSolvSetups(void *ida_mem, long int *nlinsetupsS)
1684 {
1685   IDAMem IDA_mem;
1686 
1687   if (ida_mem==NULL) {
1688     IDAProcessError(NULL, IDA_MEM_NULL, "IDAS", "IDAGetSensNumLinSolvSetups", MSG_NO_MEM);
1689     return(IDA_MEM_NULL);
1690   }
1691 
1692   IDA_mem = (IDAMem) ida_mem;
1693 
1694   if (IDA_mem->ida_sensi==SUNFALSE) {
1695     IDAProcessError(IDA_mem, IDA_NO_SENS, "IDAS", "IDAGetSensNumLinSolvSetups", MSG_NO_SENSI);
1696     return(IDA_NO_SENS);
1697   }
1698 
1699   *nlinsetupsS = IDA_mem->ida_nsetupsS;
1700 
1701   return(IDA_SUCCESS);
1702 }
1703 
1704 /*-----------------------------------------------------------------*/
1705 
IDAGetSensErrWeights(void * ida_mem,N_Vector_S eSweight)1706 int IDAGetSensErrWeights(void *ida_mem, N_Vector_S eSweight)
1707 {
1708   IDAMem IDA_mem;
1709   int is;
1710 
1711   if (ida_mem==NULL) {
1712     IDAProcessError(NULL, IDA_MEM_NULL, "IDAS", "IDAGetSensErrWeights", MSG_NO_MEM);
1713     return(IDA_MEM_NULL);
1714   }
1715 
1716   IDA_mem = (IDAMem) ida_mem;
1717 
1718   if (IDA_mem->ida_sensi==SUNFALSE) {
1719     IDAProcessError(IDA_mem, IDA_NO_SENS, "IDAS", "IDAGetSensErrWeights", MSG_NO_SENSI);
1720     return(IDA_NO_SENS);
1721   }
1722 
1723   for (is=0; is<IDA_mem->ida_Ns; is++)
1724     N_VScale(ONE, IDA_mem->ida_ewtS[is], eSweight[is]);
1725 
1726   return(IDA_SUCCESS);
1727 }
1728 
1729 /*-----------------------------------------------------------------*/
1730 
IDAGetSensStats(void * ida_mem,long int * nrSevals,long int * nrevalsS,long int * nSetfails,long int * nlinsetupsS)1731 int IDAGetSensStats(void *ida_mem, long int *nrSevals, long int *nrevalsS,
1732                       long int *nSetfails, long int *nlinsetupsS)
1733 {
1734   IDAMem IDA_mem;
1735 
1736   if (ida_mem==NULL) {
1737     IDAProcessError(NULL, IDA_MEM_NULL, "IDAS", "IDAGetSensStats", MSG_NO_MEM);
1738     return(IDA_MEM_NULL);
1739   }
1740 
1741   IDA_mem = (IDAMem) ida_mem;
1742 
1743   if (IDA_mem->ida_sensi==SUNFALSE) {
1744     IDAProcessError(IDA_mem, IDA_NO_SENS, "IDAS", "IDAGetSensStats", MSG_NO_SENSI);
1745     return(IDA_NO_SENS);
1746   }
1747 
1748   *nrSevals = IDA_mem->ida_nrSe;
1749   *nrevalsS = IDA_mem->ida_nreS;
1750   *nSetfails = IDA_mem->ida_netfS;
1751   *nlinsetupsS = IDA_mem->ida_nsetupsS;
1752 
1753   return(IDA_SUCCESS);
1754 }
1755 
1756 /*-----------------------------------------------------------------*/
1757 
IDAGetSensNumNonlinSolvIters(void * ida_mem,long int * nSniters)1758 int IDAGetSensNumNonlinSolvIters(void *ida_mem, long int *nSniters)
1759 {
1760   IDAMem IDA_mem;
1761 
1762   if (ida_mem==NULL) {
1763     IDAProcessError(NULL, IDA_MEM_NULL, "IDAS",
1764                     "IDAGetSensNumNonlinSolvIters", MSG_NO_MEM);
1765     return(IDA_MEM_NULL);
1766   }
1767 
1768   IDA_mem = (IDAMem) ida_mem;
1769 
1770   if (IDA_mem->ida_sensi==SUNFALSE) {
1771     IDAProcessError(IDA_mem, IDA_NO_SENS, "IDAS",
1772                     "IDAGetSensNumNonlinSolvIters", MSG_NO_SENSI);
1773     return(IDA_NO_SENS);
1774   }
1775 
1776   *nSniters = IDA_mem->ida_nniS;
1777 
1778   return(IDA_SUCCESS);
1779 }
1780 
1781 /*-----------------------------------------------------------------*/
1782 
IDAGetSensNumNonlinSolvConvFails(void * ida_mem,long int * nSncfails)1783 int IDAGetSensNumNonlinSolvConvFails(void *ida_mem, long int *nSncfails)
1784 {
1785   IDAMem IDA_mem;
1786 
1787   if (ida_mem==NULL) {
1788     IDAProcessError(NULL, IDA_MEM_NULL, "IDAS", "IDAGetSensNumNonlinSolvConvFails", MSG_NO_MEM);
1789     return(IDA_MEM_NULL);
1790   }
1791 
1792   IDA_mem = (IDAMem) ida_mem;
1793 
1794   if (IDA_mem->ida_sensi==SUNFALSE) {
1795     IDAProcessError(IDA_mem, IDA_NO_SENS, "IDAS", "IDAGetSensNumNonlinSolvConvFails", MSG_NO_SENSI);
1796     return(IDA_NO_SENS);
1797   }
1798 
1799   *nSncfails = IDA_mem->ida_ncfnS;
1800 
1801   return(IDA_SUCCESS);
1802 }
1803 
1804 /*-----------------------------------------------------------------*/
1805 
IDAGetSensNonlinSolvStats(void * ida_mem,long int * nSniters,long int * nSncfails)1806 int IDAGetSensNonlinSolvStats(void *ida_mem, long int *nSniters, long int *nSncfails)
1807 {
1808   IDAMem IDA_mem;
1809 
1810   if (ida_mem==NULL) {
1811     IDAProcessError(NULL, IDA_MEM_NULL, "IDAS",
1812                     "IDAGetSensNonlinSolvstats", MSG_NO_MEM);
1813     return(IDA_MEM_NULL);
1814   }
1815 
1816   IDA_mem = (IDAMem) ida_mem;
1817 
1818   if (IDA_mem->ida_sensi==SUNFALSE) {
1819     IDAProcessError(IDA_mem, IDA_NO_SENS, "IDAS",
1820                     "IDAGetSensNonlinSolvStats", MSG_NO_SENSI);
1821     return(IDA_NO_SENS);
1822   }
1823 
1824   *nSniters  = IDA_mem->ida_nniS;
1825   *nSncfails = IDA_mem->ida_ncfnS;
1826 
1827   return(IDA_SUCCESS);
1828 }
1829 
1830 /*
1831  * =================================================================
1832  * IDAGetReturnFlagName
1833  * =================================================================
1834  */
1835 
1836 
IDAGetReturnFlagName(long int flag)1837 char *IDAGetReturnFlagName(long int flag)
1838 {
1839   char *name;
1840 
1841   name = (char *)malloc(24*sizeof(char));
1842 
1843   switch(flag) {
1844   case IDA_SUCCESS:
1845     sprintf(name,"IDA_SUCCESS");
1846     break;
1847   case IDA_TSTOP_RETURN:
1848     sprintf(name,"IDA_TSTOP_RETURN");
1849     break;
1850   case IDA_ROOT_RETURN:
1851     sprintf(name,"IDA_ROOT_RETURN");
1852     break;
1853   case IDA_TOO_MUCH_WORK:
1854     sprintf(name,"IDA_TOO_MUCH_WORK");
1855     break;
1856   case IDA_TOO_MUCH_ACC:
1857     sprintf(name,"IDA_TOO_MUCH_ACC");
1858     break;
1859   case IDA_ERR_FAIL:
1860     sprintf(name,"IDA_ERR_FAIL");
1861     break;
1862   case IDA_CONV_FAIL:
1863     sprintf(name,"IDA_CONV_FAIL");
1864     break;
1865   case IDA_LINIT_FAIL:
1866     sprintf(name,"IDA_LINIT_FAIL");
1867     break;
1868   case IDA_LSETUP_FAIL:
1869     sprintf(name,"IDA_LSETUP_FAIL");
1870     break;
1871   case IDA_LSOLVE_FAIL:
1872     sprintf(name,"IDA_LSOLVE_FAIL");
1873     break;
1874   case IDA_CONSTR_FAIL:
1875     sprintf(name,"IDA_CONSTR_FAIL");
1876     break;
1877   case IDA_RES_FAIL:
1878     sprintf(name,"IDA_RES_FAIL");
1879     break;
1880   case IDA_FIRST_RES_FAIL:
1881     sprintf(name,"IDA_FIRST_RES_FAIL");
1882     break;
1883   case IDA_REP_RES_ERR:
1884     sprintf(name,"IDA_REP_RES_ERR");
1885     break;
1886   case IDA_RTFUNC_FAIL:
1887     sprintf(name,"IDA_RTFUNC_FAIL");
1888     break;
1889   case IDA_MEM_FAIL:
1890     sprintf(name,"IDA_MEM_FAIL");
1891     break;
1892   case IDA_MEM_NULL:
1893     sprintf(name,"IDA_MEM_NULL");
1894     break;
1895   case IDA_ILL_INPUT:
1896     sprintf(name,"IDA_ILL_INPUT");
1897     break;
1898   case IDA_NO_MALLOC:
1899     sprintf(name,"IDA_NO_MALLOC");
1900     break;
1901   case IDA_BAD_T:
1902     sprintf(name,"IDA_BAD_T");
1903     break;
1904   case IDA_BAD_K:
1905     sprintf(name,"IDA_BAD_K");
1906     break;
1907   case IDA_BAD_DKY:
1908     sprintf(name,"IDA_BAD_DKY");
1909     break;
1910   case IDA_BAD_EWT:
1911     sprintf(name,"IDA_BAD_EWT");
1912     break;
1913   case IDA_NO_RECOVERY:
1914     sprintf(name,"IDA_NO_RECOVERY");
1915     break;
1916   case IDA_LINESEARCH_FAIL:
1917     sprintf(name,"IDA_LINESEARCH_FAIL");
1918     break;
1919   case IDA_NO_SENS:
1920     sprintf(name,"IDA_NO_SENS");
1921     break;
1922   case IDA_SRES_FAIL:
1923     sprintf(name, "IDA_SRES_FAIL");
1924     break;
1925   case IDA_REP_SRES_ERR:
1926     sprintf(name, "IDA_REP_SRES_ERR");
1927     break;
1928   case IDA_BAD_IS:
1929     sprintf(name,"IDA_BAD_IS");
1930     break;
1931   case IDA_NO_QUAD:
1932     sprintf(name,"IDA_NO_QUAD");
1933     break;
1934   case IDA_NO_QUADSENS:
1935     sprintf(name, "IDA_NO_QUADSENS");
1936     break;
1937   case IDA_QRHS_FAIL:
1938     sprintf(name,"IDA_QRHS_FAIL");
1939     break;
1940   case IDA_REP_QRHS_ERR:
1941     sprintf(name,"IDA_REP_QRHS_ERR");
1942     break;
1943   case IDA_QSRHS_FAIL:
1944     sprintf(name, "IDA_QSRHS_FAIL");
1945     break;
1946   case IDA_REP_QSRHS_ERR:
1947     sprintf(name,"IDA_REP_QSRHS_ERR");
1948     break;
1949 
1950     /* IDAA flags follow below. */
1951   case IDA_NO_ADJ:
1952     sprintf(name, "IDA_NO_ADJ");
1953     break;
1954   case IDA_BAD_TB0:
1955     sprintf(name, "IDA_BAD_TB0");
1956     break;
1957   case IDA_REIFWD_FAIL:
1958     sprintf(name, "IDA_REIFWD_FAIL");
1959     break;
1960   case IDA_FWD_FAIL:
1961     sprintf(name, "IDA_FWD_FAIL");
1962     break;
1963   case IDA_GETY_BADT:
1964     sprintf(name, "IDA_GETY_BADT");
1965     break;
1966   case IDA_NO_BCK:
1967     sprintf(name, "IDA_NO_BCK");
1968     break;
1969   case IDA_NO_FWD:
1970     sprintf(name,"IDA_NO_FWD");
1971     break;
1972   case IDA_NLS_SETUP_FAIL:
1973     sprintf(name,"IDA_NLS_SETUP_FAIL");
1974     break;
1975   case IDA_NLS_FAIL:
1976     sprintf(name,"IDA_NLS_FAIL");
1977     break;
1978   default:
1979     sprintf(name,"NONE");
1980   }
1981 
1982   return(name);
1983 }
1984