9 #include <dolfinx/common/MPI.h>
10 #include <dolfinx/la/PETScKrylovSolver.h>
22 class PETScKrylovSolver;
45 void setF(
const std::function<
void(
const Vec, Vec)>& F, Vec b);
51 void setJ(
const std::function<
void(
const Vec, Mat)>& J, Mat Jmat);
56 void setP(
const std::function<
void(
const Vec, Mat)>& P, Mat Pmat);
62 void set_form(
const std::function<
void(Vec x)>& form);
75 const Vec dx, Vec x)>& update);
82 std::pair<int, bool>
solve(Vec x);
131 std::function<void(
const Vec x, Vec b)> _fnF;
136 std::function<void(
const Vec x, Mat J)> _fnJ;
141 std::function<void(
const Vec x, Mat P)> _fnP;
145 std::function<void(
const Vec x)> _system;
151 Mat _matJ =
nullptr, _matP =
nullptr;
163 int _krylov_iterations;
169 double _residual, _residual0;
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