DOLFINx
DOLFINx C++ interface
Loading...
Searching...
No Matches
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/petsc.h>
11#include <functional>
12#include <memory>
13#include <petscmat.h>
14#include <petscvec.h>
15#include <utility>
16
17namespace dolfinx
18{
19
20namespace la::petsc
21{
22class KrylovSolver;
23} // namespace la::petsc
24
25namespace nls::petsc
26{
27
30
32{
33public:
36 explicit NewtonSolver(MPI_Comm comm);
37
39 NewtonSolver(NewtonSolver&& solver) = default;
40
41 // Copy constructor (deleted)
42 NewtonSolver(const NewtonSolver& solver) = delete;
43
45 NewtonSolver& operator=(NewtonSolver&& solver) = default;
46
47 // Assignment operator (deleted)
48 NewtonSolver& operator=(const NewtonSolver& solver) = delete;
49
52
57 void setF(std::function<void(const Vec, Vec)> F, Vec b);
58
63 void setJ(std::function<void(const Vec, Mat)> J, Mat Jmat);
64
69 void setP(std::function<void(const Vec, Mat)> P, Mat Pmat);
70
78
86
91 void set_form(std::function<void(Vec)> form);
92
97 std::function<std::pair<double, bool>(const NewtonSolver&, const Vec)> c);
98
102 void set_update(
103 std::function<void(const NewtonSolver& solver, const Vec, Vec)> update);
104
110 std::pair<int, bool> solve(Vec x);
111
115 int iteration() const;
116
120 int krylov_iterations() const;
121
124 double residual() const;
125
128 double residual0() const;
129
131 MPI_Comm comm() const;
132
134 int max_it = 50;
135
137 double rtol = 1e-9;
138
140 double atol = 1e-10;
141
142 // FIXME: change to string to enum
144 std::string convergence_criterion = "residual";
145
147 bool report = true;
148
151
154
155private:
156 // Function for computing the residual vector. The first argument is
157 // the latest solution vector x and the second argument is the
158 // residual vector.
159 std::function<void(const Vec x, Vec b)> _fnF;
160
161 // Function for computing the Jacobian matrix operator. The first
162 // argument is the latest solution vector x and the second argument is
163 // the matrix operator.
164 std::function<void(const Vec x, Mat J)> _fnJ;
165
166 // Function for computing the preconditioner matrix operator. The
167 // first argument is the latest solution vector x and the second
168 // argument is the matrix operator.
169 std::function<void(const Vec x, Mat P)> _fnP;
170
171 // Function called before the residual and Jacobian function at each
172 // iteration.
173 std::function<void(const Vec x)> _system;
174
175 // Residual vector
176 Vec _b = nullptr;
177
178 // Jacobian matrix and preconditioner matrix
179 Mat _matJ = nullptr, _matP = nullptr;
180
181 // Function to check for convergence
182 std::function<std::pair<double, bool>(const NewtonSolver& solver,
183 const Vec r)>
184 _converged;
185
186 // Function to update the solution once convergence is reached
187 std::function<void(const NewtonSolver& solver, const Vec dx, Vec x)>
188 _update_solution;
189
190 // Accumulated number of Krylov iterations since solve began
191 int _krylov_iterations;
192
193 // Number of iterations
194 int _iteration;
195
196 // Most recent residual and initial residual
197 double _residual, _residual0;
198
199 // Linear solver
201
202 // Solution vector
203 Vec _dx = nullptr;
204
205 // MPI communicator
206 dolfinx::MPI::Comm _comm;
207};
208} // namespace nls::petsc
209} // namespace dolfinx
A duplicate MPI communicator and manage lifetime of the communicator.
Definition MPI.h:43
This class implements Krylov methods for linear systems of the form Ax = b. It is a wrapper for the K...
Definition petsc.h:451
This class defines a Newton solver for nonlinear systems of equations of the form .
Definition NewtonSolver.h:32
double relaxation_parameter
Relaxation parameter.
Definition NewtonSolver.h:153
int krylov_iterations() const
Get number of Krylov iterations elapsed since solve started.
Definition NewtonSolver.cpp:274
void set_form(std::function< void(Vec)> form)
Set the function that is called before the residual or Jacobian are computed. It is commonly used to ...
Definition NewtonSolver.cpp:126
void setF(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:91
NewtonSolver(NewtonSolver &&solver)=default
Move constructor.
double rtol
Relative tolerance.
Definition NewtonSolver.h:137
~NewtonSolver()
Destructor.
Definition NewtonSolver.cpp:79
NewtonSolver & operator=(NewtonSolver &&solver)=default
Move assignment constructor.
double atol
Absolute tolerance.
Definition NewtonSolver.h:140
double residual() const
Get current residual.
Definition NewtonSolver.cpp:281
void setP(std::function< void(const Vec, Mat)> P, Mat Pmat)
Set the function for computing the preconditioner matrix (optional).
Definition NewtonSolver.cpp:107
double residual0() const
Return initial residual.
Definition NewtonSolver.cpp:283
std::string convergence_criterion
Convergence criterion.
Definition NewtonSolver.h:144
int iteration() const
The number of Newton iterations. It can can called by functions that check for convergence during a s...
Definition NewtonSolver.cpp:279
void set_update(std::function< void(const NewtonSolver &solver, const Vec, Vec)> update)
Set function that is called after each Newton iteration to update the solution.
Definition NewtonSolver.cpp:137
void setJ(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:99
bool report
Monitor convergence.
Definition NewtonSolver.h:147
void set_convergence_check(std::function< std::pair< double, bool >(const NewtonSolver &, const Vec)> c)
Set function that is called at the end of each Newton iteration to test for convergence.
Definition NewtonSolver.cpp:131
bool error_on_nonconvergence
Throw error if solver fails to converge.
Definition NewtonSolver.h:150
MPI_Comm comm() const
Get MPI communicator.
Definition NewtonSolver.cpp:285
const la::petsc::KrylovSolver & get_krylov_solver() const
Get the internal Krylov solver used to solve for the Newton updates (const version).
Definition NewtonSolver.cpp:116
std::pair< int, bool > solve(Vec x)
Solve the nonlinear problem for given and Jacobian .
Definition NewtonSolver.cpp:143
int max_it
Maximum number of iterations.
Definition NewtonSolver.h:134
Top-level namespace.
Definition defines.h:12