% Lorenz attractor

function f = lorfun(t,s)

sigma = 10;
rho = 28;
beta = 8/3;

x = s(1,1);
y = s(2,1);
z = s(3,1);

f(1,1) = -beta*x + y*z;
f(2,1) = sigma*(z-y);
f(3,1) = -x*y + rho*y - z;