Reference documentation for deal.II version 9.5.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\}}\)
Loading...
Searching...
No Matches
kinsol.h
Go to the documentation of this file.
1// ---------------------------------------------------------------------
2//
3// Copyright (C) 2017 - 2023 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
17#ifndef dealii_sundials_kinsol_h
18#define dealii_sundials_kinsol_h
19
20#include <deal.II/base/config.h>
21
22#ifdef DEAL_II_WITH_SUNDIALS
23
24
30
31# include <deal.II/lac/vector.h>
33
34# include <boost/signals2.hpp>
35
36# include <kinsol/kinsol.h>
37# include <nvector/nvector_serial.h>
38# include <sundials/sundials_math.h>
39# include <sundials/sundials_types.h>
40
41# include <exception>
42# include <memory>
43
44
46
47// Shorthand notation for KINSOL error codes.
48# define AssertKINSOL(code) Assert(code >= 0, ExcKINSOLError(code))
49
50namespace SUNDIALS
51{
183 template <typename VectorType = Vector<double>>
184 class KINSOL
185 {
186 public:
191 {
192 public:
198 {
202 newton = KIN_NONE,
206 linesearch = KIN_LINESEARCH,
210 fixed_point = KIN_FP,
214 picard = KIN_PICARD,
215 };
216
246 const unsigned int maximum_non_linear_iterations = 200,
247 const double function_tolerance = 0.0,
248 const double step_tolerance = 0.0,
249 const bool no_init_setup = false,
250 const unsigned int maximum_setup_calls = 0,
251 const double maximum_newton_step = 0.0,
252 const double dq_relative_error = 0.0,
253 const unsigned int maximum_beta_failures = 0,
254 const unsigned int anderson_subspace_size = 0);
255
294 void
296
305
310
318
326
336
346
353
362
369
377 };
378
389
397 KINSOL(const AdditionalData &data, const MPI_Comm mpi_comm);
398
402 ~KINSOL();
403
409 unsigned int
410 solve(VectorType &initial_guess_and_solution);
411
426 std::function<void(VectorType &)> reinit_vector;
427
441 std::function<void(const VectorType &src, VectorType &dst)> residual;
442
456 std::function<void(const VectorType &src, VectorType &dst)>
458
504 std::function<void(const VectorType &current_u,
505 const VectorType &current_f)>
507
574 std::function<void(const VectorType &ycur,
575 const VectorType &fcur,
576 const VectorType &rhs,
577 VectorType & dst)>
579
619 std::function<
620 void(const VectorType &rhs, VectorType &dst, const double tolerance)>
622
668 std::function<VectorType &()> get_solution_scaling;
669
691 std::function<VectorType &()> get_function_scaling;
692
697 int,
698 << "One of the SUNDIALS KINSOL internal functions "
699 << "returned a negative error code: " << arg1
700 << ". Please consult SUNDIALS manual.");
701
702 private:
708 std::string,
709 << "Please provide an implementation for the function \""
710 << arg1 << "\"");
711
717 void
719
724
729
734
735# if DEAL_II_SUNDIALS_VERSION_GTE(6, 0, 0)
740# endif
741
742
747
753 mutable std::exception_ptr pending_exception;
754 };
755
756} // namespace SUNDIALS
757
758
760
761#endif
762
763#endif
SolutionStrategy strategy
Definition kinsol.h:304
unsigned int maximum_beta_failures
Definition kinsol.h:368
void add_parameters(ParameterHandler &prm)
Definition kinsol.cc:84
unsigned int maximum_non_linear_iterations
Definition kinsol.h:309
unsigned int anderson_subspace_size
Definition kinsol.h:376
std::function< VectorType &()> get_solution_scaling
Definition kinsol.h:668
std::function< void(const VectorType &current_u, const VectorType &current_f)> setup_jacobian
Definition kinsol.h:506
void set_functions_to_trigger_an_assert()
Definition kinsol.cc:547
MPI_Comm mpi_communicator
Definition kinsol.h:728
unsigned int solve(VectorType &initial_guess_and_solution)
Definition kinsol.cc:187
SUNContext kinsol_ctx
Definition kinsol.h:739
std::function< void(const VectorType &src, VectorType &dst)> residual
Definition kinsol.h:441
std::function< void(const VectorType &src, VectorType &dst)> iteration_function
Definition kinsol.h:457
GrowingVectorMemory< VectorType > mem
Definition kinsol.h:746
std::exception_ptr pending_exception
Definition kinsol.h:753
std::function< void(const VectorType &rhs, VectorType &dst, const double tolerance)> solve_with_jacobian
Definition kinsol.h:621
std::function< VectorType &()> get_function_scaling
Definition kinsol.h:691
std::function< void(VectorType &)> reinit_vector
Definition kinsol.h:426
std::function< void(const VectorType &ycur, const VectorType &fcur, const VectorType &rhs, VectorType &dst)> solve_jacobian_system
Definition kinsol.h:578
AdditionalData data
Definition kinsol.h:723
void * kinsol_mem
Definition kinsol.h:733
#define DEAL_II_DEPRECATED
Definition config.h:172
#define DEAL_II_NAMESPACE_OPEN
Definition config.h:472
#define DEAL_II_NAMESPACE_CLOSE
Definition config.h:473
static ::ExceptionBase & ExcKINSOLError(int arg1)
static ::ExceptionBase & ExcFunctionNotProvided(std::string arg1)
#define DeclException1(Exception1, type1, outsequence)
Definition exceptions.h:510