LightKrylov_NewtonKrylov Module



Interfaces

public interface newton

Implements the simple Newton-Krylov method for finding roots (fixed points) of a nonlinear vector-valued function , i.e. solutions such that starting from an initial guess via successive solution increments based on local linearization (the Jacobian) of the nonlinear function in the vicinity of the current solution.

Algorthmic Features

  • At iteration , the standard Newton step is computed as the solution of the linear system

where is the residual of the nonlinear function. The new guess for the fixed point is then given by:

where parametrizes the step length. The standard Newton algorithm sets .

  • The Jacobian is never assembled and the linear system is solved using one of the available iterative solvers.
  • When the residual norm does not decrease during iteration indicating that the linearization is not a very accurate model of the function's behaviour, which often happens during the initial iterations, a 1D step bisection method based on the golden ratio is implemented to dampen the step and improve convergence of the method.
  • The implementation allows for dynamic tolerances (also known as inexact Newton), where the approximation for the residual and the linear system can be solved with relaxed tolerances to reduce overall time to solution.
  • The method is suitable to both fixed points and periodic orbits via the choice of residual and corresponding Jacobian matrix. In the case of unforced periodic orbits, the period is itself an unknown that must be included in the iteration.

Advantages

  • The iterative solution of the linear systems has a comparatively low storage footprint.
  • If the Newton iteration converges, the convergence is formally asymptotically of second order. Using dynamic tolerances and line searches slightly reduce this convergence rate in exchange for a larger convergence region.

Limitations

  • The method is not guaranteed to converge if the initial guess is too far from the fixed point. If the Newton iteration diverges even with step bisection, the best suggestion is to find a better initial guess. If this is not feasible, some alternatives to improve the convergence of the Newton iteration are possible (but not implemented to date), including various line search algorithms and trust region methods (doglog, double dogleg, hookstep, ...).

References

  • Sánchez, J., Net, M., Garcıa-Archilla, B., & Simó, C. (2004). "Newton–Krylov continuation of periodic orbits for Navier–Stokes flows". Journal of Computational Physics, 201(1), 13-33.
  • Viswanath, D. (2007). "Recurrent motions within plane Couette turbulence". Journal of Fluid Mechanics, 580, 339-358.
  • Duguet, Y., Pringle, C. C. T., Kerswell, R. R. (2008). "Relative periodic orbits in transitional pipe flow". Physics of Fluids, 20(11), 114102.
  • Frantz, R. A., Loiseau, J. C., & Robinet, J. C. (2023). "Krylov methods for large-scale dynamical systems: Application in fluid dynamics". Applied Mechanics Reviews, 75(3), 030802.
  • private subroutine newton_rsp(sys, X, solver, info, tolerance, options, linear_solver_options, preconditioner, scheduler)

    Arguments

    Type IntentOptional Attributes Name
    class(abstract_system_rsp), intent(inout) :: sys

    Dynamical system for which we wish to compute a fixed point

    class(abstract_vector_rsp), intent(inout) :: X

    Initial guess for the fixed point, will be overwritten with solution

    procedure(abstract_linear_solver_rsp) :: solver

    Linear solver to be used to find Newton step

    integer, intent(out) :: info

    Information flag

    real(kind=sp), intent(in), optional :: tolerance
    type(newton_sp_opts), intent(in), optional :: options
    class(abstract_opts), intent(in), optional :: linear_solver_options

    Options for the linear solver

    class(abstract_precond_rsp), intent(in), optional :: preconditioner

    Preconditioner for the linear solver

    procedure(abstract_scheduler_sp), optional :: scheduler
  • private subroutine newton_rdp(sys, X, solver, info, tolerance, options, linear_solver_options, preconditioner, scheduler)

    Arguments

    Type IntentOptional Attributes Name
    class(abstract_system_rdp), intent(inout) :: sys

    Dynamical system for which we wish to compute a fixed point

    class(abstract_vector_rdp), intent(inout) :: X

    Initial guess for the fixed point, will be overwritten with solution

    procedure(abstract_linear_solver_rdp) :: solver

    Linear solver to be used to find Newton step

    integer, intent(out) :: info

    Information flag

    real(kind=dp), intent(in), optional :: tolerance
    type(newton_dp_opts), intent(in), optional :: options
    class(abstract_opts), intent(in), optional :: linear_solver_options

    Options for the linear solver

    class(abstract_precond_rdp), intent(in), optional :: preconditioner

    Preconditioner for the linear solver

    procedure(abstract_scheduler_dp), optional :: scheduler
  • private subroutine newton_csp(sys, X, solver, info, tolerance, options, linear_solver_options, preconditioner, scheduler)

    Arguments

    Type IntentOptional Attributes Name
    class(abstract_system_csp), intent(inout) :: sys

    Dynamical system for which we wish to compute a fixed point

    class(abstract_vector_csp), intent(inout) :: X

    Initial guess for the fixed point, will be overwritten with solution

    procedure(abstract_linear_solver_csp) :: solver

    Linear solver to be used to find Newton step

    integer, intent(out) :: info

    Information flag

    real(kind=sp), intent(in), optional :: tolerance
    type(newton_sp_opts), intent(in), optional :: options
    class(abstract_opts), intent(in), optional :: linear_solver_options

    Options for the linear solver

    class(abstract_precond_csp), intent(in), optional :: preconditioner

    Preconditioner for the linear solver

    procedure(abstract_scheduler_sp), optional :: scheduler
  • private subroutine newton_cdp(sys, X, solver, info, tolerance, options, linear_solver_options, preconditioner, scheduler)

    Arguments

    Type IntentOptional Attributes Name
    class(abstract_system_cdp), intent(inout) :: sys

    Dynamical system for which we wish to compute a fixed point

    class(abstract_vector_cdp), intent(inout) :: X

    Initial guess for the fixed point, will be overwritten with solution

    procedure(abstract_linear_solver_cdp) :: solver

    Linear solver to be used to find Newton step

    integer, intent(out) :: info

    Information flag

    real(kind=dp), intent(in), optional :: tolerance
    type(newton_dp_opts), intent(in), optional :: options
    class(abstract_opts), intent(in), optional :: linear_solver_options

    Options for the linear solver

    class(abstract_precond_cdp), intent(in), optional :: preconditioner

    Preconditioner for the linear solver

    procedure(abstract_scheduler_dp), optional :: scheduler

Subroutines

public subroutine constant_atol_dp(tol, target_tol, rnorm, iter, info)

Abstract interface to define tolerance scheduler for the Newton iteration

Arguments

Type IntentOptional Attributes Name
real(kind=dp), intent(out) :: tol

Tolerance to be used

real(kind=dp), intent(in) :: target_tol

Target tolerance

real(kind=dp), intent(in) :: rnorm

Norm of the residual of the current iterate

integer, intent(in) :: iter

Newton iteration count

integer, intent(out) :: info

Information flag

public subroutine constant_atol_sp(tol, target_tol, rnorm, iter, info)

Abstract interface to define tolerance scheduler for the Newton iteration

Arguments

Type IntentOptional Attributes Name
real(kind=sp), intent(out) :: tol

Tolerance to be used

real(kind=sp), intent(in) :: target_tol

Target tolerance

real(kind=sp), intent(in) :: rnorm

Norm of the residual of the current iterate

integer, intent(in) :: iter

Newton iteration count

integer, intent(out) :: info

Information flag

public subroutine dynamic_tol_dp(tol, target_tol, rnorm, iter, info)

Abstract interface to define tolerance scheduler for the Newton iteration

Arguments

Type IntentOptional Attributes Name
real(kind=dp), intent(out) :: tol

Tolerance to be used

real(kind=dp), intent(in) :: target_tol

Target tolerance

real(kind=dp), intent(in) :: rnorm

Norm of the residual of the current iterate

integer, intent(in) :: iter

Newton iteration count

integer, intent(out) :: info

Information flag

public subroutine dynamic_tol_sp(tol, target_tol, rnorm, iter, info)

Abstract interface to define tolerance scheduler for the Newton iteration

Arguments

Type IntentOptional Attributes Name
real(kind=sp), intent(out) :: tol

Tolerance to be used

real(kind=sp), intent(in) :: target_tol

Target tolerance

real(kind=sp), intent(in) :: rnorm

Norm of the residual of the current iterate

integer, intent(in) :: iter

Newton iteration count

integer, intent(out) :: info

Information flag