function xn = babylonian(num)
% MEC E 390 Chapter 2
% Example 2.2
% Use the Babylonian method to find the square root of a number
% Inputs:
%   num = the number to approximate the square root of
%
% Outputs:
%   xn = the estimate of the square root
% 
% Constants:
%   es = 0.001
%   x0 = 1;

clc  % clear the command window

% Initial Estimate
x0 = 1;

% Stopping criteria
es = 0.001;

% Max Iterations
maxiter = 100; % To prevent an infinite loop

% Display the top labeling line
fprintf('\n  #      Value      Error\n');

% Display the first iteration
fprintf('%3g %10.5g %10.5s\n',0,x0,'-');

% Determine second estimate
xold=x0;
xnew=(xold+num/xold)/2;

er = (xnew-xold)/xnew;
fprintf('%3g %10.5g %10.5g\n',1,xnew,er);

i = 1; % Initializing counter

% Loop to run the Babylonian method
while abs(er)>es && i<maxiter,  % test absolute relative approx. error and max iterations
    i = i+1;                    % increment the counter
    xold = xnew;                % assign the current estimate to xold
    xnew = (xold+num/xold)/2;   % calculate the new estimate
    er = (xnew-xold)/xnew;      % calculate the relative approx. error
   
    fprintf('%3g %10.5g %10.5g\n',i,xnew,er);
end

% Display the final result
disp(' ');
disp(['Total number of iterations: ' num2str(i)]);

fprintf('%5s %5.5g %2s %5.4f %5s\n','Square root of ',num,' = ',xnew,'(4 decimal places)');    % With 4 decimal places
fprintf('%5s %5.5g %2s %5.4g %5s\n','Square root of ',num,' = ',xnew,'(4 sig. digits)');    % With 4 decimal places

% Return the estimate
xn = xnew;
    


