SciML / NonlinearSolve.jl

High-performance and differentiation-enabled nonlinear solvers (Newton methods), bracketed rootfinding (bisection, Falsi), with sparsity and Newton-Krylov support.
https://docs.sciml.ai/NonlinearSolve/stable/
MIT License
216 stars 39 forks source link

Muller's method #403

Open fgittins opened 2 months ago

fgittins commented 2 months ago

What kind of problems is it mostly used for? Please describe.

Determining the roots of univariate complex functions. Muller's method is commonly used to determine the (quasi-normal) oscillation modes of neutron stars [1, 2].

Note: Given that this method only works for univariate functions, it might be too restrictive for NonlinearSolve.jl, but I thought I would suggest it anyway, since I frequently use it. Also, the algorithm is very simple so it should be easy to implement.

Describe the algorithm you’d like

The algorithm is quite simple and described well in Press et al. [3]. (See also Wikipedia [4].)

To summarise: Muller's method is a generalisation of the secant method, with the key difference being that it uses quadratic interpolation across three points (as opposed to linear interpolation among two). Solving for the roots of a quadratic is trivial and allows the method to work for complex roots.

Other implementations to know about

There is an existing Julia implementation in Roots.jl [5] and a Fortran routine in IMSL [6].

References

[1] Kokkotas & Schutz (1992) [2] Kruger, PhD Thesis (2016) [3] Press et al. (2007), Sec. 9.5.2, p. 466 [4] Wikipedia, Muller's method [5] Roots.muller [6] IMSL Math/Library Users Manual, Zanly

avik-pal commented 2 months ago

Given that this method only works for univariate functions, it might be too restrictive for NonlinearSolve.jl

This would be great to have in SimpleNonlinearSolve

fgittins commented 2 months ago

Thanks for the feedback. I had some time so I put together an initial implementation at this algorithm with some simple tests. If it looks reasonable, I can make a pull request.

One thing to note is that differently to other non-linear solvers, Muller's method requires three initial guesses. This would need to be made clear in subsequent documentation.

avik-pal commented 2 months ago

yes open the pr