The Problem
I have the following system for which I would like to determine its dynamic response, taken from these lecture notes pg.20:
where the forcing functions are $u_1 = 1-U(t - t_0)$, $u_2=0$, and $u_3=0$, $U(\cdot)$ is the unit step function and $t_0 = 2\pi/\omega_1$ where $\omega_1$ is the first undamped natural frequency.
My Attempt
I know the time domain representation of the system to be of the typical format as shown in Eq.(1), and therefore the Fourier transform is as shown in Eq.(2):
\begin{gather} \mathbf{M}\ddot{\mathbf{x}} + \mathbf{C}\dot{\mathbf{x}} + \mathbf{Kx} = \mathbf{u}(t) \tag{1} \\ % \Rightarrow [-\omega^2\mathbf{M} + i\omega\mathbf{C} + \mathbf{K}]\mathbf{X}(i\omega) = \mathbf{U}(i\omega) \tag{2} \end{gather}
where since we know that $\mathcal{L}\{1 - U(t - t_0)\} = 1 - \frac{\exp{(-st_0})}{s}$, then the Fourier transform of the right hand side of Eq.(2) is:
$$\mathbf{U}(i\omega) = \left[\left(1 - \frac{e^{-i\omega t_0}}{i\omega}\right),\, 0,\, 0,\right]^T$$
As a result, the dynamic response would be:
$$\mathbf{X}(i\omega) = [-\omega^2\mathbf{M} + i\omega\mathbf{C} + \mathbf{K}]^{-1}\mathbf{U}(i\omega) \tag{3}$$
Since the response by which I am comparing my answer is given in generalised modal coordinates $\mathbf{q}$, then the modal transformation of Eq.(3) follows where $\mathbf{V}$ is the eigenvector matrix thus we know that the modal mass, stiffness, damping matrices should be $\bar{\mathbf{M}} = \mathbf{I}$ (i.e. the identity matrix), $\bar{\mathbf{K}} = \mathbf{\Omega}^2$ (i.e. a diagonal marix containing the square of the natural frequencies), and $\bar{\mathbf{C}} = \mathbf{V}^T\mathbf{C}\mathbf{V}$ respectively. So Eq.(3) becomes:
$$\mathbf{q}(i\omega) = \left[-\omega^2\mathbf{I} + i\omega\bar{\mathbf{C}} + \mathbf{\Omega}^2 \right]^{-1}\left[\mathbf{V}^T\mathbf{U}(i\omega)\right] \tag{4}$$
Thefore, the dynamic response would be the magnitude of the imaginary numbers obtained from Eq.(4) above, for all positive $\omega$.
However, when I follow the above steps and solve the problem in MATLAB
, the answer is far from the one provided in the excercise. Admitedly, the code is very simple as well (see below) so I struggle to see where I could have gone wrong:
clear ; clc
% Input
n = 3 ; % Number of DOFs
m = 1 ; % Mass [kg]
k = 1 ; % Stiffnes [N/m]
c = 0.2 ; % Viscous damping [Ns/m]
nv = [1 0 0]' ; % Load application vector (load applied to 1st DOF)
% Matrices
M = diag(m*ones(n, 1)) ;
K = diag(2*k*ones(n,1)) - diag(k*ones(n-1,1),-1) - diag(k*ones(n-1,1),1) ;
C = diag(c*ones(n,1)) ;
% Undamped Natural Frequencies & Modes
[V, Omsq] = eig(K, M) ;
Om = diag(sqrt(abs(Omsq))) ; % Natural Frequencies
% Frequency Aanalysis
t0 = 2*pi/Om(1) ; % Time where unit step is applied
omega = linspace(0, 3, 100) ; % Frequency vector
for i = 1:numel(omega)
s = 1i*omega(i) ;
Us = V'*((1 - exp(-s*t0)/s)*nv) ;
FRF(:,i) = (s^2*eye(n) + s*(V'*C*V) + Omsq)\Us ;
end
% Plotting
figure ; plot(omega, abs(FRF)) ; grid on ; set(gca, 'YScale', 'log') ; ylim([0 10]) ;