NewtonRaphson from NonlinearSolve.jl is not throwing an error when convergence is not achieved

I think I found a bug in the NonlinearSolve.jl package. Consider the problem

using NonlinearSolve, GLMakie
f(x, p) = tanh(x / 1e-3) + x - 0.5
x0 = 1.0

prob = NonlinearProblem(f, x0)
sol = solve(prob, NewtonRaphson(), abstol=1e-3)

which returns

julia> sol.u

julia> f(sol.u, nothing)

This is clearly not a converged solution, since it is far from a zero of the function. However, no error is thrown by the solver even though the Newton iterations are not converging. This non-convergence is due to the shape of this nonlinear function

in the case where the zero is in the steep part of the graph. For an initial condition on the regions of the function dominated by x, you just end up iterating between two different values and do not achieve convergence.

I understand that this is an expected behavior of the Newton Raphson method for this function, but I think NonlinearSolve.jl should throw and error instead of letting you think that convergence was achieved.

But sol.retcode?

sol.retcode gives

julia> sol.retcode
ReturnCode.MaxIters = 4

which I don’t think is correct, right? Since if you look for sol.stats gives

julia> sol.stats
SciMLBase.NLStats(1001, 1000, 1000, 1000, 1000)

and it seems that the maximum of 1000 iterations has been reached (I assumed this number is the one by default on the package.

I was expecting an error similar to how it is handled when solving and ODE with OrdinaryDiffEq, where it is explicitly stated when the integration is not converging within tolerances for adaptive timesteps.

It stopped at max iterations and exited early as a failure. That’s what you’re seeing?

OrdinaryDiffEq also returns a retcode that describes the solution, the same retcode actually. Indeed we should make NonlinearSolve.jl throw a warning for this case, but that’s a feature request for the repo (or a PR implementation) and not really a Discourse discussion.

In case it’s not clear, this is the 4th kind of return code, it doesn’t mean there were only 4 iterations.

Oh, where could I see the meaning of these error codes?

That would be nice, if I have time next week I may try to implement it myself. The problem is that even though the no convergence information was in the retcode, it did not seem obvious to me and I assumed that the nonlinear solver succeeded. It was not until I plotted the results that I noticed that truly there was no convergence, which I think it could be misleading when using the package.

Even better than a warning, I think that if the full description of the solution object was printed in the REPL would already be okay to see at a glance if convergence happened. Now I was only getting sol.u printed in the REPL.

Anyway, thank you a lot for the answers.

1 Like

Yes. If you see the ODESolution in SciMLBase you will find the show overload. We should do the same to NonlinearSolution