Reference documentation for deal.II version GIT relicensing-446-ge11b936273 2024-04-20 12:30:02+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.h
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#ifndef dealii_sundials_arkode_h
17#define dealii_sundials_arkode_h
18
19#include <deal.II/base/config.h>
20
21
22#ifdef DEAL_II_WITH_SUNDIALS
23
29
30# ifdef DEAL_II_WITH_PETSC
33# endif
34# include <deal.II/lac/vector.h>
36
37# include <arkode/arkode.h>
38# include <nvector/nvector_serial.h>
39# ifdef DEAL_II_WITH_MPI
40# include <nvector/nvector_parallel.h>
41# endif
43
47
48# include <boost/signals2.hpp>
49
50# include <sundials/sundials_linearsolver.h>
51# include <sundials/sundials_math.h>
52
53# include <exception>
54# include <memory>
55
56
58
59
60// Shorthand notation for ARKODE error codes.
61# define AssertARKode(code) Assert(code >= 0, ExcARKodeError(code))
62
66namespace SUNDIALS
67{
328 template <typename VectorType = Vector<double>>
329 class ARKode
330 {
331 public:
336 {
337 public:
369 // Initial parameters
370 const double initial_time = 0.0,
371 const double final_time = 1.0,
372 const double initial_step_size = 1e-2,
373 const double output_period = 1e-1,
374 // Running parameters
375 const double minimum_step_size = 1e-6,
376 const unsigned int maximum_order = 5,
377 const unsigned int maximum_non_linear_iterations = 10,
378 const bool implicit_function_is_linear = false,
379 const bool implicit_function_is_time_independent = false,
380 const bool mass_is_time_independent = false,
381 const int anderson_acceleration_subspace = 3,
382 // Error parameters
383 const double absolute_tolerance = 1e-6,
384 const double relative_tolerance = 1e-5);
385
400 void
402
407
412
417
422
427
432
436 unsigned int maximum_order;
437
443
449
454
460
466
472 };
473
484
492 ARKode(const AdditionalData &data, const MPI_Comm mpi_comm);
493
497 ~ARKode();
498
506 unsigned int
507 solve_ode(VectorType &solution);
508
535 unsigned int
536 solve_ode_incrementally(VectorType &solution,
537 const double intermediate_time,
538 const bool reset_solver = false);
539
554 void
555 reset(const double t, const double h, const VectorType &y);
556
573 void *
574 get_arkode_memory() const;
575
592 std::function<
593 void(const double t, const VectorType &y, VectorType &explicit_f)>
595
612 std::function<void(const double t, const VectorType &y, VectorType &res)>
614
632 std::function<void(const double t, const VectorType &v, VectorType &Mv)>
634
669 std::function<void(const double t)> mass_times_setup;
670
697 std::function<void(const VectorType &v,
698 VectorType &Jv,
699 const double t,
700 const VectorType &y,
701 const VectorType &fy)>
703
738 std::function<
739 void(const double t, const VectorType &y, const VectorType &fy)>
741
769
793
825 std::function<void(const double t,
826 const VectorType &y,
827 const VectorType &fy,
828 const VectorType &r,
829 VectorType &z,
830 const double gamma,
831 const double tol,
832 const int lr)>
834
877 std::function<void(const double t,
878 const VectorType &y,
879 const VectorType &fy,
880 const int jok,
881 int &jcur,
882 const double gamma)>
884
912 std::function<void(const double t,
913 const VectorType &r,
914 VectorType &z,
915 const double tol,
916 const int lr)>
918
943 std::function<void(const double t)> mass_preconditioner_setup;
944
959 std::function<void(const double t,
960 const VectorType &sol,
961 const unsigned int step_number)>
963
981 std::function<bool(const double t, VectorType &sol)> solver_should_restart;
982
989 std::function<VectorType &()> get_local_tolerances;
990
1015 std::function<void(void *arkode_mem)> custom_setup;
1016
1017 private:
1023 std::string,
1024 << "Please provide an implementation for the function \""
1025 << arg1 << "\"");
1026
1030 int
1031 do_evolve_time(VectorType &solution,
1032 ::DiscreteTime &time,
1033 const bool do_reset);
1034
1041 void
1042 setup_system_solver(const VectorType &solution);
1043
1050 void
1051 setup_mass_solver(const VectorType &solution);
1052
1058 void
1060
1065
1070
1071# if DEAL_II_SUNDIALS_VERSION_GTE(6, 0, 0)
1076# endif
1077
1083
1088
1089 std::unique_ptr<internal::LinearSolverWrapper<VectorType>> linear_solver;
1090 std::unique_ptr<internal::LinearSolverWrapper<VectorType>> mass_solver;
1091
1097 mutable std::exception_ptr pending_exception;
1098
1099# ifdef DEAL_II_WITH_PETSC
1100# ifdef PETSC_USE_COMPLEX
1101 static_assert(!std::is_same_v<VectorType, PETScWrappers::MPI::Vector>,
1102 "Sundials does not support complex scalar types, "
1103 "but PETSc is configured to use a complex scalar type!");
1104
1105 static_assert(!std::is_same_v<VectorType, PETScWrappers::MPI::BlockVector>,
1106 "Sundials does not support complex scalar types, "
1107 "but PETSc is configured to use a complex scalar type!");
1108# endif // PETSC_USE_COMPLEX
1109# endif // DEAL_II_WITH_PETSC
1110 };
1111
1112
1117 int,
1118 << "One of the SUNDIALS ARKode internal functions "
1119 << " returned a negative error code: " << arg1
1120 << ". Please consult SUNDIALS manual.");
1121
1122
1123 template <typename VectorType>
1125 // Initial parameters
1126 const double initial_time,
1127 const double final_time,
1128 const double initial_step_size,
1129 const double output_period,
1130 // Running parameters
1131 const double minimum_step_size,
1132 const unsigned int maximum_order,
1133 const unsigned int maximum_non_linear_iterations,
1134 const bool implicit_function_is_linear,
1135 const bool implicit_function_is_time_independent,
1136 const bool mass_is_time_independent,
1137 const int anderson_acceleration_subspace,
1138 // Error parameters
1139 const double absolute_tolerance,
1140 const double relative_tolerance)
1141 : initial_time(initial_time)
1142 , final_time(final_time)
1143 , initial_step_size(initial_step_size)
1144 , minimum_step_size(minimum_step_size)
1145 , absolute_tolerance(absolute_tolerance)
1146 , relative_tolerance(relative_tolerance)
1147 , maximum_order(maximum_order)
1148 , output_period(output_period)
1149 , maximum_non_linear_iterations(maximum_non_linear_iterations)
1150 , implicit_function_is_linear(implicit_function_is_linear)
1151 , implicit_function_is_time_independent(
1152 implicit_function_is_time_independent)
1153 , mass_is_time_independent(mass_is_time_independent)
1154 , anderson_acceleration_subspace(anderson_acceleration_subspace)
1155 {}
1156
1157
1158
1159 template <typename VectorType>
1160 void
1162 {
1163 prm.add_parameter("Initial time", initial_time);
1164 prm.add_parameter("Final time", final_time);
1165 prm.add_parameter("Time interval between each output", output_period);
1166 prm.enter_subsection("Running parameters");
1167 prm.add_parameter("Initial step size", initial_step_size);
1168 prm.add_parameter("Minimum step size", minimum_step_size);
1169 prm.add_parameter("Maximum order of ARK", maximum_order);
1170 prm.add_parameter("Maximum number of nonlinear iterations",
1171 maximum_non_linear_iterations);
1172 prm.add_parameter("Implicit function is linear",
1173 implicit_function_is_linear);
1174 prm.add_parameter("Implicit function is time independent",
1175 implicit_function_is_time_independent);
1176 prm.add_parameter("Mass is time independent", mass_is_time_independent);
1177 prm.add_parameter("Anderson-acceleration subspace",
1178 anderson_acceleration_subspace);
1179 prm.leave_subsection();
1180 prm.enter_subsection("Error control");
1181 prm.add_parameter("Absolute error tolerance", absolute_tolerance);
1182 prm.add_parameter("Relative error tolerance", relative_tolerance);
1183 prm.leave_subsection();
1184 }
1185
1186} // namespace SUNDIALS
1187
1189#endif
1190
1191
1192#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, const bool create_path_if_needed=true)
void add_parameters(ParameterHandler &prm)
Definition arkode.h:1161
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:1124
unsigned int maximum_non_linear_iterations
Definition arkode.h:448
void setup_mass_solver(const VectorType &solution)
Definition arkode.cc:596
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
void * arkode_mem
Definition arkode.h:1069
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
void * get_arkode_memory() const
Definition arkode.cc:747
std::function< void(const double t, const VectorType &sol, const unsigned int step_number)> output_step
Definition arkode.h:962
std::function< VectorType &()> get_local_tolerances
Definition arkode.h:989
AdditionalData data
Definition arkode.h:1064
LinearSolveFunction< VectorType > solve_linearized_system
Definition arkode.h:768
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:146
MPI_Comm mpi_communicator
Definition arkode.h:1082
void reset(const double t, const double h, const VectorType &y)
Definition arkode.cc:241
SUNContext arkode_ctx
Definition arkode.h:1075
std::unique_ptr< internal::LinearSolverWrapper< VectorType > > mass_solver
Definition arkode.h:1090
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< bool(const double t, VectorType &sol)> solver_should_restart
Definition arkode.h:981
std::function< void(const double t)> mass_times_setup
Definition arkode.h:669
void set_functions_to_trigger_an_assert()
Definition arkode.cc:736
std::unique_ptr< internal::LinearSolverWrapper< VectorType > > linear_solver
Definition arkode.h:1089
std::function< void(const double t, const VectorType &v, VectorType &Mv)> mass_times_vector
Definition arkode.h:633
LinearSolveFunction< VectorType > solve_mass
Definition arkode.h:792
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:114
void setup_system_solver(const VectorType &solution)
Definition arkode.cc:386
unsigned int solve_ode_incrementally(VectorType &solution, const double intermediate_time, const bool reset_solver=false)
Definition arkode.cc:127
double last_end_time
Definition arkode.h:1087
std::function< void(void *arkode_mem)> custom_setup
Definition arkode.h:1015
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_NAMESPACE_CLOSE
Definition config.h:503
static ::ExceptionBase & ExcARKodeError(int arg1)
static ::ExceptionBase & ExcFunctionNotProvided(std::string arg1)
#define DeclException1(Exception1, type1, outsequence)
Definition exceptions.h:516
std::function< void(SundialsOperator< VectorType > &op, SundialsPreconditioner< VectorType > &prec, VectorType &x, const VectorType &b, double tol)> LinearSolveFunction