19 #include "libmesh/diff_system.h" 20 #include "libmesh/dof_map.h" 21 #include "libmesh/libmesh_logging.h" 22 #include "libmesh/linear_solver.h" 23 #include "libmesh/newton_solver.h" 24 #include "libmesh/numeric_vector.h" 25 #include "libmesh/sparse_matrix.h" 40 Real & current_residual,
46 if ((current_residual < last_residual) ||
61 (current_residual > last_residual &&
68 substepdivision = std::min(0.5, last_residual/current_residual);
69 substepdivision = std::max(substepdivision, tol*2.);
72 substepdivision = 0.5;
74 newton_iterate.
add (bx * (1.-substepdivision),
76 newton_iterate.
close();
77 bx *= substepdivision;
89 current_residual = rhs.
l2_norm();
92 << current_residual << std::endl;
96 (current_residual > last_residual)))
103 libmesh_convergence_failure();
123 Real x = bx, w = bx, v = bx;
126 Real fx = current_residual,
127 fw = current_residual,
128 fv = current_residual;
131 const unsigned int max_i = 20;
136 for (
unsigned int i=1; i <= max_i; i++)
138 Real xm = (ax+cx)*0.5;
140 Real tol2 = 2.0 * tol1;
143 if (
std::abs(x-xm) <= (tol2 - 0.5 * (cx - ax)))
151 Real r = (x-w)*(fx-fv);
152 Real q = (x-v)*(fx-fw);
153 Real p = (x-v)*q-(x-w)*r;
164 e = x >= xm ? ax-x : cx-x;
165 d = golden_ratio * e;
171 if (x+d-ax < tol2 || cx-(x+d) < tol2)
172 d =
SIGN(tol1, xm - x);
178 e = x >= xm ? ax-x : cx-x;
179 d = golden_ratio * e;
185 newton_iterate.
add (bx - u, linear_solution);
186 newton_iterate.
close();
209 fv = fw; fw = fx; fx = fu;
217 if (fu <= fw || w == x)
222 else if (fu <= fv || v == x || v == w)
231 libMesh::out <<
"Warning! Too many iterations used in Brent line search!" 239 require_residual_reduction(true),
240 require_finite_residual(true),
241 brent_line_search(true),
242 track_linear_convergence(false),
244 linear_tolerance_multiplier(1e-3),
282 LOG_SCOPE(
"solve()",
"NewtonSolver");
289 std::unique_ptr<NumericVector<Number>> linear_solution_ptr = newton_iterate.
zero_clone();
293 newton_iterate.
close();
294 linear_solution.
close();
297 #ifdef LIBMESH_ENABLE_CONSTRAINTS 328 <<
" with residual Not-a-Number" 330 libmesh_convergence_failure();
339 <<
" meets absolute tolerance " 348 if (current_residual == 0)
362 Real last_residual = current_residual;
374 << current_residual << std::endl;
377 current_linear_tolerance =
378 double(std::min (current_linear_tolerance,
388 current_linear_tolerance =
389 double(std::max(current_linear_tolerance,
398 linear_solution.
zero();
402 << current_linear_tolerance << std::endl;
405 const std::pair<unsigned int, Real> rval =
407 linear_solution, rhs, current_linear_tolerance,
415 if (linear_c_reason < 0)
420 libMesh::out <<
"Linear solver failed during Newton step, dropping out." 429 #ifdef LIBMESH_ENABLE_CONSTRAINTS 432 (
_system, &linear_solution,
true);
435 const unsigned int linear_steps = rval.first;
439 const bool linear_solve_finished =
443 libMesh::out <<
"Linear solve finished, step " << linear_steps
444 <<
", residual " << rval.second
453 newton_iterate.
add (-1., linear_solution);
454 newton_iterate.
close();
459 norm_total = newton_iterate.
l2_norm();
462 newton_iterate, norm_total,
478 current_residual = rhs.l2_norm();
481 << current_residual << std::endl;
485 linear_solve_finished &&
486 current_residual <= last_residual))
490 norm_delta, linear_solve_finished &&
491 current_residual <= last_residual);
500 last_residual, current_residual,
501 newton_iterate, linear_solution);
502 norm_delta *= steplength;
514 libMesh::out <<
" Nonlinear solver reached maximum step " 516 << current_residual << std::endl;
524 libmesh_convergence_failure();
530 norm_total = newton_iterate.
l2_norm();
538 << norm_delta / norm_total
539 <<
", |du| = " << norm_delta
545 linear_solve_finished))
549 norm_delta / steplength,
550 linear_solve_finished);
557 #ifdef LIBMESH_ENABLE_CONSTRAINTS 569 !max_nonlinear_iterations);
578 bool linear_solve_finished)
581 bool has_converged =
false;
587 has_converged =
true;
595 has_converged =
true;
599 if (!linear_solve_finished)
601 return has_converged;
608 has_converged =
true;
616 has_converged =
true;
619 return has_converged;
624 Real current_residual,
626 bool linear_solve_finished)
631 libMesh::out <<
" Nonlinear solver converged, step " << step_num
632 <<
", residual " << current_residual
639 << current_residual <<
" > " 647 libMesh::out <<
" Nonlinear solver converged, step " << step_num
648 <<
", residual reduction " 663 if (!linear_solve_finished)
669 libMesh::out <<
" Nonlinear solver converged, step " << step_num
670 <<
", absolute step size " 688 libMesh::out <<
" Nonlinear solver converged, step " << step_num
689 <<
", relative step size " bool continue_after_max_iterations
Defaults to true, telling the DiffSolver to continue rather than exit when a solve has reached its ma...
virtual void init() override
The initialization function.
NewtonSolver(sys_type &system)
Constructor.
std::unique_ptr< LinearSolver< Number > > _linear_solver
The LinearSolver defines the interface used to solve the linear_implicit system.
Real absolute_step_tolerance
The DiffSolver should exit after the full nonlinear step norm is reduced to either less than absolute...
bool quiet
The DiffSolver should not print anything to libMesh::out unless quiet is set to false; default is tru...
static constexpr Real TOLERANCE
virtual ~NewtonSolver()
Destructor.
virtual std::unique_ptr< NumericVector< T > > zero_clone() const =0
unsigned int max_nonlinear_iterations
The DiffSolver should exit in failure if max_nonlinear_iterations is exceeded and continue_after_max_...
Real absolute_residual_tolerance
The DiffSolver should exit after the residual is reduced to either less than absolute_residual_tolera...
Real max_solution_norm
The largest solution norm which the DiffSolver has yet seen will be stored here, to be used for stopp...
Real max_residual_norm
The largest nonlinear residual which the DiffSolver has yet seen will be stored here, to be used for stopping criteria based on relative_residual_tolerance.
bool _exact_constraint_enforcement
Whether we should enforce exact constraints globally during a solve.
virtual unsigned int solve() override
This method performs a solve, using an inexact Newton-Krylov method with line search.
The DiffSolver achieved the desired absolute step size tolerance.
NumericVector< Number > * rhs
The system matrix.
LinearConvergenceReason
Linear solver convergence flags (taken from the PETSc flags).
virtual void reinit()
The reinitialization function.
ADRealEigenVector< T, D, asd > sqrt(const ADRealEigenVector< T, D, asd > &)
The libMesh namespace provides an interface to certain functionality in the library.
The linear solver used by the DiffSolver failed to find a solution.
const SparseMatrix< Number > * request_matrix(std::string_view mat_name) const
unsigned int _solve_result
Initialized to zero.
unsigned int _inner_iterations
The number of inner iterations used by the last solve.
double linear_tolerance_multiplier
The tolerance for linear solves is kept below this multiplier (which defaults to 1e-3) times the norm...
virtual void zero()=0
Set all entries to zero.
ADRealEigenVector< T, D, asd > abs(const ADRealEigenVector< T, D, asd > &)
virtual void init()
The initialization function.
This is a generic class that defines a solver to handle ImplicitSystem classes, including NonlinearIm...
The DiffSolver reached the maximum allowed number of nonlinear iterations before satisfying any conve...
This base class can be inherited from to provide interfaces to linear solvers from different packages...
The DiffSolver achieved the desired relative step size tolerance.
virtual void assembly(bool, bool, bool=false, bool=false)
Assembles a residual in rhs and/or a jacobian in matrix, as requested.
sys_type & _system
A reference to the system we are solving.
virtual Real l2_norm() const =0
unsigned int max_linear_iterations
Each linear solver step should exit after max_linear_iterations is exceeded.
bool require_finite_residual
If this is set to true, the solver is forced to test the residual after each Newton step...
std::unique_ptr< NumericVector< Number > > solution
Data structure to hold solution values.
The DiffSolver achieved the desired relative residual tolerance.
bool track_linear_convergence
If set to true, check for convergence of the linear solve.
virtual void close()=0
Calls the NumericVector's internal assembly routines, ensuring that the values are consistent across ...
double minimum_linear_tolerance
The tolerance for linear solves is kept above this minimum.
bool verbose
The DiffSolver may print a lot more to libMesh::out if verbose is set to true; default is false...
virtual void update()
Update the local values to reflect the solution on neighboring processors.
virtual void close()=0
Calls the SparseMatrix's internal assembly routines, ensuring that the values are consistent across p...
DIE A HORRIBLE DEATH HERE typedef LIBMESH_DEFAULT_SCALAR_TYPE Real
bool continue_after_backtrack_failure
Defaults to false, telling the DiffSolver to throw an error when the backtracking scheme fails to fin...
The DiffSolver achieved the desired absolute residual tolerance.
Real minsteplength
If the quasi-Newton step length must be reduced to below this factor to give a residual reduction...
virtual void reinit() override
The reinitialization function.
unsigned int _outer_iterations
The number of outer iterations used by the last solve.
SparseMatrix< Number > * matrix
The system matrix.
bool test_convergence(Real current_residual, Real step_norm, bool linear_solve_finished)
bool brent_line_search
If require_residual_reduction is true, the solver may reduce step lengths when required.
The DiffSolver failed to find a descent direction by backtracking (See newton_solver.C)
bool require_residual_reduction
If this is set to true, the solver is forced to test the residual after each Newton step...
bool on_command_line(std::string arg)
const std::string & name() const
A default or invalid solve result.
virtual void add(const numeric_index_type i, const T value)=0
Adds value to the vector entry specified by i.
Real line_search(Real tol, Real last_residual, Real ¤t_residual, NumericVector< Number > &newton_iterate, const NumericVector< Number > &linear_solution)
This does a line search in the direction opposite linear_solution to try and minimize the residual of...
void print_convergence(unsigned int step_num, Real current_residual, Real step_norm, bool linear_solve_finished)
This prints output for the convergence criteria based on by the given residual and step size...
double initial_linear_tolerance
Any required linear solves will at first be done with this tolerance; the DiffSolver may tighten the ...
const DofMap & get_dof_map() const
Real relative_residual_tolerance
Real relative_step_tolerance
std::unique_ptr< LinearSolutionMonitor > linear_solution_monitor
Pointer to functor which is called right after each linear solve.
void enforce_constraints_exactly(const System &system, NumericVector< Number > *v=nullptr, bool homogeneous=false) const
Constrains the numeric vector v, which represents a solution defined on the mesh. ...
Manages consistently variables, degrees of freedom, coefficient vectors, and matrices for implicit sy...