%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% Compute position and velocity of damped pendulum  with sinusoidal forcing
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

theta0		= 0.5;			% initial position
theta_dot0	= -1;			% initial velocity
tspan		= [0 20];		% time boundaries
x0		= [theta0;theta_dot0];

global alpha w02 f_input
alpha	= 1;
w02	= 20;

f_input = inline('sawtooth(5*t)');

[t,x]	= ode45('f4',tspan, x0);

x1	= x(:,1);
x2	= x(:,2);

subplot(3,1,1)
plot(t,x1);
xlabel('time in s'); ylabel('position (rad)'); 
title('simple damped pendulum with custom forcing');
subplot(3,1,2)
plot(t,x2);
xlabel('time in s'); ylabel('velocity (rad/s)'); 
subplot(3,1,3)
plot(t,f_input(t));
xlabel('time in s'); ylabel('input');