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.h
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#ifndef dealii_sundials_arkode_h
18#define dealii_sundials_arkode_h
19
20#include <deal.II/base/config.h>
21
22
23#ifdef DEAL_II_WITH_SUNDIALS
24
30
31# ifdef DEAL_II_WITH_PETSC
34# endif
35# include <deal.II/lac/vector.h>
37
38# include <arkode/arkode.h>
39# include <nvector/nvector_serial.h>
40# ifdef DEAL_II_WITH_MPI
41# include <nvector/nvector_parallel.h>
42# endif
44
47
48# include <boost/signals2.hpp>
49
50# include <sundials/sundials_linearsolver.h>
51# include <sundials/sundials_math.h>
52# include <sundials/sundials_types.h>
53
54# include <exception>
55# include <memory>
56
57
59
60
61// Shorthand notation for ARKODE error codes.
62# define AssertARKode(code) Assert(code >= 0, ExcARKodeError(code))
63
67namespace SUNDIALS
68{
329 template <typename VectorType = Vector<double>>
330 class ARKode
331 {
332 public:
337 {
338 public:
370 // Initial parameters
371 const double initial_time = 0.0,
372 const double final_time = 1.0,
373 const double initial_step_size = 1e-2,
374 const double output_period = 1e-1,
375 // Running parameters
376 const double minimum_step_size = 1e-6,
377 const unsigned int maximum_order = 5,
378 const unsigned int maximum_non_linear_iterations = 10,
379 const bool implicit_function_is_linear = false,
380 const bool implicit_function_is_time_independent = false,
381 const bool mass_is_time_independent = false,
382 const int anderson_acceleration_subspace = 3,
383 // Error parameters
384 const double absolute_tolerance = 1e-6,
385 const double relative_tolerance = 1e-5);
386
401 void
403
408
413
418
423
428
433
437 unsigned int maximum_order;
438
444
450
455
461
467
473 };
474
485
493 ARKode(const AdditionalData &data, const MPI_Comm mpi_comm);
494
498 ~ARKode();
499
507 unsigned int
508 solve_ode(VectorType &solution);
509
536 unsigned int
537 solve_ode_incrementally(VectorType & solution,
538 const double intermediate_time,
539 const bool reset_solver = false);
540
555 void
556 reset(const double t, const double h, const VectorType &y);
557
574 void *
575 get_arkode_memory() const;
576
593 std::function<
594 void(const double t, const VectorType &y, VectorType &explicit_f)>
596
613 std::function<void(const double t, const VectorType &y, VectorType &res)>
615
633 std::function<void(const double t, const VectorType &v, VectorType &Mv)>
635
670 std::function<void(const double t)> mass_times_setup;
671
698 std::function<void(const VectorType &v,
699 VectorType & Jv,
700 const double t,
701 const VectorType &y,
702 const VectorType &fy)>
704
739 std::function<
740 void(const double t, const VectorType &y, const VectorType &fy)>
742
770
794
826 std::function<void(const double t,
827 const VectorType &y,
828 const VectorType &fy,
829 const VectorType &r,
830 VectorType & z,
831 const double gamma,
832 const double tol,
833 const int lr)>
835
878 std::function<void(const double t,
879 const VectorType &y,
880 const VectorType &fy,
881 const int jok,
882 int & jcur,
883 const double gamma)>
885
913 std::function<void(const double t,
914 const VectorType &r,
915 VectorType & z,
916 const double tol,
917 const int lr)>
919
944 std::function<void(const double t)> mass_preconditioner_setup;
945
960 std::function<void(const double t,
961 const VectorType & sol,
962 const unsigned int step_number)>
964
982 std::function<bool(const double t, VectorType &sol)> solver_should_restart;
983
990 std::function<VectorType &()> get_local_tolerances;
991
1016 std::function<void(void *arkode_mem)> custom_setup;
1017
1018 private:
1024 std::string,
1025 << "Please provide an implementation for the function \""
1026 << arg1 << "\"");
1027
1031 int
1032 do_evolve_time(VectorType & solution,
1033 ::DiscreteTime &time,
1034 const bool do_reset);
1035
1042 void
1043 setup_system_solver(const VectorType &solution);
1044
1051 void
1052 setup_mass_solver(const VectorType &solution);
1053
1059 void
1061
1066
1071
1072# if DEAL_II_SUNDIALS_VERSION_GTE(6, 0, 0)
1077# endif
1078
1084
1089
1090 std::unique_ptr<internal::LinearSolverWrapper<VectorType>> linear_solver;
1091 std::unique_ptr<internal::LinearSolverWrapper<VectorType>> mass_solver;
1092
1098 mutable std::exception_ptr pending_exception;
1099
1100# ifdef DEAL_II_WITH_PETSC
1101# ifdef PETSC_USE_COMPLEX
1102 static_assert(!std::is_same<VectorType, PETScWrappers::MPI::Vector>::value,
1103 "Sundials does not support complex scalar types, "
1104 "but PETSc is configured to use a complex scalar type!");
1105
1106 static_assert(
1107 !std::is_same<VectorType, PETScWrappers::MPI::BlockVector>::value,
1108 "Sundials does not support complex scalar types, "
1109 "but PETSc is configured to use a complex scalar type!");
1110# endif // PETSC_USE_COMPLEX
1111# endif // DEAL_II_WITH_PETSC
1112 };
1113
1114
1119 int,
1120 << "One of the SUNDIALS ARKode internal functions "
1121 << " returned a negative error code: " << arg1
1122 << ". Please consult SUNDIALS manual.");
1123
1124
1125 template <typename VectorType>
1127 // Initial parameters
1128 const double initial_time,
1129 const double final_time,
1130 const double initial_step_size,
1131 const double output_period,
1132 // Running parameters
1133 const double minimum_step_size,
1134 const unsigned int maximum_order,
1135 const unsigned int maximum_non_linear_iterations,
1136 const bool implicit_function_is_linear,
1137 const bool implicit_function_is_time_independent,
1138 const bool mass_is_time_independent,
1139 const int anderson_acceleration_subspace,
1140 // Error parameters
1141 const double absolute_tolerance,
1142 const double relative_tolerance)
1143 : initial_time(initial_time)
1144 , final_time(final_time)
1145 , initial_step_size(initial_step_size)
1146 , minimum_step_size(minimum_step_size)
1147 , absolute_tolerance(absolute_tolerance)
1148 , relative_tolerance(relative_tolerance)
1149 , maximum_order(maximum_order)
1150 , output_period(output_period)
1151 , maximum_non_linear_iterations(maximum_non_linear_iterations)
1152 , implicit_function_is_linear(implicit_function_is_linear)
1153 , implicit_function_is_time_independent(
1154 implicit_function_is_time_independent)
1155 , mass_is_time_independent(mass_is_time_independent)
1156 , anderson_acceleration_subspace(anderson_acceleration_subspace)
1157 {}
1158
1159
1160
1161 template <typename VectorType>
1162 void
1164 {
1165 prm.add_parameter("Initial time", initial_time);
1166 prm.add_parameter("Final time", final_time);
1167 prm.add_parameter("Time interval between each output", output_period);
1168 prm.enter_subsection("Running parameters");
1169 prm.add_parameter("Initial step size", initial_step_size);
1170 prm.add_parameter("Minimum step size", minimum_step_size);
1171 prm.add_parameter("Maximum order of ARK", maximum_order);
1172 prm.add_parameter("Maximum number of nonlinear iterations",
1173 maximum_non_linear_iterations);
1174 prm.add_parameter("Implicit function is linear",
1175 implicit_function_is_linear);
1176 prm.add_parameter("Implicit function is time independent",
1177 implicit_function_is_time_independent);
1178 prm.add_parameter("Mass is time independent", mass_is_time_independent);
1179 prm.add_parameter("Anderson-acceleration subspace",
1180 anderson_acceleration_subspace);
1181 prm.leave_subsection();
1182 prm.enter_subsection("Error control");
1183 prm.add_parameter("Absolute error tolerance", absolute_tolerance);
1184 prm.add_parameter("Relative error tolerance", relative_tolerance);
1185 prm.leave_subsection();
1186 }
1187
1188} // namespace SUNDIALS
1189
1191#endif
1192
1193
1194#endif
void add_parameter(const std::string &entry, ParameterType &parameter, const std::string &documentation="", const Patterns::PatternBase &pattern= *Patterns::Tools::Convert< ParameterType >::to_pattern(), const bool has_to_be_set=false)
void enter_subsection(const std::string &subsection)
void add_parameters(ParameterHandler &prm)
Definition arkode.h:1163
AdditionalData(const double initial_time=0.0, const double final_time=1.0, const double initial_step_size=1e-2, const double output_period=1e-1, const double minimum_step_size=1e-6, const unsigned int maximum_order=5, const unsigned int maximum_non_linear_iterations=10, const bool implicit_function_is_linear=false, const bool implicit_function_is_time_independent=false, const bool mass_is_time_independent=false, const int anderson_acceleration_subspace=3, const double absolute_tolerance=1e-6, const double relative_tolerance=1e-5)
Definition arkode.h:1126
unsigned int maximum_non_linear_iterations
Definition arkode.h:449
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
void * arkode_mem
Definition arkode.h:1070
std::exception_ptr pending_exception
Definition arkode.h:1098
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< VectorType &()> get_local_tolerances
Definition arkode.h:990
AdditionalData data
Definition arkode.h:1065
std::function< void(const VectorType &v, VectorType &Jv, const double t, const VectorType &y, const VectorType &fy)> jacobian_times_vector
Definition arkode.h:703
LinearSolveFunction< VectorType > solve_linearized_system
Definition arkode.h:769
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::unique_ptr< internal::LinearSolverWrapper< VectorType > > mass_solver
Definition arkode.h:1091
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< bool(const double t, VectorType &sol)> solver_should_restart
Definition arkode.h:982
std::function< void(const double t)> mass_times_setup
Definition arkode.h:670
std::function< void(const double t, const VectorType &sol, const unsigned int step_number)> output_step
Definition arkode.h:963
void set_functions_to_trigger_an_assert()
Definition arkode.cc:708
std::unique_ptr< internal::LinearSolverWrapper< VectorType > > linear_solver
Definition arkode.h:1090
std::function< void(const double t, const VectorType &v, VectorType &Mv)> mass_times_vector
Definition arkode.h:634
LinearSolveFunction< VectorType > solve_mass
Definition arkode.h:793
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
double last_end_time
Definition arkode.h:1088
std::function< void(void *arkode_mem)> custom_setup
Definition arkode.h:1016
#define DEAL_II_NAMESPACE_OPEN
Definition config.h:472
#define DEAL_II_NAMESPACE_CLOSE
Definition config.h:473
static ::ExceptionBase & ExcARKodeError(int arg1)
static ::ExceptionBase & ExcFunctionNotProvided(std::string arg1)
#define DeclException1(Exception1, type1, outsequence)
Definition exceptions.h:510
std::function< void(SundialsOperator< VectorType > &op, SundialsPreconditioner< VectorType > &prec, VectorType &x, const VectorType &b, double tol)> LinearSolveFunction