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