16#ifndef dealii_sundials_kinsol_h
17#define dealii_sundials_kinsol_h
21#ifdef DEAL_II_WITH_SUNDIALS
33# include <boost/signals2.hpp>
35# include <kinsol/kinsol.h>
36# if DEAL_II_SUNDIALS_VERSION_LT(4, 1, 0)
37# include <kinsol/kinsol_impl.h>
39# include <nvector/nvector_serial.h>
40# include <sundials/sundials_math.h>
41# include <sundials/sundials_types.h>
49# define AssertKINSOL(code) Assert(code >= 0, ExcKINSOLError(code))
204 template <
typename VectorType = Vector<
double>>
409 const MPI_Comm & mpi_comm = MPI_COMM_WORLD);
422 solve(VectorType &initial_guess_and_solution);
446 std::function<
int(
const VectorType &src, VectorType &dst)>
residual;
462 std::function<
int(
const VectorType &src, VectorType &dst)>
510 std::function<
int(
const VectorType ¤t_u,
const VectorType ¤t_f)>
571 std::function<
int(
const VectorType &ycur,
572 const VectorType &fcur,
573 const VectorType &rhs,
617 int(
const VectorType &rhs, VectorType &dst,
const double tolerance)>
681 <<
"One of the SUNDIALS KINSOL internal functions "
682 <<
"returned a negative error code: " << arg1
683 <<
". Please consult SUNDIALS manual.");
693 <<
"Please provide an implementation for the function \""
unsigned int maximum_setup_calls
double maximum_newton_step
SolutionStrategy strategy
unsigned int maximum_beta_failures
AdditionalData(const SolutionStrategy &strategy=linesearch, const unsigned int maximum_non_linear_iterations=200, const double function_tolerance=0.0, const double step_tolerance=0.0, const bool no_init_setup=false, const unsigned int maximum_setup_calls=0, const double maximum_newton_step=0.0, const double dq_relative_error=0.0, const unsigned int maximum_beta_failures=0, const unsigned int anderson_subspace_size=0)
void add_parameters(ParameterHandler &prm)
unsigned int maximum_non_linear_iterations
unsigned int anderson_subspace_size
double function_tolerance
std::function< VectorType &()> get_solution_scaling
std::function< int(const VectorType &rhs, VectorType &dst, const double tolerance)> solve_with_jacobian
std::function< int(const VectorType ¤t_u, const VectorType ¤t_f)> setup_jacobian
void set_functions_to_trigger_an_assert()
KINSOL(const AdditionalData &data=AdditionalData(), const MPI_Comm &mpi_comm=MPI_COMM_WORLD)
std::function< int(const VectorType &ycur, const VectorType &fcur, const VectorType &rhs, VectorType &dst)> solve_jacobian_system
unsigned int solve(VectorType &initial_guess_and_solution)
GrowingVectorMemory< VectorType > mem
std::function< int(const VectorType &src, VectorType &dst)> residual
std::function< VectorType &()> get_function_scaling
std::function< void(VectorType &)> reinit_vector
std::function< int(const VectorType &src, VectorType &dst)> iteration_function
#define DEAL_II_NAMESPACE_OPEN
#define DEAL_II_DEPRECATED_EARLY
#define DEAL_II_NAMESPACE_CLOSE
static ::ExceptionBase & ExcKINSOLError(int arg1)
static ::ExceptionBase & ExcFunctionNotProvided(std::string arg1)
DeclException1(ExcARKodeError, int,<< "One of the SUNDIALS ARKode internal functions "<< " returned a negative error code: "<< arg1<< ". Please consult SUNDIALS manual.")