Automatic Differentiation
 
Loading...
Searching...
No Matches
check_flag_sundials.hpp
Go to the documentation of this file.
1#ifndef STAN_MATH_PRIM_ERR_CHECK_FLAG_SUNDIALS_HPP
2#define STAN_MATH_PRIM_ERR_CHECK_FLAG_SUNDIALS_HPP
3
6#include <kinsol/kinsol.h>
7#include <cvodes/cvodes.h>
8#include <array>
9
10namespace stan {
11namespace math {
12
13#define CHECK_CVODES_CALL(call) cvodes_check(call, #call)
14#define CHECK_IDAS_CALL(call) idas_check(call, #call)
15#define CHECK_KINSOL_CALL(call) kinsol_check(call, #call)
16
28inline std::array<std::string, 2> cvodes_flag_msg(int flag) {
29 std::array<std::string, 2> msg;
30 switch (flag) {
31 case -1:
32 msg = {"CV_TOO_MUCH_WORK",
33 "The solver took mxstep internal steps but could not reach tout"};
34 break; // NOLINT
35 case -2:
36 msg = {"CV_TOO_MUCH_ACC",
37 "The solver could not satisfy the accuracy demanded by the user "
38 "for some internal step"};
39 break; // NOLINT
40 case -3:
41 msg = {"CV_ERR_FAILURE",
42 "Error test failures occurred too many times during one internal "
43 "time step or minimum step size was reached"};
44 break; // NOLINT
45 case -4:
46 msg = {"CV_CONV_FAILURE",
47 "Convergence test failures occurred too many times during one "
48 "internal time step or minimum step size was reached"};
49 break; // NOLINT
50 case -8:
51 msg = {"CV_RHSFUNC_FAIL",
52 "The right-hand side function failed in an unrecoverable manner"};
53 break; // NOLINT
54 case -9:
55 msg = {"CV_FIRST_RHSFUNC_ERR",
56 "The right-hand side function failed at the first call"};
57 break; // NOLINT
58 case -10:
59 msg = {"CV_REPTD_RHSFUNC_ERR",
60 "The right-hand side function had repetead recoverable errors"};
61 break; // NOLINT
62 case -11:
63 msg = {"CV_UNREC_RHSFUNC_ERR",
64 "The right-hand side function had a recoverable error, but no "
65 "recovery is possible"};
66 break; // NOLINT
67 case -27:
68 msg = {"CV_TOO_CLOSE",
69 "The output and initial times are too close to each other"};
70 break; // NOLINT
71 default:
72 switch (flag) {
73 case -5:
74 msg = {"CV_LINIT_FAIL",
75 "The linear solver's initialization function failed"};
76 break; // NOLINT
77 case -6:
78 msg = {"CV_LSETUP_FAIL",
79 "The linear solver's setup function failed in an "
80 "unrecoverable manner"};
81 break; // NOLINT
82 case -7:
83 msg = {"CV_LSOLVE_FAIL",
84 "The linear solver's solve function failed in an "
85 "unrecoverable manner"};
86 break; // NOLINT
87 case -20:
88 msg = {"CV_MEM_FAIL", "A memory allocation failed"};
89 break; // NOLINT
90 case -21:
91 msg = {"CV_MEM_NULL", "The cvode_mem argument was NULL"};
92 break; // NOLINT
93 case -22:
94 msg = {"CV_ILL_INPUT", "One of the function inputs is illegal"};
95 break; // NOLINT
96 case -23:
97 msg = {"CV_NO_MALLOC",
98 "The CVODE memory block was not allocated by a call to "
99 "CVodeMalloc"};
100 break; // NOLINT
101 case -24:
102 msg = {"CV_BAD_K",
103 "The derivative order k is larger than the order used"};
104 break; // NOLINT
105 case -25:
106 msg = {"CV_BAD_T", "The time t s outside the last step taken"};
107 break; // NOLINT
108 case -26:
109 msg = {"CV_BAD_DKY", "The output derivative vector is NULL"};
110 break; // NOLINT
111 case -40:
112 msg = {"CV_BAD_IS",
113 "The sensitivity index is larger than the number of "
114 "sensitivities computed"};
115 break; // NOLINT
116 case -41:
117 msg = {"CV_NO_SENS",
118 "Forward sensitivity integration was not activated"};
119 break; // NOLINT
120 case -42:
121 msg = {"CV_SRHSFUNC_FAIL",
122 "The sensitivity right-hand side function failed in an "
123 "unrecoverable manner"};
124 break; // NOLINT
125 case -43:
126 msg = {"CV_FIRST_SRHSFUNC_ER",
127 "The sensitivity right-hand side function failed at the first "
128 "call"};
129 break; // NOLINT
130 case -44:
131 msg = {"CV_REPTD_SRHSFUNC_ER",
132 "The sensitivity ight-hand side function had repetead "
133 "recoverable errors"};
134 break; // NOLINT
135 case -45:
136 msg = {"CV_UNREC_SRHSFUNC_ER",
137 "The sensitivity right-hand side function had a recoverable "
138 "error, but no recovery is possible"};
139 break; // NOLINT
140 case -101:
141 msg = {"CV_ADJMEM_NULL", "The cvadj_mem argument was NULL"};
142 break; // NOLINT
143 case -103:
144 msg = {"CV_BAD_TB0",
145 "The final time for the adjoint problem is outside the "
146 "interval over which the forward problem was solved"};
147 break; // NOLINT
148 case -104:
149 msg = {"CV_BCKMEM_NULL",
150 "The cvodes memory for the backward problem was not created"};
151 break; // NOLINT
152 case -105:
153 msg = {"CV_REIFWD_FAIL",
154 "Reinitialization of the forward problem failed at the first "
155 "checkpoint"};
156 break; // NOLINT
157 case -106:
158 msg = {
159 "CV_FWD_FAIL",
160 "An error occured during the integration of the forward problem"};
161 break; // NOLINT
162 case -107:
163 msg = {"CV_BAD_ITASK", "Wrong task for backward integration"};
164 break; // NOLINT
165 case -108:
166 msg = {"CV_BAD_TBOUT",
167 "The desired output time is outside the interval over which "
168 "the forward problem was solved"};
169 break; // NOLINT
170 case -109:
171 msg = {"CV_GETY_BADT", "Wrong time in interpolation function"};
172 break; // NOLINT
173 }
174 }
175 return msg;
176}
177
186inline void cvodes_check(int flag, const char* func_name) {
187 if (flag < 0) {
188 std::ostringstream ss;
189 ss << func_name << " failed with error flag " << flag << ": \n"
190 << cvodes_flag_msg(flag).at(1) << ".";
191 throw std::domain_error(ss.str());
192 }
193}
194
195inline std::array<std::string, 2> idas_flag_msg(int flag) {
196 std::array<std::string, 2> msg;
197 switch (flag) {
198 case -1:
199 msg = {"IDA_TOO_MUCH_WORK",
200 "The solver took mxstep internal steps but could not reach tout."};
201 break; // NOLINT
202 case -2:
203 msg = {"IDA_TOO_MUCH_ACC",
204 "The solver could not satisfy the accuracy demanded by the user "
205 "for some internal step."};
206 break; // NOLINT
207 case -3:
208 msg = {"IDA_ERR_FAIL",
209 "Error test failures occurred too many times during one internal "
210 "time step or minimum step size was reached."};
211 break; // NOLINT
212 case -4:
213 msg = {"IDA_CONV_FAIL",
214 "Convergence test failures occurred too many times during one "
215 "internal time step or minimum step size was reached."};
216 break; // NOLINT
217 case -5:
218 msg = {"IDA_LINIT_FAIL",
219 "The linear solver’s initialization function failed."};
220 break; // NOLINT
221 case -6:
222 msg = {"IDA_LSETUP_FAIL",
223 "The linear solver’s setup function failed in an unrecoverable "
224 "manner."};
225 break; // NOLINT
226 case -7:
227 msg = {"IDA_LSOLVE_FAIL",
228 "The linear solver’s solve function failed in an unrecoverable "
229 "manner."};
230 break; // NOLINT
231 case -8:
232 msg = {"IDA_RES_FAIL",
233 "The user-provided residual function failed in an unrecoverable "
234 "manner."};
235 break; // NOLINT
236 case -9:
237 msg = {"IDA_REP_RES_FAIL",
238 "The user-provided residual function repeatedly returned a "
239 "recoverable error flag, but the solver was unable to recover."};
240 break; // NOLINT
241 case -10:
242 msg = {"IDA_RTFUNC_FAIL",
243 "The rootfinding function failed in an unrecoverable manner."};
244 break; // NOLINT
245 case -11:
246 msg = {"IDA_CONSTR_FAIL",
247 "The inequality constraints were violated and the solver was "
248 "unable to recover."};
249 break; // NOLINT
250 case -12:
251 msg = {"IDA_FIRST_RES_FAIL",
252 "The user-provided residual function failed recoverably on the "
253 "first call."};
254 break; // NOLINT
255 case -13:
256 msg = {"IDA_LINESEARCH_FAIL", "The line search failed."};
257 break; // NOLINT
258 default:
259 switch (flag) {
260 case -14:
261 msg = {"IDA_NO_RECOVERY",
262 "The residual function, linear solver setup function, or "
263 "linear solver solve function had a recoverable failure, but "
264 "IDACalcIC could not recover."};
265 break; // NOLINT
266 case -15:
267 msg = {"IDA_NLS_INIT_FAIL",
268 "The nonlinear solver’s init routine failed."};
269 break; // NOLINT
270 case -16:
271 msg = {"IDA_NLS_SETUP_FAIL",
272 "The nonlinear solver’s setup routine failed."};
273 break; // NOLINT
274 case -20:
275 msg = {"IDA_MEM_NULL", "The ida mem argument was NULL."};
276 break; // NOLINT
277 case -21:
278 msg = {"IDA_MEM_FAIL", "A memory allocation failed."};
279 break; // NOLINT
280 case -22:
281 msg = {"IDA_ILL_INPUT", "One of the function inputs is illegal."};
282 break; // NOLINT
283 case -23:
284 msg = {"IDA_NO_MALLOC",
285 "The idas memory was not allocated by a call to IDAInit."};
286 break; // NOLINT
287 case -24:
288 msg = {"IDA_BAD_EWT", "Zero value of some error weight component."};
289 break; // NOLINT
290 case -25:
291 msg = {"IDA_BAD_K", "The k-th derivative is not available."};
292 break; // NOLINT
293 case -26:
294 msg = {"IDA_BAD_T", "The time t is outside the last step taken."};
295 break; // NOLINT
296 case -27:
297 msg = {
298 "IDA_BAD_DKY",
299 "The vector argument where derivative should be stored is NULL."};
300 break; // NOLINT
301 case -30:
302 msg = {"IDA_NO_QUAD", "Quadratures were not initialized."};
303 break; // NOLINT
304 case -31:
305 msg = {"IDA_QRHS_FAIL",
306 "The user-provided right-hand side function for quadratures "
307 "failed in an unrecoverable manner."};
308 break; // NOLINT
309 case -32:
310 msg = {"IDA_FIRST_QRHS_ERR",
311 "The user-provided right-hand side function for quadratures "
312 "failed -in an unrecoverable manner on the first call."};
313 break; // NOLINT
314 case -33:
315 msg = {"IDA_REP_QRHS_ERR",
316 "The user-provided right-hand side repeatedly returned a re- "
317 "coverable error flag, but the solver was unable to recover."};
318 break; // NOLINT
319 case -40:
320 msg = {"IDA_NO_SENS", "Sensitivities were not initialized."};
321 break; // NOLINT
322 case -41:
323 msg = {"IDA_SRES_FAIL",
324 "The user-provided sensitivity residual function failed in an "
325 "unrecoverable manner."};
326 break; // NOLINT
327 case -42:
328 msg = {"IDA_REP_SRES_ERR",
329 "The user-provided sensitivity residual function repeatedly "
330 "re- turned a recoverable error flag, but the solver was "
331 "unable to recover."};
332 break; // NOLINT
333 case -43:
334 msg = {"IDA_BAD_IS", "The sensitivity identifier is not valid."};
335 break; // NOLINT
336 case -50:
337 msg = {"IDA_NO_QUADSENS",
338 "Sensitivity-dependent quadratures were not initialized."};
339 break; // NOLINT
340 case -51:
341 msg = {"IDA_QSRHS_FAIL",
342 "The user-provided sensitivity-dependent quadrature right- "
343 "hand side function failed in an unrecoverable manner."};
344 break; // NOLINT
345 case -52:
346 msg = {"IDA_FIRST_QSRHS_ERR",
347 "The user-provided sensitivity-dependent quadrature right- "
348 "hand side function failed in an unrecoverable manner on the "
349 "first call."};
350 break; // NOLINT
351 case -53:
352 msg = {"IDA_REP_QSRHS_ERR",
353 "The user-provided sensitivity-dependent quadrature right- "
354 "hand side repeatedly returned a recoverable error flag, but "
355 "the solver was unable to recover."};
356 break; // NOLINT
357 }
358 }
359 return msg;
360}
361
362inline void idas_check(int flag, const char* func_name) {
363 if (flag < 0) {
364 std::ostringstream ss;
365 ss << func_name << " failed with error flag " << flag << ": \n"
366 << idas_flag_msg(flag).at(1);
367 throw std::domain_error(ss.str());
368 }
369}
370
380inline void kinsol_check(int flag, const char* func_name) {
381 if (flag < 0) {
382 std::ostringstream ss;
383 ss << "algebra_solver failed with error flag " << flag << ".";
384 throw std::runtime_error(ss.str());
385 }
386}
387
401inline void kinsol_check(int flag, const char* func_name,
402 long int max_num_steps) { // NOLINT(runtime/int)
403 std::ostringstream ss;
404 if (flag == -6) {
405 domain_error("algebra_solver", "maximum number of iterations",
406 max_num_steps, "(", ") was exceeded in the solve.");
407 } else if (flag == -11) {
408 ss << "The linear solver’s setup function failed "
409 << "in an unrecoverable manner.";
410 throw std::runtime_error(ss.str());
411 } else if (flag < 0) {
412 ss << "algebra_solver failed with error flag " << flag << ".";
413 throw std::runtime_error(ss.str());
414 }
415}
416
417} // namespace math
418} // namespace stan
419#endif
std::array< std::string, 2 > idas_flag_msg(int flag)
std::array< std::string, 2 > cvodes_flag_msg(int flag)
Map cvodes error flag to acutally error msg.
void idas_check(int flag, const char *func_name)
void cvodes_check(int flag, const char *func_name)
Throws a std::domain_error exception when a Sundial function fails (i.e.
void kinsol_check(int flag, const char *func_name)
Throws an exception message when the functions in KINSOL fails.
void domain_error(const char *function, const char *name, const T &y, const char *msg1, const char *msg2)
The lgamma implementation in stan-math is based on either the reentrant safe lgamma_r implementation ...