% step3s.m
% solve the linearized system
% simplified version 

beta, x0, t1, u0   % these data must be given before the computation starts

arr=arrhen(Xbar(2)); 
arrprime=arr*gamma*Xbar(2)^(-2);
A=[-1-arr -Xbar(1)*arrprime; alpha*arr -1-beta+alpha*Xbar(1)*arrprime] 
lambda=eig(A)   % check the eigenvalues
B=[0;beta]; 
Ubar=Xbar(2)-(1-Xbar(2)+alpha*Xbar(1)*arr)/beta 

[t,x]=ode45('linear',[0 t1],x0);  % solve equations 

figure(3), clf % clear figure
plot(t,x)                         % solution curves

figure(4), clf % clear figure
plot(x(:,1),x(:,2))               % phase portrait 

