clear all;
close all;

% Parameters of the system:
a = 1;
b = 5;
h = 0.1;
rInterval = a : h : b;
Y0 = [0;1]; 

% Solving the system by means of MATLAB ODE solver
option = odeset('AbsTol',10^(-10),'RelTol',10^(-10));

omega = 0.7630790233;
[r,y] = ode23(@system,rInterval,Y0,option,omega);
figure(1); 
plot(r,y(:,1)); xlabel('r'); ylabel('first eigenvector'); 
error = abs(y(length(r),1));
fprintf('Error of the first eigenvalue is %10.9f\n',error);

omega = 1.5556700046;
[r,y] = ode23(@system,rInterval,Y0,option,omega);
figure(2); 
plot(r,y(:,1)); xlabel('r'); ylabel('second eigenvector'); 
error = abs(y(length(r),1));
fprintf('Error of the second eigenvalue is %10.9f\n',error);

omega = 2.3412414511;
[r,y] = ode23(@system,rInterval,Y0,option,omega);
figure(3); 
plot(r,y(:,1)); xlabel('r'); ylabel('third eigenvector'); 
error = abs(y(length(r),1));
fprintf('Error of the third eigenvalue is %10.9f\n',error);

omega = 3.1214853270;
[r,y] = ode23(@system,rInterval,Y0,option,omega);
figure(4); 
plot(r,y(:,1)); xlabel('r'); ylabel('fourtheigenvector'); 
error = abs(y(length(r),1));
fprintf('Error of the fourth eigenvalue is %10.9f\n',error);

omega = 3.8961027211;
[r,y] = ode23(@system,rInterval,Y0,option,omega);
figure(5); 
plot(r,y(:,1)); xlabel('r'); ylabel('fifth eigenvector'); 
error = abs(y(length(r),1));
fprintf('Error of the fifth eigenvalue is %10.9f\n',error);
