Nonlinear Solver Iterator Interface

There is an iterator form of the nonlinear solver which somewhat mirrors the DiffEq integrator interface:

import NonlinearSolve as NLS
import NonlinearSolveBase as NLSB

f(u, p) = u .* u .- 2.0
u0 = 1.5
probB = NLS.NonlinearProblem(f, u0)

nlcache = NLS.init(probB, NLS.NewtonRaphson())
GeneralizedFirstOrderAlgorithmCache(
    alg = NewtonRaphson(
        descent = NewtonDescent(),
        jacobian_reuse = JacobianReuse(
            max_age = -1,
            max_residual_ratio = 0.1
        ),
        autodiff = AutoForwardDiff(),
        vjp_autodiff = AutoFiniteDiff(
            fdtype = Val{:forward}(),
            fdjtype = Val{:forward}(),
            fdhtype = Val{:hcentral}(),
            dir = true
        ),
        jvp_autodiff = AutoForwardDiff(),
        concrete_jac = Val{false}()
    ),
    u = 1.5,
    residual = 0.25,
    inf-norm(residual) = 0.25,
    nsteps = 0,
    retcode = Default
)

init takes the same keyword arguments as solve, but it returns a cache object that satisfies typeof(nlcache) <: AbstractNonlinearSolveCache. There are two cache forms:

  • Native iterative algorithms return a stepping cache. Call step! to advance it, or solve! to run it to completion.
  • Algorithms without a SciMLBase.__init method, such as the SimpleNonlinearSolve algorithms, return NonlinearSolveNoInitCache. This cache stores the problem and options but no iteration state, so call solve! directly; step!, get_fu, and get_nsteps are unavailable.

The iterator interface supports:

CommonSolve.step! — Method
step!(cache::AbstractNonlinearSolveCache, args...; kwargs...)

Perform one step of a nonlinear solver and mutate cache in place.

The public wrapper first checks whether the cache is still active, then calls NonlinearSolveBase.InternalAPI.step!, updates the step counters, and enforces a time limit when one is configured. It returns the value produced by the algorithm-specific implementation, which is commonly nothing.

Arguments

  • cache::AbstractNonlinearSolveCache: the initialized stepping cache to advance.
  • args...: positional arguments forwarded to the algorithm-specific implementation.

Keywords

  • recompute_jacobian::Union{Nothing, Bool} = nothing: whether to recompute a Jacobian when the algorithm uses one. nothing delegates the choice to the algorithm. This keyword is ignored or rejected by algorithms according to their own interface.
  • evaluate_residual::Bool = true: a hint that the algorithm may skip the residual evaluation at the newly accepted iterate. It is honored only when supports_deferred_residual returns true; call refresh_residual! before reading a deferred residual.
  • kwargs...: additional algorithm-specific keyword arguments.

Returns

The value returned by InternalAPI.step!. The cache is mutated in place. Calling step! on a terminated cache does nothing and returns nothing.

Extension Rules

Solver packages implement InternalAPI.step!, not this CommonSolve.step! method. The implementation must leave the cache state consistent with get_u, get_fu, and the termination cache. The wrapper owns the top-level nsteps and stats.nsteps increments, so an implementation should not increment those counters for the same step.

Examples

import NonlinearSolve

prob = NonlinearSolve.NonlinearProblem((u, p) -> u^2 - p, 1.0, 2.0)
cache = NonlinearSolve.init(prob, NonlinearSolve.NewtonRaphson())
NonlinearSolve.step!(cache)
source

We can perform 10 steps of the Newton-Raphson solver with the following:

for i in 1:10
    NLS.step!(nlcache)
end

Code that accepts any nonlinear algorithm can detect the second form and choose the complete solve path. The stepping branch can use solve_cache! when it needs the allocation-sensitive cache result:

function solve_any(prob, alg)
    cache = NLS.init(prob, alg)
    if cache isa NLSB.NonlinearSolveNoInitCache
        return NLS.solve!(cache)
    end
    NLSB.solve_cache!(cache)
    return NLS.solve!(cache)
end

simple_sol = solve_any(probB, NLS.SimpleNewtonRaphson())
retcode: Success
u: 1.414213562373095

We currently don't implement a Base.iterate interface but that will be added in the future.