clear all
% numerical intergration of Ex 6.3.9 using Euler scheme
% for y>1 a small time step needs to used (cubic nonlinearity)

%timestep
dt = 0.0001;
%final time
T = 2.5;
%no of steps
tsteps = ceil(T/dt);
% initial condition
x0 = 1;
y0=  0.2;

time =0;
xn = x0;
yn = y0;

% save the discrete time and solution values
tdata = time;
xdata = x0;
ydata = y0;

%the Euler method
for n=1:tsteps,
    fnx = yn^3-4*xn;
    fny = yn^3-yn-3*xn;
    
    xn = xn + dt*fnx;
    yn = yn + dt*fny;
    
    time = dt*n;

	% save data
	tdata = [tdata time];
	xdata = [xdata xn];
    ydata = [ydata yn];
    
end;	

plot(xdata,ydata,'k')
xlabel(['\bf x'],'fontsize',12)
ylabel(['\bf y'],'fontsize',12)
    