Untitled

clear all
close all

% Define the frequency ratio values (between 0 and 3 as per assignment
FreqRatio=linspace(0,3,5000); 
% Define values for damping ratio (zeta)
zeta=[0.01 0.1 0.2 0.7 1]; 

% Plot deformation response factor (Rd) for different zeta values
for k = 1:length(zeta)

    Rd = abs(1./(1-FreqRatio.^2+i*2*zeta(k)*FreqRatio));  

    % Plot Rd for each zeta
    figure(1)
    plot(FreqRatio,Rd)  
    axis([0 3 0 5])  
    hold on

end

xlabel('Frequency ratio \omega / \omega_{n}')
ylabel('Deformation response factor, \it{R_d = u_0 / (u_{st})_0}') 
xline(1, 'HandleVisibility','off') % Vertical line at FreqRatio = 1
yline(1, 'HandleVisibility','off') % Horizontal line at Rd = 1
legend('\zeta=0.01', '\zeta=0.1', '\zeta=0.2', '\zeta=0.7', '\zeta=1')


% Plot phase angle (phi) for different zeta values
for j=1:length(zeta)

    phi=angle((1-FreqRatio.^2+i*2*zeta(j)*FreqRatio))*180/pi; % Calculate phi in degrees to get the right scale
    
    % Plot phi for each zeta
    figure(2)
    plot(FreqRatio,phi)
    axis([0 3 0 180])  
 
    hold on

end

xlabel('Frequency ratio \omega / \omega_{n}')
ylabel('Phase angle \phi') 
xline(1, 'HandleVisibility','off') % Vertical line at FreqRatio = 1
yline(90, 'HandleVisibility','off') % Horizontal line at phi = 90 degrees
legend('\zeta=0.01', '\zeta=0.1', '\zeta=0.2', '\zeta=0.7', '\zeta=1')

% Determine vibration amplitude

zeta=0.05;
omega=10*2*pi;
e=0.002; % From lecture 2
k=2.25e6; % From lecture 2
mr=250; % From lecture 2
m=814; % From lecture 2

omega_n=sqrt(k/m); % Natural circular frequency

p0=mr*e*omega^2; 

% Calculate dynamic amplitude u0
u0=(p0/k)*abs(1/(1-(omega/omega_n)^2+i*2*zeta*(omega/omega_n))) % Dynamic amplitudide