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