Hacker News Board Hacker News, browsed like a classic message board
Live data from Hacker News

Nonlinearsolve.jl: Fast and Robust Solvers for Nonlinear Equations in Julia

arxiv.org

Post reply on HN
Print view View original

21–30 of 37 posts

Re: Nonlinearsolve.jl: Fast and Robust Solvers for Nonlinear Equations in Julia

#21

Post by ubj »

Bimos wrote: post #7
FreakLegion wrote: post #6

Earlier quoted context omitted.

> all nonlinear equation systems are nonconvex Maybe you have something more particular in mind when you say "systems", but not all nonlinear functions are non-convex. Least squares, for example, is nonlinear and convex. Also note that IPOPT, while wonderful, is a local solver. It may not be limited to convex problems, but those are the only ones it's guaranteed to solve to optimality.

The feasible region {x | f(x) = 0} is nonconvex no matter whether f is convex.

I'm assuming you're referring to nonlinear f(x) because this statement is trivially false for linear f(x).

But consider the function f(x) = max(0, g(x) - c) where the following holds:

* g(x) is nonlinear, positive definite in x, and convex.

* c > 0

Then f(x) is nonlinear and convex (it's the pointwise maximum of convex functions), and the set {x : f(x) = 0} is a convex set.

Re: Nonlinearsolve.jl: Fast and Robust Solvers for Nonlinear Equations in Julia

#22

Post by adgjlsfhk1 »

wenc wrote: post #20
ChrisRackauckas wrote: post #14

Earlier quoted context omitted.

Definitely an interesting response I didn't see coming. This is one of the reasons to share preprints! IPOPT is a nonlinear optimization tool, not a nonlinear solver. Generally from what I've seen it's not a good idea to solve nonlinear systems with a nonlinear optimizer, but you can phrase it as a "constraint satisfaction problem", i.e. a nonlinear optimization with a trivial `maximize 0` loss function but with equa…

Yes tribal knowledge says NLP solvers aren’t tailored to NLE solvers as NLE solvers are. There have been papers in this past showing this but in my experience, it’s not obvious in practice. Nonlinear systems can often surprise us. The initial point (and randomness) has more influence on the solution time than almost anything else. I forget the exact internal formulation of the problem in IPOPT but for an NLE it’s lik…

isn't MA57 just a standard good sparse linear solver? For sparse problems, NonlinearSolve.jl (via LinearSolve.jl) will either use KLU or Cholmod both of which are also quite fast. Also, even if the linear solver is the bottleneck, that doesn't mean different algorithms won't make a differenc.e Different nonlinear solve algorithms can lead to needing fewer iterations (e.g. fewer linear solves).

Re: Nonlinearsolve.jl: Fast and Robust Solvers for Nonlinear Equations in Julia

#23

Post by wenc »

adgjlsfhk1 wrote: post #22
wenc wrote: post #20

Earlier quoted context omitted.

Yes tribal knowledge says NLP solvers aren’t tailored to NLE solvers as NLE solvers are. There have been papers in this past showing this but in my experience, it’s not obvious in practice. Nonlinear systems can often surprise us. The initial point (and randomness) has more influence on the solution time than almost anything else. I forget the exact internal formulation of the problem in IPOPT but for an NLE it’s lik…

isn't MA57 just a standard good sparse linear solver? For sparse problems, NonlinearSolve.jl (via LinearSolve.jl) will either use KLU or Cholmod both of which are also quite fast. Also, even if the linear solver is the bottleneck, that doesn't mean different algorithms won't make a differenc.e Different nonlinear solve algorithms can lead to needing fewer iterations (e.g. fewer linear solves).

I don’t know. I’m curious is all.

MA57 seems generic but I’ve tried different linear solvers with IPOPT even ones that were supposed to be better, but MA57 still performed the best.

Nonlinear numerical problems involve a bit of art. Theory doesn’t always pan out in practice.

Re: Nonlinearsolve.jl: Fast and Robust Solvers for Nonlinear Equations in Julia

#24

Post by selimthegrim »

ChrisRackauckas wrote: post #13
Bimos wrote: post #8

On GPU, it can only solve small instances which one single GPU thread can handle, and can only benefit from GPU by solving multiple instances, right?

No, there are two types of GPU usage: the ensemble form and the array form. In the array form, you simply make your state variable a GPU and can run all internal operations on the GPU. This is used for GPU acceleration of `f` functions where you have a large state and a relatively "regular" function in `f`, such as neural networks or PDE discretizations. Thus you can see this in action for example on things like the…

This sounds fascinating and I’d love to help.

Re: Nonlinearsolve.jl: Fast and Robust Solvers for Nonlinear Equations in Julia

#25

Post by kxyvr »

ChrisRackauckas wrote: post #19
sobellian wrote: post #17

Earlier quoted context omitted.

You can also put it as a minimization problem to minimize ||g(u)||, but I assume that does not magically make anything more efficient.

That makes your Newton method require calculating the Hessian of the objective function, which is a second derivative. A normal Newton method for nonlinear systems uses the Jacobian, which is the first derivative (the connection is that the Hessian is the Jacobian of the gradient). But indeed, formulating it so that you have discontinuities (absolute values) and require second derivatives instead of the first for the…

This is not entirely true. A typical trick for accomplishing this is to simply use the Hessian approximation `g'(x)*g'(x)` (star denotes adjoint), thus cutting off the second-derivative evaluation of g, and then use a modified CG method like Newton-Krylov for line-search methods or Steihaug-Toint for trust-region methods. It's very robust, matrix-free, and doesn't require an excessive amount of code to implement. Further, the method can be preconditioned with a positive definite preconditioner.

And, to be clear, there's no guarantee that a zero residual be found, only that `g'(x)*g(x)=0`, but that tends to be the nature of nonlinear equations.

Re: Nonlinearsolve.jl: Fast and Robust Solvers for Nonlinear Equations in Julia

#26

Post by kxyvr »

wenc wrote: post #23
adgjlsfhk1 wrote: post #22

Earlier quoted context omitted.

isn't MA57 just a standard good sparse linear solver? For sparse problems, NonlinearSolve.jl (via LinearSolve.jl) will either use KLU or Cholmod both of which are also quite fast. Also, even if the linear solver is the bottleneck, that doesn't mean different algorithms won't make a differenc.e Different nonlinear solve algorithms can lead to needing fewer iterations (e.g. fewer linear solves).

I don’t know. I’m curious is all. MA57 seems generic but I’ve tried different linear solvers with IPOPT even ones that were supposed to be better, but MA57 still performed the best. Nonlinear numerical problems involve a bit of art. Theory doesn’t always pan out in practice.

MA57 is a sparse symmetric solver. IPOPT uses this because it solves a modified KKT system that includes both the Hessian and constraint information simultaneously. This differs from composite step methods, which split the optimization step into a step for feasibility (quasi-normal step) and a step for optimality (tangential step). An algorithm like NITRO uses a composite step method and as a result uses a different kind of linear system solver.

Re: Nonlinearsolve.jl: Fast and Robust Solvers for Nonlinear Equations in Julia

#27

Post by ChrisRackauckas »

kxyvr wrote: post #25
ChrisRackauckas wrote: post #19

Earlier quoted context omitted.

That makes your Newton method require calculating the Hessian of the objective function, which is a second derivative. A normal Newton method for nonlinear systems uses the Jacobian, which is the first derivative (the connection is that the Hessian is the Jacobian of the gradient). But indeed, formulating it so that you have discontinuities (absolute values) and require second derivatives instead of the first for the…

This is not entirely true. A typical trick for accomplishing this is to simply use the Hessian approximation `g'(x)*g'(x)` (star denotes adjoint), thus cutting off the second-derivative evaluation of g, and then use a modified CG method like Newton-Krylov for line-search methods or Steihaug-Toint for trust-region methods. It's very robust, matrix-free, and doesn't require an excessive amount of code to implement. Fur…

> It's very robust, matrix-free, and doesn't require an excessive amount of code to implement. Further, the method can be preconditioned with a positive definite preconditioner.

It is not robust, in fact it's well-known as a numerically unstable method. g'(x)*g'(x)` is numerically unstable because it squares the condition number of the matrix. If you have a condition number of 1e-10, which is true for many examples in scientific problems (like the DFN battery model in the paper), the condition number of J'J is 1e-20, which effectively means a linear solve is unable to retrieve any digits of accuracy with 64-bit floating point numbers. So while I agree that you can use symmetric linear solvers as a small performance improvement in some cases, you have to be very careful when you do this as it is a very numerically unstable method. For some details, see https://math.stackexchange.com/a/2874902

Also, if you look at Figure 5 in the paper, you will see that there is a way to choose operator conditioning https://docs.sciml.ai/LinearSolve/stable/basics/OperatorAssu..., and if you tell it to assume the operator is well-conditioned it will actually use this trick. But we do not default to that because of the ill-conditioning.

Re: Nonlinearsolve.jl: Fast and Robust Solvers for Nonlinear Equations in Julia

#28

Post by avikpal1410 »

kxyvr wrote: post #25
ChrisRackauckas wrote: post #19

Earlier quoted context omitted.

That makes your Newton method require calculating the Hessian of the objective function, which is a second derivative. A normal Newton method for nonlinear systems uses the Jacobian, which is the first derivative (the connection is that the Hessian is the Jacobian of the gradient). But indeed, formulating it so that you have discontinuities (absolute values) and require second derivatives instead of the first for the…

This is not entirely true. A typical trick for accomplishing this is to simply use the Hessian approximation `g'(x)*g'(x)` (star denotes adjoint), thus cutting off the second-derivative evaluation of g, and then use a modified CG method like Newton-Krylov for line-search methods or Steihaug-Toint for trust-region methods. It's very robust, matrix-free, and doesn't require an excessive amount of code to implement. Fur…

That is precisely how trust region methods work for nonlinear least squares problems and nonlinear systems. MINPACK, which both SciPy and Matlab use, defaults to this kind of TR scheme (without the matrix-free part). We do have TR comparisons in the paper.

Note that we don't however solve the normal form J'J system as that is generally bad (unless user opts in to doing so) and instead use least squares formulation which is more numeraically stable

Re: Nonlinearsolve.jl: Fast and Robust Solvers for Nonlinear Equations in Julia

#29

Post by avikpal1410 »

akutlay wrote: post #12
wenc wrote: post #2

Would be interested to see this compared to IPOPT. https://coin-or.github.io/Ipopt/ IPOPT isn't tailored to NLEs specifically, but it does solve NLPs well. It also uses some amazing Fortran linear algebra routines from Harwell (MA57).

I dont know much about this but the fact that they did not mention IPOPT or Highs (which I would think them as two of the most popular software used for these problems) made me question this too.

IPOPT as mentioned in other posts solves Nonlinear Optimization Problems. Working around it to make it solve nonlinear equations or nonlinear least squares problems typically doesn't go well. We will add benchmarks for the NLLS part (which we don't talk about in the paper) and that will contain IPOPT comparisons.

For some comparison (this is not a benchmark me or anyone on this paper have written), see https://juliapackagecomparisons.github.io/comparisons/math/n.... Specialized solvers for NLLS pretty much always beat general optimization solvers.

Re: Nonlinearsolve.jl: Fast and Robust Solvers for Nonlinear Equations in Julia

#30

Post by kxyvr »

ChrisRackauckas wrote: post #27
kxyvr wrote: post #25

Earlier quoted context omitted.

This is not entirely true. A typical trick for accomplishing this is to simply use the Hessian approximation `g'(x)*g'(x)` (star denotes adjoint), thus cutting off the second-derivative evaluation of g, and then use a modified CG method like Newton-Krylov for line-search methods or Steihaug-Toint for trust-region methods. It's very robust, matrix-free, and doesn't require an excessive amount of code to implement. Fur…

> It's very robust, matrix-free, and doesn't require an excessive amount of code to implement. Further, the method can be preconditioned with a positive definite preconditioner. It is not robust, in fact it's well-known as a numerically unstable method. g'(x)*g'(x)` is numerically unstable because it squares the condition number of the matrix. If you have a condition number of 1e-10, which is true for many examples i…

Yes, care must be taken to stabilize the algorithm. However, I do not believe it to be true that no digits of accuracy can be obtained. The algorithm falls back to the Cauchy (steepest-descent) point on the first iteration and this will give reduction. That said, I'm willing to solve this particular problem and see. Where's your code for the DFN battery model? The page here leads to a dead github link:

https://help.juliahub.com/batteries/stable/api/#JuliaSimBatt...

Post reply on HN