|
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\}}\)
Go to the documentation of this file.
16 #ifndef dealii_sundials_ida_h
17 #define dealii_sundials_ida_h
22 #ifdef DEAL_II_WITH_SUNDIALS
28 # ifdef DEAL_II_WITH_PETSC
35 # ifdef DEAL_II_SUNDIALS_WITH_IDAS
36 # include <idas/idas.h>
41 # include <sundials/sundials_config.h>
42 # if DEAL_II_SUNDIALS_VERSION_LT(3, 0, 0)
43 # include <ida/ida_spbcgs.h>
44 # include <ida/ida_spgmr.h>
45 # include <ida/ida_sptfqmr.h>
47 # include <boost/signals2.hpp>
49 # include <nvector/nvector_serial.h>
50 # include <sundials/sundials_math.h>
51 # include <sundials/sundials_types.h>
59 # define AssertIDA(code) Assert(code >= 0, ExcIDAError(code))
234 template <
typename VectorType = Vector<
double>>
398 "Ignore algebraic terms for error computations",
400 "Indicate whether or not to suppress algebraic variables "
401 "in the local error test.");
405 static std::string ic_type_str =
"use_y_diff";
407 "Correction type at initial time",
409 "This is one of the following three options for the "
410 "initial condition calculation. \n"
411 " none: do not try to make initial conditions consistent. \n"
412 " use_y_diff: compute the algebraic components of y and differential\n"
413 " components of y_dot, given the differential components of y. \n"
414 " This option requires that the user specifies differential and \n"
415 " algebraic components in the function get_differential_components.\n"
416 " use_y_dot: compute all components of y, given y_dot.",
418 prm.
add_action(
"Correction type at initial time",
419 [&](
const std::string &
value) {
420 if (
value ==
"use_y_diff")
422 else if (
value ==
"use_y_dot")
424 else if (
value ==
"none")
430 static std::string reset_type_str =
"use_y_diff";
432 "Correction type after restart",
434 "This is one of the following three options for the "
435 "initial condition calculation. \n"
436 " none: do not try to make initial conditions consistent. \n"
437 " use_y_diff: compute the algebraic components of y and differential\n"
438 " components of y_dot, given the differential components of y. \n"
439 " This option requires that the user specifies differential and \n"
440 " algebraic components in the function get_differential_components.\n"
441 " use_y_dot: compute all components of y, given y_dot.",
443 prm.
add_action(
"Correction type after restart",
444 [&](
const std::string &
value) {
445 if (
value ==
"use_y_diff")
447 else if (
value ==
"use_y_dot")
449 else if (
value ==
"none")
585 const MPI_Comm mpi_comm = MPI_COMM_WORLD);
638 std::function<
int(
const double t,
674 std::function<
int(
const double t,
725 std::function<void(
const double t,
728 const unsigned int step_number)>
773 <<
"One of the SUNDIALS IDA internal functions "
774 <<
" returned a negative error code: " << arg1
775 <<
". Please consult SUNDIALS manual.");
785 <<
"Please provide an implementation for the function \""
838 # ifdef DEAL_II_WITH_PETSC
839 # ifdef PETSC_USE_COMPLEX
841 "Sundials does not support complex scalar types, "
842 "but PETSc is configured to use a complex scalar type!");
846 "Sundials does not support complex scalar types, "
847 "but PETSc is configured to use a complex scalar type!");
848 # endif // PETSC_USE_COMPLEX
849 # endif // DEAL_II_WITH_PETSC
856 #endif // DEAL_II_WITH_SUNDIALS
unsigned int maximum_order
InitialConditionCorrection ic_type
InitialConditionCorrection reset_type
std::function< void(VectorType &)> reinit_vector
std::function< int(const VectorType &rhs, VectorType &dst)> solve_jacobian_system
double relative_tolerance
std::function< void(const double t, const VectorType &sol, const VectorType &sol_dot, const unsigned int step_number)> output_step
unsigned maximum_non_linear_iterations_ic
SymmetricTensor< 2, dim, Number > e(const Tensor< 2, dim, Number > &F)
std::function< bool(const double t, VectorType &sol, VectorType &sol_dot)> solver_should_restart
GrowingVectorMemory< VectorType > mem
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 double absolute_tolerance=1e-6, const double relative_tolerance=1e-5, const bool ignore_algebraic_terms_for_errors=true, const InitialConditionCorrection &ic_type=use_y_diff, const InitialConditionCorrection &reset_type=use_y_diff, const unsigned int maximum_non_linear_iterations_ic=5)
std::function< IndexSet()> differential_components
#define DEAL_II_NAMESPACE_OPEN
void set_functions_to_trigger_an_assert()
std::function< int(const double t, const VectorType &y, const VectorType &y_dot, const double alpha)> setup_jacobian
bool ignore_algebraic_terms_for_errors
DeclException1(ExcIDAError, int,<< "One of the SUNDIALS IDA internal functions "<< " returned a negative error code: "<< arg1<< ". Please consult SUNDIALS manual.")
void add_action(const std::string &entry, const std::function< void(const std::string &value)> &action)
InitialConditionCorrection
unsigned int maximum_non_linear_iterations
static ::ExceptionBase & ExcInternalError()
unsigned int solve_dae(VectorType &solution, VectorType &solution_dot)
void reset(const double t, const double h, VectorType &y, VectorType &yp)
void add_parameter(const std::string &entry, ParameterType ¶meter, const std::string &documentation="", const Patterns::PatternBase &pattern= *Patterns::Tools::Convert< ParameterType >::to_pattern(), const bool has_to_be_set=false)
IDA(const AdditionalData &data=AdditionalData(), const MPI_Comm mpi_comm=MPI_COMM_WORLD)
std::function< int(const double t, const VectorType &y, const VectorType &y_dot, VectorType &res)> residual
void enter_subsection(const std::string &subsection)
double absolute_tolerance
void add_parameters(ParameterHandler &prm)
std::function< VectorType &()> get_local_tolerances
#define DEAL_II_NAMESPACE_CLOSE
#define AssertThrow(cond, exc)