solving a non linear equation

I am trying to solve this equation (Colebrook). I have created a problem to solve it directly. However, that is not the method I can use for this exercise. This is program I solved it initially.

function [f]= moody(ed,Re)
if Re<0
  error(sprintf('Reynolds number = %f cannot be negative',Re));
elseif Re<2000
  f = 64/Re;  return      %  laminar flow
end
if ed>0.05
  warning(sprintf('epsilon/diameter ratio = %f is not on Moody chart',ed));
end
findf = inline('1.0/sqrt(f) + 2.0*log10( ed/3.7 + 2.51/( Re*sqrt(f)) )','f','ed','Re');
fi = 1/(1.8*log10(6.9/Re + (ed/3.7)^1.11))^2;   %  initial guess at f

I need to use this bisection method to solve the Colebrook equation. which this is the bisection method.

function bisection(f,a,b,tolerance,maxiter)
%bisection method for finding root of f(x)=0
%calculate function at the initial interval and mid-point
a0=a; b0=b;%initial bracket
iter=0;
error=abs(b0-a0);
f_a0=feval(f,a0);
f_b0=feval(f,b0);
if f_a0*f_b0>0 
         disp ('method does not work, change the initial bracket ')
         return
    else
%repeat the calculation until error<tolerance
        while  error>tolerance
            f_a=feval(f,a);
            f_b=feval(f,b);
                c=(a+b)*0.5;
            f_c=feval(f,c);
                if f_c*f_a<0
                    b=c;
                elseif f_b*f_c<0
                    a=c;
                end
            error=abs(b-a);
            iter=iter+1;
                if iter>maxiter
                 disp('no solution within maxiter. increase maxiter')
                    break
                end
         end
end
%
if iter<=maxiter
 disp(['iter=  ' num2str(iter)])
 disp(['root is   ' num2str(c)])
end
fplot(f,[a0,b0])
xlabel('x'),ylabel('f(x)');grid on
dfTol = 5e-6;
f = fzero(findf,fi,optimset('TolX',dfTol,'Display','off'),ed,Re);

Here is where I am at now.

function [f]=friction_factor(f)
Re=10e4;
ed=0.03;

f=(1.0/((sqrt(f)+2.0*log10((ed/3.7)+(2.51/(Re*sqrt(f)))))))^2;

end

I think there is something small that I am not seeing. I just went blank after writing the much. I would appreciate any help to get to the next steps.

thank you

[Roger Stafford]

It seems a shame that you should be messing around with iterative methods when you can obtain an explicit solution using the 'lambertw' function of matlab.

Answer 1

[Brand]

% Bisection Method in MATLAB
a=input('Enter function with right hand side zero:','s');
f=inline(a);

xl=input('Enter the first value of guess interval:') ;
xu=input('Enter the end value of guess interval:');
tol=input('Enter the allowed error:');

if f(xu)*f(xl)<0
else
    fprintf('The guess is incorrect! Enter new guesses\n');
    xl=input('Enter the first value of guess interval:\n') ;
    xu=input('Enter the end value of guess interval:\n');
end

fori=2:1000
xr=(xu+xl)/2;
if f(xu)*f(xr)<0
    xl=xr;
else
    xu=xr;
end

if f(xl)*f(xr)<0
    xu=xr;
else
    xl=xr;
end

xnew(1)=0;
xnew(i)=xr;
if abs((xnew(i)-xnew(i-1))/xnew(i))<tol,break,end
end
str = ['The required root of the equation is: ', num2str(xr), '']

Answer 2

[Brand]

% Bisection Method in MATLAB
a=input('Enter function with right hand side zero:','s');
f=inline(a);

xl=input('Enter the first value of guess interval:') ;
xu=input('Enter the end value of guess interval:');
tol=input('Enter the allowed error:');

if f(xu)*f(xl)<0
else
    fprintf('The guess is incorrect! Enter new guesses\n');
    xl=input('Enter the first value of guess interval:\n') ;
    xu=input('Enter the end value of guess interval:\n');
end

fori=2:1000
xr=(xu+xl)/2;
if f(xu)*f(xr)<0
    xl=xr;
else
    xu=xr;
end

if f(xl)*f(xr)<0
    xu=xr;
else
    xl=xr;
end

xnew(1)=0;
xnew(i)=xr;
if abs((xnew(i)-xnew(i-1))/xnew(i))<tol,break,end
end
str = ['The required root of the equation is: ', num2str(xr), '']

[Walter Roberson]

What is the difference between your two answers?