Consider the system $\text{x” + x = sin3t}$ , with initial conditions $\text{x(0) = 0}$ and $\text{x’(0) = 0}$  

Take the Laplace Transform

${\text{s}}^2\text{X(s)}+\text{X(s)}=3/({\text{s}}^2+3^2)$ 

Solve for X(s) in the frequency domain

$\text{X(s)}=3/\left[\right({\text{s}}^2+1\left)\right({\text{s}}^2+9\left)\right]=(3/8)[1/({\text{s}}^2+1\left)\right]-(1/8)[3/({\text{s}}^2+9\left)\right]$ 

Inverse transform to get back to the time domain

$\text{X(t)}=(3/8)\text{sint}-(1/8)\text{sin3t}$ 

The response will be made up of two bits. The first, will be a sinusoidal response at the natural frequency $[{\text{ω}}_{\text{n}}={\text{(k/m)}}^{\text{1/2}}={\text{(1/1)}}^{\text{1/2}}=1]$ and the other will be a sinusoidal response at frequency $\text{ω = 3}$, from the frequency of the external input.

We can treat the system just described as a static deflection problem, where the weight of each mass will be balanced by the force in each spring.  The weight of each mass can be expressed as a column matrix [F] of three force terms found by multiplying the acceleration of gravity (g = 9.8 m/sec²) times each mass.

Then the equation

$\left[\text{k}\right]\left[\text{x}\right]=\left[\text{F}\right]$ is 

$\left[\begin{array}{ccc}\mathrm{80}& -30& 0\\ -30& 70& -40\\ 0& -40& 60\end{array}\right]\left[\begin{array}{c}{x}_1\\ x_2\\ x_3\end{array}\right]=\left[\begin{array}{c}9.8\\ 29.4\\ 98\end{array}\right]$ 

The deflection that the system will experience under the effect of the weight [X] can be found using matrix operations involving the inverse of [K]:  $\left[\text{X}\right]={\left[\text{K}\right]}^{-1}\left[\text{F}\right]$

This can be done easily using a computer mathematics routine such as MATLAB by defining the two matrices, inverting the K matrix and multiplying it by the F matrix

>> K = [80 -30 0; -30 70 -40; 0 -40 60]

K =

    80   -30     0

   -30    70   -40

     0   -40    60

>> F = [9.8; 29.4; 98]

F =

    9.8000

   29.4000

   98.0000

>> X = inv(K)*F

X =

    1.2727

    3.0673

    3.6782

So mass 1 will deflect 1.2727, mass 2 will deflect 3.0673 and mass 3 will deflect 3.6782 under these static weights.

We can also view this as a dynamic problem, where the system is displace and released.  This will cause all three bits to oscillate in simple harmonic motion, i.e. it will look like a cosine wave. 

We can assume a solution that looks like a displacement of

${\text{x}}_{\text{n}}={\text{A}}_{\text{n}}\text{cos}{\text{ ω}}_{\text{n}}\text{t}$ 

So the first derivative of displacement (velocity) is

${\text{x}}_{\text{n}}\text{’}={\text{A}}_{\text{n}}{{\text{ω}}_{\text{n}}}^{\mathrm{}}\text{cos}{\text{ ω}}_{\text{n}}\text{t}$ 

And the second derivative of displacement (acceleration) is

${\text{x}}_{\text{n}}\text{’’}=-{\text{A}}_{\text{n}}{{\text{ω}}_{\text{n}}}^2\text{cos}{\text{ ω}}_{\text{n}}\text{t}$ 

Thus, we can also write 

${\text{x}}_{\text{n}}\text{’’}=-{{\text{ω}}_{\text{n}}}^2{\text{x}}_{\text{n}}$ 

Since the matrix expression is

$\left[\text{M}\right]\left[\text{X’’}\right]+\left[\text{K}\right]\left[\text{X}\right]=\left[0\right]$ 

The above allows us to express the matrix equation as

$-{{\text{ω}}_{\text{n}}}^2\left[\text{M}\right]\left[\text{X}\right]+\left[\text{K}\right]\left[\text{X}\right]=\left[0\right]$ 

Let $\text{λ}={{\text{ω}}_{\text{n}}}^2$ and define a new matrix $\left[\text{A}\right]={\left[\text{M}\right]}^{-1}\left[\text{K}\right]$  and this develops into an eigenvalue form $\left\{{\text{λ}}_{\text{n}}\right[\text{I}]-[\text{A}\left]\right\}\left[\text{X}\right]=0$ 

So taking the determinant $\text{det{λ[I] – [A]} = 0}$ yields the system’s characteristic equation

Plugging in the values we had in the example problem

$\left[\begin{array}{c}\text{M}\end{array}\right]=\left[\begin{array}{ccc}1& 0& 0\\ 0& 3& 0\\ 0& 0& 10\end{array}\right]\left[\begin{array}{c}\text{K}\end{array}\right]=\left[\begin{array}{ccc}80& -30& 0\\ -30& 70& -40\\ 0& -40& 60\end{array}\right]$

We can use MATLAB to find the system’s eigenvalues via the following commands

>> M = [1 0 0; 0 3 0; 0 0 10]

M =

1   0   0

0   3   0

0   0   10

>>K = [80 -30 0; -30 70 -40; 0 -40 60]

K =

80   -30    0

-30   70   -40

0    -40    60

>>A = inv(M)*K

A =

80.0000    -30.0000          0

-10.0000     23.3333    -13.3333

0          -4.0000       6.0000

>>lambda = eig(A)

lambda =

84.9248

21.6116

2.7969

And the results tell us that our system eigenvalues are 2.8, 21.6 and 84.9

Recalling that $\text{λ}$ is related to our system natural frequency ${\text{ω}}_{\text{n}}$ by $\text{λ}={{\text{ω}}_{\text{n}}}^2$ or ${\text{ω}}_{\text{n}}={\left(\text{λ}\right)}^{1/2}$ 

We can find the natural frequencies:

${\text{ω}}_1={\left(2.8\right)}^{1/2}=1.67\text{ rad/sec}$ 

${\text{ω}}_2={\left(21.6\right)}^{1/2}=4.65\text{ rad/sec}$ 

${\text{ω}}_3={\left(84.9\right)}^{1/2}=9.2\text{ rad/sec}$ 

Now recall that each eigenvalue has an associated eigenvector. These eigenvectors correspond to the mode shapes of the vibratory system. In this case, they tell us that at the first natural frequency, the associated mode shape, from the first eigenvector will indicate the relative motions of the masses.

We can use MATLAB to get the eigenvalues and eigenvectors using

[V,D] = eig(A)

V =

-0.9868    -0.4456    0.2360

0.1620    -0.8672    0.6074

-0.0082     0.2222    0.7585

D =

84.9248         0            0

0               21.6116      0

0                    0        2.7969

Normalising these eigenvectors yields the vibratory mode shapes

$\left[X_1\right]=\left(1/0.7585\right)\left[\begin{array}{c}0.2360\\ 0.6074\\ 0.7585\end{array}\right]=\left[\begin{array}{c}0.3\\ 0.8\\ 1.0\end{array}\right]$ 

$\left[X_2\right]=\left(1/-.8672\right)\left[\begin{array}{c}-0.4456\\ -0.8672\\ 0.2222\end{array}\right]=\left[\begin{array}{c}0.5\\ 1.0\\ -0.3\end{array}\right]$ 

$\left[X_3\right]=\left(1/-.9868\right)\left[\begin{array}{c}-0.9868\\ 0.1620\\ -0.0082\end{array}\right]=\left[\begin{array}{c}1.0\\ -0.2\\ 0.0\end{array}\right]$ 

The first eigenvector corresponds to the vibrational mode shape occurring at the first natural frequency of 1.67 rad/sec (16 hz), which means the motion of the masses is as shown

Vibratory modes for 2^nd and 3^rd natural frequencies can be found similarly.  A negative sign means vibratory motion is out-of-phase (motion is up for one mass, down for the other).

Consider the racecar on a transporter trailer as shown, with vehicle mass (m = 10kg) and a spring rate (k = 100 N/m) and damping (b = 20N-sec/m) created by the attachments of the car to the trailer.  Let u and y be the deflections of the trailer and the car respectively.

We will use solution Approach 2 as developed in this topic.

Assume that the system is at rest until t = 0, at which time the trailer begins moving with a constant velocity u’ = 1 m/sec.

Let us apply Newton’s Second Law (F = ma) to construct the system equation

$\text{my” + b(y’ – u’) + k(y – u) = 0}$ 

Rearrange

$\text{my” + by’ + ky = bu’ + ku}$

Divide by m

$\text{y” + (b/m)y’ + (k/m)y = (b/m)u’ + (k/m)u}$

Rewrite as

$\text{y”}+{\text{a}}_1\text{y’}+{\text{a}}_2\text{y}={\text{b}}_0\text{y”}+{\text{b}}_1\text{u’}+{\text{b}}_2\text{u}$

And we find that this is in one of the forms previously described in solution Approach 2

Where:

${\text{a}}_1=\text{b/m}$ 

${\text{a}}_2=\text{k/m}$ 

${\text{b}}_0=\text{0}$ 

${\text{b}}_1=\text{b/m}$ 

${\text{b}}_2=\text{k/m}$ 

${\text{β}}_0={\text{b}}_0=0$ 

${\text{β}}_1={\text{b}}_1-{\text{a}}_1{\text{b}}_0=\text{b/m}$ 

${\text{β}}_2={\text{b}}_2-{\text{a}}_1{\text{β}}_1-{\text{a}}_2{\text{β}}_0=\text{(k/m)}-{\text{(b/m)}}^2$

We previously developed:

So

We previously said that if we could define [A], [B], [C] and [D] then we could get a solution.  And the easiest way to do so is in MATLAB using the ss command to tell it that we are giving it a state space form, with commands as follows:

>>A = [0 1; -10 -2]

>>B = [2; 6]

>>C = [1  0]

>>D = 0

>>sys = ss(A,B,C,D)

These commands define our four state space form matrices and then create a system response (which we named sys) using the matrices we just defined

We can also create a plot of this response, from time 0 to time 4 in steps of 0.01 using the commands

>> t = 0:0.01:4

>> u = t

>>lsim(sys,u,t)

Where the blue line shows the motion of the trailer, u.  And the red line shows the motion of the racecar, y, which begins by lagging, then oscillating, and finally stabilises to moving in unison with the trailer.

MATLAB also lets us go back and forth between the state space and the transfer fucntion expressions using the commands tf2ss (transfer function to state space) and ss2tf (state space to transfer function). 

 

Let us use the state space matrices fromt eh previous example and change to the transfer fucntion using the command ss2tf(A,B,C,D,iu) where iu is the number of inputs to the system (1 in our case).

 

>> A = [0 1 0 0; 0 0 1 0; 0 0 0 1; -10 -15 -12 -6]

A =

     0     1     0     0

     0     0     1     0

     0     0     0     1

   -10   -15   -12    -6

>> B = [0; 0; 0; 1]

B =

     0

     0

     0

     1

>> C = [10 15 5 0]

C =

    10    15     5     0

>> D = 0

D =

     0

>> [num,den] = ss2tf(A,B,C,D,1)

num =

     0     0     5    15    10

den =

    1.0000    6.0000   12.0000   15.0000   10.0000

These results would mean that the transfer function was

$\text{Y(s)}/\text{U(s)}=(0{\text{s}}^4+0{\text{s}}^3+5{\text{s}}^2+15\text{s}+10(/(1{\text{s}}^4+6{\text{s}}^3+12{\text{s}}^2+15\text{s}+10)$ 

Which can be seen to be the same transfer function that was generated in the previous example.

 

We can reverse this process if we define the numerator of the transfer function as num, define the denominator of the transfer function as den, and change to state space using the command  [A,B,C,D] = tf2ss(num,den)

>> [A,B,C,D] = tf2ss(num,den)

A =

   -6.0000  -12.0000  -15.0000  -10.0000

    1.0000         0         0         0

         0    1.0000         0         0

         0         0    1.0000         0

B =

     1

     0

     0

     0

C =

     0     5    15    10

D =

     0

 

And we are back where we started.

Find the natural frequencies of the system of three masses shown in the figure, with a spring and damper between each and force f₁t applied to mass 1.

Since we are after a dynamic solution, we will use solution Approach 1 as developed in this topic.

Write the equations of motion:

${\text{m}}_1{\text{x}}_1\text{”}+{\text{k}}_1({\text{x}}_1-{\text{x}}_2)+{\text{c}}_1({\text{x}}_1\text{”}-{\text{x}}_2\text{”})=-{\text{f}}_1\text{(t)}$ 

${\text{m}}_2{\text{x}}_2\text{”}+{\text{k}}_1({\text{x}}_2-{\text{x}}_1){}_{\mathrm{}}+{\text{k}}_2({\text{x}}_2-{\text{x}}_3)+{\text{c}}_1({\text{x}}_2\text{’}-{\text{x}}_1\text{’})+{\text{c}}_2({\text{x}}_2\text{’}-{\text{x}}_3\text{’})=0$ 

${\text{m}}_3{\text{x}}_3\text{”}+{\text{k}}_2({\text{x}}_3-{\text{x}}_2){}_{\mathrm{}}+{\text{k}}_3{\text{x}}_3+{\text{c}}_2({\text{x}}_3\text{’}-{\text{x}}_2\text{’}){\text{x}}_2+{\text{k}}_3{\text{x}}_3\text{’}=0$  

Arrange as a matrix equation

Develop a state space form by defining

And

The inverse of the mass matrix is

Recall from what we did before that the state space matrix A should achieve [Y’] = [A][Y] and this can be achieved from

As discussed in solution Approach 1, we could now do an eigen-solution using matrix A and generate the eigenvalues (which as previously discussed lead to the natural frequencies) and the eigenvectors (which as previously discussed lead to the vibration mode shapes).