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-2021, 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 input and output
19  * functions for the adjoint module in the IDAS solver.
20  * -----------------------------------------------------------------
21  */
22 
23 
24 /*
25  * =================================================================
26  * IMPORTED HEADER FILES
27  * =================================================================
28  */
29 
30 #include <stdio.h>
31 #include <stdlib.h>
32 
33 #include "idas_impl.h"
34 #include <sundials/sundials_types.h>
35 
36 /*
37  * =================================================================
38  * IDAA PRIVATE CONSTANTS
39  * =================================================================
40  */
41 
42 #define ONE         RCONST(1.0)
43 
44 /*
45  * -----------------------------------------------------------------
46  * Optional input functions for ASA
47  * -----------------------------------------------------------------
48  */
49 
50 /*
51  * -----------------------------------------------------------------
52  * IDAAdjSetNoSensi
53  * -----------------------------------------------------------------
54  * Disables the forward sensitivity analysis in IDASolveF.
55  * -----------------------------------------------------------------
56  */
57 
IDAAdjSetNoSensi(void * ida_mem)58 int IDAAdjSetNoSensi(void *ida_mem)
59 {
60   IDAMem IDA_mem;
61   IDAadjMem IDAADJ_mem;
62 
63   /* Is ida_mem valid? */
64   if (ida_mem == NULL) {
65     IDAProcessError(NULL, IDA_MEM_NULL, "IDAA", "IDAAdjSetNoSensi", MSGAM_NULL_IDAMEM);
66     return IDA_MEM_NULL;
67   }
68   IDA_mem = (IDAMem) ida_mem;
69 
70   /* Is ASA initialized? */
71   if (IDA_mem->ida_adjMallocDone == SUNFALSE) {
72     IDAProcessError(IDA_mem, IDA_NO_ADJ, "IDAA", "IDAAdjSetNoSensi",  MSGAM_NO_ADJ);
73     return(IDA_NO_ADJ);
74   }
75   IDAADJ_mem = IDA_mem->ida_adj_mem;
76 
77   IDAADJ_mem->ia_storeSensi = SUNFALSE;
78 
79   return(IDA_SUCCESS);
80 }
81 
82 /*
83  * -----------------------------------------------------------------
84  * Optional input functions for backward integration
85  * -----------------------------------------------------------------
86  */
87 
IDASetNonlinearSolverB(void * ida_mem,int which,SUNNonlinearSolver NLS)88 int IDASetNonlinearSolverB(void *ida_mem, int which, SUNNonlinearSolver NLS)
89 {
90   IDAMem IDA_mem;
91   IDAadjMem IDAADJ_mem;
92   IDABMem IDAB_mem;
93   void *ida_memB;
94 
95   /* Check if ida_mem exists */
96   if (ida_mem == NULL) {
97     IDAProcessError(NULL, IDA_MEM_NULL, "IDAA",
98                     "IDASetNonlinearSolverB", MSGAM_NULL_IDAMEM);
99     return(IDA_MEM_NULL);
100   }
101   IDA_mem = (IDAMem) ida_mem;
102 
103   /* Was ASA initialized? */
104   if (IDA_mem->ida_adjMallocDone == SUNFALSE) {
105     IDAProcessError(IDA_mem, IDA_NO_ADJ, "IDAA",
106                     "IDASetNonlinearSolverB", MSGAM_NO_ADJ);
107     return(IDA_NO_ADJ);
108   }
109   IDAADJ_mem = IDA_mem->ida_adj_mem;
110 
111   /* Check the value of which */
112   if ( which >= IDAADJ_mem->ia_nbckpbs ) {
113     IDAProcessError(IDA_mem, IDA_ILL_INPUT, "IDAA",
114                     "IDASetNonlinearSolverB", MSGAM_BAD_WHICH);
115     return(IDA_ILL_INPUT);
116   }
117 
118   /* Find the IDABMem entry in the linked list corresponding to 'which' */
119   IDAB_mem = IDAADJ_mem->IDAB_mem;
120   while (IDAB_mem != NULL) {
121     if ( which == IDAB_mem->ida_index ) break;
122     /* advance */
123     IDAB_mem = IDAB_mem->ida_next;
124   }
125 
126   ida_memB = (void *) (IDAB_mem->IDA_mem);
127 
128   return(IDASetNonlinearSolver(ida_memB, NLS));
129 }
130 
IDASetUserDataB(void * ida_mem,int which,void * user_dataB)131 int IDASetUserDataB(void *ida_mem, int which, void *user_dataB)
132 {
133   IDAMem IDA_mem;
134   IDAadjMem IDAADJ_mem;
135   IDABMem IDAB_mem;
136 
137   /* Is ida_mem valid? */
138   if (ida_mem == NULL) {
139     IDAProcessError(NULL, IDA_MEM_NULL, "IDAA", "IDASetUserDataB", MSGAM_NULL_IDAMEM);
140     return IDA_MEM_NULL;
141   }
142   IDA_mem = (IDAMem) ida_mem;
143 
144   /* Is ASA initialized? */
145   if (IDA_mem->ida_adjMallocDone == SUNFALSE) {
146     IDAProcessError(IDA_mem, IDA_NO_ADJ, "IDAA", "IDASetUserDataB",  MSGAM_NO_ADJ);
147     return(IDA_NO_ADJ);
148   }
149   IDAADJ_mem = IDA_mem->ida_adj_mem;
150 
151   /* Check the value of which */
152   if ( which >= IDAADJ_mem->ia_nbckpbs ) {
153     IDAProcessError(IDA_mem, IDA_ILL_INPUT, "IDAA", "IDASetUserDataB", MSGAM_BAD_WHICH);
154     return(IDA_ILL_INPUT);
155   }
156 
157   /* Find the IDABMem entry in the linked list corresponding to 'which'. */
158   IDAB_mem = IDAADJ_mem->IDAB_mem;
159   while (IDAB_mem != NULL) {
160     if( which == IDAB_mem->ida_index ) break;
161     /* advance */
162     IDAB_mem = IDAB_mem->ida_next;
163   }
164 
165   /* Set user data for this backward problem. */
166   IDAB_mem->ida_user_data = user_dataB;
167 
168   return(IDA_SUCCESS);
169 }
170 
IDASetMaxOrdB(void * ida_mem,int which,int maxordB)171 int IDASetMaxOrdB(void *ida_mem, int which, int maxordB)
172 {
173   IDAMem IDA_mem;
174   IDAadjMem IDAADJ_mem;
175   IDABMem IDAB_mem;
176   void *ida_memB;
177 
178   /* Is ida_mem valid? */
179   if (ida_mem == NULL) {
180     IDAProcessError(NULL, IDA_MEM_NULL, "IDAA", "IDASetMaxOrdB", MSGAM_NULL_IDAMEM);
181     return IDA_MEM_NULL;
182   }
183   IDA_mem = (IDAMem) ida_mem;
184 
185   /* Is ASA initialized? */
186   if (IDA_mem->ida_adjMallocDone == SUNFALSE) {
187     IDAProcessError(IDA_mem, IDA_NO_ADJ, "IDAA", "IDASetMaxOrdB",  MSGAM_NO_ADJ);
188     return(IDA_NO_ADJ);
189   }
190   IDAADJ_mem = IDA_mem->ida_adj_mem;
191 
192   /* Check the value of which */
193   if ( which >= IDAADJ_mem->ia_nbckpbs ) {
194     IDAProcessError(IDA_mem, IDA_ILL_INPUT, "IDAA", "IDASetMaxOrdB", MSGAM_BAD_WHICH);
195     return(IDA_ILL_INPUT);
196   }
197 
198   /* Find the IDABMem entry in the linked list corresponding to 'which'. */
199   IDAB_mem = IDAADJ_mem->IDAB_mem;
200   while (IDAB_mem != NULL) {
201     if( which == IDAB_mem->ida_index ) break;
202     /* advance */
203     IDAB_mem = IDAB_mem->ida_next;
204   }
205   ida_memB = (void *) IDAB_mem->IDA_mem;
206 
207   return IDASetMaxOrd(ida_memB, maxordB);
208 }
209 
IDASetMaxNumStepsB(void * ida_mem,int which,long int mxstepsB)210 int IDASetMaxNumStepsB(void *ida_mem, int which, long int mxstepsB)
211 {
212   IDAMem IDA_mem;
213   IDAadjMem IDAADJ_mem;
214   IDABMem IDAB_mem;
215   void *ida_memB;
216 
217   /* Is ida_mem valid? */
218   if (ida_mem == NULL) {
219     IDAProcessError(NULL, IDA_MEM_NULL, "IDAA", "IDASetMaxNumStepsB", MSGAM_NULL_IDAMEM);
220     return IDA_MEM_NULL;
221   }
222   IDA_mem = (IDAMem) ida_mem;
223 
224   /* Is ASA initialized? */
225   if (IDA_mem->ida_adjMallocDone == SUNFALSE) {
226     IDAProcessError(IDA_mem, IDA_NO_ADJ, "IDAA", "IDASetMaxNumStepsB",  MSGAM_NO_ADJ);
227     return(IDA_NO_ADJ);
228   }
229   IDAADJ_mem = IDA_mem->ida_adj_mem;
230 
231   /* Check the value of which */
232   if ( which >= IDAADJ_mem->ia_nbckpbs ) {
233     IDAProcessError(IDA_mem, IDA_ILL_INPUT, "IDAA", "IDASetMaxNumStepsB", MSGAM_BAD_WHICH);
234     return(IDA_ILL_INPUT);
235   }
236 
237   /* Find the IDABMem entry in the linked list corresponding to 'which'. */
238   IDAB_mem = IDAADJ_mem->IDAB_mem;
239   while (IDAB_mem != NULL) {
240     if( which == IDAB_mem->ida_index ) break;
241     /* advance */
242     IDAB_mem = IDAB_mem->ida_next;
243   }
244   ida_memB = (void *) IDAB_mem->IDA_mem;
245 
246   return IDASetMaxNumSteps(ida_memB, mxstepsB);
247 }
248 
IDASetInitStepB(void * ida_mem,int which,realtype hinB)249 int IDASetInitStepB(void *ida_mem, int which, realtype hinB)
250 {
251   IDAMem IDA_mem;
252   IDAadjMem IDAADJ_mem;
253   IDABMem IDAB_mem;
254   void *ida_memB;
255 
256   /* Is ida_mem valid? */
257   if (ida_mem == NULL) {
258     IDAProcessError(NULL, IDA_MEM_NULL, "IDAA", "IDASetInitStepB", MSGAM_NULL_IDAMEM);
259     return IDA_MEM_NULL;
260   }
261   IDA_mem = (IDAMem) ida_mem;
262 
263   /* Is ASA initialized? */
264   if (IDA_mem->ida_adjMallocDone == SUNFALSE) {
265     IDAProcessError(IDA_mem, IDA_NO_ADJ, "IDAA", "IDASetInitStepB",  MSGAM_NO_ADJ);
266     return(IDA_NO_ADJ);
267   }
268   IDAADJ_mem = IDA_mem->ida_adj_mem;
269 
270   /* Check the value of which */
271   if ( which >= IDAADJ_mem->ia_nbckpbs ) {
272     IDAProcessError(IDA_mem, IDA_ILL_INPUT, "IDAA", "IDASetInitStepB", MSGAM_BAD_WHICH);
273     return(IDA_ILL_INPUT);
274   }
275 
276   /* Find the IDABMem entry in the linked list corresponding to 'which'. */
277   IDAB_mem = IDAADJ_mem->IDAB_mem;
278   while (IDAB_mem != NULL) {
279     if( which == IDAB_mem->ida_index ) break;
280     /* advance */
281     IDAB_mem = IDAB_mem->ida_next;
282   }
283   ida_memB = (void *) IDAB_mem->IDA_mem;
284 
285   return IDASetInitStep(ida_memB, hinB);
286 }
287 
IDASetMaxStepB(void * ida_mem,int which,realtype hmaxB)288 int IDASetMaxStepB(void *ida_mem, int which, realtype hmaxB)
289 {
290   IDAMem IDA_mem;
291   IDAadjMem IDAADJ_mem;
292   IDABMem IDAB_mem;
293   void *ida_memB;
294 
295   /* Is ida_mem valid? */
296   if (ida_mem == NULL) {
297     IDAProcessError(NULL, IDA_MEM_NULL, "IDAA", "IDASetMaxStepB", MSGAM_NULL_IDAMEM);
298     return IDA_MEM_NULL;
299   }
300   IDA_mem = (IDAMem) ida_mem;
301 
302   /* Is ASA initialized? */
303   if (IDA_mem->ida_adjMallocDone == SUNFALSE) {
304     IDAProcessError(IDA_mem, IDA_NO_ADJ, "IDAA", "IDASetMaxStepB",  MSGAM_NO_ADJ);
305     return(IDA_NO_ADJ);
306   }
307   IDAADJ_mem = IDA_mem->ida_adj_mem;
308 
309   /* Check the value of which */
310   if ( which >= IDAADJ_mem->ia_nbckpbs ) {
311     IDAProcessError(IDA_mem, IDA_ILL_INPUT, "IDAA", "IDASetMaxStepB", MSGAM_BAD_WHICH);
312     return(IDA_ILL_INPUT);
313   }
314 
315   /* Find the IDABMem entry in the linked list corresponding to 'which'. */
316   IDAB_mem = IDAADJ_mem->IDAB_mem;
317   while (IDAB_mem != NULL) {
318     if( which == IDAB_mem->ida_index ) break;
319     /* advance */
320     IDAB_mem = IDAB_mem->ida_next;
321   }
322   ida_memB = (void *) IDAB_mem->IDA_mem;
323 
324   return IDASetMaxStep(ida_memB, hmaxB);
325 }
326 
IDASetSuppressAlgB(void * ida_mem,int which,booleantype suppressalgB)327 int IDASetSuppressAlgB(void *ida_mem, int which,  booleantype suppressalgB)
328 {
329   IDAMem IDA_mem;
330   IDAadjMem IDAADJ_mem;
331   IDABMem IDAB_mem;
332   void *ida_memB;
333 
334   /* Is ida_mem valid? */
335   if (ida_mem == NULL) {
336     IDAProcessError(NULL, IDA_MEM_NULL, "IDAA", "IDASetSuppressAlgB", MSGAM_NULL_IDAMEM);
337     return IDA_MEM_NULL;
338   }
339   IDA_mem = (IDAMem) ida_mem;
340 
341   /* Is ASA initialized? */
342   if (IDA_mem->ida_adjMallocDone == SUNFALSE) {
343     IDAProcessError(IDA_mem, IDA_NO_ADJ, "IDAA", "IDASetSuppressAlgB",  MSGAM_NO_ADJ);
344     return(IDA_NO_ADJ);
345   }
346   IDAADJ_mem = IDA_mem->ida_adj_mem;
347 
348   /* Check the value of which */
349   if ( which >= IDAADJ_mem->ia_nbckpbs ) {
350     IDAProcessError(IDA_mem, IDA_ILL_INPUT, "IDAA", "IDASetSuppressAlgB", MSGAM_BAD_WHICH);
351     return(IDA_ILL_INPUT);
352   }
353 
354   /* Find the IDABMem entry in the linked list corresponding to 'which'. */
355   IDAB_mem = IDAADJ_mem->IDAB_mem;
356   while (IDAB_mem != NULL) {
357     if( which == IDAB_mem->ida_index ) break;
358     /* advance */
359     IDAB_mem = IDAB_mem->ida_next;
360   }
361   ida_memB = (void *) IDAB_mem->IDA_mem;
362 
363   return IDASetSuppressAlg(ida_memB, suppressalgB);
364 }
365 
IDASetIdB(void * ida_mem,int which,N_Vector idB)366 int IDASetIdB(void *ida_mem, int which,  N_Vector idB)
367 {
368   IDAMem IDA_mem;
369   IDAadjMem IDAADJ_mem;
370   IDABMem IDAB_mem;
371   void *ida_memB;
372 
373   /* Is ida_mem valid? */
374   if (ida_mem == NULL) {
375     IDAProcessError(NULL, IDA_MEM_NULL, "IDAA", "IDASetIdB", MSGAM_NULL_IDAMEM);
376     return IDA_MEM_NULL;
377   }
378   IDA_mem = (IDAMem) ida_mem;
379 
380   /* Is ASA initialized? */
381   if (IDA_mem->ida_adjMallocDone == SUNFALSE) {
382     IDAProcessError(IDA_mem, IDA_NO_ADJ, "IDAA", "IDASetIdB",  MSGAM_NO_ADJ);
383     return(IDA_NO_ADJ);
384   }
385   IDAADJ_mem = IDA_mem->ida_adj_mem;
386 
387   /* Check the value of which */
388   if ( which >= IDAADJ_mem->ia_nbckpbs ) {
389     IDAProcessError(IDA_mem, IDA_ILL_INPUT, "IDAA", "IDASetIdB", MSGAM_BAD_WHICH);
390     return(IDA_ILL_INPUT);
391   }
392 
393   /* Find the IDABMem entry in the linked list corresponding to 'which'. */
394   IDAB_mem = IDAADJ_mem->IDAB_mem;
395   while (IDAB_mem != NULL) {
396     if( which == IDAB_mem->ida_index ) break;
397     /* advance */
398     IDAB_mem = IDAB_mem->ida_next;
399   }
400   ida_memB = (void *) IDAB_mem->IDA_mem;
401 
402   return IDASetId(ida_memB, idB);
403 }
404 
IDASetConstraintsB(void * ida_mem,int which,N_Vector constraintsB)405 int IDASetConstraintsB(void *ida_mem, int which,  N_Vector constraintsB)
406 {
407   IDAMem IDA_mem;
408   IDAadjMem IDAADJ_mem;
409   IDABMem IDAB_mem;
410   void *ida_memB;
411 
412   /* Is ida_mem valid? */
413   if (ida_mem == NULL) {
414     IDAProcessError(NULL, IDA_MEM_NULL, "IDAA", "IDASetConstraintsB", MSGAM_NULL_IDAMEM);
415     return IDA_MEM_NULL;
416   }
417   IDA_mem = (IDAMem) ida_mem;
418 
419   /* Is ASA initialized? */
420   if (IDA_mem->ida_adjMallocDone == SUNFALSE) {
421     IDAProcessError(IDA_mem, IDA_NO_ADJ, "IDAA", "IDASetConstraintsB",  MSGAM_NO_ADJ);
422     return(IDA_NO_ADJ);
423   }
424   IDAADJ_mem = IDA_mem->ida_adj_mem;
425 
426   /* Check the value of which */
427   if ( which >= IDAADJ_mem->ia_nbckpbs ) {
428     IDAProcessError(IDA_mem, IDA_ILL_INPUT, "IDAA", "IDASetConstraintsB", MSGAM_BAD_WHICH);
429     return(IDA_ILL_INPUT);
430   }
431 
432   /* Find the IDABMem entry in the linked list corresponding to 'which'. */
433   IDAB_mem = IDAADJ_mem->IDAB_mem;
434   while (IDAB_mem != NULL) {
435     if( which == IDAB_mem->ida_index ) break;
436     /* advance */
437     IDAB_mem = IDAB_mem->ida_next;
438   }
439   ida_memB = (void *) IDAB_mem->IDA_mem;
440 
441   return IDASetConstraints(ida_memB, constraintsB);
442 }
443 /*
444  * ----------------------------------------------------------------
445  * Input quadrature functions for ASA
446  * ----------------------------------------------------------------
447  */
448 
IDASetQuadErrConB(void * ida_mem,int which,int errconQB)449 int IDASetQuadErrConB(void *ida_mem, int which, int errconQB)
450 {
451   IDAMem IDA_mem;
452   IDAadjMem IDAADJ_mem;
453   IDABMem IDAB_mem;
454   void *ida_memB;
455 
456   /* Is ida_mem valid? */
457   if (ida_mem == NULL) {
458     IDAProcessError(NULL, IDA_MEM_NULL, "IDAA", "IDASetQuadErrConB", MSGAM_NULL_IDAMEM);
459     return IDA_MEM_NULL;
460   }
461   IDA_mem = (IDAMem) ida_mem;
462 
463   /* Is ASA initialized? */
464   if (IDA_mem->ida_adjMallocDone == SUNFALSE) {
465     IDAProcessError(IDA_mem, IDA_NO_ADJ, "IDAA", "IDASetQuadErrConB",  MSGAM_NO_ADJ);
466     return(IDA_NO_ADJ);
467   }
468   IDAADJ_mem = IDA_mem->ida_adj_mem;
469 
470   /* Check the value of which */
471   if ( which >= IDAADJ_mem->ia_nbckpbs ) {
472     IDAProcessError(IDA_mem, IDA_ILL_INPUT, "IDAA", "IDASetQuadErrConB", MSGAM_BAD_WHICH);
473     return(IDA_ILL_INPUT);
474   }
475 
476   /* Find the IDABMem entry in the linked list corresponding to 'which'. */
477   IDAB_mem = IDAADJ_mem->IDAB_mem;
478   while (IDAB_mem != NULL) {
479     if( which == IDAB_mem->ida_index ) break;
480     /* advance */
481     IDAB_mem = IDAB_mem->ida_next;
482   }
483   ida_memB = (void *) IDAB_mem->IDA_mem;
484 
485   return IDASetQuadErrCon(ida_memB, errconQB);
486 
487 }
488 
489 /*
490  * -----------------------------------------------------------------
491  * Optional output functions for backward integration
492  * -----------------------------------------------------------------
493  */
494 
495 /*
496  * IDAGetAdjIDABmem
497  *
498  * This function returns a (void *) pointer to the IDAS
499  * memory allocated for the backward problem. This pointer can
500  * then be used to call any of the IDAGet* IDAS routines to
501  * extract optional output for the backward integration phase.
502  */
503 
IDAGetAdjIDABmem(void * ida_mem,int which)504 void *IDAGetAdjIDABmem(void *ida_mem, int which)
505 {
506   IDAMem IDA_mem;
507   IDAadjMem IDAADJ_mem;
508   IDABMem IDAB_mem;
509   void *ida_memB;
510 
511   /* Is ida_mem valid? */
512   if (ida_mem == NULL) {
513     IDAProcessError(NULL, 0, "IDAA", "IDAGetAdjIDABmem", MSGAM_NULL_IDAMEM);
514     return(NULL);
515   }
516   IDA_mem = (IDAMem) ida_mem;
517 
518   /* Is ASA initialized? */
519   if (IDA_mem->ida_adjMallocDone == SUNFALSE) {
520     IDAProcessError(IDA_mem, 0, "IDAA", "IDAGetAdjIDABmem",  MSGAM_NO_ADJ);
521     return(NULL);
522   }
523   IDAADJ_mem = IDA_mem->ida_adj_mem;
524 
525   /* Check the value of which */
526   if ( which >= IDAADJ_mem->ia_nbckpbs ) {
527     IDAProcessError(IDA_mem, 0, "IDAA", "IDAGetAdjIDABmem", MSGAM_BAD_WHICH);
528     return(NULL);
529   }
530 
531   /* Find the IDABMem entry in the linked list corresponding to 'which'. */
532   IDAB_mem = IDAADJ_mem->IDAB_mem;
533   while (IDAB_mem != NULL) {
534     if( which == IDAB_mem->ida_index ) break;
535     /* advance */
536     IDAB_mem = IDAB_mem->ida_next;
537   }
538   ida_memB = (void *) IDAB_mem->IDA_mem;
539 
540   return(ida_memB);
541 }
542 
543 /*
544  * IDAGetAdjCheckPointsInfo
545  *
546  * Loads an array of nckpnts structures of type IDAadjCheckPointRec
547  * defined below.
548  *
549  * The user must allocate space for ckpnt (ncheck+1).
550  */
551 
IDAGetAdjCheckPointsInfo(void * ida_mem,IDAadjCheckPointRec * ckpnt)552 int IDAGetAdjCheckPointsInfo(void *ida_mem, IDAadjCheckPointRec *ckpnt)
553 {
554   IDAMem IDA_mem;
555   IDAadjMem IDAADJ_mem;
556   CkpntMem ck_mem;
557   int i;
558 
559   /* Is ida_mem valid? */
560   if (ida_mem == NULL) {
561     IDAProcessError(NULL, IDA_MEM_NULL, "IDAA", "IDAGetAdjCheckPointsInfo", MSGAM_NULL_IDAMEM);
562     return(IDA_MEM_NULL);
563   }
564   IDA_mem = (IDAMem) ida_mem;
565 
566   /* Is ASA initialized? */
567   if (IDA_mem->ida_adjMallocDone == SUNFALSE) {
568     IDAProcessError(IDA_mem, IDA_NO_ADJ, "IDAA", "IDAGetAdjCheckPointsInfo",  MSGAM_NO_ADJ);
569     return(IDA_NO_ADJ);
570   }
571   IDAADJ_mem = IDA_mem->ida_adj_mem;
572 
573   i=0;
574   ck_mem = IDAADJ_mem->ck_mem;
575   while (ck_mem != NULL) {
576 
577     ckpnt[i].my_addr = (void *) ck_mem;
578     ckpnt[i].next_addr = (void *) ck_mem->ck_next;
579     ckpnt[i].t0 = ck_mem->ck_t0;
580     ckpnt[i].t1 = ck_mem->ck_t1;
581     ckpnt[i].nstep = ck_mem->ck_nst;
582     ckpnt[i].order = ck_mem->ck_kk;
583     ckpnt[i].step = ck_mem->ck_hh;
584 
585     ck_mem = ck_mem->ck_next;
586     i++;
587   }
588 
589   return(IDA_SUCCESS);
590 }
591 
592 /* IDAGetConsistentICB
593  *
594  * Returns the consistent initial conditions computed by IDACalcICB or
595  * IDACalcICBS
596  *
597  * It must be preceded by a successful call to IDACalcICB or IDACalcICBS
598  * for 'which' backward problem.
599  */
600 
IDAGetConsistentICB(void * ida_mem,int which,N_Vector yyB0_mod,N_Vector ypB0_mod)601 int IDAGetConsistentICB(void *ida_mem, int which, N_Vector yyB0_mod, N_Vector ypB0_mod)
602 {
603   IDAMem IDA_mem;
604   IDAadjMem IDAADJ_mem;
605   IDABMem IDAB_mem;
606   void *ida_memB;
607   int flag;
608 
609   /* Is ida_mem valid? */
610   if (ida_mem == NULL) {
611     IDAProcessError(NULL, IDA_MEM_NULL, "IDAA", "IDAGetConsistentICB", MSGAM_NULL_IDAMEM);
612     return IDA_MEM_NULL;
613   }
614   IDA_mem = (IDAMem) ida_mem;
615 
616   /* Is ASA initialized? */
617   if (IDA_mem->ida_adjMallocDone == SUNFALSE) {
618     IDAProcessError(IDA_mem, IDA_NO_ADJ, "IDAA", "IDAGetConsistentICB",  MSGAM_NO_ADJ);
619     return(IDA_NO_ADJ);
620   }
621   IDAADJ_mem = IDA_mem->ida_adj_mem;
622 
623   /* Check the value of which */
624   if ( which >= IDAADJ_mem->ia_nbckpbs ) {
625     IDAProcessError(IDA_mem, IDA_ILL_INPUT, "IDAA", "IDAGetConsistentICB", MSGAM_BAD_WHICH);
626     return(IDA_ILL_INPUT);
627   }
628 
629   /* Find the IDABMem entry in the linked list corresponding to 'which'. */
630   IDAB_mem = IDAADJ_mem->IDAB_mem;
631   while (IDAB_mem != NULL) {
632     if( which == IDAB_mem->ida_index ) break;
633     /* advance */
634     IDAB_mem = IDAB_mem->ida_next;
635   }
636   ida_memB = (void *) IDAB_mem->IDA_mem;
637 
638   flag = IDAGetConsistentIC(ida_memB, yyB0_mod, ypB0_mod);
639 
640   return(flag);
641 }
642 
643 
644 /*
645  * -----------------------------------------------------------------
646  * Undocumented development user-callable functions
647  * -----------------------------------------------------------------
648  */
649 
650 /*
651  * -----------------------------------------------------------------
652  * IDAGetAdjDataPointHermite
653  * -----------------------------------------------------------------
654  * Returns the 2 vectors stored for cubic Hermite interpolation at
655  * the data point 'which'. The user must allocate space for yy and
656  * yd.
657  *
658  * Returns IDA_MEM_NULL if ida_mem is NULL, IDA_ILL_INPUT if the
659  * interpolation type previously specified is not IDA_HERMITE or
660  * IDA_SUCCESS otherwise.
661  *
662  */
IDAGetAdjDataPointHermite(void * ida_mem,int which,realtype * t,N_Vector yy,N_Vector yd)663 int IDAGetAdjDataPointHermite(void *ida_mem, int which,
664                               realtype *t, N_Vector yy, N_Vector yd)
665 
666 {
667   IDAMem IDA_mem;
668   IDAadjMem IDAADJ_mem;
669   DtpntMem *dt_mem;
670   HermiteDataMem content;
671 
672     /* Is ida_mem valid? */
673   if (ida_mem == NULL) {
674     IDAProcessError(NULL, IDA_MEM_NULL, "IDAA", "IDAGetAdjDataPointHermite", MSGAM_NULL_IDAMEM);
675     return(IDA_MEM_NULL);
676   }
677   IDA_mem = (IDAMem) ida_mem;
678 
679   /* Is ASA initialized? */
680   if (IDA_mem->ida_adjMallocDone == SUNFALSE) {
681     IDAProcessError(IDA_mem, IDA_NO_ADJ, "IDAA", "IDAGetAdjDataPointHermite",  MSGAM_NO_ADJ);
682     return(IDA_NO_ADJ);
683   }
684   IDAADJ_mem = IDA_mem->ida_adj_mem;
685 
686   dt_mem = IDAADJ_mem->dt_mem;
687 
688   if (IDAADJ_mem->ia_interpType != IDA_HERMITE) {
689     IDAProcessError(IDA_mem, IDA_ILL_INPUT, "IDAA", "IDAGetAdjDataPointHermite",  MSGAM_WRONG_INTERP);
690     return(IDA_ILL_INPUT);
691   }
692 
693   *t = dt_mem[which]->t;
694   content = (HermiteDataMem) dt_mem[which]->content;
695 
696   if (yy != NULL) N_VScale(ONE, content->y, yy);
697   if (yd != NULL) N_VScale(ONE, content->yd, yd);
698 
699   return(IDA_SUCCESS);
700 }
701 
702 /*
703  * IDAGetAdjDataPointPolynomial
704  *
705  * Returns the vector stored for polynomial interpolation at the
706  * data point 'which'. The user must allocate space for y.
707  *
708  * Returns IDA_MEM_NULL if ida_mem is NULL, IDA_ILL_INPUT if the
709  * interpolation type previously specified is not IDA_POLYNOMIAL or
710  * IDA_SUCCESS otherwise.
711  */
712 
713 
IDAGetAdjDataPointPolynomial(void * ida_mem,int which,realtype * t,int * order,N_Vector y)714 int IDAGetAdjDataPointPolynomial(void *ida_mem, int which,
715                                  realtype *t, int *order, N_Vector y)
716 {
717   IDAMem IDA_mem;
718   IDAadjMem IDAADJ_mem;
719   DtpntMem *dt_mem;
720   PolynomialDataMem content;
721   /* Is ida_mem valid? */
722   if (ida_mem == NULL) {
723     IDAProcessError(NULL, IDA_MEM_NULL, "IDAA", "IDAGetAdjDataPointPolynomial", MSGAM_NULL_IDAMEM);
724     return(IDA_MEM_NULL);
725   }
726   IDA_mem = (IDAMem) ida_mem;
727 
728   /* Is ASA initialized? */
729   if (IDA_mem->ida_adjMallocDone == SUNFALSE) {
730     IDAProcessError(IDA_mem, IDA_NO_ADJ, "IDAA", "IDAGetAdjDataPointPolynomial", MSGAM_NO_ADJ);
731     return(IDA_NO_ADJ);
732   }
733   IDAADJ_mem = IDA_mem->ida_adj_mem;
734 
735   dt_mem = IDAADJ_mem->dt_mem;
736 
737   if (IDAADJ_mem->ia_interpType != IDA_POLYNOMIAL) {
738     IDAProcessError(IDA_mem, IDA_ILL_INPUT, "IDAA", "IDAGetAdjDataPointPolynomial", MSGAM_WRONG_INTERP);
739     return(IDA_ILL_INPUT);
740   }
741 
742   *t = dt_mem[which]->t;
743   content = (PolynomialDataMem) dt_mem[which]->content;
744 
745   if (y != NULL) N_VScale(ONE, content->y, y);
746 
747   *order = content->order;
748 
749   return(IDA_SUCCESS);
750 }
751 
752 /*
753  * IDAGetAdjCurrentCheckPoint
754  *
755  * Returns the address of the 'active' check point.
756  */
757 
IDAGetAdjCurrentCheckPoint(void * ida_mem,void ** addr)758 int IDAGetAdjCurrentCheckPoint(void *ida_mem, void **addr)
759 {
760   IDAMem IDA_mem;
761   IDAadjMem IDAADJ_mem;
762 
763   if (ida_mem == NULL) {
764     IDAProcessError(NULL, IDA_MEM_NULL, "IDAA", "IDAGetAdjCurrentCheckPoint", MSGAM_NULL_IDAMEM);
765     return(IDA_MEM_NULL);
766   }
767   IDA_mem = (IDAMem) ida_mem;
768 
769   /* Is ASA initialized? */
770   if (IDA_mem->ida_adjMallocDone == SUNFALSE) {
771     IDAProcessError(IDA_mem, IDA_NO_ADJ, "IDAA", "IDAGetAdjCurrentCheckPoint",  MSGAM_NO_ADJ);
772     return(IDA_NO_ADJ);
773   }
774   IDAADJ_mem = IDA_mem->ida_adj_mem;
775 
776   *addr = (void *) IDAADJ_mem->ia_ckpntData;
777 
778   return(IDA_SUCCESS);
779 }
780 
781 
782 
783