function [] = harm_osc (x0,v0,time)

% solve the system x'= v; v'=-k/m*x; 
% plot position x(t) and velocity v(t);
% plot trajectory of equi-energy points (x(t),v(t)) in the plane as 
%                    a function of time.

options = odeset('RelTol',5e-3,'AbsTol',[5e-3 5e-3]);
[T,Y] = ode45('osc',[0 time],[x0 v0],options);
figure(1);
plot(T,Y(:,1),'b-',T,Y(:,2),'ro')
title('x ---  v ..., Tolerance 0.005')
figure(2);
plot(Y(:,1),Y(:,2),'-.')
title('Phase plot (x,v), Tolerance 0.005')

options = odeset('RelTol',5e-2,'AbsTol',[5e-2 5e-2]);
[T,Y] = ode45('osc',[0 time],[x0 v0],options);
figure(3);
plot(T,Y(:,1),'b-',T,Y(:,2),'ro')
title('x ---  v ..., Tolerance 0.05')
figure(4);
plot(Y(:,1),Y(:,2),'-.')
title('Phase plot (x,v),  Tolerance 0.05')


%in connection with slides: molec_dyna.ppt