9. MATLAB Code: F-14 Wing Sweep
% ME380R - Robot Mechanism Design - Final Project
% Michael Estridge
% 4/27/2020
clc;
clear all;
close all;
global f
f = 1;
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%%
% F-14 Wing Sweep Analysis:
% Linear Actuator Gear Train Analysis
tf = 60; % end time, s
dt = 1; % time step, s
t = 0:dt:tf; % time vector
Win = 7400*sind(3*t); % motor input velocity, rpm
N_in1 = 1; % equivalent worm gear teeth
N_out1 = 20; % equivalent output gear teeth
TPI = .125; % lead screw thread per inch
GR =N_in1/N_out1*TPI; % gear ratio
Win2 = Win/60; % motor input velocity, rot/s
a_dot = Win2*GR; % linear actuation velocity, in/s
a_ddot(1) = 0; % initial actuation acceleration, in/s^
a(1) = 3.75*12; % initial actuator link length, in
for i=1:numel(t)-1
a(i+1) = a(i)+a_dot(i+1)*dt; % linear actuator link length, in
a_ddot(i+1) = (a_dot(i+1)-a_dot(i))/dt; % linear actuation acceleration, in/s^2
end
% 3 Link Mechanism Analysis
b = 3*12; % wing link, in
c = 6*12; % ground link, in
t0 = 170; % mechanism offset angle WRT o2(deg)
p = 8*12; % wing COM position relative to wing-ground pin
tp = -65; % wing COM angle relative to b
% t = 0:1:numel(a)-1; % time, s
% Position:
[Ax,Ay,t2,t3] = wing_sweep_plot(a,b,c,t0,p,tp);
figure(f)
plot(t,Win)
title('Motor Speed vs Time')
xlabel('Time (s)')
ylabel('Motor Speed (rpm)')
f = f+1;
figure(f)
plot(t,a)
title('Actuator Link Length vs Time')
xlabel('Time (s)')
ylabel('Input Link Length (in)')
f = f+1;
figure(f)
plot(t,a_dot)
title('Actuator Speed vs Time')
xlabel('Time (s)')
ylabel('Linear Actuator Speed (in/s)')
f = f+1;
figure(f)
plot(a,t2)
title('Actuator Angle vs Length')
xlabel('Input Link Length (in)')
ylabel('Angle (deg)')
f = f+1;
figure(f)
plot(t,t2)
title('Actuator Angle vs Time')
xlabel('Time (s)')
ylabel('Angle (deg)')
f = f+1;
figure(f)
plot(a,t3)
title('Output Angle vs Actuator Length')
xlabel('Input Link Length (in)')
ylabel('Angle (deg)')
f = f+1;
figure(f)
plot(t,t3)
title('Output Angle vs Time')
xlabel('Time (s)')
ylabel('Angle (deg)')
f = f+1;
figure(f)
plot(a,t3+109)
title('Sweep Angle vs Actuator Length')
xlabel('Input Link Length (in)')
ylabel('Angle (deg)')
f = f+1;
figure(f)
plot(t,t3+109)
title('Sweep Angle vs Time')
xlabel('Time (s)')
ylabel('Angle (deg)')
f = f+1;
% Velocity:
[Wa,Wb,V_A,V_COM] = wing_sweep_vel(a,b,p,t2,t3,a_dot);
figure(f)
plot(a,Wa)
title('Input Link Rotational Velocity vs Actuator Length')
xlabel('Input Link Length (in)')
ylabel('Rotational Velocity (rad/s)')
f = f+1;
figure(f)
plot(t,Wa)
title('Input Link Rotational Velocity vs Time')
xlabel('Time (s)')
ylabel('Rotational Velocity (rad/s)')
f = f+1;
figure(f)
plot(a,Wb)
title('Output Link Rotational Velocity vs Actuator Length')
xlabel('Input Link Length (in)')
ylabel('Rotational Velocity (rad/s)')
f = f+1;
figure(f)
plot(t,Wb)
title('Output Link Rotational Velocity vs Time')
xlabel('Time (s)')
ylabel('Rotational Velocity (rad/s)')
f = f+1;
figure(f)
plot(a,V_A)
title('Point A Velocity vs Actuator Length')
xlabel('Input Link Length (in)')
ylabel('Linear Velocity (in/s)')
f = f+1;
figure(f)
plot(t,V_A)
title('Point A Velocity vs Time')
xlabel('Time (s)')
ylabel('Linear Velocity (in/s)')
f = f+1;
figure(f)
plot(a,V_COM)
title('Wing COM Velocity vs Actuator Length')
xlabel('Input Link Length (in)')
ylabel('Linear Velocity (in/s)')
f = f+1;
figure(f)
plot(t,V_COM)
title('Wing COM Velocity vs Time')
xlabel('Time (s)')
ylabel('Linear Velocity (in/s)')
f = f+1;
% Mechanical Advantage
for i = 1:numel(a)
Y(i) = 180+t3(i)-(t2(i)-180);
MA(i) = sind(Y(i))*GR^-1;
end
figure(f)
plot(a,MA)
title('Mechanical Advantage vs Actuator Length')
xlabel('Input Link Length (in)')
ylabel('MA')
f = f+1;
figure(f)
plot(t,MA)
title('Mechanical Advantage vs Time')
xlabel('Time (s)')
ylabel('MA')
f = f+1;
% Acceleration
[Aa,Ab,A_A,A_COM] = wing_sweep_accel(a,b,p,t2,t3,a_dot,Wa,Wb,a_ddot);
figure(f)
plot(a,a_ddot)
title('Input Link Linear Acceleration vs Actuator Length')
xlabel('Input Link Length (in)')
ylabel('Linear Acceleration (in/s^2)')
f = f+1;
figure(f)
plot(t,a_ddot)
title('Input Link Linear Acceleration vs Time')
xlabel('Time (s)')
ylabel('Linear Acceleration (in/s^2)')
f = f+1;
figure(f)
plot(a,Aa)
title('Input Link Rotational Acceleration vs Actuator Length')
xlabel('Input Link Length (in)')
ylabel('Rotational Acceleration (rad/s^2)')
f = f+1;
figure(f)
plot(t,Aa)
title('Input Link Rotational Acceleration vs Time')
xlabel('Time (s)')
ylabel('Rotational Acceleration (rad/s^2)')
f = f+1;
figure(f)
plot(a,Ab)
title('Output Link Rotational Acceleration vs Actuator Length')
xlabel('Input Link Length (in)')
ylabel('Rotational Acceleration (rad/s^2)')
f = f+1;
figure(f)
plot(t,Ab)
title('Output Link Rotational Acceleration vs Time')
xlabel('Time (s)')
ylabel('Rotational Acceleration (rad/s^2)')
f = f+1;
figure(f)
plot(a,A_A)
title('Point A Acceleration vs Actuator Length')
xlabel('Input Link Length (in)')
ylabel('Linear Acceleration (in/s^2)')
f = f+1;
figure(f)
plot(t,A_A)
title('Point A Acceleration vs Time')
xlabel('Time (s)')
ylabel('Linear Acceleration (in/s^2)')
f = f+1;
figure(f)
plot(a,A_COM)
title('Wing COM Acceleration vs Actuator Length')
xlabel('Input Link Length (in)')
ylabel('Linear Acceleration (in/s^2)')
f = f+1;
figure(f)
plot(t,A_COM)
title('Wing COM Acceleration vs Time')
xlabel('Time (s)')
ylabel('Linear Acceleration (in/s^2)')
f = f+1;
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%%
% Rotation Function:
function [Xr,Yr] = rotate_Z(X,Y,t0)
Rz = [cosd(t0) -sind(t0);sind(t0) cosd(t0)];
R = [X;Y];
R = Rz*R;
Xr = R(1);
Yr = R(2);
end
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%%
% Wing Sweep Plot Function:
function [Ax,Ay,t2,t3] = wing_sweep_plot(a,b,c,t0,p,tp)
global f
div = 10;
B = linspace(0,b,div);
[o4x,o4y] = rotate_Z(c,0,t0);
for i=1:numel(a)
[t2(i),t3(i)] = wing_sweep_pos(a(i),b,c,t0);
A_x(i) = a(i)*cosd(t2(i));
A_y(i) = a(i)*sind(t2(i));
P_x(i) = p*cosd(tp+t3(i))+o4x;
P_y(i) = p*sind(tp+t3(i))+o4y;
end
X_all = [0 A_x o4x P_x];
Y_all = [0 A_y o4y P_y];
for i=1:numel(a)
A = linspace(0,a(i),div);
Ax(i) = A_x(i);
Ay(i) = A_y(i);
Px(i) = P_x(i);
Py(i) = P_y(i);
link1x(i,:) = A.*cosd(t2(i));
link1y(i,:) = A.*sind(t2(i));
link3x(i,:) = o4x+B.*cosd(t3(i));
link3y(i,:) = o4y+B.*sind(t3(i));
Px1(i,:) = linspace(Ax(i),Px(i),div);
Py1(i,:) = linspace(Ay(i),Py(i),div);
Px2(i,:) = linspace(o4x,Px(i),div);
Py2(i,:) = linspace(o4y,Py(i),div);
figure(f)
plot(Ax(i),Ay(i),'ob')
hold on
plot(Px(i),Py(i),'og')
plot(link1x(i,:),link1y(i,:),'k',link3x(i,:),link3y(i,:),'k',...
Px1(i,:),Py1(i,:),'k',Px2(i,:),Py2(i,:),'k')
plot(0,0,'ok',o4x,o4y,'ok')
plot(Ax,Ay,'b')
plot(Px,Py,'g')
grid on
title('Wing Sweep Position Analysis')
xlim([1.2*(min(X_all)-1) 1.2*(max(X_all)+1)])
ylim([1.2*(min(Y_all)-1) 1.2*(max(Y_all)+1)])
axis equal;
hold off
end
figure(f)
legend('A','COM');
axis equal;
f = f+1;
end
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%%
% Wing Sweep Acceleration Function:
function [Aa,Ab,A_A,A_COM] = wing_sweep_accel(a,b,p,t2,t3,a_dot,Wa,Wb,a_ddot)
for i=1:numel(a)
M = [-a(i)*sind(t2(i)) b*sind(t3(i));a(i)*cosd(t2(i)) -b*cosd(t3(i))];
N(1) = -a_ddot(i)*cosd(t2(i))+2*a_dot(i)*Wa(i)*sind(t2(i))+...
Wa(i)^2*a(i)*cosd(t2(i))-Wb(i)^2*b*cosd(t3(i));
N(2) = -a_ddot(i)*sind(t2(i))+2*a_dot(i)*Wa(i)*cosd(t2(i))+...
Wa(i)^2*a(i)*sind(t2(i))-Wb(i)^2*b*sind(t3(i));
[X] = M\(N.');
Aa(i) = X(1);
Ab(i) = X(2);
A_Ax = -Ab(i)*b*sind(t3(i))-Wb(i)^2*b*cosd(t3(i));
A_Ay = Ab(i)*b*cosd(t3(i))-Wb(i)^2*b*sind(t3(i));
A_A(i) = sqrt(A_Ax^2+A_Ay^2);
A_COMx = -Ab(i)*p*sind(t3(i))-Wb(i)^2*p*cosd(t3(i));
A_COMy = Ab(i)*p*cosd(t3(i))-Wb(i)^2*p*sind(t3(i));
A_COM(i) = sqrt(A_COMx^2+A_COMy^2);
end
end
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%%
% Wing Sweep Velocity Function:
function [Wa,Wb,V_A,V_COM] = wing_sweep_vel(a,b,p,t2,t3,a_dot)
for i=1:numel(a)
M = [-a(i)*sind(t2(i)) b*sind(t3(i));a(i)*cosd(t2(i)) -b*cosd(t3(i))];
N = [-a_dot(i)*cosd(t2(i));-a_dot(i)*sind(t2(i))];
[X] = M\N;
Wa(i) = X(1);
Wb(i) = X(2);
V_A(i) = b*Wb(i);
V_COM(i) = p*Wb(i);
end
end
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%%
% Wing Sweep Position Function:
function [t2,t3] = wing_sweep_pos(a,b,c,t0)
A = acosd((a^2-b^2-c^2)/(-2*b*c));
B = acosd((b^2-a^2-c^2)/(-2*a*c));
C = acosd((c^2-a^2-b^2)/(-2*a*b));
t2 = B+t0;
t3 = t0-180-A;
end
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%