%MAE 82 
%Project 4 solution

clear all
close all
%% PROBLEM #1
% 
% %Part 1
% 
% t1 = 0:2*pi/100:6*pi;
% 
% y1 = -0.5.*t1.*sin(t1);
% 
% figure(1)
% plot(t1,y1)
% xlabel('Time [s]')
% ylabel('Deflection y(t)')
% title('Problem #1 - Part 1')
% 
% t2 = 0:2*pi/100:6*pi;
% 
% y2 = (1/3).*cos(t2)-(1/3).*cos(2.*t2);
% 
% figure(2)
% plot(t2,y2)
% xlabel('Time [s]')
% ylabel('Deflection y(t)')
% title('Problem #1 - Part 2')

%% PROBLEM #2

% %Part a
% b = 1;
% a = 4;
% alpha = 1;
% 
% %phase 1
% ta1= 0:pi/100:pi;
% ya1 = (1+4/15).*sin(ta1)-(1/15).*sin(4.*ta1);
% %phase 2
% ta2 = pi:pi/100:3*pi/2;
% ya2 = -(9/15).*sin(2.*ta2)-(1/12).*sin(4.*ta2);
% 
% %find velocity at beginning of phase 3
% NewVelocity = -(18/15).*cos(2*3*pi/2)-(4/12).*cos(4*3*pi/2)
% figure(3)
% plot(ta1,ya1)
% xlabel('Time [s]')
% ylabel('Deflection y(t)')
% title('Problem #2 - Part a')
% hold on
% plot(ta2,ya2)
% hold off
% %Part b
% b = 64;
% a = 4;
% %phase 1
% tb1= 0:pi/1000:pi/7.5;
% yb1 = (11/96).*sin(8.*tb1)+(1/48).*sin(4.*tb1);
% NewICb1 = (88/96)*cos(8*pi/7.5) + (4/48)*cos(4*pi/7.5)
% %phase 2
% t = pi/7.5;
% A = [cos(2*t) sin(2*t);
%     -2*sin(2*t) 2*cos(2*t)];
% d = [(1/12)*sin(4*t);(1/3)*cos(4*t)+NewICb1];
% x = zeros(2,1);
% x = A\d;
% c1 = x(1)
% c2 = x(2)
% tt = 2.143;
% tb2 = t:pi/10000:tt;
% yb2 = c1.*cos(2.*tb2) + c2.*sin(2*tb2)-(1/12).*sin(4*tb2);
% 
% %find velocity at beginning of phase 3
% NewVelocity = -2*c1*sin(2*tt)+2*c2*cos(2*tt)-(1/3)*cos(4*tt)
% figure(4)
% plot(tb1,yb1)
% xlabel('Time [s]')
% ylabel('Deflection y(t)')
% title('Problem #2 - Part b')
% hold on
% plot(tb2,yb2)
% %Part c
% b = 36;
% a = 25;
% %phase 1
% tc1= 0:pi/10000:.571;
% yc1 = (2/15).*sin(6.*tc1)+(1/20).*sin(4.*tc1);
% t=.571;
% NewICc1 = (12/15)*cos(6*t) + (1/5)*cos(4*t)
% %phase 2
% A = [cos(5*t) sin(5*t);
%     -5*sin(5*t) 5*cos(5*t)];
% d = [-(1/9)*sin(4*t);-(4/9)*cos(4*t)+NewICc1];
% x = zeros(2,1);
% x = A\d;
% c1 = x(1)
% c2 = x(2)
% tt = 1.241;
% tc2 = t:pi/10000:tt;
% yc2 = c1.*cos(5.*tc2) + c2.*sin(5*tc2)+(1/9).*sin(4*tc2);
% % %find velocity at beginning of phase 3
% NewVelocity = -5*c1*sin(5*tt)+5*c2*cos(5*tt)-(4/9)*cos(4*tt)
% figure(5)
% plot(tc1,yc1)
% xlabel('Time [s]')
% ylabel('Deflection y(t)')
% title('Problem #2 - Part c')
% hold on
% plot(tc2,yc2)


%% PROBLEM #3

%First Phase
c = .1;
M = abs(sqrt((c^2/4)-4));
N = 4*c^2 + 36;
A = -c/N;
B = -3/N;

c1 = -A;
c2 = (1 - c*A/2 - 4*B)/M;

t = 0:pi/100:pi/2;

y = exp(-c.*t./2).*(c1.*cos(M.*t) + c2.*sin(M.*t)) + A.*cos(4.*t) + B.*sin(4.*t);

newIC = exp(-c*pi/4)*(c1*c*M/2 +c2*M);
figure(5)
plot(t,y)
xlabel('Time [s]')
ylabel('Deflection y(t)')
title('Problem #3 - c = 0.1')
hold on

%Second Phase

M = abs(sqrt((c^2/4)-1));
N = 16*c^2 + 180;
A = -5*c/N;
B = -15/N;

c2 = -(A/M)*exp(c*pi/4);
c1 = 4*B*exp(c*pi/4)-c*c2/2/M - newIC*exp(c*pi/4)/M;
tt=3*pi/2;
t = pi/2:pi/100:tt;

y = -exp(-c.*t./2).*(c1.*cos(M.*t) + c2.*sin(M.*t)) + A.*cos(4.*t) + B.*sin(4.*t);
t1=tt-.01;
y1 = -exp(-c.*t1./2).*(c1.*cos(M.*t1) + c2.*sin(M.*t1)) + A.*cos(4.*t1) + B.*sin(4.*t1);
y2 = -exp(-c.*tt./2).*(c1.*cos(M.*tt) + c2.*sin(M.*tt)) + A.*cos(4.*tt) + B.*sin(4.*tt);
NewVelocity = (y2-y1)/(tt-t1)
plot(t,y)