OGS
NumLib::NonlinearSolver< NonlinearSolverTag::Newton > Class Referencefinal

## Detailed Description

Find a solution to a nonlinear equation using the Newton-Raphson method.

Definition at line 75 of file NonlinearSolver.h.

#include <NonlinearSolver.h>

Inheritance diagram for NumLib::NonlinearSolver< NonlinearSolverTag::Newton >:
[legend]
Collaboration diagram for NumLib::NonlinearSolver< NonlinearSolverTag::Newton >:
[legend]

## Public Types

using System = NonlinearSystem< NonlinearSolverTag::Newton >
Type of the nonlinear equation system to be solved. More...

## Public Member Functions

NonlinearSolver (GlobalLinearSolver &linear_solver, int const maxiter, double const damping=1.0)

~NonlinearSolver ()

void setEquationSystem (System &eq, ConvergenceCriterion &conv_crit)

void calculateNonEquilibriumInitialResiduum (std::vector< GlobalVector * > const &x, std::vector< GlobalVector * > const &x_prev, int const process_id) override

NonlinearSolverStatus solve (std::vector< GlobalVector * > &x, std::vector< GlobalVector * > const &x_prev, std::function< void(int, std::vector< GlobalVector * > const &)> const &postIterationCallback, int const process_id) override

void compensateNonEquilibriumInitialResiduum (bool const value)

Public Member Functions inherited from NumLib::NonlinearSolverBase
virtual ~NonlinearSolverBase ()=default

## Private Attributes

GlobalLinearSolver_linear_solver

System_equation_system = nullptr

ConvergenceCriterion_convergence_criterion = nullptr

int const _maxiter
maximum number of iterations More...

double const _damping

GlobalVector_r_neq = nullptr
non-equilibrium initial residuum. More...

std::size_t _res_id = 0u
ID of the residual vector. More...

std::size_t _J_id = 0u
ID of the Jacobian matrix. More...

std::size_t _minus_delta_x_id = 0u
ID of the $$-\Delta x$$ vector. More...

std::size_t _x_new_id
ID of the vector storing $$x - (-\Delta x)$$. More...

std::size_t _r_neq_id = 0u

bool _compensate_non_equilibrium_initial_residuum = false

## ◆ System

Type of the nonlinear equation system to be solved.

Definition at line 80 of file NonlinearSolver.h.

## ◆ NonlinearSolver()

 NumLib::NonlinearSolver< NonlinearSolverTag::Newton >::NonlinearSolver ( GlobalLinearSolver & linear_solver, int const maxiter, double const damping = 1.0 )
inlineexplicit

Constructs a new instance.

Parameters
 linear_solver the linear solver used by this nonlinear solver. maxiter the maximum number of iterations used to solve the equation. damping A positive damping factor.
_damping

Definition at line 90 of file NonlinearSolver.h.

93  : _linear_solver(linear_solver), _maxiter(maxiter), _damping(damping)
94  {
95  }
int const _maxiter
maximum number of iterations

## ◆ ~NonlinearSolver()

 NumLib::NonlinearSolver< NonlinearSolverTag::Newton >::~NonlinearSolver ( )

Definition at line 465 of file NonlinearSolver.cpp.

466 {
467  if (_r_neq != nullptr)
468  {
470  }
471 }
GlobalVector * _r_neq
non-equilibrium initial residuum.
virtual void releaseVector(GlobalVector const &x)=0
static NUMLIB_EXPORT VectorProvider & provider

## ◆ calculateNonEquilibriumInitialResiduum()

 void NumLib::NonlinearSolver< NonlinearSolverTag::Newton >::calculateNonEquilibriumInitialResiduum ( std::vector< GlobalVector * > const & x, std::vector< GlobalVector * > const & x_prev, int const process_id )
overridevirtual

Implements NumLib::NonlinearSolverBase.

Definition at line 214 of file NonlinearSolver.cpp.

218 {
220  {
221  return;
222  }
223
224  _equation_system->assemble(x, x_prev, process_id);
226  _equation_system->getResidual(*x[process_id], *x_prev[process_id], *_r_neq);
227 }
virtual void assemble(std::vector< GlobalVector * > const &x, std::vector< GlobalVector * > const &x_prev, int const process_id)=0
virtual void getResidual(GlobalVector const &x, GlobalVector const &x_prev, GlobalVector &res) const =0
virtual GlobalVector & getVector(std::size_t &id)=0
Get an uninitialized vector with the given id.

## ◆ compensateNonEquilibriumInitialResiduum()

 void NumLib::NonlinearSolver< NonlinearSolverTag::Newton >::compensateNonEquilibriumInitialResiduum ( bool const value )
inline

Definition at line 119 of file NonlinearSolver.h.

120  {
122  }

## ◆ setEquationSystem()

 void NumLib::NonlinearSolver< NonlinearSolverTag::Newton >::setEquationSystem ( System & eq, ConvergenceCriterion & conv_crit )
inline

Set the nonlinear equation system that will be solved. TODO doc

Definition at line 101 of file NonlinearSolver.h.

102  {
103  _equation_system = &eq;
104  _convergence_criterion = &conv_crit;
105  }

## ◆ solve()

 NonlinearSolverStatus NumLib::NonlinearSolver< NonlinearSolverTag::Newton >::solve ( std::vector< GlobalVector * > & x, std::vector< GlobalVector * > const & x_prev, std::function< void(int, std::vector< GlobalVector * > const &)> const & postIterationCallback, int const process_id )
overridevirtual

Assemble and solve the equation system.

Parameters
 x in: the initial guess, out: the solution. x_prev previous time step solution. postIterationCallback called after each iteration if set. process_id usually used in staggered schemes.
Return values
 true if the equation system could be solved false otherwise

Implements NumLib::NonlinearSolverBase.

Definition at line 229 of file NonlinearSolver.cpp.

235 {
236  namespace LinAlg = MathLib::LinAlg;
237  auto& sys = *_equation_system;
238
240  auto& minus_delta_x =
243
244  bool error_norms_met = false;
245
246  // TODO be more efficient
247  // init minus_delta_x to the right size
248  LinAlg::copy(*x[process_id], minus_delta_x);
249
251
252  int iteration = 1;
253  for (; iteration <= _maxiter; ++iteration, _convergence_criterion->reset())
254  {
255  BaseLib::RunTime timer_dirichlet;
256  double time_dirichlet = 0.0;
257
258  BaseLib::RunTime time_iteration;
259  time_iteration.start();
260
261  timer_dirichlet.start();
262  sys.computeKnownSolutions(*x[process_id], process_id);
263  sys.applyKnownSolutions(*x[process_id]);
264  time_dirichlet += timer_dirichlet.elapsed();
265
266  sys.preIteration(iteration, *x[process_id]);
267
268  BaseLib::RunTime time_assembly;
269  time_assembly.start();
270  try
271  {
272  sys.assemble(x, x_prev, process_id);
273  }
274  catch (AssemblyException const& e)
275  {
276  ERR("Abort nonlinear iteration. Repeating timestep. Reason: {:s}",
277  e.what());
278  error_norms_met = false;
279  iteration = _maxiter;
280  break;
281  }
282  sys.getResidual(*x[process_id], *x_prev[process_id], res);
283  sys.getJacobian(J);
284  INFO("[time] Assembly took {:g} s.", time_assembly.elapsed());
285
286  // Subtract non-equilibrium initial residuum if set
287  if (_r_neq != nullptr)
288  LinAlg::axpy(res, -1, *_r_neq);
289
290  minus_delta_x.setZero();
291
292  timer_dirichlet.start();
293  sys.applyKnownSolutionsNewton(J, res, minus_delta_x);
294  time_dirichlet += timer_dirichlet.elapsed();
295  INFO("[time] Applying Dirichlet BCs took {:g} s.", time_dirichlet);
296
297  if (!sys.isLinear() && _convergence_criterion->hasResidualCheck())
298  {
300  }
301
302  BaseLib::RunTime time_linear_solver;
303  time_linear_solver.start();
304  bool iteration_succeeded = _linear_solver.solve(J, res, minus_delta_x);
305  INFO("[time] Linear solver took {:g} s.", time_linear_solver.elapsed());
306
307  if (!iteration_succeeded)
308  {
309  ERR("Newton: The linear solver failed.");
310  }
311  else
312  {
313  // TODO could be solved in a better way
314  // cf.
315  // https://www.petsc.org/release/docs/manualpages/Vec/VecWAXPY.html
316
317  // Copy pointers, replace the one for the given process id.
318  std::vector<GlobalVector*> x_new{x};
319  x_new[process_id] =
321  *x[process_id], _x_new_id);
322  LinAlg::axpy(*x_new[process_id], -_damping, minus_delta_x);
323
324  if (postIterationCallback)
325  {
326  postIterationCallback(iteration, x_new);
327  }
328
329  switch (sys.postIteration(*x_new[process_id]))
330  {
332  break;
334  ERR("Newton: The postIteration() hook reported a "
335  "non-recoverable error.");
336  iteration_succeeded = false;
337  break;
339  INFO(
340  "Newton: The postIteration() hook decided that this "
341  "iteration has to be repeated.");
342  // TODO introduce some onDestroy hook.
344  *x_new[process_id]);
345  continue; // That throws the iteration result away.
346  }
347
348  LinAlg::copy(*x_new[process_id],
349  *x[process_id]); // copy new solution to x
351  *x_new[process_id]);
352  }
353
354  if (!iteration_succeeded)
355  {
356  // Don't compute further error norms, but break here.
357  error_norms_met = false;
358  break;
359  }
360
361  if (sys.isLinear())
362  {
363  error_norms_met = true;
364  }
365  else
366  {
368  {
369  // Note: x contains the new solution!
370  _convergence_criterion->checkDeltaX(minus_delta_x,
371  *x[process_id]);
372  }
373
374  error_norms_met = _convergence_criterion->isSatisfied();
375  }
376
377  INFO("[time] Iteration #{:d} took {:g} s.", iteration,
378  time_iteration.elapsed());
379
380  if (error_norms_met)
381  {
382  break;
383  }
384
385  // Avoid increment of the 'iteration' if the error norms are not met,
386  // but maximum number of iterations is reached.
387  if (iteration >= _maxiter)
388  {
389  break;
390  }
391  }
392
393  if (iteration > _maxiter)
394  {
395  ERR("Newton: Could not solve the given nonlinear system within {:d} "
396  "iterations",
397  _maxiter);
398  }
399
403
404  return {error_norms_met, iteration};
405 }
void INFO(char const *fmt, Args const &... args)
Definition: Logging.h:32
void ERR(char const *fmt, Args const &... args)
Definition: Logging.h:42
Count the running time.
Definition: RunTime.h:29
double elapsed() const
Get the elapsed time in seconds.
Definition: RunTime.h:42
void start()
Start the timer.
Definition: RunTime.h:32
bool solve(EigenMatrix &A, EigenVector &b, EigenVector &x)
virtual void checkResidual(GlobalVector const &residual)=0
Check if the residual satisfies the convergence criterion.
virtual bool hasResidualCheck() const =0
virtual bool isSatisfied() const
Tell if the convergence criterion is satisfied.
virtual void checkDeltaX(GlobalVector const &minus_delta_x, GlobalVector const &x)=0
virtual bool hasDeltaXCheck() const =0
virtual GlobalMatrix & getMatrix(std::size_t &id)=0
Get an uninitialized matrix with the given id.
virtual void releaseMatrix(GlobalMatrix const &A)=0
std::size_t _J_id
ID of the Jacobian matrix.
std::size_t _x_new_id
ID of the vector storing .
std::size_t _res_id
ID of the residual vector.
void axpy(PETScVector &y, double const a, PETScVector const &x)
Definition: LinAlg.cpp:57
void copy(PETScVector const &x, PETScVector &y)
Definition: LinAlg.cpp:37
static NUMLIB_EXPORT MatrixProvider & provider

## ◆ _compensate_non_equilibrium_initial_residuum

 bool NumLib::NonlinearSolver< NonlinearSolverTag::Newton >::_compensate_non_equilibrium_initial_residuum = false
private

Enables computation of the non-equilibrium initial residuum $$r_{\rm neq}$$ before the first time step. The forces are zero if the external forces are in equilibrium with the initial state/initial conditions. During the simulation the new residuum reads $$\tilde r = r - r_{\rm neq}$$.

Definition at line 152 of file NonlinearSolver.h.

## ◆ _convergence_criterion

 ConvergenceCriterion* NumLib::NonlinearSolver< NonlinearSolverTag::Newton >::_convergence_criterion = nullptr
private

Definition at line 129 of file NonlinearSolver.h.

## ◆ _damping

 double const NumLib::NonlinearSolver< NonlinearSolverTag::Newton >::_damping
private

A positive damping factor. The default value 1.0 gives a non-damped Newton method. Common values are in the range 0.5 to 0.7 for somewhat conservative method and seldom become smaller than 0.2 for very conservative approach.

Definition at line 136 of file NonlinearSolver.h.

## ◆ _equation_system

 System* NumLib::NonlinearSolver< NonlinearSolverTag::Newton >::_equation_system = nullptr
private

Definition at line 126 of file NonlinearSolver.h.

## ◆ _J_id

 std::size_t NumLib::NonlinearSolver< NonlinearSolverTag::Newton >::_J_id = 0u
private

ID of the Jacobian matrix.

Definition at line 140 of file NonlinearSolver.h.

## ◆ _linear_solver

 GlobalLinearSolver& NumLib::NonlinearSolver< NonlinearSolverTag::Newton >::_linear_solver
private

Definition at line 125 of file NonlinearSolver.h.

## ◆ _maxiter

 int const NumLib::NonlinearSolver< NonlinearSolverTag::Newton >::_maxiter
private

maximum number of iterations

Definition at line 130 of file NonlinearSolver.h.

## ◆ _minus_delta_x_id

 std::size_t NumLib::NonlinearSolver< NonlinearSolverTag::Newton >::_minus_delta_x_id = 0u
private

ID of the $$-\Delta x$$ vector.

Definition at line 141 of file NonlinearSolver.h.

## ◆ _r_neq

 GlobalVector* NumLib::NonlinearSolver< NonlinearSolverTag::Newton >::_r_neq = nullptr
private

non-equilibrium initial residuum.

Definition at line 138 of file NonlinearSolver.h.

## ◆ _r_neq_id

 std::size_t NumLib::NonlinearSolver< NonlinearSolverTag::Newton >::_r_neq_id = 0u
private

ID of the non-equilibrium initial residuum vector.

Definition at line 144 of file NonlinearSolver.h.

## ◆ _res_id

 std::size_t NumLib::NonlinearSolver< NonlinearSolverTag::Newton >::_res_id = 0u
private

ID of the residual vector.

Definition at line 139 of file NonlinearSolver.h.

## ◆ _x_new_id

 std::size_t NumLib::NonlinearSolver< NonlinearSolverTag::Newton >::_x_new_id
private
Initial value:
=
0u

ID of the vector storing $$x - (-\Delta x)$$.

Definition at line 142 of file NonlinearSolver.h.

The documentation for this class was generated from the following files: