Reference documentation for deal.II version 9.2.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\}}\)
arkode.h
Go to the documentation of this file.
1 //-----------------------------------------------------------
2 //
3 // Copyright (C) 2017 - 2020 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 #ifndef dealii_sundials_arkode_h
17 #define dealii_sundials_arkode_h
18 
19 #include <deal.II/base/config.h>
20 
21 #include <deal.II/base/mpi.h>
22 
23 #ifdef DEAL_II_WITH_SUNDIALS
24 
26 # include <deal.II/base/exceptions.h>
27 # include <deal.II/base/logstream.h>
29 # ifdef DEAL_II_WITH_PETSC
32 # endif
33 # include <deal.II/lac/vector.h>
35 
36 # include <arkode/arkode.h>
37 # include <arkode/arkode_impl.h>
38 # include <nvector/nvector_serial.h>
39 # ifdef DEAL_II_WITH_MPI
40 # include <nvector/nvector_parallel.h>
41 # endif
42 # include <boost/signals2.hpp>
43 
44 # include <sundials/sundials_math.h>
45 # include <sundials/sundials_types.h>
46 
47 # include <memory>
48 
49 
51 
52 // Shorthand notation for ARKODE error codes.
53 # define AssertARKode(code) Assert(code >= 0, ExcARKodeError(code))
54 
58 namespace SUNDIALS
59 {
312  template <typename VectorType = Vector<double>>
313  class ARKode
314  {
315  public:
320  {
321  public:
349  // Initial parameters
350  const double initial_time = 0.0,
351  const double final_time = 1.0,
352  const double initial_step_size = 1e-2,
353  const double output_period = 1e-1,
354  // Running parameters
355  const double minimum_step_size = 1e-6,
356  const unsigned int maximum_order = 5,
357  const unsigned int maximum_non_linear_iterations = 10,
358  const bool implicit_function_is_linear = false,
359  const bool implicit_function_is_time_independent = false,
360  // Error parameters
361  const double absolute_tolerance = 1e-6,
362  const double relative_tolerance = 1e-5)
375  {}
376 
419  void
421  {
422  prm.add_parameter("Initial time", initial_time);
423  prm.add_parameter("Final time", final_time);
424  prm.add_parameter("Time interval between each output", output_period);
425 
426  prm.enter_subsection("Running parameters");
427  prm.add_parameter("Initial step size", initial_step_size);
428  prm.add_parameter("Minimum step size", minimum_step_size);
429  prm.add_parameter("Maximum order of ARK", maximum_order);
430  prm.add_parameter("Maximum number of nonlinear iterations",
432  prm.add_parameter("Implicit function is linear",
434  prm.add_parameter("Implicit function is time independent",
436  prm.leave_subsection();
437 
438  prm.enter_subsection("Error control");
439  prm.add_parameter("Absolute error tolerance", absolute_tolerance);
440  prm.add_parameter("Relative error tolerance", relative_tolerance);
441  prm.leave_subsection();
442  }
443 
447  double initial_time;
448 
452  double final_time;
453 
458 
463 
468 
473 
477  unsigned int maximum_order;
478 
483 
489 
494 
500  };
501 
514  const MPI_Comm mpi_comm = MPI_COMM_WORLD);
515 
519  ~ARKode();
520 
525  unsigned int
526  solve_ode(VectorType &solution);
527 
542  void
543  reset(const double t, const double h, const VectorType &y);
544 
549  std::function<void(VectorType &)> reinit_vector;
550 
567  std::function<
568  int(const double t, const VectorType &y, VectorType &explicit_f)>
570 
587  std::function<int(const double t, const VectorType &y, VectorType &res)>
589 
674  std::function<int(const int convfail,
675  const double t,
676  const double gamma,
677  const VectorType &ypred,
678  const VectorType &fpred,
679  bool & j_is_current)>
681 
726  std::function<int(const double t,
727  const double gamma,
728  const VectorType &ycur,
729  const VectorType &fcur,
730  const VectorType &rhs,
731  VectorType & dst)>
733 
734 
768  std::function<int(const double t)> setup_mass;
769 
788  std::function<int(const VectorType &rhs, VectorType &dst)>
790 
805  std::function<void(const double t,
806  const VectorType & sol,
807  const unsigned int step_number)>
809 
827  std::function<bool(const double t, VectorType &sol)> solver_should_restart;
828 
835  std::function<VectorType &()> get_local_tolerances;
836 
840  DeclException1(ExcARKodeError,
841  int,
842  << "One of the SUNDIALS ARKode internal functions "
843  << " returned a negative error code: " << arg1
844  << ". Please consult SUNDIALS manual.");
845 
846 
847  private:
852  DeclException1(ExcFunctionNotProvided,
853  std::string,
854  << "Please provide an implementation for the function \""
855  << arg1 << "\"");
856 
862  void
864 
869 
873  void *arkode_mem;
874 
878  N_Vector yy;
879 
883  N_Vector abs_tolls;
884 
891 
896 
897 # ifdef DEAL_II_WITH_PETSC
898 # ifdef PETSC_USE_COMPLEX
900  "Sundials does not support complex scalar types, "
901  "but PETSc is configured to use a complex scalar type!");
902 
903  static_assert(
905  "Sundials does not support complex scalar types, "
906  "but PETSc is configured to use a complex scalar type!");
907 # endif // PETSC_USE_COMPLEX
908 # endif // DEAL_II_WITH_PETSC
909  };
910 
911 } // namespace SUNDIALS
912 
914 #endif
915 
916 
917 #endif
internal::QGaussLobatto::gamma
long double gamma(const unsigned int n)
Definition: quadrature_lib.cc:96
SUNDIALS::ARKode::AdditionalData::initial_time
double initial_time
Definition: arkode.h:447
SUNDIALS::ARKode
Definition: arkode.h:313
SUNDIALS::ARKode::AdditionalData::maximum_non_linear_iterations
unsigned int maximum_non_linear_iterations
Definition: arkode.h:488
SUNDIALS::ARKode::abs_tolls
N_Vector abs_tolls
Definition: arkode.h:883
SUNDIALS::ARKode::reset
void reset(const double t, const double h, const VectorType &y)
Definition: arkode.cc:353
SUNDIALS::ARKode::AdditionalData::absolute_tolerance
double absolute_tolerance
Definition: arkode.h:467
VectorType
MPI_Comm
Physics::Elasticity::Kinematics::e
SymmetricTensor< 2, dim, Number > e(const Tensor< 2, dim, Number > &F)
SUNDIALS::ARKode::AdditionalData::maximum_order
unsigned int maximum_order
Definition: arkode.h:477
SUNDIALS::ARKode::setup_jacobian
std::function< int(const int convfail, const double t, const double gamma, const VectorType &ypred, const VectorType &fpred, bool &j_is_current)> setup_jacobian
Definition: arkode.h:680
SUNDIALS
Definition: arkode.h:58
SUNDIALS::ARKode::~ARKode
~ARKode()
Definition: arkode.cc:253
GrowingVectorMemory
Definition: vector_memory.h:320
bool
conditional_ostream.h
SUNDIALS::ARKode::AdditionalData::AdditionalData
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 double absolute_tolerance=1e-6, const double relative_tolerance=1e-5)
Definition: arkode.h:348
SUNDIALS::ARKode::solve_ode
unsigned int solve_ode(VectorType &solution)
Definition: arkode.cc:271
SUNDIALS::ARKode::AdditionalData
Definition: arkode.h:319
SUNDIALS::ARKode::AdditionalData::relative_tolerance
double relative_tolerance
Definition: arkode.h:472
mpi.h
SUNDIALS::ARKode::AdditionalData::output_period
double output_period
Definition: arkode.h:482
SUNDIALS::ARKode::output_step
std::function< void(const double t, const VectorType &sol, const unsigned int step_number)> output_step
Definition: arkode.h:808
DEAL_II_NAMESPACE_OPEN
#define DEAL_II_NAMESPACE_OPEN
Definition: config.h:358
SUNDIALS::ARKode::DeclException1
DeclException1(ExcARKodeError, int,<< "One of the SUNDIALS ARKode internal functions "<< " returned a negative error code: "<< arg1<< ". Please consult SUNDIALS manual.")
SUNDIALS::ARKode::get_local_tolerances
std::function< VectorType &()> get_local_tolerances
Definition: arkode.h:835
petsc_block_vector.h
SUNDIALS::ARKode::AdditionalData::add_parameters
void add_parameters(ParameterHandler &prm)
Definition: arkode.h:420
SUNDIALS::ARKode::AdditionalData::initial_step_size
double initial_step_size
Definition: arkode.h:457
SUNDIALS::ARKode::AdditionalData::implicit_function_is_time_independent
bool implicit_function_is_time_independent
Definition: arkode.h:499
parameter_handler.h
SUNDIALS::ARKode::AdditionalData::minimum_step_size
double minimum_step_size
Definition: arkode.h:462
SUNDIALS::ARKode::solve_jacobian_system
std::function< int(const double t, const double gamma, const VectorType &ycur, const VectorType &fcur, const VectorType &rhs, VectorType &dst)> solve_jacobian_system
Definition: arkode.h:732
SUNDIALS::ARKode::explicit_function
std::function< int(const double t, const VectorType &y, VectorType &explicit_f)> explicit_function
Definition: arkode.h:569
exceptions.h
SUNDIALS::ARKode::arkode_mem
void * arkode_mem
Definition: arkode.h:873
SUNDIALS::ARKode::reinit_vector
std::function< void(VectorType &)> reinit_vector
Definition: arkode.h:549
value
static const bool value
Definition: dof_tools_constraints.cc:433
SUNDIALS::ARKode::mem
GrowingVectorMemory< VectorType > mem
Definition: arkode.h:895
SUNDIALS::ARKode::AdditionalData::implicit_function_is_linear
bool implicit_function_is_linear
Definition: arkode.h:493
SUNDIALS::ARKode::ARKode
ARKode(const AdditionalData &data=AdditionalData(), const MPI_Comm mpi_comm=MPI_COMM_WORLD)
Definition: arkode.cc:239
vector.h
SUNDIALS::ARKode::set_functions_to_trigger_an_assert
void set_functions_to_trigger_an_assert()
Definition: arkode.cc:494
SUNDIALS::ARKode::solve_mass_system
std::function< int(const VectorType &rhs, VectorType &dst)> solve_mass_system
Definition: arkode.h:789
SUNDIALS::ARKode::communicator
MPI_Comm communicator
Definition: arkode.h:890
SUNDIALS::ARKode::yy
N_Vector yy
Definition: arkode.h:878
ParameterHandler::add_parameter
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)
Definition: parameter_handler.h:2328
petsc_vector.h
vector_memory.h
ParameterHandler
Definition: parameter_handler.h:845
config.h
ParameterHandler::enter_subsection
void enter_subsection(const std::string &subsection)
Definition: parameter_handler.cc:927
SUNDIALS::ARKode::setup_mass
std::function< int(const double t)> setup_mass
Definition: arkode.h:768
SUNDIALS::ARKode::AdditionalData::final_time
double final_time
Definition: arkode.h:452
SUNDIALS::ARKode::solver_should_restart
std::function< bool(const double t, VectorType &sol)> solver_should_restart
Definition: arkode.h:827
DEAL_II_NAMESPACE_CLOSE
#define DEAL_II_NAMESPACE_CLOSE
Definition: config.h:359
ParameterHandler::leave_subsection
void leave_subsection()
Definition: parameter_handler.cc:941
logstream.h
SUNDIALS::ARKode::implicit_function
std::function< int(const double t, const VectorType &y, VectorType &res)> implicit_function
Definition: arkode.h:588
SUNDIALS::ARKode::data
AdditionalData data
Definition: arkode.h:868
int