function r = linearFractInterpolationF(fcn,a,b,c,tol,itMax)
% Root finding using linear fractional interpolation
% a, b, & c are 3 initial guesses to the root.
% See p.235 of Heath, example 5.14
% K. Ming Leung, 03/13/03

%clear all;fcn='x2m4sinx'; format long;
fa=feval(fcn,a);
fb=feval(fcn,b);
fc=feval(fcn,c);
it=1; h=10; 
r=[];       % to store intermediate results
while it <= itMax & abs(h) > tol;
    ac=a-c; bc=b-c;
    h = fc*ac*bc*(fa-fb)/...
        (fa*ac*(fc-fb)-fb*bc*(fc-fa));
    a=b; fa=fb;     % a replaced by old b
    b=c; fb=fc;     % b replaced by old c
    c=c+h; fc=feval(fcn,c); % new c
    r=[r; [it c fc]]; it=it+1;
end;


