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