function dxdt = catalysis(t,x)
% USAGE: dxdt = catalysis(t,x)
global k1 km1 k2 e0
dxdt = zeros(size(x));

s = x(1);
e = x(2);
c = x(3);
p = x(4);

dxdt(1) = -k1*e*s + km1*c;
dxdt(2) = -k1*e*s + km1*c + k2*c;
dxdt(3) =  k1*e*s - km1*c - k2*c;
dxdt(4) =  k2*c;

return