# Ridders' method in Julia

numerical-analysis root-finding julia

Ridders' method is a root-finding method based on the regula falsi method that uses an exponential function to fit a given function bracketed between $x_0$ and $x_1$.

In the algorithm, in each iteration first the function is evaluated at a third point $x_2 = (x_0+x_1)/2$, then applies the unique exponential function

$f(x_0) -2f(x_2)e^Q +f(x_1)e^{2Q} = 0$

to transform the function into a straight line.

Solving the quadratic equation for $e^Q$ gives

$e^Q = {f(x_2) + \mathrm{sign}(f(x_1)) \sqrt{f(x_2)^2-f(x_0)f(x_1)} \over f(x_1)} \ ,$

and applying the false position method to $f(x_0), f(x_2)e^Q$ and $f(x_1)e^{2Q}$ gives the new guess

$x = x_2 + (x_2-x_0) {\mathrm{sign}(f(x_0)-x(x_1)) f(x_2) \over \sqrt{f(x_2)^2-f(x_0)f(x_1)}} \ .$

This equation converges quadratically, therefore the number of significant digits approximately doubles with each iteration (two function evaluations).

After each iteration, $x_0$ and $x_1$ will be redefined, being the new guess one of them and the other one defined according to the sign of the evaluations of $f$ at $x_0$, $x_1$ and $x_2$.

In Julia:

function ridders(f::Function, x0::Number, x1::Number, args::Tuple=();
xtol::AbstractFloat=1e-5, ytol=2eps(Float64),
maxiter::Integer=50)
y1 = f(x1,args...)
y0 = f(x0,args...)
for _ in 1:maxiter
x2 = mean([x0,x1])
y2 = f(x2,args...)
x = x2 + (x2-x0) * sign(y0-y1)*y2/sqrt(y2^2-y0*y1)
# x-tolerance.
if min(abs(x-x0),abs(x-x1)) < xtol
return x
end
y = f(x,args...)
# y-tolerance.
if abs(y) < ytol
return x
end
if sign(y2) != sign(y)
x0 = x2
y0 = y2
x1 = x
y1 = y
elseif sign(y1) != sign(y)
x0 = x
y0 = y
else
x1 = x
y1 = y
end
end
error("Max iteration exceeded")
end

The code follows the same structure as the false position method implementation, supports complex numbers and defined accuracy in $x$ and $y$ dimensions.

As an example, for the equation $f(x)=x^2/12+x-4 \quad (1 \leq x \leq 5)$, the solution is $x = \sqrt{84}-6$:

julia> (f(x)=x^2/12+x-4; x0=1; x1=5; ridders(f,x0,x1))
3.1651513899117836
julia> sqrt(84)-6
3.1651513899116797