%% By Kittiphan Techakittiroj (13 Feb 1996)

%% Step input
timeend = 20;
figure(1);

h = 0.01;
[output1, time1] = rk4('input1', 'hw2_f', 'hw2_g', [0 0], h, timeend);

h = 0.1;
[output2, time2] = rk4('input1', 'hw2_f', 'hw2_g', [0 0], h, timeend);

h = 0.5;
[output3, time3] = rk4('input1', 'hw2_f', 'hw2_g', [0 0], h, timeend);

h = 1;
[output4, time4] = rk4('input1', 'hw2_f', 'hw2_g', [0 0], h, timeend);

h = 2;
[output5, time5] = rk4('input1', 'hw2_f', 'hw2_g', [0 0], h, timeend);

plot(time1, input1(time1), ...
     'y', time1, output1, 'r', time2, output2, 'g', ...
     time3, output3, 'b', time4, output4, 'm', ...
     time5, output5, 'c');
title('Input is u(t), different step size (0.01, 0.1, 0.5, 1, 2)');
xlabel('time(s)');  ylabel('displacement(m)');  axis;


%% Sin(0.1*t)*u(t)
timeend = 100;
figure(2);

h = 0.01;
[output1, time1] = rk4('input2', 'hw2_f', 'hw2_g', [0 0], h, timeend);

h = 0.1;
[output2, time2] = rk4('input2', 'hw2_f', 'hw2_g', [0 0], h, timeend);

h = 0.5;
[output3, time3] = rk4('input2', 'hw2_f', 'hw2_g', [0 0], h, timeend);

h = 1;
[output4, time4] = rk4('input2', 'hw2_f', 'hw2_g', [0 0], h, timeend);

h = 2;
[output5, time5] = rk4('input2', 'hw2_f', 'hw2_g', [0 0], h, timeend);

plot(time1, input2(time1), ...
     'y', time1, output1, 'r', time2, output2, 'g', ...
     time3, output3, 'b', time4, output4, 'm', ...
     time5, output5, 'c');
title('Input is sin(0.1*t), different step size (0.01, 0.1, 0.5, 1, 2)');
xlabel('time(s)');  ylabel('displacement(m)');  axis;


%% Sin(1*t)*u(t)
timeend = 10;
figure(3);

h = 0.01;
[output1, time1] = rk4('input3', 'hw2_f', 'hw2_g', [0 0], h, timeend);

h = 0.1;
[output2, time2] = rk4('input3', 'hw2_f', 'hw2_g', [0 0], h, timeend);

h = 0.5;
[output3, time3] = rk4('input3', 'hw2_f', 'hw2_g', [0 0], h, timeend);

h = 1;
[output4, time4] = rk4('input3', 'hw2_f', 'hw2_g', [0 0], h, timeend);

h = 2;
[output5, time5] = rk4('input3', 'hw2_f', 'hw2_g', [0 0], h, timeend);

plot(time1, input3(time1), ...
     'y', time1, output1, 'r', time2, output2, 'g', ...
     time3, output3, 'b', time4, output4, 'm', ...
     time5, output5, 'c');
title('Input is sin(1*t), different step size (0.01, 0.1, 0.5, 1, 2)');
xlabel('time(s)');  ylabel('displacement(m)');  axis;

%% Sin(10*t)*u(t)
figure(4);

h = 0.01;
[output1, time1] = rk4('input4', 'hw2_f', 'hw2_g', [0 0], h, timeend);

h = 0.1;
[output2, time2] = rk4('input4', 'hw2_f', 'hw2_g', [0 0], h, timeend);

h = 0.5;
[output3, time3] = rk4('input4', 'hw2_f', 'hw2_g', [0 0], h, timeend);

h = 1;
[output4, time4] = rk4('input4', 'hw2_f', 'hw2_g', [0 0], h, timeend);

h = 2;
[output5, time5] = rk4('input4', 'hw2_f', 'hw2_g', [0 0], h, timeend);

plot(time1, input4(time1), ...
     'y', time1, output1, 'r', time2, output2, 'g', ...
     time3, output3, 'b', time4, output4, 'm', ...
     time5, output5, 'c');
title('Input is sin(10*t), different step size (0.01, 0.1, 0.5, 1, 2)');
xlabel('time(s)');  ylabel('displacement(m)');  axis;


figure(5);

timeend = 20;
h = 0.01;
[output1, time1] = rk4('input2', 'hw2_f', 'hw2_g', [0 0], h, timeend);
[output2, time2] = rk4('input3', 'hw2_f', 'hw2_g', [0 0], h, timeend);
[output3, time3] = rk4('input4', 'hw2_f', 'hw2_g', [0 0], h, timeend);

plot(time1, input2(time1), 'y.', time1, output1, 'r-', ...
     time2, input3(time2), 'm.', time2, output2, 'g-', ...
     time3, input4(time3), 'c.', time3, output3, 'b-');
title('Input is sin(w*t), different frequency (0.1, 1, 10)');
xlabel('time(s)');  ylabel('displacement(m)');  axis;
