function yprime = ozone(t,y);
% ozone(t,y) returns rate of change of chemistry of atmosphere
%	equation.  ozone-dr in ozone folder

k1=1.6E-14; k2=2.0E-11; k3=3.6E-13; k4=1.3E-14; k5=1.4E-10;

yp1 = -k1*y(1)*y(2) ...
      +k2*y(3)*y(5) ...
      -k4*y(1)*y(6) ...
      -k5*y(1)*y(8);
      
yp2 = -k1*y(1)*y(2); 

yp3 = -k2*y(3)*y(5) ...
      +k5*y(1)*y(8); 
      
yp4 =  k1*y(1)*y(2) ...
       +k2*y(3)*y(5);
       
yp5 =  k1*y(1)*y(2) ...
      -k2*y(3)*y(5) ...
      -k3*y(5)*y(6); 
      
yp6 = -k3*y(5)*y(6) ...
      -k4*y(1)*y(6) ...
      +k5*y(1)*y(8);

yp7 =  k3*y(5)*y(6);

yp8 =  k4*y(1)*y(6) ...
      -k5*y(1)*y(8);
      
% return a column vector
yprime = [ yp1; yp2; yp3; yp4; yp5; yp6; yp7; yp8];