function lorenza

% Needs R to be defined as a global variable in the command window.

global A B R

A = 10; B = 8/3; 

close all  % close all figures.

[t,u] = ode45(@lore,[0,100],30*rand(1,3));

% Time plots
figh(1) = figure;
dfigpos = get(figh(1),'pos');
height = dfigpos(4);
width = dfigpos(3);
left = 300;
bottom = 400;
pos = [left,bottom,width, height];
set(figh(1),'pos',pos);
plot(t,u)
xlabel('t')
ylabel('x, y, & z')
title('Time plots')
legend('x','y','z')

% Phase space plot
figh(2) = figure;
left = left + 50;
bottom = bottom - 50;
pos = [left,bottom,width, height];
set(figh(2),'pos',pos);
plot3(u(:,1),u(:,2),u(:,3))
xlabel('x')
ylabel('y')
zlabel('z')
title('Phase space plot')

% Eliminate transients to show the limiting behavior.
figh(3) = figure;
left = left + 50;
bottom = bottom - 50;
pos = [left,bottom,width, height];
set(figh(3),'pos',pos);
N = find(t>50);
plot3(u(N,1),u(N,2),u(N,3))
xlabel('x')
ylabel('y')
zlabel('z')
title('Limiting behavior (t > 50)')


function upr = lore(t,u)

global A B R

x = u(1);
y = u(2);
z = u(3);

xpr = A*(y-x);
ypr = R*x - y - x*z;
zpr = -B*z + x*y;

upr = [xpr;ypr;zpr];

