Well, functions like `x -> exp(-x^3)-1` or `x -> (x^3 + 1) - 1` have
obvious roots. The root of the test function above is not that obvious. To
tell the rest of the story: I once tried to symbolically solve it with
Mathematica
Mathematica> Solve[Log[x] + x*x/(2*Exp[1]) - 2*x/Sqrt[Exp[1]] + 1 == 0,
x]
Solve:: This system cannot be solved with the methods available to Solve
and, being unsuccessful, asked on their mailing list. People from Wolfram
took it up and promised to tweaked the solver for a class of expressions
this function belongs to. Up to now (Version 9.0) this has not happened.
You were right about the stop condition, function `ridders` used up all
1000 loops. Correcting that, timings are comparable:
julia> @time [ridders(sin, 3.0, 4.0) for i in 1:10000];
elapsed time: 0.167489776 seconds (40471888 bytes allocated)
@time [bisect_root(sin, 3.0, 4.0) for i in 1:10000];
elapsed time: 0.127175753 seconds (26311888 bytes allocated)
both finding the root within 4.44e-16 accuracy.
As a side remark, `bisect_root` as is should check whether the sign changes
between `x1` and `x2`, otherwise it returns a point on the boundary:
julia> bisect_root(sin, 1.0, 2.0)
(1.9999999999999998,2.0)
And another small suggestion: In my version of `bisect` I check in every
loop whether `sign(fn(xm)) == 0`. This helps to stop early in some cases
such as when 0.0 is a root.
It could also be helpful to return the number of iterations (or function
calls).
On Tuesday, February 25, 2014 5:01:12 AM UTC+1, Jason Merrill wrote:
>
> Because I can't leave well enough alone...
>
> julia> function ridder_loop()
> accum = 0.0
> for i = 1:10000
> accum += ridders(sin, 2.0, 4.0)
> end
> accum
> end
>
> julia> function bisect_loop()
> accum = 0.0
> for i = 1:10000
> accum += bisect_root(sin, 2.0, 4.0)[1]
> end
> accum
> end
>
> julia> @time ridder_loop()
> elapsed time: 10.430907316 seconds (1596537204 bytes allocated)
> 31415.92653589205
>
> julia> @time ridder_loop()
> elapsed time: 10.440722588 seconds (1596480048 bytes allocated)
> 31415.92653589205
>
> julia> @time bisect_loop()
> elapsed time: 0.152900673 seconds (25984456 bytes allocated)
> 31415.92653589205
>
> julia> @time bisect_loop()
> elapsed time: 0.144590689 seconds (25920064 bytes allocated)
> 31415.92653589205
>
> So for this example, admittedly with an overly simplistic objective
> function, bisection is winning by a factor of about 70.
>
>