% 16.07, Fall 2004
% Undamped spring-mass example
% Last updated September 7th, 2004
%
% State vector and state derivative vector:
% X' and X are vectors in --column format--!!!!!!
% X  = [ x
%        x_dot ]
% X' = [ x_dot
%        x_dot_dot ]

% set some tolerances for the MATLAB integrator ode45
options = odeset('RelTol', 1e-4, 'AbsTol', 1e-4, 'MaxStep', .05);

% initial conditions
X0 = [ 0.0;             % mass position [m]
       0.0 ];           % mass velocity [m/s]
   
% how long do we run for?
tspan = [0, 400];       % starting and ending time [s]


% run the sim and store the solution
[t, X] = ode45(@getDerivs, tspan, X0, options);

% plot mass position vs. time
figure;
plot(t, X(:,1));
xlabel('time [s]');
ylabel('position [m]');
title('Mass Position vs. Time');
grid on;

% plot mass velocity vs. time
figure;
plot(t, X(:,2));
xlabel('time [s]');
ylabel('velocity [m/s]');
title('Mass Velocity vs. Time');
grid on;


% replicate forcing function vs. time
A = 1;                              % from the state derivative function
forcing = zeros(1, length(t));
for count = 1:length(t)
    if(t(count) < 100)
        omega = .1;
    elseif(t(count) < 150)
        omega = 1;
    else
        omega = 1 + .1*(t(count) - 150);
    end;    
    
    forcing(count) = A*sin(omega*t(count));
end

% plot forcing function vs. time
figure;
plot(t, forcing)
xlabel('time [s]');
ylabel('force [N]');
title('Forcing Function vs. Time');
grid on;