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