NLsolve.jl icon indicating copy to clipboard operation
NLsolve.jl copied to clipboard

Add info about the function in the README

Open dpsanders opened this issue 6 years ago • 9 comments

The function in the README has precisely 80 roots in (-5..5)^2 (assuming IntervalRootFinding is working correctly...):

julia> using IntervalArithmetic, IntervalRootFinding, StaticArrays

julia> f(x, y) = SVector( (x+3)*(y^3-7)+18, sin(y*exp(x)-1) )
f (generic function with 1 method)

julia> f(xx) = f(xx...)
f (generic function with 2 methods)

julia> rts = roots(f, IntervalBox(-5..5, 2), Newton, 1e-20)
80-element Array{IntervalRootFinding.Root{IntervalArithmetic.IntervalBox{2,Float64}},1}:
 Root([4.99884, 4.99885] × [1.68094, 1.68095], :unique)
 Root([4.98639, 4.9864] × [1.68053, 1.68054], :unique)
 Root([4.97379, 4.9738] × [1.68011, 1.68012], :unique)
 Root([4.96103, 4.96104] × [1.67968, 1.67969], :unique)
 Root([4.9481, 4.94811] × [1.67925, 1.67926], :unique)
 Root([4.935, 4.93501] × [1.67881, 1.67882], :unique)
 Root([4.92172, 4.92173] × [1.67836, 1.67837], :unique)
 Root([4.90826, 4.90827] × [1.6779, 1.67791], :unique)
 Root([4.89461, 4.89462] × [1.67743, 1.67744], :unique)
 Root([4.88077, 4.88078] × [1.67696, 1.67697], :unique)
 ⋮
 Root([2.91809, 2.9181] × [1.58188, 1.58189], :unique)
 Root([2.80939, 2.8094] × [1.57427, 1.57428], :unique)
 Root([2.68705, 2.68706] × [1.56525, 1.56526], :unique)
 Root([2.54714, 2.54715] × [1.55431, 1.55432], :unique)
 Root([2.3837, 2.38371] × [1.5406, 1.54061], :unique)
 Root([2.18717, 2.18718] × [1.5226, 1.52261], :unique)
 Root([1.94053, 1.94054] × [1.49727, 1.49728], :unique)
 Root([1.60901, 1.60902] × [1.45725, 1.45726], :unique)
 Root([1.10116, 1.10117] × [1.377, 1.37701], :unique)
 Root([0, 0] × [1, 1], :unique)

julia> all([root.status for root in rts] .== :unique)
true

dpsanders avatar Mar 12 '18 19:03 dpsanders

(This is with the latest PR to IntervalRootFinding, changing bisection to not bisect exactly in the centre of an interval.)

dpsanders avatar Mar 12 '18 19:03 dpsanders

Thanks! I'll try to write something intelligent about the problem :)

pkofod avatar Mar 12 '18 19:03 pkofod

By the way, where did that function come from? It's pretty nasty ;)

dpsanders avatar Mar 12 '18 20:03 dpsanders

Do you know any other software with which we can check if that number of roots is correct? Mathematica doesn't seem to be able to do it...

dpsanders avatar Mar 12 '18 20:03 dpsanders

It's pretty nasty ;)

It's very nasty indeed! And above all, solving it is a very bad example of what this package claims it can do :) It's old... You'd have to ask svillemot, he added it in Dec 24, 2013 :) https://github.com/JuliaNLSolvers/NLsolve.jl/commit/c434bd80b69bc3a159d681162ee88b9d886a116b#diff-04c6e90faac2675aa89e2176d2eec7d8

pkofod avatar Mar 12 '18 20:03 pkofod

cc @svillemot ^

dpsanders avatar Mar 12 '18 20:03 dpsanders

I found this:

https://www.sciencedirect.com/science/article/pii/S0377042700007111

which effectively seems to do what IntervalRootFinding.jl does, only without intervals.

It has some nice looking hard test cases.

dpsanders avatar Mar 12 '18 20:03 dpsanders

I actually don't remember why I picked this particular example, sorry. But clearly I was not aware that it had 80 roots! I guess it should be replaced by a better-behaving one…

svillemot avatar Mar 12 '18 22:03 svillemot

I actually don't remember why I picked this particular example, sorry. But clearly I was not aware that it had 80 roots! I guess it should be replaced by a better-behaving one…

It does serve as a warning though, we might just want to make it explicit :)

pkofod avatar Mar 13 '18 08:03 pkofod