Loading [MathJax]/extensions/TeX/newcommand.js
 Reference documentation for deal.II version 9.3.3
\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\}}
All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Modules Pages
kinsol.h
Go to the documentation of this file.
1//-----------------------------------------------------------
2//
3// Copyright (C) 2017 - 2021 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_kinsol_h
17#define dealii_sundials_kinsol_h
18
19#include <deal.II/base/config.h>
20
21#ifdef DEAL_II_WITH_SUNDIALS
22
23
27# include <deal.II/base/mpi.h>
29
30# include <deal.II/lac/vector.h>
32
33# include <boost/signals2.hpp>
34
35# include <kinsol/kinsol.h>
36# if DEAL_II_SUNDIALS_VERSION_LT(4, 1, 0)
37# include <kinsol/kinsol_impl.h>
38# endif
39# include <nvector/nvector_serial.h>
40# include <sundials/sundials_math.h>
41# include <sundials/sundials_types.h>
42
43# include <memory>
44
45
47
48// Shorthand notation for KINSOL error codes.
49# define AssertKINSOL(code) Assert(code >= 0, ExcKINSOLError(code))
50
51namespace SUNDIALS
52{
204 template <typename VectorType = Vector<double>>
205 class KINSOL
206 {
207 public:
212 {
213 public:
219 {
223 newton = KIN_NONE,
227 linesearch = KIN_LINESEARCH,
231 fixed_point = KIN_FP,
235 picard = KIN_PICARD,
236 };
237
267 const unsigned int maximum_non_linear_iterations = 200,
268 const double function_tolerance = 0.0,
269 const double step_tolerance = 0.0,
270 const bool no_init_setup = false,
271 const unsigned int maximum_setup_calls = 0,
272 const double maximum_newton_step = 0.0,
273 const double dq_relative_error = 0.0,
274 const unsigned int maximum_beta_failures = 0,
275 const unsigned int anderson_subspace_size = 0);
276
315 void
317
326
331
339
347
357
367
374
383
390
398 };
399
409 const MPI_Comm & mpi_comm = MPI_COMM_WORLD);
410
414 ~KINSOL();
415
421 unsigned int
422 solve(VectorType &initial_guess_and_solution);
423
431 std::function<void(VectorType &)> reinit_vector;
432
446 std::function<int(const VectorType &src, VectorType &dst)> residual;
447
462 std::function<int(const VectorType &src, VectorType &dst)>
464
510 std::function<int(const VectorType &current_u, const VectorType &current_f)>
512
571 std::function<int(const VectorType &ycur,
572 const VectorType &fcur,
573 const VectorType &rhs,
574 VectorType & dst)>
576
616 std::function<
617 int(const VectorType &rhs, VectorType &dst, const double tolerance)>
619
658 std::function<VectorType &()> get_solution_scaling;
659
674 std::function<VectorType &()> get_function_scaling;
675
680 int,
681 << "One of the SUNDIALS KINSOL internal functions "
682 << "returned a negative error code: " << arg1
683 << ". Please consult SUNDIALS manual.");
684
685
686 private:
692 std::string,
693 << "Please provide an implementation for the function \""
694 << arg1 << "\"");
695
701 void
703
708
713
717 N_Vector solution;
718
722 N_Vector u_scale;
723
727 N_Vector f_scale;
728
733
738 };
739
740} // namespace SUNDIALS
741
742
744
745#endif
746
747#endif
unsigned int maximum_setup_calls
Definition: kinsol.h:366
SolutionStrategy strategy
Definition: kinsol.h:325
unsigned int maximum_beta_failures
Definition: kinsol.h:389
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)
Definition: kinsol.cc:57
void add_parameters(ParameterHandler &prm)
Definition: kinsol.cc:84
unsigned int maximum_non_linear_iterations
Definition: kinsol.h:330
unsigned int anderson_subspace_size
Definition: kinsol.h:397
std::function< VectorType &()> get_solution_scaling
Definition: kinsol.h:658
std::function< int(const VectorType &rhs, VectorType &dst, const double tolerance)> solve_with_jacobian
Definition: kinsol.h:618
std::function< int(const VectorType &current_u, const VectorType &current_f)> setup_jacobian
Definition: kinsol.h:511
void set_functions_to_trigger_an_assert()
Definition: kinsol.cc:607
KINSOL(const AdditionalData &data=AdditionalData(), const MPI_Comm &mpi_comm=MPI_COMM_WORLD)
Definition: kinsol.cc:340
N_Vector solution
Definition: kinsol.h:717
std::function< int(const VectorType &ycur, const VectorType &fcur, const VectorType &rhs, VectorType &dst)> solve_jacobian_system
Definition: kinsol.h:575
unsigned int solve(VectorType &initial_guess_and_solution)
Definition: kinsol.cc:376
GrowingVectorMemory< VectorType > mem
Definition: kinsol.h:737
MPI_Comm communicator
Definition: kinsol.h:732
std::function< int(const VectorType &src, VectorType &dst)> residual
Definition: kinsol.h:446
N_Vector f_scale
Definition: kinsol.h:727
std::function< VectorType &()> get_function_scaling
Definition: kinsol.h:674
std::function< void(VectorType &)> reinit_vector
Definition: kinsol.h:431
AdditionalData data
Definition: kinsol.h:707
std::function< int(const VectorType &src, VectorType &dst)> iteration_function
Definition: kinsol.h:463
N_Vector u_scale
Definition: kinsol.h:722
void * kinsol_mem
Definition: kinsol.h:712
#define DEAL_II_NAMESPACE_OPEN
Definition: config.h:402
#define DEAL_II_DEPRECATED_EARLY
Definition: config.h:170
#define DEAL_II_NAMESPACE_CLOSE
Definition: config.h:403
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.")