function xdot=reac4(t,x,e0,s0,c0,k,r,d)
xdot(1)=k*(e0-x(1))*x(2)-(r+d)*x(1);
xdot(2)=-k*(e0-x(1))*x(2)+r*x(1);

