Recreating Blended Bisection and False Position algorithm from reference

Question

I'm trying to recreate a blended bisection algorithm (Algorithm 3) from the website below (link takes you to exact section of the algorithm I'm referencing)

https://www.mdpi.com/2227-7390/7/11/1118/htm#sec3-mathematics-07-01118

I'm not quite sure if what I've typed out currently is correct and I'm stuck on line 29 of the algorithm from the website where I'm not sure what it means especially with the intersection symbol.

Code so far

/* Math function to test on */
function fn(x) {
    //return x * x - x - 2; /* root for this is x = 2 */
    return x*x*x-2; /* root for this is x = (2)^(1/3) */
}

function blendedMethod(a, b, eps, maxIterations, fn) {
    let k = 0,
        r, fa, fb, ba, bb, eps_a;

    do {
        let m = (a + b) * .5;
        let eps_m = Math.abs(fn(m));

        let fn_a = fn(a),
            fn_r;
        let s = a - ((fn_a * (b - a)) / (fn(b) - fn_a));
        let eps_s = Math.abs(fn(s));

        if (eps_m < eps_s) {
            r = m;
            fn_r = fn(r);
            eps_a = eps_m;
            if (fn_a * fn_r < 0) {
                ba = a;
                bb = r;
            } else {
                ba = r;
                bb = b;
            }

        } else {
            r = s;
            fn_r = fn(r)
            eps_a = eps_s;
            if (fn_a * fn_r < 0) {
                fa = a;
                fb = r;
            } else {
                fa = r;
                fb = b;
            }

            /* line 29 here! */
            /* [a, b] = [ba, bb] ∩ [fa, fb] */
            /* either fa,fb or ba,bb haven't yet been defined */
            /* so this will fail with NaN */
            a = Math.max(ba, fa);
            b = Math.min(bb, fb);
        }

        r = r;
        eps_a = Math.abs(fn_r)
        k = k + 1;

    } while (Math.abs(fn(r)) > eps || k < maxIterations)

    /* think this is the root to return*/
    return r;
}

console.log(blendedMethod(1,4,0.00001,1000,fn));

EDIT: Fixed some errors, only problem is that this algorithm defines either fa,fb or ba,bb inside the conditional statements without defining the other two. So by the time it comes to these calculations below, it fail with NaN and messes up for the next iterations. a = Math.max(ba,fa); b = Math.min(bb,fb);

∩ generally denotes intersection. So it seems to be saying you need the intersection of [ba, bb] and [fa,fb]. The paper implies that every [a,b] is an interval. Does it makes sense to say [ba, bb] ∩ [fa, fb] is equivalent to `[Math.max(ba, fa), Math.min(bb, fb)]` — rfestag, Jul 18 '20 at 15:39
I think you're correct, only problem is that this algorithm defines either fa,fb or ba,bb inside the conditional statements without defining the other two. So by the time it comes to these calculations below, it fail with NaN and messes up for the next iterations. a = Math.max(ba,fa); b = Math.min(bb,fb); — C9C, Jul 18 '20 at 16:23
The paper is re-inventing the wheel. The regula falsi method can be improved using weighted means like in the Illinois variant. A blend of secant and bisection was developed by Dekker some 50 years ago and became famous as the `fzeroin` algorithm in matlab and other libraries. — Lutz Lehmann, Jul 18 '20 at 16:50
@LutzLehmann Interesting! Have you any links to the source code of both algorithms either in JS, Java, Python, or C#. I would be curious to check them out! — C9C, Jul 18 '20 at 16:55
@C9C - True, although you can use the `min`, `max` conditionally (i,e, only if you have ba, bb, fa, and fb). Otherwise, use whichever pair you have (since you should have at least one of the pairs). — rfestag, Jul 18 '20 at 16:56
Some discussion with links to original sources in https://scicomp.stackexchange.com/questions/25005/dekkers-method-and-fixed-further-border/29711#29711, for the Illinois regula falsi against the basic method https://stackoverflow.com/questions/22273751/regula-falsi-algorithm/22284632#22284632 — Lutz Lehmann, Jul 18 '20 at 17:34

score 1 · Accepted Answer · answered Jul 20 '20 at 17:56

You are right in that this intersection makes no sense. There is in every step only one sub-interval defined. As all intervals are successive nested subsets, the stale, old values of the interval that was not set in the current loop is still a superset of the new interval. The new interval could be directly set in each branch. Or the method selection branch could be totally separated from the interval selection branch.

The implementation is not very economic as 6 or more function values are computed where only 2 evaluations are needed. The idea being that the dominating factor in the time complexity are the function evaluations, so that a valid metric for any root finding algorithm is the number of function evaluations. To that end, always keep points and function value as pair, generate them as a pair, assign them as a pair.

let fn_a =f(a), fn_b=f(b)
do {
    let m = (a + b) * .5;
    let fm = f(m);
    let s = a - (fn_a * (b - a)) / (fn_b - fn_a)
    let fn_s = f(s);
    let c,fn_c; 
    // method selection
    if(Math.abs(fn_m) < Math.abs(fn_s)) {
        c = m; fn_c = fn_m;
    } else {
        c = s; fn_c = fn_s;
    }
    // sub-interval selection
    if( fn_a*fn_c > 0 ) {
        a = c; fn_a = fn_c;
    } else {
        b = c; fn_b = fn_c; 
    }
while( Math.abs(b-a) > eps );

It is also not clear in what way the blended method avoids or alleviates the shortcomings of the basis algorithms. To avoid the stalling (deviation from a secant step) of the regula falsi method it would be better to introduce a stalling counter and apply a bisection step after 1 or 2 stalled steps. Or just simply alternate the false position and bisection steps. Both variants ensure the reduction of the bracketing interval.

Known effective modifications of the regula falsi method are on one hand the variations like the Illinois variant that add a weight factor to the function values, thus shifting the root approximation towards the repeated, stalled interval bound. On the other hand there are more general algorithms that combine the ideas of the bracketing interval and reverse interpolation like the Dekker and Brent methods.

Recreating Blended Bisection and False Position algorithm from reference

1 Answers1