libMesh
newton_solver.C
Go to the documentation of this file.
1 // The libMesh Finite Element Library.
2 // Copyright (C) 2002-2017 Benjamin S. Kirk, John W. Peterson, Roy H. Stogner
3 
4 // This library is free software; you can redistribute it and/or
5 // modify it under the terms of the GNU Lesser General Public
6 // License as published by the Free Software Foundation; either
7 // version 2.1 of the License, or (at your option) any later version.
8 
9 // This library is distributed in the hope that it will be useful,
10 // but WITHOUT ANY WARRANTY; without even the implied warranty of
11 // MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
12 // Lesser General Public License for more details.
13 
14 // You should have received a copy of the GNU Lesser General Public
15 // License along with this library; if not, write to the Free Software
16 // Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA
17 
18 
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"
26 
27 namespace libMesh
28 {
29 
30 // SIGN from Numerical Recipes
31 template <typename T>
32 inline
33 T SIGN(T a, T b)
34 {
35  return b >= 0 ? std::abs(a) : -std::abs(a);
36 }
37 
39  Real last_residual,
40  Real & current_residual,
41  NumericVector<Number> & newton_iterate,
42  const NumericVector<Number> & linear_solution)
43 {
44  // Take a full step if we got a residual reduction or if we
45  // aren't substepping
46  if ((current_residual < last_residual) ||
48  (!require_finite_residual || !libmesh_isnan(current_residual))))
49  return 1.;
50 
51  // The residual vector
53 
54  Real ax = 0.; // First abscissa, don't take negative steps
55  Real cx = 1.; // Second abscissa, don't extrapolate steps
56 
57  // Find bx, a step length that gives lower residual than ax or cx
58  Real bx = 1.;
59 
60  while (libmesh_isnan(current_residual) ||
61  (current_residual > last_residual &&
63  {
64  // Reduce step size to 1/2, 1/4, etc.
65  Real substepdivision;
66  if (brent_line_search && !libmesh_isnan(current_residual))
67  {
68  substepdivision = std::min(0.5, last_residual/current_residual);
69  substepdivision = std::max(substepdivision, tol*2.);
70  }
71  else
72  substepdivision = 0.5;
73 
74  newton_iterate.add (bx * (1.-substepdivision),
75  linear_solution);
76  newton_iterate.close();
77  bx *= substepdivision;
78  if (verbose)
79  libMesh::out << " Shrinking Newton step to "
80  << bx << std::endl;
81 
82  // Check residual with fractional Newton step
83  _system.assembly (true, false);
84 
85  rhs.close();
86  current_residual = rhs.l2_norm();
87  if (verbose)
88  libMesh::out << " Current Residual: "
89  << current_residual << std::endl;
90 
91  if (bx/2. < minsteplength &&
92  (libmesh_isnan(current_residual) ||
93  (current_residual > last_residual)))
94  {
95  libMesh::out << "Inexact Newton step FAILED at step "
96  << _outer_iterations << std::endl;
97 
99  {
100  libmesh_convergence_failure();
101  }
102  else
103  {
104  libMesh::out << "Continuing anyway ..." << std::endl;
106  return bx;
107  }
108  }
109  } // end while (current_residual > last_residual)
110 
111  // Now return that reduced-residual step, or use Brent's method to
112  // find a more optimal step.
113 
114  if (!brent_line_search)
115  return bx;
116 
117  // Brent's method adapted from Numerical Recipes in C, ch. 10.2
118  Real e = 0.;
119 
120  Real x = bx, w = bx, v = bx;
121 
122  // Residuals at bx
123  Real fx = current_residual,
124  fw = current_residual,
125  fv = current_residual;
126 
127  // Max iterations for Brent's method loop
128  const unsigned int max_i = 20;
129 
130  // for golden ratio steps
131  const Real golden_ratio = 1.-(std::sqrt(5.)-1.)/2.;
132 
133  for (unsigned int i=1; i <= max_i; i++)
134  {
135  Real xm = (ax+cx)*0.5;
136  Real tol1 = tol * std::abs(x) + tol*tol;
137  Real tol2 = 2.0 * tol1;
138 
139  // Test if we're done
140  if (std::abs(x-xm) <= (tol2 - 0.5 * (cx - ax)))
141  return x;
142 
143  Real d;
144 
145  // Construct a parabolic fit
146  if (std::abs(e) > tol1)
147  {
148  Real r = (x-w)*(fx-fv);
149  Real q = (x-v)*(fx-fw);
150  Real p = (x-v)*q-(x-w)*r;
151  q = 2. * (q-r);
152  if (q > 0.)
153  p = -p;
154  else
155  q = std::abs(q);
156  if (std::abs(p) >= std::abs(0.5*q*e) ||
157  p <= q * (ax-x) ||
158  p >= q * (cx-x))
159  {
160  // Take a golden section step
161  e = x >= xm ? ax-x : cx-x;
162  d = golden_ratio * e;
163  }
164  else
165  {
166  // Take a parabolic fit step
167  d = p/q;
168  if (x+d-ax < tol2 || cx-(x+d) < tol2)
169  d = SIGN(tol1, xm - x);
170  }
171  }
172  else
173  {
174  // Take a golden section step
175  e = x >= xm ? ax-x : cx-x;
176  d = golden_ratio * e;
177  }
178 
179  Real u = std::abs(d) >= tol1 ? x+d : x + SIGN(tol1,d);
180 
181  // Assemble the residual at the new steplength u
182  newton_iterate.add (bx - u, linear_solution);
183  newton_iterate.close();
184  bx = u;
185  if (verbose)
186  libMesh::out << " Shrinking Newton step to "
187  << bx << std::endl;
188 
189  _system.assembly (true, false);
190 
191  rhs.close();
192  Real fu = current_residual = rhs.l2_norm();
193  if (verbose)
194  libMesh::out << " Current Residual: "
195  << fu << std::endl;
196 
197  if (fu <= fx)
198  {
199  if (u >= x)
200  ax = x;
201  else
202  cx = x;
203  v = w; w = x; x = u;
204  fv = fw; fw = fx; fx = fu;
205  }
206  else
207  {
208  if (u < x)
209  ax = u;
210  else
211  cx = u;
212  if (fu <= fw || w == x)
213  {
214  v = w; w = u;
215  fv = fw; fw = fu;
216  }
217  else if (fu <= fv || v == x || v == w)
218  {
219  v = u;
220  fv = fu;
221  }
222  }
223  }
224 
225  if (!quiet)
226  libMesh::out << "Warning! Too many iterations used in Brent line search!"
227  << std::endl;
228  return bx;
229 }
230 
231 
233  : Parent(s),
236  brent_line_search(true),
238  minsteplength(1e-5),
241 {
242 }
243 
244 
245 
247 {
248 }
249 
250 
251 
253 {
254  Parent::init();
255 
256  if (libMesh::on_command_line("--solver_system_names"))
257  _linear_solver->init((_system.name()+"_").c_str());
258  else
259  _linear_solver->init();
260 
261  _linear_solver->init_names(_system);
262 }
263 
264 
265 
267 {
268  Parent::reinit();
269 
270  _linear_solver->clear();
271 
272  _linear_solver->init_names(_system);
273 }
274 
275 
276 
277 unsigned int NewtonSolver::solve()
278 {
279  LOG_SCOPE("solve()", "NewtonSolver");
280 
281  // Reset any prior solve result
283 
284  NumericVector<Number> & newton_iterate = *(_system.solution);
285 
286  UniquePtr<NumericVector<Number>> linear_solution_ptr = newton_iterate.zero_clone();
287  NumericVector<Number> & linear_solution = *linear_solution_ptr;
288  NumericVector<Number> & rhs = *(_system.rhs);
289 
290  newton_iterate.close();
291  linear_solution.close();
292  rhs.close();
293 
294 #ifdef LIBMESH_ENABLE_CONSTRAINTS
296 #endif
297 
298  SparseMatrix<Number> & matrix = *(_system.matrix);
299 
300  // Set starting linear tolerance
301  Real current_linear_tolerance = initial_linear_tolerance;
302 
303  // Start counting our linear solver steps
304  _inner_iterations = 0;
305 
306  // Now we begin the nonlinear loop
309  {
310  if (verbose)
311  libMesh::out << "Assembling the System" << std::endl;
312 
313  _system.assembly(true, true);
314  rhs.close();
315  Real current_residual = rhs.l2_norm();
316 
317  if (libmesh_isnan(current_residual))
318  {
319  libMesh::out << " Nonlinear solver DIVERGED at step "
321  << " with residual Not-a-Number"
322  << std::endl;
323  libmesh_convergence_failure();
324  continue;
325  }
326 
327  if (current_residual <= absolute_residual_tolerance)
328  {
329  if (verbose)
330  libMesh::out << "Linear solve unnecessary; residual "
331  << current_residual
332  << " meets absolute tolerance "
334  << std::endl;
335 
336  // We're not doing a solve, but other code may reuse this
337  // matrix.
338  matrix.close();
339 
341  if (current_residual == 0)
342  {
345  if (absolute_step_tolerance > 0)
347  if (relative_step_tolerance > 0)
349  }
350 
351  break;
352  }
353 
354  // Prepare to take incomplete steps
355  Real last_residual = current_residual;
356 
357  max_residual_norm = std::max (current_residual,
359 
360  // Compute the l2 norm of the whole solution
361  Real norm_total = newton_iterate.l2_norm();
362 
364 
365  if (verbose)
366  libMesh::out << "Nonlinear Residual: "
367  << current_residual << std::endl;
368 
369  // Make sure our linear tolerance is low enough
370  current_linear_tolerance = std::min (current_linear_tolerance,
371  current_residual * linear_tolerance_multiplier);
372 
373  // But don't let it be too small
374  if (current_linear_tolerance < minimum_linear_tolerance)
375  {
376  current_linear_tolerance = minimum_linear_tolerance;
377  }
378 
379  // If starting the nonlinear solve with a really good initial guess, we dont want to set an absurd linear tolerance
380  current_linear_tolerance = std::max(current_linear_tolerance, absolute_residual_tolerance / current_residual / 10.0);
381 
382  // At this point newton_iterate is the current guess, and
383  // linear_solution is now about to become the NEGATIVE of the next
384  // Newton step.
385 
386  // Our best initial guess for the linear_solution is zero!
387  linear_solution.zero();
388 
389  if (verbose)
390  libMesh::out << "Linear solve starting, tolerance "
391  << current_linear_tolerance << std::endl;
392 
393  // Solve the linear system.
394  const std::pair<unsigned int, Real> rval =
395  _linear_solver->solve (matrix, _system.request_matrix("Preconditioner"),
396  linear_solution, rhs, current_linear_tolerance,
398 
400  {
401  LinearConvergenceReason linear_c_reason = _linear_solver->get_converged_reason();
402 
403  // Check if something went wrong during the linear solve
404  if (linear_c_reason < 0)
405  {
406  // The linear solver failed somehow
408  // Print a message
409  libMesh::out << "Linear solver failed during Newton step, dropping out."
410  << std::endl;
411  break;
412  }
413  }
414 
415  // We may need to localize a parallel solution
416  _system.update ();
417  // The linear solver may not have fit our constraints exactly
418 #ifdef LIBMESH_ENABLE_CONSTRAINTS
420  (_system, &linear_solution, /* homogeneous = */ true);
421 #endif
422 
423  const unsigned int linear_steps = rval.first;
424  libmesh_assert_less_equal (linear_steps, max_linear_iterations);
425  _inner_iterations += linear_steps;
426 
427  const bool linear_solve_finished =
428  !(linear_steps == max_linear_iterations);
429 
430  if (verbose)
431  libMesh::out << "Linear solve finished, step " << linear_steps
432  << ", residual " << rval.second
433  << std::endl;
434 
435  // Compute the l2 norm of the nonlinear update
436  Real norm_delta = linear_solution.l2_norm();
437 
438  if (verbose)
439  libMesh::out << "Trying full Newton step" << std::endl;
440  // Take a full Newton step
441  newton_iterate.add (-1., linear_solution);
442  newton_iterate.close();
443 
444  if (this->linear_solution_monitor.get())
445  {
446  // Compute the l2 norm of the whole solution
447  norm_total = newton_iterate.l2_norm();
448  rhs.close();
449  (*this->linear_solution_monitor)(linear_solution, norm_delta,
450  newton_iterate, norm_total,
451  rhs, rhs.l2_norm(), _outer_iterations);
452  }
453 
454  // Check residual with full Newton step, if that's useful for determining
455  // whether to line search, whether to quit early, or whether to die after
456  // hitting our max iteration count
457  if (this->require_residual_reduction ||
458  this->require_finite_residual ||
459  _outer_iterations+1 < max_nonlinear_iterations ||
461  {
462  _system.assembly(true, false);
463 
464  rhs.close();
465  current_residual = rhs.l2_norm();
466  if (verbose)
467  libMesh::out << " Current Residual: "
468  << current_residual << std::endl;
469 
470  // don't fiddle around if we've already converged
471  if (test_convergence(current_residual, norm_delta,
472  linear_solve_finished &&
473  current_residual <= last_residual))
474  {
475  if (!quiet)
476  print_convergence(_outer_iterations, current_residual,
477  norm_delta, linear_solve_finished &&
478  current_residual <= last_residual);
480  break; // out of _outer_iterations for loop
481  }
482  }
483 
484  // since we're not converged, backtrack if necessary
485  Real steplength =
486  this->line_search(std::sqrt(TOLERANCE),
487  last_residual, current_residual,
488  newton_iterate, linear_solution);
489  norm_delta *= steplength;
490 
491  // Check to see if backtracking failed,
492  // and break out of the nonlinear loop if so...
494  {
496  break; // out of _outer_iterations for loop
497  }
498 
499  if (_outer_iterations + 1 >= max_nonlinear_iterations)
500  {
501  libMesh::out << " Nonlinear solver reached maximum step "
502  << max_nonlinear_iterations << ", latest evaluated residual "
503  << current_residual << std::endl;
505  {
507  libMesh::out << " Continuing..." << std::endl;
508  }
509  else
510  {
511  libmesh_convergence_failure();
512  }
513  continue;
514  }
515 
516  // Compute the l2 norm of the whole solution
517  norm_total = newton_iterate.l2_norm();
518 
520 
521  // Print out information for the
522  // nonlinear iterations.
523  if (verbose)
524  libMesh::out << " Nonlinear step: |du|/|u| = "
525  << norm_delta / norm_total
526  << ", |du| = " << norm_delta
527  << std::endl;
528 
529  // Terminate the solution iteration if the difference between
530  // this iteration and the last is sufficiently small.
531  if (test_convergence(current_residual, norm_delta / steplength,
532  linear_solve_finished))
533  {
534  if (!quiet)
535  print_convergence(_outer_iterations, current_residual,
536  norm_delta / steplength,
537  linear_solve_finished);
539  break; // out of _outer_iterations for loop
540  }
541  } // end nonlinear loop
542 
543  // The linear solver may not have fit our constraints exactly
544 #ifdef LIBMESH_ENABLE_CONSTRAINTS
546 #endif
547 
548  // We may need to localize a parallel solution
549  _system.update ();
550 
551  // Make sure we are returning something sensible as the
552  // _solve_result, except in the edge case where we weren't really asked to
553  // solve.
555  !max_nonlinear_iterations);
556 
557  return _solve_result;
558 }
559 
560 
561 
562 bool NewtonSolver::test_convergence(Real current_residual,
563  Real step_norm,
564  bool linear_solve_finished)
565 {
566  // We haven't converged unless we pass a convergence test
567  bool has_converged = false;
568 
569  // Is our absolute residual low enough?
570  if (current_residual < absolute_residual_tolerance)
571  {
573  has_converged = true;
574  }
575 
576  // Is our relative residual low enough?
577  if ((current_residual / max_residual_norm) <
579  {
581  has_converged = true;
582  }
583 
584  // For incomplete linear solves, it's not safe to test step sizes
585  if (!linear_solve_finished)
586  {
587  return has_converged;
588  }
589 
590  // Is our absolute Newton step size small enough?
591  if (step_norm < absolute_step_tolerance)
592  {
594  has_converged = true;
595  }
596 
597  // Is our relative Newton step size small enough?
598  if (step_norm / max_solution_norm <
600  {
602  has_converged = true;
603  }
604 
605  return has_converged;
606 }
607 
608 
609 void NewtonSolver::print_convergence(unsigned int step_num,
610  Real current_residual,
611  Real step_norm,
612  bool linear_solve_finished)
613 {
614  // Is our absolute residual low enough?
615  if (current_residual < absolute_residual_tolerance)
616  {
617  libMesh::out << " Nonlinear solver converged, step " << step_num
618  << ", residual " << current_residual
619  << std::endl;
620  }
622  {
623  if (verbose)
624  libMesh::out << " Nonlinear solver current_residual "
625  << current_residual << " > "
626  << (absolute_residual_tolerance) << std::endl;
627  }
628 
629  // Is our relative residual low enough?
630  if ((current_residual / max_residual_norm) <
632  {
633  libMesh::out << " Nonlinear solver converged, step " << step_num
634  << ", residual reduction "
635  << current_residual / max_residual_norm
636  << " < " << relative_residual_tolerance
637  << std::endl;
638  }
640  {
641  if (verbose)
642  libMesh::out << " Nonlinear solver relative residual "
643  << (current_residual / max_residual_norm)
644  << " > " << relative_residual_tolerance
645  << std::endl;
646  }
647 
648  // For incomplete linear solves, it's not safe to test step sizes
649  if (!linear_solve_finished)
650  return;
651 
652  // Is our absolute Newton step size small enough?
653  if (step_norm < absolute_step_tolerance)
654  {
655  libMesh::out << " Nonlinear solver converged, step " << step_num
656  << ", absolute step size "
657  << step_norm
658  << " < " << absolute_step_tolerance
659  << std::endl;
660  }
661  else if (absolute_step_tolerance)
662  {
663  if (verbose)
664  libMesh::out << " Nonlinear solver absolute step size "
665  << step_norm
666  << " > " << absolute_step_tolerance
667  << std::endl;
668  }
669 
670  // Is our relative Newton step size small enough?
671  if (step_norm / max_solution_norm <
673  {
674  libMesh::out << " Nonlinear solver converged, step " << step_num
675  << ", relative step size "
676  << (step_norm / max_solution_norm)
677  << " < " << relative_step_tolerance
678  << std::endl;
679  }
680  else if (relative_step_tolerance)
681  {
682  if (verbose)
683  libMesh::out << " Nonlinear solver relative step size "
684  << (step_norm / max_solution_norm)
685  << " > " << relative_step_tolerance
686  << std::endl;
687  }
688 }
689 
690 } // namespace libMesh
bool continue_after_max_iterations
Defaults to true, telling the DiffSolver to continue rather than exit when a solve has reached its ma...
Definition: diff_solver.h:174
NewtonSolver(sys_type &system)
Constructor.
double abs(double a)
Real minimum_linear_tolerance
The tolerance for linear solves is kept above this minimum.
Definition: diff_solver.h:215
Real absolute_step_tolerance
The DiffSolver should exit after the full nonlinear step norm is reduced to either less than absolute...
Definition: diff_solver.h:203
bool quiet
The DiffSolver should not print anything to libMesh::out unless quiet is set to false; default is tru...
Definition: diff_solver.h:162
virtual Real l2_norm() const =0
virtual void init() libmesh_override
The initialization function.
virtual ~NewtonSolver()
Destructor.
unsigned int max_nonlinear_iterations
The DiffSolver should exit in failure if max_nonlinear_iterations is exceeded and continue_after_max_...
Definition: diff_solver.h:156
Real absolute_residual_tolerance
The DiffSolver should exit after the residual is reduced to either less than absolute_residual_tolera...
Definition: diff_solver.h:191
static UniquePtr< DiffSolver > build(sys_type &s)
Factory method.
Definition: diff_solver.C:53
Real max_solution_norm
The largest solution norm which the DiffSolver has yet seen will be stored here, to be used for stopp...
Definition: diff_solver.h:296
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.
Definition: diff_solver.h:303
The DiffSolver achieved the desired absolute step size tolerance.
Definition: diff_solver.h:251
NumericVector< Number > * rhs
The system matrix.
virtual void reinit()
The reinitialization function.
Definition: diff_solver.C:60
static const Real TOLERANCE
The libMesh namespace provides an interface to certain functionality in the library.
Real initial_linear_tolerance
Any required linear solves will at first be done with this tolerance; the DiffSolver may tighten the ...
Definition: diff_solver.h:210
The linear solver used by the DiffSolver failed to find a solution.
Definition: diff_solver.h:282
long double max(long double a, double b)
unsigned int _solve_result
Initialized to zero.
Definition: diff_solver.h:326
const std::string & name() const
Definition: system.h:1998
unsigned int _inner_iterations
The number of inner iterations used by the last solve.
Definition: diff_solver.h:313
virtual void zero()=0
Set all entries to zero.
libmesh_assert(j)
std::unique_ptr< T > UniquePtr
Definition: auto_ptr.h:46
virtual unsigned int solve() libmesh_override
This method performs a solve, using an inexact Newton-Krylov method with line search.
virtual void reinit() libmesh_override
The reinitialization function.
virtual void init()
The initialization function.
Definition: diff_solver.C:69
This is a generic class that defines a solver to handle ImplicitSystem classes, including NonlinearIm...
Definition: diff_solver.h:70
The DiffSolver reached the maximum allowed number of nonlinear iterations before satisfying any conve...
Definition: diff_solver.h:270
This base class can be inherited from to provide interfaces to linear solvers from different packages...
Definition: linear_solver.h:58
The DiffSolver achieved the desired relative step size tolerance.
Definition: diff_solver.h:257
T SIGN(T a, T b)
Definition: newton_solver.C:33
virtual void assembly(bool, bool, bool=false, bool=false)
Assembles a residual in rhs and/or a jacobian in matrix, as requested.
PetscErrorCode Vec x
const DofMap & get_dof_map() const
Definition: system.h:2030
sys_type & _system
A reference to the system we are solving.
Definition: diff_solver.h:318
unsigned int max_linear_iterations
Each linear solver step should exit after max_linear_iterations is exceeded.
Definition: diff_solver.h:148
bool require_finite_residual
If this is set to true, the solver is forced to test the residual after each Newton step...
The DiffSolver achieved the desired relative residual tolerance.
Definition: diff_solver.h:245
bool track_linear_convergence
If set to true, check for convergence of the linear solve.
UniquePtr< NumericVector< Number > > solution
Data structure to hold solution values.
Definition: system.h:1523
void enforce_constraints_exactly(const System &system, NumericVector< Number > *v=libmesh_nullptr, bool homogeneous=false) const
Constrains the numeric vector v, which represents a solution defined on the mesh. ...
Definition: dof_map.h:1816
Real linear_tolerance_multiplier
The tolerance for linear solves is kept below this multiplier (which defaults to 1e-3) times the norm...
virtual void close()=0
Calls the NumericVector&#39;s internal assembly routines, ensuring that the values are consistent across ...
const SparseMatrix< Number > * request_matrix(const std::string &mat_name) const
bool verbose
The DiffSolver may print a lot more to libMesh::out if verbose is set to true; default is false...
Definition: diff_solver.h:168
virtual void update()
Update the local values to reflect the solution on neighboring processors.
Definition: system.C:425
virtual void close()=0
Calls the SparseMatrix&#39;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...
Definition: diff_solver.h:180
The DiffSolver achieved the desired absolute residual tolerance.
Definition: diff_solver.h:239
bool libmesh_isnan(float a)
Real minsteplength
If the quasi-Newton step length must be reduced to below this factor to give a residual reduction...
bool on_command_line(const std::string &arg)
Definition: libmesh.C:921
virtual UniquePtr< NumericVector< T > > zero_clone() const =0
unsigned int _outer_iterations
The number of outer iterations used by the last solve.
Definition: diff_solver.h:308
const Parallel::Communicator & comm() const
OStreamProxy out
SparseMatrix< Number > * matrix
The system matrix.
UniquePtr< LinearSolver< Number > > _linear_solver
The LinearSolver defines the interface used to solve the linear_implicit system.
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.
LinearConvergenceReason
Linear solver convergence flags (taken from the PETSc flags)
UniquePtr< LinearSolutionMonitor > linear_solution_monitor
Pointer to functor which is called right after each linear solve.
Definition: diff_solver.h:288
The DiffSolver failed to find a descent direction by backtracking (See newton_solver.C)
Definition: diff_solver.h:276
bool require_residual_reduction
If this is set to true, the solver is forced to test the residual after each Newton step...
Definition: newton_solver.h:98
A default or invalid solve result.
Definition: diff_solver.h:227
virtual void add(const numeric_index_type i, const T value)=0
Adds value to each entry of the vector.
Real line_search(Real tol, Real last_residual, Real &current_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...
Definition: newton_solver.C:38
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...
long double min(long double a, double b)
Real relative_residual_tolerance
Definition: diff_solver.h:192
This class provides a specific system class.