function prob_2_20
% 
% Written by:
% -- 
% John L. Weatherwax                2005-05-07
% 
% email: wax@alum.mit.edu
% 
% Please send comments and especially bug reports to the
% above email address.
% 
%-----
  
close all;   
  
addpath( fullfile( pwd, 'Prob_2_20_Files' ) );

global tcond cond;

% x_1 = y(1), x_2 = y(2), y = y(3)
options = odeset('AbsTol',1e-20);
[t,y] = mode15s(@odes,[0 8e5],[0; 1; 0],options);
plot(t,y(:,1:2))
figure
semilogx(t,y(:,3))

figure;
loglog(tcond,cond,'*-');
title( 'Time V.S. Condition Number' );

figure; 
loglog(t(1:end-1),diff(t));

%===================================================
function dydt = odes(t,y)
k = [8.4303270e-10 2.9002673e+11 2.4603642e+10 8.7600580e-06];
dydt = zeros(3,1);
dydt(1) = -k(1)*y(1) + k(2)*y(3);
dydt(2) = -k(4)*y(2) + k(3)*y(3);
dydt(3) =  k(1)*y(1) + k(4)*y(2) - (k(2) + k(3))*y(3);