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
petsc_snes.h
Go to the documentation of this file.
1//-----------------------------------------------------------
2//
3// Copyright (C) 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#ifndef dealii_petsc_snes_h
17#define dealii_petsc_snes_h
18
19#include <deal.II/base/config.h>
20
21#ifdef DEAL_II_WITH_PETSC
22# include <deal.II/base/mpi.h>
25
29
30# include <petscsnes.h>
31
32# include <exception>
33# if defined(DEAL_II_HAVE_CXX20)
34# include <concepts>
35# endif
36
37
39
40namespace PETScWrappers
41{
46 {
47 public:
53 using real_type = PetscReal;
54
75 // Running parameters
76 const std::string &options_prefix = "",
77 const std::string &snes_type = "",
78 const std::string &snes_linesearch_type = "",
81 const real_type step_tolerance = 0,
82 const int maximum_non_linear_iterations = -1,
83 const int max_n_function_evaluations = -1)
92 {}
93
97 void
99
103 std::string options_prefix;
104
112 std::string snes_type;
113
121
128
135
142
149
156 };
157
246 template <typename VectorType = PETScWrappers::VectorBase,
247 typename PMatrixType = PETScWrappers::MatrixBase,
248 typename AMatrixType = PMatrixType>
251 std::constructible_from<
252 VectorType,
254 std::constructible_from<
255 PMatrixType,
257 std::constructible_from<AMatrixType, Mat>))
259 {
260 public:
266 using real_type = PetscReal;
267
272 const MPI_Comm mpi_comm = PETSC_COMM_WORLD);
273
278
284 operator SNES() const;
285
289 SNES
291
297
301 void
303
307 void
309
321 void
322 set_matrix(PMatrixType &P);
323
328 void
329 set_matrices(AMatrixType &A, PMatrixType &P);
330
338 unsigned int
339 solve(VectorType &x);
340
350 unsigned int
351 solve(VectorType &x, PMatrixType &P);
352
363 unsigned int
364 solve(VectorType &x, AMatrixType &A, PMatrixType &P);
365
374 std::function<void(const VectorType &x, VectorType &res)> residual;
375
385 std::function<void(const VectorType &x, AMatrixType &A, PMatrixType &P)>
387
400 std::function<void(const VectorType & x,
401 const unsigned int step_number,
402 const real_type f_norm)>
404
418 std::function<void(const VectorType &x)> setup_jacobian;
419
431 std::function<void(const VectorType &src, VectorType &dst)>
433
451 std::function<void(const VectorType &x, real_type &energy_value)> energy;
452
453 private:
457 SNES snes;
458
464
469
475 mutable std::exception_ptr pending_exception;
476 };
477
478} // namespace PETScWrappers
479
481
482#endif // DEAL_II_WITH_PETSC
483
484#endif
void add_parameters(ParameterHandler &prm)
std::function< void(const VectorType &src, VectorType &dst)> solve_with_jacobian
Definition petsc_snes.h:432
void set_matrices(AMatrixType &A, PMatrixType &P)
std::function< void(const VectorType &x, real_type &energy_value)> energy
Definition petsc_snes.h:451
SmartPointer< PMatrixType, NonlinearSolver > P
Definition petsc_snes.h:463
unsigned int solve(VectorType &x, AMatrixType &A, PMatrixType &P)
std::function< void(const VectorType &x, VectorType &res)> residual
Definition petsc_snes.h:374
unsigned int solve(VectorType &x)
std::exception_ptr pending_exception
Definition petsc_snes.h:475
NonlinearSolver(const NonlinearSolverData &data=NonlinearSolverData(), const MPI_Comm mpi_comm=PETSC_COMM_WORLD)
MPI_Comm get_mpi_communicator() const
void set_matrix(PMatrixType &P)
std::function< void(const VectorType &x, const unsigned int step_number, const real_type f_norm)> monitor
Definition petsc_snes.h:403
unsigned int solve(VectorType &x, PMatrixType &P)
std::function< void(const VectorType &x, AMatrixType &A, PMatrixType &P)> jacobian
Definition petsc_snes.h:386
std::function< void(const VectorType &x)> setup_jacobian
Definition petsc_snes.h:418
SmartPointer< AMatrixType, NonlinearSolver > A
Definition petsc_snes.h:462
void reinit(const NonlinearSolverData &data)
#define DEAL_II_NAMESPACE_OPEN
Definition config.h:472
#define DEAL_II_CXX20_REQUIRES(condition)
Definition config.h:160
#define DEAL_II_NAMESPACE_CLOSE
Definition config.h:473
NonlinearSolverData(const std::string &options_prefix="", const std::string &snes_type="", const std::string &snes_linesearch_type="", const real_type absolute_tolerance=0, const real_type relative_tolerance=0, const real_type step_tolerance=0, const int maximum_non_linear_iterations=-1, const int max_n_function_evaluations=-1)
Definition petsc_snes.h:74