Reference documentation for deal.II version GIT relicensing-487-ge9eb5ab491 2024-04-25 07:20:02+00:00
\(\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// SPDX-License-Identifier: LGPL-2.1-or-later
4// Copyright (C) 2023 by the deal.II authors
5//
6// This file is part of the deal.II library.
7//
8// Part of the source code is dual licensed under Apache-2.0 WITH
9// LLVM-exception OR LGPL-2.1-or-later. Detailed license information
10// governing the source code and code contributions can be found in
11// LICENSE.md and CONTRIBUTING.md at the top level directory of deal.II.
12//
13// ------------------------------------------------------------------------
14
15#ifndef dealii_petsc_snes_h
16#define dealii_petsc_snes_h
17
18#include <deal.II/base/config.h>
19
20#ifdef DEAL_II_WITH_PETSC
21# include <deal.II/base/mpi.h>
24
28
29# include <petscsnes.h>
30
31# include <exception>
32
33# if defined(DEAL_II_HAVE_CXX20)
34# include <concepts>
35# endif
36
37
39
40namespace PETScWrappers
41{
157
257 template <typename VectorType = PETScWrappers::VectorBase,
258 typename PMatrixType = PETScWrappers::MatrixBase,
259 typename AMatrixType = PMatrixType>
262 std::constructible_from<
263 VectorType,
265 std::constructible_from<
266 PMatrixType,
268 std::constructible_from<AMatrixType, Mat>))
269 class NonlinearSolver
270 {
271 public:
277 using real_type = PetscReal;
278
282 NonlinearSolver(const NonlinearSolverData &data = NonlinearSolverData(),
283 const MPI_Comm mpi_comm = PETSC_COMM_WORLD);
284
288 ~NonlinearSolver();
289
295 operator SNES() const;
296
300 SNES
301 petsc_snes();
302
307 get_mpi_communicator() const;
308
312 void
313 reinit();
314
319 void
320 reinit(const NonlinearSolverData &data);
321
333 void
334 set_matrix(PMatrixType &P);
335
341 void
342 set_matrices(AMatrixType &A, PMatrixType &P);
343
351 unsigned int
352 solve(VectorType &x);
353
363 unsigned int
364 solve(VectorType &x, PMatrixType &P);
365
376 unsigned int
377 solve(VectorType &x, AMatrixType &A, PMatrixType &P);
378
387 std::function<void(const VectorType &x, VectorType &res)> residual;
388
398 std::function<void(const VectorType &x, AMatrixType &A, PMatrixType &P)>
399 jacobian;
400
413 std::function<void(const VectorType &x,
414 const unsigned int step_number,
415 const real_type f_norm)>
416 monitor;
417
431 std::function<void(const VectorType &x)> setup_jacobian;
432
444 std::function<void(const VectorType &src, VectorType &dst)>
445 solve_with_jacobian;
446
464 std::function<void(const VectorType &x, real_type &energy_value)> energy;
465
466 private:
470 SNES snes;
471
477
481 bool need_dummy_assemble;
482
488 mutable std::exception_ptr pending_exception;
489 };
490
491} // namespace PETScWrappers
492
494
495#endif // DEAL_II_WITH_PETSC
496
497#endif
void add_parameters(ParameterHandler &prm)
#define DEAL_II_NAMESPACE_OPEN
Definition config.h:502
#define DEAL_II_CXX20_REQUIRES(condition)
Definition config.h:177
#define DEAL_II_NAMESPACE_CLOSE
Definition config.h:503
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