`roots` returns NaNs for roots with complicated function

[alecloudenback]

I am trying to build a more robust internal rate of return solver for ActuaryUtilities.jl and couldn't get IntervalRootFinding to pass some tests finding the roots. The following case failed for me.

cfs = [t % 10 == 0 ? -10 : 1.5 for t in 0:99]
f(i) =  sum(cfs .* [1/(1+i[1])^t for t in 0:length(cfs)-1])

This worked:

julia> roots(f,0..2)
1-element Vector{Root{Interval{Float64}}}:
 Root([0.0646316, 0.0646317], :unique)

But the following takes considerable time and returns a bunch of roots that are actually `N⚓

 julia> roots(f,-2..2)
19662-element Vector{Root{Interval{Float64}}}:
 Root([-1.00046, -1.00045], :unknown)
 Root([-0.999951, -0.99995], :unknown)
 Root([-1.00047, -1.00046], :unknown)
 Root([-1.00066, -1.00065], :unknown)
 Root([-0.999964, -0.999963], :unknown)
 Root([-1.00029, -1.00028], :unknown)
 Root([-1.0004, -1.00039], :unknown)
 ⋮
 Root([-0.999918, -0.999917], :unknown)
 Root([-0.999854, -0.999853], :unknown)
 Root([-0.999729, -0.999728], :unknown)
 Root([-0.999867, -0.999866], :unknown)
 Root([-1.00032, -1.00031], :unknown)
 Root([-1.00051, -1.0005], :unknown)

f shoots up towards infinity around -0.5, so maybe related to #161?

[dpsanders]

What a nasty function!

You can write f(i) = sum(cfs .* [1/(1+i)^t for t in 0:length(cfs)-1])

without the index [1], by the way.

The function diverges when i == -1, so I don't think the :unknown are unexpected. You can speed up the calculation using

roots(x -> 1/f(x), -2..2, Bisection, 1e-3)

although this will not prove the uniqueness of the root. You can also try using more precision using e.g

f(big(-1.00066.. -1.00065))

to try to exclude some of those :unknown, but when there is a jump from -infinity to +infinity I don't think you can ever exclude it using interval methods alone. I think it will be difficult to have any automated way of detecting the divergence to infinity in a function like this.

If you have prior knowledge about the function then you should use that to exclude these badly-behaved places from consideration.

[alecloudenback]

My goal is to return the nearest root to zero between -2..2 or nothing if no root within that range. I could in practice justify limiting the range to -1..1, but that doesn't cut at the root of the problem.

Unfortunately, I don't know in advance how well behaved a series of cashflows will be. Do you know of techniques where I could efficiently search outwards from zero?

[dpsanders]

I wonder if Roots.jl would be better behaved for this particular problem?

[alecloudenback]

Thanks, Roots.jl does seem to work better in this situation. Do you want me to leave the issue open?

[dpsanders]

OK I'm glad to hear it. Yes let's leave it open since it's a useful test case to think about, thanks!

[Kolaru]

This is related to #59, as you probably want to limit the number of unknown intervals to significant less than 19K.

Then filtering the resulting intervals to find the root closest to 0 should not be a problem.