function [] = ch00_fun ()
%%%%% ch01_fun.m %%%%%
% Run function and study the effects of each line.

disp(sprintf('**** Exercises for practicing MatLab ****\n'));
disp(sprintf('**** FUNCTIONS ****\n'));

%M-files can be either scripts or functions. 
%Scripts are simply files containing a sequence of MATLAB statements. 
%Functions make use of their own local variables and accept input arguments.

%A line at the top of a function M-file contains the syntax definition. 
%The name of a function, as defined in the first line of the M-file, should
%be the same as the name of the file without the .m extension. 

%For example, the existence of a file on disk called stat.m with

%    function [mean,stdev] = stat(x)
%    n = length(x);
%    mean = sum(x)/n;
%    stdev = sqrt(sum((x-mean).^2/n));

%defines a new function called stat that calculates the mean and standard 
%deviation of a vector x. The variables within the body of the function are 
%all local variables.
%Subfunctions can be placed also in a single file for a local use.

disp('   n         fact     Stirling      abs_err      rel_err');
ch01_3stir_err(1,10);

disp(sprintf('\n%%Also COMPUTER PROBLEM 1.1'))
disp('%The absolute error grows with n.')
disp('%The relative error shrinks with n.')
return;

function [] = ch01_3stir_err (min, max)

% calculate and plot the absolute and relative error of Stirling's 
% approximation for an integer interval [min, max]

x=[min:max];
y_stir_aerr=[min:max];
y_stir_rerr=[min:max];

for i=min:max 
   [y_stir_aerr(i-min+1),y_stir_rerr(i-min+1)]= calc_err (i); 
end

semilogy(x,y_stir_aerr,'r-',x,y_stir_rerr,'b-');
title(['Absolute and Relative error of of Stirling''s approximation'])
xlabel('integer interval')
ylabel('error')
text(6,10^3,'absolute error')
text(7.5,10^-1,'relative error')
return;


function [stir_aerr,stir_rerr] =  calc_err (n)
% calculate and return absolute and relative error of Stirling's approximation

fact=1;
for i=1:n fact=i*fact; end

stir= sqrt(2*pi*n)*(n/exp(1))^n;

stir_aerr= abs(stir-fact);
stir_rerr= stir_aerr/fact;

disp(sprintf('%4g %12g   %10.3e   %10.3e   %10.3e', n, fact, stir, stir_aerr, stir_rerr'));
return;