function dydt=Lorenz(t,y)
%page 334, Lorenz equation
s=10; r=28; b=8/3;
dydt=[-s*y(1)+s*y(2);
      -y(1)*y(3)+r*y(1)-y(2);
      y(1)*y(2)-b*y(3)];
  