function fig312

nt=100;
tmax=120;
t=linspace(0,tmax,nt);

% get(gcf)
%set(gcf,'Position', [831 504 515 280]);
set(gcf,'Position', [1665 1089 638 257]);

ep=0.1;
for it=1:nt
	k=2+cos(ep*t(it));
	ya(it)=1/sqrt(3*k);
end;

%  initial values

y10=0; y20=1; 
y0=[y10 y20];

%  calculate solution using a MATLAB routine
[tt,y] = ode45(@rhs,[0 tmax],y0);  titleX='ode45';
%[t,y] = ode23s(@rhs,[0 tmax],y0);  titleX='ode23s';

hold on
box on
grid on

plot(t,ya,'--','Linewidth',1.3)
plot(tt,y(:,1),'-r','Linewidth',1)

axis([0 tmax -0.601 0.6])
%loc='NorthWest';
loc='SouthEast';

xlabel('t-axis','FontSize',14,'FontWeight','bold')
ylabel('Solution','FontSize',14,'FontWeight','bold')

set(gca,'FontSize',14);
%legend(' Exact Solution',' y_0(t)+\epsilon y_1(t)','Location',loc);
set(findobj(gcf,'tag','legend'),'FontSize',14); 


%  define f(t,y)
function dy=rhs(t,y)
dy=zeros(2,1);
ep=0.1; 
k=2+cos(ep*t);
dy(1) = y(2);
dy(2) = -k*k*y(1);