Getting to zero allocations in NonlinearSolve.jl

I am in the process of migrating some code from NLsolve.jl to NonlinearSolve.jl. I have written a simple Newton-Raphson test, but I cannot get the allocations down to zero; both solve! and reinit! allocate. Can this be improved?

using NonlinearSolve

N = 1000
u0 = rand(N)
f!(J, u, p) = @. J = (u - 1)

problem = NonlinearProblem{true}(f!, u0)
solver = NewtonRaphson(
  autodiff = AutoForwardDiff(;
    chunksize = NonlinearSolve.pickchunksize(problem.u0),
  ),
)
cache = init(problem, solver)

solve!(cache)
reinit!(cache)

println(@allocated solve!(cache))
println(@allocated reinit!(cache))
1 Like

The caching API minimizes allocations but doesn’t make them go down to zero. The allocations in the above code would come from LU factorization IIRC.

1 Like

Zero allocations is only possible for StaticArrays but that wouldn’t be performant for such large arrays

may i ask why the LU factorization allocates? in theory, an LU factorization only requires to allocate 1 matrix and 1 Integer vector.

Julia interface for mutating LU accepts only a matrix and allocates a new permutation vector each time.

If worried about allocations, GitHub - DynareJulia/FastLapackInterface.jl allows you to pre-allocate temporaries and reuse them…

It’s quite a bit slower though in many cases.

Hi Chris, would you know by any chance if LinearSolve.jl supports pre-allocation? At least for the native julia backend?

If you use the init interface it will allocate all caches required for solve!, as described in