DOLFIN-X
DOLFIN-X C++ interface
NewtonSolver.h
1 // Copyright (C) 2005-2021 Garth N. Wells
2 //
3 // This file is part of DOLFINX (https://www.fenicsproject.org)
4 //
5 // SPDX-License-Identifier: LGPL-3.0-or-later
6 
7 #pragma once
8 
9 #include <dolfinx/common/MPI.h>
10 #include <dolfinx/la/PETScKrylovSolver.h>
11 #include <functional>
12 #include <memory>
13 #include <petscmat.h>
14 #include <petscvec.h>
15 #include <utility>
16 
17 namespace dolfinx
18 {
19 
20 namespace la
21 {
22 class PETScKrylovSolver;
23 } // namespace la
24 
25 namespace nls
26 {
27 
30 
32 {
33 public:
36  explicit NewtonSolver(MPI_Comm comm);
37 
39  virtual ~NewtonSolver();
40 
45  void setF(const std::function<void(const Vec, Vec)>& F, Vec b);
46 
51  void setJ(const std::function<void(const Vec, Mat)>& J, Mat Jmat);
52 
56  void setP(const std::function<void(const Vec, Mat)>& P, Mat Pmat);
57 
62  void set_form(const std::function<void(Vec x)>& form);
63 
67  void
68  set_convergence_check(const std::function<std::pair<double, bool>(
69  const nls::NewtonSolver& solver, const Vec r)>& c);
70 
74  void set_update(const std::function<void(const nls::NewtonSolver& solver,
75  const Vec dx, Vec x)>& update);
76 
82  std::pair<int, bool> solve(Vec x);
83 
87  int iteration() const;
88 
92  int krylov_iterations() const;
93 
96  double residual() const;
97 
100  double residual0() const;
101 
103  MPI_Comm mpi_comm() const;
104 
106  int max_it = 50;
107 
109  double rtol = 1e-9;
110 
112  double atol = 1e-10;
113 
114  // FIXME: change to string to enum
116  std::string convergence_criterion = "residual";
117 
119  bool report = true;
120 
123 
125  double relaxation_parameter = 1.0;
126 
127 private:
128  // Function for computing the residual vector. The first argument is
129  // the latest solution vector x and the second argument is the
130  // residual vector.
131  std::function<void(const Vec x, Vec b)> _fnF;
132 
133  // Function for computing the Jacobian matrix operator. The first
134  // argument is the latest solution vector x and the second argument is
135  // the matrix operator.
136  std::function<void(const Vec x, Mat J)> _fnJ;
137 
138  // Function for computing the preconditioner matrix operator. The
139  // first argument is the latest solution vector x and the second
140  // argument is the matrix operator.
141  std::function<void(const Vec x, Mat P)> _fnP;
142 
143  // Function called before the residual and Jacobian function at each
144  // iteration.
145  std::function<void(const Vec x)> _system;
146 
147  // Residual vector
148  Vec _b = nullptr;
149 
150  // Jacobian matrix and preconditioner matrix
151  Mat _matJ = nullptr, _matP = nullptr;
152 
153  // Function to check for convergence
154  std::function<std::pair<double, bool>(const nls::NewtonSolver& solver,
155  const Vec r)>
156  _converged;
157 
158  // Function to update the solution once convergence is reached
159  std::function<void(const nls::NewtonSolver& solver, const Vec dx, Vec x)>
160  _update_solution;
161 
162  // Accumulated number of Krylov iterations since solve began
163  int _krylov_iterations;
164 
165  // Number of iterations
166  int _iteration;
167 
168  // Most recent residual and initial residual
169  double _residual, _residual0;
170 
171  // Linear solver
172  la::PETScKrylovSolver _solver;
173 
174  // Solution vector
175  Vec _dx = nullptr;
176 
177  // MPI communicator
178  dolfinx::MPI::Comm _mpi_comm;
179 };
180 } // namespace nls
181 } // namespace dolfinx
A duplicate MPI communicator and manage lifetime of the communicator.
Definition: MPI.h:36
This class implements Krylov methods for linear systems of the form Ax = b. It is a wrapper for the K...
Definition: PETScKrylovSolver.h:30
This class defines a Newton solver for nonlinear systems of equations of the form .
Definition: NewtonSolver.h:32
int max_it
Maximum number of iterations.
Definition: NewtonSolver.h:106
void setF(const std::function< void(const Vec, Vec)> &F, Vec b)
Set the function for computing the residual and the vector to the assemble the residual into.
Definition: NewtonSolver.cpp:92
std::pair< int, bool > solve(Vec x)
Solve the nonlinear problem for given and Jacobian .
Definition: NewtonSolver.cpp:135
void setP(const std::function< void(const Vec, Mat)> &P, Mat Pmat)
Set the function for computing the preconditioner matrix (optional)
Definition: NewtonSolver.cpp:108
double relaxation_parameter
Relaxation parameter.
Definition: NewtonSolver.h:125
double atol
Absolute tolerance.
Definition: NewtonSolver.h:112
void set_update(const std::function< void(const nls::NewtonSolver &solver, const Vec dx, Vec x)> &update)
Set function that is called after each Newton iteration to update the solution.
Definition: NewtonSolver.cpp:128
NewtonSolver(MPI_Comm comm)
Create nonlinear solver.
Definition: NewtonSolver.cpp:65
void set_form(const std::function< void(Vec x)> &form)
Set the function that is called before the residual or Jacobian are computed. It is commonly used to ...
Definition: NewtonSolver.cpp:116
int iteration() const
The number of Newton interations. It can can called by functions that check for convergence during a ...
Definition: NewtonSolver.cpp:264
int krylov_iterations() const
Return number of Krylov iterations elapsed since solve started.
Definition: NewtonSolver.cpp:262
bool report
Monitor convergence.
Definition: NewtonSolver.h:119
double rtol
Relative tolerance.
Definition: NewtonSolver.h:109
double residual0() const
Return initial residual.
Definition: NewtonSolver.cpp:268
virtual ~NewtonSolver()
Destructor.
Definition: NewtonSolver.cpp:80
bool error_on_nonconvergence
Throw error if solver fails to converge.
Definition: NewtonSolver.h:122
double residual() const
Return current residual.
Definition: NewtonSolver.cpp:266
std::string convergence_criterion
Convergence criterion.
Definition: NewtonSolver.h:116
MPI_Comm mpi_comm() const
Return MPI communicator.
Definition: NewtonSolver.cpp:270
void set_convergence_check(const std::function< std::pair< double, bool >(const nls::NewtonSolver &solver, const Vec r)> &c)
Set function that is called at the end of each Newton iteration to test for convergence.
Definition: NewtonSolver.cpp:121
void setJ(const std::function< void(const Vec, Mat)> &J, Mat Jmat)
Set the function for computing the Jacobian (dF/dx) and the matrix to assemble the residual into.
Definition: NewtonSolver.cpp:100