9#include <dolfinx/common/MPI.h>
10#include <dolfinx/la/petsc.h>
57 void setF(std::function<
void(
const Vec, Vec)> F, Vec b);
63 void setJ(std::function<
void(
const Vec, Mat)> J, Mat Jmat);
69 void setP(std::function<
void(
const Vec, Mat)> P, Mat Pmat);
91 void set_form(std::function<
void(Vec)> form);
97 std::function<std::pair<double, bool>(
const NewtonSolver&,
const Vec)> c);
103 std::function<
void(
const NewtonSolver& solver,
const Vec, Vec)> update);
110 std::pair<int, bool>
solve(Vec x);
131 MPI_Comm
comm()
const;
159 std::function<void(
const Vec x, Vec b)> _fnF;
164 std::function<void(
const Vec x, Mat J)> _fnJ;
169 std::function<void(
const Vec x, Mat P)> _fnP;
173 std::function<void(
const Vec x)> _system;
179 Mat _matJ =
nullptr, _matP =
nullptr;
182 std::function<std::pair<double, bool>(
const NewtonSolver& solver,
187 std::function<void(
const NewtonSolver& solver,
const Vec dx, Vec x)>
191 int _krylov_iterations;
197 double _residual, _residual0;
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