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. 
>
>

Reply via email to