Reference documentation for deal.II version 9.5.0
\(\newcommand{\dealvcentcolon}{\mathrel{\mathop{:}}}\) \(\newcommand{\dealcoloneq}{\dealvcentcolon\mathrel{\mkern-1.2mu}=}\) \(\newcommand{\jump}[1]{\left[\!\left[ #1 \right]\!\right]}\) \(\newcommand{\average}[1]{\left\{\!\left\{ #1 \right\}\!\right\}}\)
Loading...
Searching...
No Matches
arkode.cc
Go to the documentation of this file.
1// ---------------------------------------------------------------------
2//
3// Copyright (C) 2017 - 2023 by the deal.II authors
4//
5// This file is part of the deal.II library.
6//
7// The deal.II library is free software; you can use it, redistribute
8// it, and/or modify it under the terms of the GNU Lesser General
9// Public License as published by the Free Software Foundation; either
10// version 2.1 of the License, or (at your option) any later version.
11// The full text of the license can be found in the file LICENSE.md at
12// the top level directory of deal.II.
13//
14// ---------------------------------------------------------------------
15
16
17
18#include <deal.II/base/config.h>
19
21
22#ifdef DEAL_II_WITH_SUNDIALS
23
26
30# ifdef DEAL_II_WITH_TRILINOS
33# endif
34# ifdef DEAL_II_WITH_PETSC
37# endif
38
42
43# include <arkode/arkode_arkstep.h>
44# include <sunlinsol/sunlinsol_spgmr.h>
45# include <sunnonlinsol/sunnonlinsol_fixedpoint.h>
46
47# include <iostream>
48
50
51namespace SUNDIALS
52{
53 template <typename VectorType>
55 : ARKode(data, MPI_COMM_SELF)
56 {}
57
58
59 template <typename VectorType>
61 const MPI_Comm mpi_comm)
62 : data(data)
63 , arkode_mem(nullptr)
65 , arkode_ctx(nullptr)
66# endif
67 , mpi_communicator(mpi_comm)
68 , last_end_time(data.initial_time)
69 , pending_exception(nullptr)
70 {
72
73# if DEAL_II_SUNDIALS_VERSION_GTE(6, 0, 0)
74 // SUNDIALS will always duplicate communicators if we provide them. This
75 // can cause problems if SUNDIALS is configured with MPI and we pass along
76 // MPI_COMM_SELF in a serial application as MPI won't be
77 // initialized. Hence, work around that by just not providing a
78 // communicator in that case.
79 const int status =
80 SUNContext_Create(mpi_communicator == MPI_COMM_SELF ? nullptr :
82 &arkode_ctx);
83 (void)status;
84 AssertARKode(status);
85# endif
86 }
87
88
89
90 template <typename VectorType>
92 {
93 ARKStepFree(&arkode_mem);
94
95# if DEAL_II_SUNDIALS_VERSION_GTE(6, 0, 0)
96 const int status = SUNContext_Free(&arkode_ctx);
97 (void)status;
98 AssertARKode(status);
99# endif
100
101 Assert(pending_exception == nullptr, ExcInternalError());
102 }
103
104
105
106 template <typename VectorType>
107 unsigned int
108 ARKode<VectorType>::solve_ode(VectorType &solution)
109 {
110 DiscreteTime time(data.initial_time,
111 data.final_time,
112 data.initial_step_size);
113
114 return do_evolve_time(solution, time, /* force_solver_restart = */ true);
115 }
116
117
118
119 template <typename VectorType>
120 unsigned int
122 const double intermediate_time,
123 const bool reset_solver)
124 {
126 intermediate_time > last_end_time,
128 "The requested intermediate time is smaller than the last requested "
129 "intermediate time."));
130
131 const bool do_reset = reset_solver || arkode_mem == nullptr;
132 DiscreteTime time(last_end_time, intermediate_time, data.initial_step_size);
133 return do_evolve_time(solution, time, do_reset);
134 }
135
136
137
138 template <typename VectorType>
139 int
141 DiscreteTime &time,
142 const bool do_reset)
143 {
144 if (do_reset)
145 {
146 reset(time.get_current_time(), time.get_next_step_size(), solution);
147 if (output_step)
148 output_step(time.get_current_time(),
149 solution,
150 time.get_step_number());
151 }
152 else
153 {
154 // If we don't do a full reset then we still need to fix the end time.
155 // In SUNDIALS 6 and later, SUNDIALS will not do timesteps if the
156 // current time is past the set end point (i.e., ARKStepEvolve will
157 // return ARK_TSTOP_RETURN).
158 const int status = ARKStepSetStopTime(arkode_mem, time.get_end_time());
159 (void)status;
160 AssertARKode(status);
161 }
162
163 auto solution_nvector = internal::make_nvector_view(solution
165 ,
166 arkode_ctx
167# endif
168 );
169
170 while (!time.is_at_end())
171 {
172 time.set_desired_next_step_size(data.output_period);
173
174 // Having set up all of the ancillary things, finally call the main
175 // ARKode function. One we return, check that what happened:
176 // - If we have a pending recoverable exception, ignore it if SUNDIAL's
177 // return code was zero -- in that case, SUNDIALS managed to indeed
178 // recover and we no longer need the exception
179 // - If we have any other exception, rethrow it
180 // - If no exception, test that SUNDIALS really did successfully return
181 Assert(pending_exception == nullptr, ExcInternalError());
182 double actual_next_time;
183 const auto status = ARKStepEvolve(arkode_mem,
184 time.get_next_time(),
185 solution_nvector,
186 &actual_next_time,
187 ARK_NORMAL);
188 if (pending_exception)
189 {
190 try
191 {
192 std::rethrow_exception(pending_exception);
193 }
194 catch (const RecoverableUserCallbackError &exc)
195 {
196 pending_exception = nullptr;
197 if (status == 0)
198 /* just eat the exception */;
199 else
200 throw;
201 }
202 catch (...)
203 {
204 pending_exception = nullptr;
205 throw;
206 }
207 }
208
209 AssertARKode(status);
210
211 // Then reflect this time advancement in our own DiscreteTime object:
212 time.set_next_step_size(actual_next_time - time.get_current_time());
213 time.advance_time();
214
215 // Finally check whether resets or output calls are desired at this
216 // time:
217 while (solver_should_restart(time.get_current_time(), solution))
218 reset(time.get_current_time(),
220 solution);
221
222 if (output_step)
223 output_step(time.get_current_time(),
224 solution,
225 time.get_step_number());
226 }
227 last_end_time = time.get_current_time();
228 return time.get_step_number();
229 }
230
231
232
233 template <typename VectorType>
234 void
235 ARKode<VectorType>::reset(const double current_time,
236 const double current_time_step,
237 const VectorType &solution)
238 {
239 last_end_time = current_time;
240 int status;
241 (void)status;
242
243# if DEAL_II_SUNDIALS_VERSION_GTE(6, 0, 0)
244 status = SUNContext_Free(&arkode_ctx);
245 AssertARKode(status);
246 // Same comment applies as in class constructor:
247 status =
248 SUNContext_Create(mpi_communicator == MPI_COMM_SELF ? nullptr :
249 &mpi_communicator,
250 &arkode_ctx);
251 AssertARKode(status);
252# endif
253
254 if (arkode_mem)
255 {
256 ARKStepFree(&arkode_mem);
257 // Initialization is version-dependent: do that in a moment
258 }
259
260 // just a view on the memory in solution, all write operations on yy by
261 // ARKODE will automatically be mirrored to solution
262 auto initial_condition_nvector = internal::make_nvector_view(solution
264 ,
265 arkode_ctx
266# endif
267 );
268
269 Assert(explicit_function || implicit_function,
270 ExcFunctionNotProvided("explicit_function || implicit_function"));
271
272 auto explicit_function_callback =
273 [](realtype tt, N_Vector yy, N_Vector yp, void *user_data) -> int {
274 Assert(user_data != nullptr, ExcInternalError());
275 ARKode<VectorType> &solver =
276 *static_cast<ARKode<VectorType> *>(user_data);
277
278 auto *src_yy = internal::unwrap_nvector_const<VectorType>(yy);
279 auto *dst_yp = internal::unwrap_nvector<VectorType>(yp);
280
282 solver.explicit_function,
283 solver.pending_exception,
284 tt,
285 *src_yy,
286 *dst_yp);
287 };
288
289
290 auto implicit_function_callback =
291 [](realtype tt, N_Vector yy, N_Vector yp, void *user_data) -> int {
292 Assert(user_data != nullptr, ExcInternalError());
293 ARKode<VectorType> &solver =
294 *static_cast<ARKode<VectorType> *>(user_data);
295
296 auto *src_yy = internal::unwrap_nvector_const<VectorType>(yy);
297 auto *dst_yp = internal::unwrap_nvector<VectorType>(yp);
298
300 solver.implicit_function,
301 solver.pending_exception,
302 tt,
303 *src_yy,
304 *dst_yp);
305 };
306
307 arkode_mem = ARKStepCreate(explicit_function ? explicit_function_callback :
308 ARKRhsFn(nullptr),
309 implicit_function ? implicit_function_callback :
310 ARKRhsFn(nullptr),
311 current_time,
312 initial_condition_nvector
313# if DEAL_II_SUNDIALS_VERSION_GTE(6, 0, 0)
314 ,
315 arkode_ctx
316# endif
317 );
318 Assert(arkode_mem != nullptr, ExcInternalError());
319
320 if (get_local_tolerances)
321 {
322 const auto abs_tols = internal::make_nvector_view(get_local_tolerances()
324 ,
325 arkode_ctx
326# endif
327 );
328 status =
329 ARKStepSVtolerances(arkode_mem, data.relative_tolerance, abs_tols);
330 AssertARKode(status);
331 }
332 else
333 {
334 status = ARKStepSStolerances(arkode_mem,
335 data.relative_tolerance,
336 data.absolute_tolerance);
337 AssertARKode(status);
338 }
339
340 // for version 4.0.0 this call must be made before solver settings
341 status = ARKStepSetUserData(arkode_mem, this);
342 AssertARKode(status);
343
344 setup_system_solver(solution);
345
346 setup_mass_solver(solution);
347
348 status = ARKStepSetInitStep(arkode_mem, current_time_step);
349 AssertARKode(status);
350
351 status = ARKStepSetStopTime(arkode_mem, data.final_time);
352 AssertARKode(status);
353
354 status = ARKStepSetOrder(arkode_mem, data.maximum_order);
355 AssertARKode(status);
356
357 if (custom_setup)
358 custom_setup(arkode_mem);
359 }
360
361
362
363 template <typename VectorType>
364 void
365 ARKode<VectorType>::setup_system_solver(const VectorType &solution)
366 {
367 // no point in setting up a solver when there is no linear system
368 if (!implicit_function)
369 return;
370
371 int status;
372 (void)status;
373 // Initialize solver
374 // currently only iterative linear solvers are supported
375 if (jacobian_times_vector)
376 {
377 SUNLinearSolver sun_linear_solver;
378 if (solve_linearized_system)
379 {
380 linear_solver =
381 std::make_unique<internal::LinearSolverWrapper<VectorType>>(
382 solve_linearized_system,
383 pending_exception
384# if DEAL_II_SUNDIALS_VERSION_GTE(6, 0, 0)
385 ,
386 arkode_ctx
387# endif
388 );
389 sun_linear_solver = *linear_solver;
390 }
391 else
392 {
393 // use default solver from SUNDIALS
394 // TODO give user options
395 auto y_template = internal::make_nvector_view(solution
397 ,
398 arkode_ctx
399# endif
400 );
401# if DEAL_II_SUNDIALS_VERSION_LT(6, 0, 0)
402 sun_linear_solver =
403 SUNLinSol_SPGMR(y_template,
404 PREC_NONE,
405 0 /*krylov subvectors, 0 uses default*/);
406# else
407 sun_linear_solver =
408 SUNLinSol_SPGMR(y_template,
409 PREC_NONE,
410 0 /*krylov subvectors, 0 uses default*/,
411 arkode_ctx);
412# endif
413 }
414 status = ARKStepSetLinearSolver(arkode_mem, sun_linear_solver, nullptr);
415 AssertARKode(status);
416
417 auto jacobian_times_vector_callback = [](N_Vector v,
418 N_Vector Jv,
419 realtype t,
420 N_Vector y,
421 N_Vector fy,
422 void * user_data,
423 N_Vector) -> int {
424 Assert(user_data != nullptr, ExcInternalError());
425 ARKode<VectorType> &solver =
426 *static_cast<ARKode<VectorType> *>(user_data);
427
428 auto *src_v = internal::unwrap_nvector_const<VectorType>(v);
429 auto *src_y = internal::unwrap_nvector_const<VectorType>(y);
430 auto *src_fy = internal::unwrap_nvector_const<VectorType>(fy);
431
432 auto *dst_Jv = internal::unwrap_nvector<VectorType>(Jv);
433
436 solver.pending_exception,
437 *src_v,
438 *dst_Jv,
439 t,
440 *src_y,
441 *src_fy);
442 };
443
444 auto jacobian_times_vector_setup_callback =
445 [](realtype t, N_Vector y, N_Vector fy, void *user_data) -> int {
446 Assert(user_data != nullptr, ExcInternalError());
447 ARKode<VectorType> &solver =
448 *static_cast<ARKode<VectorType> *>(user_data);
449
450 auto *src_y = internal::unwrap_nvector_const<VectorType>(y);
451 auto *src_fy = internal::unwrap_nvector_const<VectorType>(fy);
452
455 solver.pending_exception,
456 t,
457 *src_y,
458 *src_fy);
459 };
460 status = ARKStepSetJacTimes(arkode_mem,
461 jacobian_times_setup ?
462 jacobian_times_vector_setup_callback :
463 ARKLsJacTimesSetupFn(nullptr),
464 jacobian_times_vector_callback);
465 AssertARKode(status);
466 if (jacobian_preconditioner_solve)
467 {
468 auto solve_with_jacobian_callback = [](realtype t,
469 N_Vector y,
470 N_Vector fy,
471 N_Vector r,
472 N_Vector z,
473 realtype gamma,
474 realtype delta,
475 int lr,
476 void * user_data) -> int {
477 Assert(user_data != nullptr, ExcInternalError());
478 ARKode<VectorType> &solver =
479 *static_cast<ARKode<VectorType> *>(user_data);
480
481 auto *src_y = internal::unwrap_nvector_const<VectorType>(y);
482 auto *src_fy = internal::unwrap_nvector_const<VectorType>(fy);
483 auto *src_r = internal::unwrap_nvector_const<VectorType>(r);
484
485 auto *dst_z = internal::unwrap_nvector<VectorType>(z);
486
489 solver.pending_exception,
490 t,
491 *src_y,
492 *src_fy,
493 *src_r,
494 *dst_z,
495 gamma,
496 delta,
497 lr);
498 };
499
500 auto jacobian_solver_setup_callback = [](realtype t,
501 N_Vector y,
502 N_Vector fy,
503 booleantype jok,
504 booleantype *jcurPtr,
505 realtype gamma,
506 void *user_data) -> int {
507 Assert(user_data != nullptr, ExcInternalError());
508 ARKode<VectorType> &solver =
509 *static_cast<ARKode<VectorType> *>(user_data);
510
511 auto *src_y = internal::unwrap_nvector_const<VectorType>(y);
512 auto *src_fy = internal::unwrap_nvector_const<VectorType>(fy);
513
516 solver.pending_exception,
517 t,
518 *src_y,
519 *src_fy,
520 jok,
521 *jcurPtr,
522 gamma);
523 };
524
525 status = ARKStepSetPreconditioner(arkode_mem,
526 jacobian_preconditioner_setup ?
527 jacobian_solver_setup_callback :
528 ARKLsPrecSetupFn(nullptr),
529 solve_with_jacobian_callback);
530 AssertARKode(status);
531 }
532 if (data.implicit_function_is_linear)
533 {
534 status = ARKStepSetLinear(
535 arkode_mem, data.implicit_function_is_time_independent ? 0 : 1);
536 AssertARKode(status);
537 }
538 }
539 else
540 {
541 auto y_template = internal::make_nvector_view(solution
543 ,
544 arkode_ctx
545# endif
546 );
547
548# if DEAL_II_SUNDIALS_VERSION_LT(6, 0, 0)
549 SUNNonlinearSolver fixed_point_solver =
550 SUNNonlinSol_FixedPoint(y_template,
551 data.anderson_acceleration_subspace);
552# else
553 SUNNonlinearSolver fixed_point_solver =
554 SUNNonlinSol_FixedPoint(y_template,
555 data.anderson_acceleration_subspace,
556 arkode_ctx);
557# endif
558
559 status = ARKStepSetNonlinearSolver(arkode_mem, fixed_point_solver);
560 AssertARKode(status);
561 }
562
563 status =
564 ARKStepSetMaxNonlinIters(arkode_mem, data.maximum_non_linear_iterations);
565 AssertARKode(status);
566 }
567
568
569
570 template <typename VectorType>
571 void
572 ARKode<VectorType>::setup_mass_solver(const VectorType &solution)
573 {
574 int status;
575 (void)status;
576
577 if (mass_times_vector)
578 {
579 SUNLinearSolver sun_mass_linear_solver;
580 if (solve_mass)
581 {
582 mass_solver =
583 std::make_unique<internal::LinearSolverWrapper<VectorType>>(
584 solve_mass,
585
586 pending_exception
587# if DEAL_II_SUNDIALS_VERSION_GTE(6, 0, 0)
588 ,
589 arkode_ctx
590# endif
591 );
592 sun_mass_linear_solver = *mass_solver;
593 }
594 else
595 {
596 auto y_template = internal::make_nvector_view(solution
598 ,
599 arkode_ctx
600# endif
601 );
602# if DEAL_II_SUNDIALS_VERSION_LT(6, 0, 0)
603 sun_mass_linear_solver =
604 SUNLinSol_SPGMR(y_template,
605 PREC_NONE,
606 0 /*krylov subvectors, 0 uses default*/);
607# else
608 sun_mass_linear_solver =
609 SUNLinSol_SPGMR(y_template,
610 PREC_NONE,
611 0 /*krylov subvectors, 0 uses default*/,
612 arkode_ctx);
613# endif
614 }
615 booleantype mass_time_dependent =
616 data.mass_is_time_independent ? SUNFALSE : SUNTRUE;
617 status = ARKStepSetMassLinearSolver(arkode_mem,
618 sun_mass_linear_solver,
619 nullptr,
620 mass_time_dependent);
621 AssertARKode(status);
622
623 auto mass_matrix_times_vector_setup_callback =
624 [](realtype t, void *mtimes_data) -> int {
625 Assert(mtimes_data != nullptr, ExcInternalError());
626 ARKode<VectorType> &solver =
627 *static_cast<ARKode<VectorType> *>(mtimes_data);
628
630 solver.mass_times_setup, solver.pending_exception, t);
631 };
632
633 auto mass_matrix_times_vector_callback =
634 [](N_Vector v, N_Vector Mv, realtype t, void *mtimes_data) -> int {
635 Assert(mtimes_data != nullptr, ExcInternalError());
636 ARKode<VectorType> &solver =
637 *static_cast<ARKode<VectorType> *>(mtimes_data);
638
639 auto *src_v = internal::unwrap_nvector_const<VectorType>(v);
640 auto *dst_Mv = internal::unwrap_nvector<VectorType>(Mv);
641
643 solver.mass_times_vector,
644 solver.pending_exception,
645 t,
646 *src_v,
647 *dst_Mv);
648 };
649
650 status = ARKStepSetMassTimes(arkode_mem,
651 mass_times_setup ?
652 mass_matrix_times_vector_setup_callback :
653 ARKLsMassTimesSetupFn(nullptr),
654 mass_matrix_times_vector_callback,
655 this);
656 AssertARKode(status);
657
658 if (mass_preconditioner_solve)
659 {
660 auto mass_matrix_solver_setup_callback =
661 [](realtype t, void *user_data) -> int {
662 Assert(user_data != nullptr, ExcInternalError());
663 ARKode<VectorType> &solver =
664 *static_cast<ARKode<VectorType> *>(user_data);
665
668 };
669
670 auto solve_with_mass_matrix_callback = [](realtype t,
671 N_Vector r,
672 N_Vector z,
673 realtype delta,
674 int lr,
675 void *user_data) -> int {
676 Assert(user_data != nullptr, ExcInternalError());
677 ARKode<VectorType> &solver =
678 *static_cast<ARKode<VectorType> *>(user_data);
679
680 auto *src_r = internal::unwrap_nvector_const<VectorType>(r);
681 auto *dst_z = internal::unwrap_nvector<VectorType>(z);
682
685 solver.pending_exception,
686 t,
687 *src_r,
688 *dst_z,
689 delta,
690 lr);
691 };
692
693 status =
694 ARKStepSetMassPreconditioner(arkode_mem,
695 mass_preconditioner_setup ?
696 mass_matrix_solver_setup_callback :
697 ARKLsMassPrecSetupFn(nullptr),
698 solve_with_mass_matrix_callback);
699 AssertARKode(status);
700 }
701 }
702 }
703
704
705
706 template <typename VectorType>
707 void
709 {
710 solver_should_restart = [](const double, VectorType &) -> bool {
711 return false;
712 };
713 }
714
715
716
717 template <typename VectorType>
718 void *
720 {
721 return arkode_mem;
722 }
723
724 template class ARKode<Vector<double>>;
725 template class ARKode<BlockVector<double>>;
726
729
730# ifdef DEAL_II_WITH_MPI
731
732# ifdef DEAL_II_WITH_TRILINOS
735
736# endif // DEAL_II_WITH_TRILINOS
737
738# ifdef DEAL_II_WITH_PETSC
739# ifndef PETSC_USE_COMPLEX
742# endif // PETSC_USE_COMPLEX
743# endif // DEAL_II_WITH_PETSC
744
745# endif // DEAL_II_WITH_MPI
746
747} // namespace SUNDIALS
748
750
751#endif
#define AssertARKode(code)
Definition arkode.h:62
double get_next_step_size() const
void set_desired_next_step_size(const double time_step_size)
unsigned int get_step_number() const
bool is_at_end() const
double get_current_time() const
void set_next_step_size(const double time_step_size)
double get_next_time() const
double get_previous_step_size() const
void advance_time()
double get_end_time() const
void setup_mass_solver(const VectorType &solution)
Definition arkode.cc:572
std::function< void(const double t)> mass_preconditioner_setup
Definition arkode.h:944
std::function< void(const double t, const VectorType &y, VectorType &explicit_f)> explicit_function
Definition arkode.h:595
std::exception_ptr pending_exception
Definition arkode.h:1098
ARKode(const AdditionalData &data=AdditionalData())
Definition arkode.cc:54
std::function< void(const double t, const VectorType &y, const VectorType &fy, const int jok, int &jcur, const double gamma)> jacobian_preconditioner_setup
Definition arkode.h:884
void * get_arkode_memory() const
Definition arkode.cc:719
std::function< void(const VectorType &v, VectorType &Jv, const double t, const VectorType &y, const VectorType &fy)> jacobian_times_vector
Definition arkode.h:703
std::function< void(const double t, const VectorType &y, VectorType &res)> implicit_function
Definition arkode.h:614
int do_evolve_time(VectorType &solution, ::DiscreteTime &time, const bool do_reset)
Definition arkode.cc:140
MPI_Comm mpi_communicator
Definition arkode.h:1083
void reset(const double t, const double h, const VectorType &y)
Definition arkode.cc:235
SUNContext arkode_ctx
Definition arkode.h:1076
std::function< void(const double t, const VectorType &y, const VectorType &fy, const VectorType &r, VectorType &z, const double gamma, const double tol, const int lr)> jacobian_preconditioner_solve
Definition arkode.h:834
std::function< void(const double t)> mass_times_setup
Definition arkode.h:670
void set_functions_to_trigger_an_assert()
Definition arkode.cc:708
std::function< void(const double t, const VectorType &v, VectorType &Mv)> mass_times_vector
Definition arkode.h:634
std::function< void(const double t, const VectorType &r, VectorType &z, const double tol, const int lr)> mass_preconditioner_solve
Definition arkode.h:918
std::function< void(const double t, const VectorType &y, const VectorType &fy)> jacobian_times_setup
Definition arkode.h:741
unsigned int solve_ode(VectorType &solution)
Definition arkode.cc:108
void setup_system_solver(const VectorType &solution)
Definition arkode.cc:365
unsigned int solve_ode_incrementally(VectorType &solution, const double intermediate_time, const bool reset_solver=false)
Definition arkode.cc:121
#define DEAL_II_NAMESPACE_OPEN
Definition config.h:472
#define DEAL_II_SUNDIALS_VERSION_GTE(major, minor, patch)
Definition config.h:359
#define DEAL_II_NAMESPACE_CLOSE
Definition config.h:473
static ::ExceptionBase & ExcFunctionNotProvided(std::string arg1)
#define Assert(cond, exc)
static ::ExceptionBase & ExcInternalError()
static ::ExceptionBase & RecoverableUserCallbackError()
static ::ExceptionBase & ExcMessage(std::string arg1)
#define AssertThrow(cond, exc)
int call_and_possibly_capture_exception(const F &f, std::exception_ptr &eptr, Args &&...args)
Definition utilities.h:46