# Nonlinear systems: NLsolve alternative (non-convergence issue)

Hello Julia community

I’m having trouble with the convergence of NLsolve when solving a large system of linear equations (the transition between steady states in a heterogeneous-agent model in economics). Even in its simpler form (coarse grids, not many periods), the solver can’t seem to get past e-03/e-04 and I need at least e-07. I tried using AD but given I have many storage variables that have to be defined as arrays of DualNumbers, the speed decreased quite a bit and the precision didn’t change much. I used all the three method available in NLSolve and the accuracy was pretty much the same (though speed varied a lot, anderson being the fastest by far).

Matlab can solve the model reasonably fast using the standard fsolve with the Levenberg-Marquardt method, and I was wondering if there’s any reasonably straightforward way to use this for solving a nonlinear system in Julia. Alternatively, do you know any other solver that could be a better fit than NLsolve given my needs? I googled but didn’t find much besides JuMP (using the standard constant objective function trick), and given its solvers are installed independently I don’t wanna spend much time going through its setup if it’s not gonna be that good, so I rather ask you first.

TL;DR: I can’t get NLsolve to converge, are there good alternatives for solving large non-linear systems? (ideally suited for heterogeneous-agent economics models)

It would help to know more about the problem:

• what nonlinear solver did you try?
• what linear solver did you try (in case of newton…)?

BifurcationKit.jl can be used to solve large scale problems and study their bifurcations (see the tutorials. You can also try SIAMFANLEquations.jl or JuMP.jl

1 Like

Sorry, I don’t follow. I stated I have only tried the NLsolve.jl package with all available methods (Newton, trust region and Anderson acceleration, default specs in all) and none worked, either with or without Forward Differentiation.

Additional info: I tried solving the model in Matlab using trust region instead of Levenberg-Marquardt and indeed the accuracy decreased (e-05/e-06) but still is 100x than NLsolve’s (with the same method, trust region), so it’s not only the method that is causing the non-convergence.

Gonna check your suggestions. Do you happen to know Ipopt, by the way? It’s the only solver that supports non-linear programming in JuMP.jl besides Artelys Knitro (and I’m not gonna pay for that one, at least not unless I’m 100% sure it will solve my problem).

perhaps an homotopy would work