Using Simulink to analyze two degrees of freedom system

Spring 2009 Compiled on January 31, 2024 at 3:47am

Abstract

A two degrees of freedom system consisting of two masses connected by springs and subject to 3 diﬀerent type of input forces is analyzed and simulated using Simulink

1 Introduction and Theory

It is required to ﬁnd

x_{1} (t)

and

x_{2} (t)

analytically and then to use Matlab’s Simulink software for the analysis.

The mathematical model of the system is ﬁrst developed and the equation of motions obtained using Lagrangian formulation then the analytical solution is found by solving the resulting coupled second order diﬀerential equations for

m_{1}

and

m_{2}

. Next, a simulink model is developed to implement the diﬀerential equations and the output

x_{1} (t)

and

x_{2} (t)

from Simulink is shown and compared to the output from the analytical solution.

2 Analytical solution

Assuming positive is downwards and that

x_{2} > x_{1}

, force-balance equations for

m_{1}

results in

m_{1} {\ddot{x}}_{1} = F (t) + k_{1} (x_{2} - x_{1})

And force-balance equations for

m_{2}

results in

m_{2} {\ddot{x}}_{2} = - k_{1} (x_{2} - x_{1}) - k_{2} x_{2}

Hence the EQM for the system become

\begin{aligned} m_{1} {\ddot{x}}_{1} - k_{1} x_{2} + k_{1} x_{1} & = F (t) \\ m_{2} {\ddot{x}}_{2} + (k_{1} + k_{2}) x_{2} - k_{1} x_{1} & = 0 \end{aligned}

Or in matrix form

(\begin{matrix} m_{1} & 0 \\ 0 & m_{2} \end{matrix}) (\begin{matrix} {\ddot{x}}_{1} \\ {\ddot{x}}_{2} \end{matrix}) + (\begin{matrix} k_{1} & - k_{1} \\ - k_{1} & (k_{1} + k_{2}) \end{matrix}) (\begin{matrix} x_{1} \\ x_{2} \end{matrix}) = (\begin{matrix} F (t) \\ 0 \end{matrix})

The above can be written in matrix form as

M \ddot{x} + K x = F

Where

\ddot{x}, x, F

are 2 by 1 vectors and

M

and

K

are the mass and stiﬀness matrices. The solution to the above is

\begin{matrix} (1) & x = x_{h} + x_{p} \end{matrix}

2.1 Finding the homogenous solution

We start by ﬁnding

x_{h}

from the following

(\begin{matrix} m_{1} & 0 \\ 0 & m_{2} \end{matrix}) (\begin{matrix} {\ddot{x}}_{1} \\ {\ddot{x}}_{2} \end{matrix}) + (\begin{matrix} k_{1} & - k_{1} \\ - k_{1} & (k_{1} + k_{2}) \end{matrix}) (\begin{matrix} x_{1} \\ x_{2} \end{matrix}) = (\begin{matrix} 0 \\ 0 \end{matrix})

Assuming that

x_{1} (t) = A_{1} \cos (ω t + ϕ)

and

x_{2} (t) = A_{2} \cos (ω t + ϕ)

, hence

{\dot{x}}_{1} (t) = - ω A_{1} \sin (ω t + ϕ)

and

{\dot{x}}_{2} (t) = - ω A_{2} \sin (ω t + ϕ)

and

{\ddot{x}}_{1} (t) = - ω^{2} A_{1} \cos (ω t + ϕ)

and

{\dot{x}}_{2} (t) = - ω^{2} A_{2} \cos (ω t + ϕ)

. Substituting the above values in the above system results in

(\begin{matrix} m_{1} & 0 \\ 0 & m_{2} \end{matrix}) (\begin{matrix} - ω^{2} A_{1} \cos (ω t + ϕ) \\ - ω^{2} A_{2} \cos (ω t + ϕ) \end{matrix}) + (\begin{matrix} k_{1} & - k_{1} \\ - k_{1} & (k_{1} + k_{2}) \end{matrix}) (\begin{matrix} A_{1} \cos (ω t + ϕ) \\ A_{2} \cos (ω t + ϕ) \end{matrix}) = (\begin{matrix} 0 \\ 0 \end{matrix})

Divide by

\cos (ω t + ϕ)

since not zero (else no solution exist) we obtain

(\begin{matrix} m_{1} & 0 \\ 0 & m_{2} \end{matrix}) (\begin{matrix} - ω^{2} A_{1} \\ - ω^{2} A_{2} \end{matrix}) + (\begin{matrix} k_{1} & - k_{1} \\ - k_{1} & (k_{1} + k_{2}) \end{matrix}) (\begin{matrix} A_{1} \\ A_{2} \end{matrix}) = (\begin{matrix} 0 \\ 0 \end{matrix})

Rewrite the above as

\begin{aligned} (\begin{array}{c} - ω^{2} m_{1} A_{1} \\ - ω^{2} m_{2} A_{2} \end{array}) + (\begin{array}{c} k_{1} A_{1} - k_{1} A_{2} \\ - k_{1} A_{1} + (k_{1} + k_{2}) A_{2} \end{array}) & = (\begin{array}{c} 0 \\ 0 \end{array}) \\ (\begin{array}{c} - ω^{2} m_{1} A_{1} + k_{1} A_{1} - k_{1} A_{2} \\ - ω^{2} m_{2} A_{2} - k_{1} A_{1} + (k_{1} + k_{2}) A_{2} \end{array}) & = (\begin{array}{c} 0 \\ 0 \end{array}) \\ (\begin{array}{c} (- ω^{2} m_{1} + k_{1}) A_{1} - k_{1} A_{2} \\ - ω^{2} m_{2} A_{2} + (k_{1} + k_{2}) A_{2} - k_{1} A_{1} \end{array}) & = (\begin{array}{c} 0 \\ 0 \end{array}) \\ (2) & (\begin{array}{c} - ω^{2} m_{1} + k_{1} & - k_{1} \\ - k_{1} & - ω^{2} m_{2} + (k_{1} + k_{2}) \end{array}) (\begin{array}{c} A_{1} \\ A_{2} \end{array}) & = (\begin{array}{c} 0 \\ 0 \end{array}) \end{aligned}

From the last equation above, we see that to obtain a solution we must have

| \begin{matrix} - ω^{2} m_{1} + k_{1} & - k_{1} \\ - k_{1} & - ω^{2} m_{2} + (k_{1} + k_{2}) \end{matrix} | = 0

since if we had

(\begin{matrix} A_{1} \\ A_{2} \end{matrix}) = 0

then no solution will exist. Therefore, taking the determinant and setting it to zero results in

\begin{aligned} (- ω^{2} m_{1} + k_{1}) (- ω^{2} m_{2} + (k_{1} + k_{2})) - k_{1}^{2} & = 0 \\ k_{1}^{2} + k_{1} k_{2} - ω^{2} k_{1} m_{1} - ω^{2} k_{1} m_{2} - ω^{2} k_{2} m_{1} + ω^{4} m_{1} m_{2} - k_{1}^{2} & = 0 \\ ω^{4} m_{1} m_{2} + ω^{2} (- k_{1} m_{1} - k_{2} m_{1} - k_{1} m_{2}) + (k_{1}^{2} + k_{1} k_{2} - k_{1}^{2}) & = 0 \end{aligned}

Let

ω^{4} = λ^{2}

, hence the above becomes

λ^{2} m_{1} m_{2} + λ (- k_{1} m_{1} - k_{2} m_{1} - k_{1} m_{2}) + (k_{1}^{2} + k_{1} k_{2} - k_{1}^{2}) = 0

Solving for

λ

gives

\begin{matrix} (3) & \begin{matrix} λ_{1} = \frac{1}{2 m_{1} m_{2}} (k_{1} m_{1} + k_{1} m_{2} + k_{2} m_{1} + \sqrt{k_{1}^{2} m_{1}^{2} + k_{1}^{2} m_{2}^{2} + k_{2}^{2} m_{1}^{2} + 2 k_{1} k_{2} m_{1}^{2} + 2 k_{1}^{2} m_{1} m_{2} - 2 k_{1} k_{2} m_{1} m_{2}}) \\ λ_{2} = \frac{1}{2 m_{1} m_{2}} (k_{1} m_{1} + k_{1} m_{2} + k_{2} m_{1} - \sqrt{k_{1}^{2} m_{1}^{2} + k_{1}^{2} m_{2}^{2} + k_{2}^{2} m_{1}^{2} + 2 k_{1} k_{2} m_{1}^{2} + 2 k_{1}^{2} m_{1} m_{2} - 2 k_{1} k_{2} m_{1} m_{2}}) \end{matrix} \end{matrix}

For each of the above solutions, we obtain a diﬀerent

(\begin{matrix} A_{1} \\ A_{2} \end{matrix})

from equation (2) as follows

For

λ_{1}

, (2) becomes

\begin{aligned} (\begin{array}{c} - λ_{1} m_{1} + k_{1} & - k_{1} \\ - k_{1} & - λ_{1} m_{2} + (k_{1} + k_{2}) \end{array}) (\begin{array}{c} A_{1} \\ A_{2} \end{array}) & = (\begin{array}{c} 0 \\ 0 \end{array}) \\ (\begin{array}{c} (- λ_{1} m_{1} + k_{1}) A_{1} - k_{1} A_{2} \\ - k_{1} A_{1} - λ_{1} m_{2} A_{2} + (k_{1} + k_{2}) A_{2} \end{array}) & = (\begin{array}{c} 0 \\ 0 \end{array}) \end{aligned}

From the ﬁrst equation above, we have

\begin{aligned} (- λ_{1} m_{1} + k_{1}) A_{1} - k_{1} A_{2} & = 0 \\ \frac{- λ_{1} m_{1} + k_{1}}{k_{1}} & = {(\frac{A_{2}}{A_{1}})}^{(1)} \end{aligned}

Similarly for

λ_{2}

\frac{- λ_{2} m_{1} + k_{1}}{k_{1}} = {(\frac{A_{2}}{A_{1}})}^{(2)}

Let

\begin{aligned} (4) & r_{1} & = \frac{- λ_{1} m_{1} + k_{1}}{k_{1}} = {(\frac{A_{2}}{A_{1}})}^{(1)} \\ r_{2} & = \frac{- λ_{2} m_{1} + k_{1}}{k_{1}} = {(\frac{A_{2}}{A_{1}})}^{(2)} \end{aligned}

Hence now

x_{h}

can be written as

{(\begin{matrix} x_{1} \\ x_{2} \end{matrix})}_{h} = (\begin{matrix} A_{1}^{(1)} \cos (ω_{1} t + ϕ_{1}) + A_{1}^{(2)} \cos (ω_{2} t + ϕ_{2}) \\ A_{2}^{(1)} \cos (ω_{1} t + ϕ_{1}) + A_{2}^{(2)} \cos (ω_{2} t + ϕ_{2}) \end{matrix})

But

A_{2}^{(1)} = r_{1} A_{1}^{(1)}

and

A_{2}^{(2)} = r_{2} A_{1}^{(2)}

, hence the above becomes

\begin{matrix} (5) & {(\begin{matrix} x_{1} \\ x_{2} \end{matrix})}_{h} = (\begin{matrix} A_{1}^{(1)} \cos (ω_{1} t + ϕ_{1}) + A_{1}^{(2)} \cos (ω_{2} t + ϕ_{2}) \\ r_{1} A_{1}^{(1)} \cos (ω_{1} t + ϕ_{1}) + r_{2} A_{1}^{(2)} \cos (ω_{2} t + ϕ_{2}) \end{matrix}) \end{matrix}

Now, given numerical values for

k_{1}, k_{2}, m_{1}, m_{2}

we can ﬁnd

ω_{1}, ω_{2}

from (3) above, and next ﬁnd

r_{1}, r_{2}

from (4). Hence (5) contains 4 unknowns,

A_{1}^{(1)}, A_{1}^{(2)}, ϕ_{1}, ϕ_{2}

which now can be found from initial conditions (after we ﬁnd the particular solution) which we will now proceed to do.

2.2 Finding particular solutions

For each of the above, we ﬁnd

x_{p}

and then add it to

x_{h}

found above in (5) to obtain (1).

2.2.1 Finding the particular solution for unit impulse input

Using the standard response for a unit impulse which for a single degree of freedom system is

x (t) = \frac{1}{m ω_{n}} \sin ω_{n} t

, then we write

x_{p}

x_{p} = {(\begin{matrix} x_{1} \\ x_{2} \end{matrix})}_{p} = (\begin{matrix} \frac{1}{m ω_{1}} \sin ω_{1} t + \frac{1}{m ω_{2}} \sin ω_{2} t \\ 0 \end{matrix})

Hence, the general solution becomes

\begin{matrix} (6) & (\begin{matrix} x_{1} \\ x_{2} \end{matrix}) = (\begin{matrix} A_{1}^{(1)} \cos (ω_{1} t + ϕ_{1}) + A_{1}^{(2)} \cos (ω_{2} t + ϕ_{2}) \\ r_{1} A_{1}^{(1)} \cos (ω_{1} t + ϕ_{1}) + r_{2} A_{1}^{(2)} \cos (ω_{2} t + ϕ_{2}) \end{matrix}) + (\begin{matrix} \frac{1}{m ω_{1}} \sin ω_{1} t + \frac{1}{m ω_{2}} \sin ω_{2} t \\ 0 \end{matrix}) \end{matrix}

2.2.2 Finding the particular solution for unit step input

Since unit step is

1

for

t > 0

, then, using convolution we write

\begin{aligned} x_{p} (t) & = \int_{0}^{t} f (τ) h (t - τ) d τ \\ = \int_{0}^{t} h (t - τ) d τ \\ = \int_{0}^{t} \frac{1}{m ω_{n}} \sin ω_{n} (t - τ) d τ \\ = \frac{1}{m ω_{n}} {[\frac{- \cos (ω_{n} (t - τ))}{- ω_{n}}]}_{0}^{t} \\ = \frac{1}{m ω_{n}^{2}} {[\cos (ω_{n} (t - τ))]}_{0}^{t} \\ = \frac{1}{m ω_{n}^{2}} [1 - \cos (ω_{n} t)] \end{aligned}

Then, since now we have 2 natural frequencies, we can write

x_{p}

x_{p} = {(\begin{matrix} x_{1} \\ x_{2} \end{matrix})}_{p} = (\begin{matrix} \frac{1}{m ω_{1}^{2}} [1 - \cos (ω_{1} t)] + \frac{1}{m ω_{2}^{2}} [1 - \cos (ω_{2} t)] \\ 0 \end{matrix})

Hence, the general solution becomes

(\begin{matrix} x_{1} \\ x_{2} \end{matrix}) = (\begin{matrix} A_{1}^{(1)} \cos (ω_{1} t + ϕ_{1}) + A_{1}^{(2)} \cos (ω_{2} t + ϕ_{2}) \\ r_{1} A_{1}^{(1)} \cos (ω_{1} t + ϕ_{1}) + r_{2} A_{1}^{(2)} \cos (ω_{2} t + ϕ_{2}) \end{matrix}) + (\begin{matrix} \frac{1}{m ω_{1}^{2}} [1 - \cos (ω_{1} t)] + \frac{1}{m ω_{2}^{2}} [1 - \cos (ω_{2} t)] \\ 0 \end{matrix})

2.2.3 Finding the particular solution for

\sin ω t

In this case, we guess that

x_{1 p} = c_{1} \cos ω t + c_{2} \sin ω t

, and since there is no forcing function being applied directly on

m_{2}

then

x_{2 p} = 0

hence

x_{p} = {(\begin{matrix} x_{1} (t) \\ x_{2} (t) \end{matrix})}_{p} = (\begin{matrix} c_{1} \cos ω t + c_{2} \sin ω t \\ 0 \end{matrix})

Then

{\dot{x}}_{1 p} (t) = - ω c_{1} \sin ω t + ω c_{2} \cos ω t

and

{\ddot{x}}_{1 p} (t) = - ω^{2} c_{1} \cos ω t + ω^{2} c_{2} \sin ω t

and now we substitute these into the original ODE for

x_{1}

which is

\begin{aligned} m_{1} {\ddot{x}}_{1} - k_{1} x_{2} + k_{1} x_{1} & = F (t) \\ m_{1} {\ddot{x}}_{1 p} - k_{1} x_{2 p} + k_{1} x_{1 p} & = \sin ω t \end{aligned}

We obtain the following

\begin{aligned} m_{1} (- ω^{2} c_{1} \cos ω t + ω^{2} c_{2} \sin ω t) - k_{1} (0) + k_{1} (c_{1} \cos ω t + c_{2} \sin ω t) & = \sin ω t \\ \cos ω t (- ω^{2} c_{1} m_{1} + k_{1} c_{1}) + \sin ω t (m_{1} ω^{2} c_{2} + k_{1} c_{2}) & = \sin ω t \end{aligned}

Hence by comparing coeﬃcients, we obtain

\begin{aligned} - ω^{2} c_{1} m_{1} + k_{1} c_{1} & = 0 \\ m_{1} ω^{2} c_{2} + k_{1} c_{2} & = 1 \end{aligned}

\begin{aligned} c_{1} (- ω^{2} m_{1} + k_{1}) & = 0 \\ c_{2} (m_{1} ω^{2} + k_{1}) & = 1 \end{aligned}

c_{1}

must be zero since

k_{1} - ω^{2} m_{1} = 0

only when

ω = ω_{n}

and we assume that this is not the case here. Hence

\begin{aligned} c_{1} & = 0 \\ c_{2} & = \frac{1}{(m_{1} ω^{2} + k_{1})} \end{aligned}

Therefore

x_{p}

becomes

x_{p} = {(\begin{matrix} x_{1} (t) \\ x_{2} (t) \end{matrix})}_{p} = (\begin{matrix} \frac{1}{(m_{1} ω^{2} + k_{1})} \sin ω t \\ 0 \end{matrix})

And the general solution becomes

\begin{matrix} (8) & (\begin{matrix} x_{1} \\ x_{2} \end{matrix}) = (\begin{matrix} A_{1}^{(1)} \cos (ω_{1} t + ϕ_{1}) + A_{1}^{(2)} \cos (ω_{2} t + ϕ_{2}) \\ r_{1} A_{1}^{(1)} \cos (ω_{1} t + ϕ_{1}) + r_{2} A_{1}^{(2)} \cos (ω_{2} t + ϕ_{2}) \end{matrix}) + (\begin{matrix} \frac{1}{(m_{1} ω^{2} + k_{1})} \sin ω t \\ 0 \end{matrix}) \end{matrix}

3 Simulink simulation and block diagrams

In simulink, we will directly solve the system from the original formulation

\begin{aligned} m_{1} {\ddot{x}}_{1} - k_{1} x_{2} + k_{1} x_{1} & = F (t) \\ m_{2} {\ddot{x}}_{2} + (k_{1} + k_{2}) x_{2} - k_{1} x_{1} & = 0 \end{aligned}

\begin{aligned} {\ddot{x}}_{1} - \frac{k_{1}}{m_{1}} x_{2} + \frac{k_{1}}{m_{1}} x_{1} & = \frac{F (t)}{m_{1}} \\ {\ddot{x}}_{2} + \frac{(k_{1} + k_{2})}{m_{2}} x_{2} - \frac{k_{1}}{m_{2}} x_{1} & = 0 \end{aligned}

Hence

\begin{aligned} {\ddot{x}}_{1} & = \frac{F (t)}{m_{1}} + \frac{k_{1}}{m_{1}} x_{2} - \frac{k_{1}}{m_{1}} x_{1} \\ {\ddot{x}}_{2} & = - \frac{(k_{1} + k_{2})}{m_{2}} x_{2} + \frac{k_{1}}{m_{2}} x_{1} \end{aligned}

3.1 Unit step simulink diagram and output

For an initial run with parameters

m_{1} = m_{2} = k_{1} = k_{2} = 1

I get this warning below

And this is the output for

x_{1} (t)

and

x_{2} (t)

for the unit step response

3.1.1 Veriﬁcation of result from Simulink by Numerically solving the diﬀerential equations

To verify the above output from Simulink, I solved the same coupled diﬀerential equations for zero initial conditions numerically (using a numerical diﬀerential equation solver) and plotted the solution for

x_{1} (t)

and

x_{2} (t)

and the result matches that shown above by simulink. Here is the code the plot as a result of this veriﬁcation

4 Unit impulse simulink diagram and output

4.0.1 Veriﬁcation of result from Simulink by Numerically solving the diﬀerential equations

To verify the above output from Simulink, The same coupled diﬀerential equations were solved numerically for zero initial conditions numerically and the solution plotted for

x_{1} (t)

and

x_{2} (t)

and the result was found to match that shown above by simulink. Here is the code used to do the veriﬁcation.

4.1

\sin (ω t)

input simulink diagram and output

For an initial run with parameters

m_{1} = m_{2} = k_{1} = k_{2} = 1

this is the output for

x_{1} (t)

and

x_{2} (t)

and showing the input signal at the same time

5 Discussion

A coupled system of two masses and springs was analyzed using Simulink. The simulation was done for one set of parameters (masses and stiﬀness). Simulink made the simulation of this system under diﬀerent loading conditions easy to do. The 2 masses response were recorded using simulink scope and the signals captured on the same plot to make it easy to compare the response of the ﬁrst mass to the second mass.

The analytical analysis was more time consuming than actually making the simulation in simulink. The ability to easily change diﬀerent sources to the system was useful as well as the ability to change the frequency of the input and immediately see the eﬀect on the response.

This was my ﬁrst project using Simulink, and I can see that this tool will be useful to learn more as it allows one to easily analyze engineering problems.