pd.m

he basic idea of the method of bisection is to start with an initial interval, [a0,b0], that is chos

%
%  CHBE 230 Computational Methods
%  Assignment 2     
%  Date: February 14th 2007
%            
%  By:   Moetamed Shariati - Saba
%  St #: 50233055     
%


function SABA

L=200 ; D=6 ;                % ft , in

Ft_To_Meter_CF= 0.3048 ; In_To_Meter_CF = 0.0254 ; % Conversion factors

h=0.01;
TOL= 0.0001;

p= 1000 ;                    % KG/(M^3)  Water density
Mu = 0.001 ; E = 0.0004 ;    % Constants
g  = 9.81 ;                  % acceleration due to gravity

D = D * In_To_Meter_CF;      % converting unit of diameter of the pipe
L = L * Ft_To_Meter_CF;      % converting unit of length of the pipe
E = E * Ft_To_Meter_CF;      % converting unit of surface roughness of the 
                             % pipe
Vl = 0.001 ; Vu = 3;         % Lower and upper velocities  ( m^3/s )

Inc_Num=10;                  % Number of points in the overall interval

B=[];                        % a Null matrix to save subintervals found by 
                             % incrememntal search method
                             
ROOTS=[];                    % a null matrix to save the roots found by secant method.

PDROP = 10600 ;              % Pressure drop Pascal  ( Pa )
                             % defining the function which we want to find
                             % its roots.
                             
Z =@(V) PDROP- p.*g.*L./(2.*D.*g).* ( V.^2).* (-1.8 .* log10( (6.9./(p.*V.*D./Mu))+((E./D)./3.7).^1.11)).^(-2);

Num_ROOT = 0 ;               % initial number of roots

U=linspace (Vl,Vu,Inc_Num);  % initiation the intervals for search incremental search method

                             % incremental search method:
                             % output is number of roots and
                             % the intervals where the roots are located
                             % in.
                             
for i= 1:(Inc_Num-1)
    if ((Z(U(i)).*Z(U(i+1)))<= (0))
        Num_ROOT = Num_ROOT +1;
        B(Num_ROOT,1)=U(i); 
        B(Num_ROOT,2)=U(i+1);
       
    end
    
end
                             % finding the exact roots, by using secant
                             % method.
for K = 1:Num_ROOT
    L=B(K,1); U=B(K,2);
    M=(L+U)/2;
    N=L;
    ROOT=M-((Z(M)*(M-N))/(Z(M)-Z(N)));
    while ( abs(ROOT-M) > TOL)
        N=M;
        M=ROOT;
        ROOT=M-((Z(M)*(M-N))/(Z(M)-Z(N)));
    end
    ROOTS(K)=ROOT;
end
ROOTS
end


-----------------------------------


ROOTS =

    1.6353