Eisenstat-Walker Newton-Krylov solver

[termi-official]

Needs a new LinearSolve.jl release first.

Checklist

• Appropriate tests were added
• Any code changes were done in a way that does not break public API
• All documentation related to code changes were updated
• The new code follows the
  contributor guidelines, in particular the SciML Style Guide and
  COLPRAC.
• Any new documentation only uses public API

Additional context

Related issue is here #319 . I left the default to not add any forcing.

[termi-official]

This MWE suggests that it speeds up things

using NonlinearSolveFirstOrder
using LinearAlgebra
using SciMLLogging
using LinearSolve
using BenchmarkTools

f(x, p) = Tridiagonal(-1.0ones(p-1), 2ones(p), -1.0ones(p-1))*x +  x .* log.(x)

N = 10000

prob = NonlinearProblem(f, rand(N) .+ 1.0, N)
sol1 = solve(prob, NewtonRaphson(; linsolve=KrylovJL_CG()); verbose=NonlinearVerbosity(SciMLLogging.All()))
sol2 = solve(prob, NewtonRaphson(; linsolve=KrylovJL_CG(), forcing = EisenstatWalkerForcing2()); verbose=NonlinearVerbosity(SciMLLogging.All()))

solver = NewtonRaphson(; linsolve=KrylovJL_CG())
@btime solve($prob, $solver) # 44ms


solver = NewtonRaphson(; linsolve=KrylovJL_CG(), forcing = EisenstatWalkerForcing2())
@btime solve($prob, $solver) # 12ms

However, for the Brussellator problem in SciMLBenchmarks this is not the case anymore, because the initial guess is too far away from the root. There is still a speedup in the order of 2x for some cases.