Modifying the Ridders method with a good first estimate

numerical methodsroots

I am using the ridder's method to approximate a numerical function.

It has been working very well for my needs but I'm feeling that it is inefficient in that it is taking more iterations than necessary.

One main constraint is that the function will only take integer values, so I have to round the answer and get the solution when (upper bound – lower bound = 1).

For my function, I already have the upper and lower bound and I also have a good first estimate of the solution that I'm trying to use to decrease the number of iterations. Here is an example of why I feel it is inefficient:

# x0 x2 x1 f(x0) f(x2) f(x1)
1 -2231 2230 -46 -1381.95 774.00 -0.17
2 -46 2230 -46 -0.17 774.00 -0.17
3 -46 1092 -46 -0.17 437.28 -0.17
4 -46 523 -46 -0.17 230.37 -0.17
5 -46 239 -46 -0.17 118.18 -0.17
6 -46 97 -46 -0.17 59.84 -0.17
7 -46 26 -46 -0.17 30.16 -0.17
8 -46 -10 -46 -0.17 15.09 -0.17
9 -46 -28 -46 -0.17 7.55 -0.17
10 -46 -37 -46 -0.17 3.71 -0.17
11 -46 -42 -46 -0.17 1.56 -0.17
12 -46 -44 -46 -0.17 0.70 -0.17
13 -46 -45 -46 -0.17 0.27 -0.17

Where x0 and x2 are the lower and upper bounds respectively, and x1 is the current estimate of the root. As you can see it converges to -46 which is the correct solution in my case (since it is the closest to 0), but I do not like the fact that it has to go through 13 iterations to get the correct solution.

The initial estimate is usually a good estimate, but it is impossible how far away from the needed solution it is.

Is there any way I can update the formula to better use my initial estimate and reduce the number of iterations? (Note: I still would like to use the ridder's method for its robustness and ease of use).

Note: the main equation that I feel could be updated is based on the false position method given below:

$$ x_{new} = x_1 + (x_1 – x_0) \times \frac{sign(f(x0))\times f(x1)}{\sqrt{f(x1)^2 – f(x0) \times f(x2)}} $$

where $x_1$ is our current estimate, $x_0$ and $x_2$ are the lower and upper bounds respectively. The $x_{new}$ is then used as either the upper or lower bound in the next iteration.

Update: Here is a matlab of the method that I am using, the main modification is the if statement at first (line 13) and the roundings to convert to integers. You can compare this to the official method here:

function xZero = rootSolve(func,xLow,xUpp, xest)

maxIter = 50;
fLow = feval(func,xLow);
fUpp = feval(func,xUpp);
xZero = [];

tol = 10*eps;

if (flow * fupp < 0)
    for i=1:maxIter
        xMid = round(0.5*(xLow+xUpp));
        if (i==1) xMid = xest % change to new estimate
        fMid = feval(func,xMid);
        s = sqrt(fMid*fMid - fLow*fUpp);
        if s==0.0, break; end
        xTmp = (xMid-xLow)*fMid/s;
        if fLow >= fUpp
            xNew = xMid + round(xTmp);
        else
            xNew = xMid - round(xTmp);
        end
        fNew = feval(func,xNew);
        if abs(fNew)<tol
            xZero = xNew;
            break
        else if (xupp - xlow) < 2
            if (abs(flow) < abs(fupp))
                xZero = xlow
            else
                xZero = xupp
            end
            break
        end

        
        %Update
        if sign(fMid) ~= sign(fNew)
            xLow = xMid;
            fLow = fMid;
            xUpp = xZero;
            fUpp = fNew;
        elseif sign(fLow) ~= sign(fNew)
            xUpp = xZero;
            fUpp = fNew;
        elseif sign(fUpp) ~= sign(fNew)
            xLow = xZero;
            fLow = fNew;
        else
            error('Something bad happened in riddersMethod!');
        end

    end
else
    if fLow == 0.0
        xZero = xLow;
    elseif fUpp == 0.0
        xZero = xUpp;
    else
        error('Root must be bracketed in Ridder''s Method!');
    end
end

Best Answer

It seems like your issue might be that you are rounding xNew into a point you have already computed rather than towards a new point. One simple fix would be to simply round towards the middle of your interval.

Another alternative option is to use tolerance appropriately. Although it is often unstated, it's usually good practice to apply a small $\Delta x$ to your root estimate so that it guarantees you will not just keep getting estimates which differ by less than your desired error (in this case an error of $1$). To handle this, one typically uses something like:

$$\hat x_3=x_1+(x_1-x_0)\frac{\operatorname{sign}(f(x_0))f(x_1)}{\sqrt{f(x_1)^2-f(x_0)f(x_2)}}$$

$$x_3=\operatorname{round}(\hat x_3+0.5\operatorname{sign}(x_1-\hat x_3))\tag{$\star$}$$

This will round $\hat x_3$ towards $x_1$, the middle of your interval. That way if you have something like $x_0=-46$ and $\hat x_3=-45.8$, you'll end up with $x_3=-45$ instead of $x_3=-46$ (which doesn't help), and in the example you showed this should give you almost immediate convergence.

$(\star):$ For non-integers, you wouldn't round the result, and instead of the $0.5$ term, you would use something like $10^{-14}$ for example.