Modeling a two-mass, spring, damper system

I'm trying to model a system with two masses, two springs, two dampers, and one applied force using transfer functions. I'll then be inputting it into Simulink.

The system looks like this but there is a force applied to the right edge of ${ m }_{ 2 }$ pointing towards the right.

I already found the two differential equations of the system. Then I found one transfer function from each differential equation and ended up with these:

$\frac { { U }_{ 1 }(s) }{ { U }_{ 2 }(s) }$ and $\frac { { U }_{ 1 }(s){ U }_{ 2 }(s) }{ F(s) }$

I'm a little unsure of how to proceed from here. The end result should be the response of ${ U }_{ 1 }$ and ${ U }_{ 2 }$.

I haven't included the exact equations because I don't really want an exact answer (I'll never learn that way!). If I've left out any information please let me know. Thanks for any help!

EDIT: I guess the general question here is, how do you model this type of system using transfer functions?

EDIT 2: Here are my differential equations. Note that my variables differ slightly from the linked image.

${ \ddot { x } }_{ 1 }=\frac { { k }_{ 12 } }{ { m }_{ 1 } } \left( { x }_{ 2 }-{ x }_{ 1 } \right) +\frac { { d }_{ 12 } }{ { m }_{ 1 } } \left( { \dot { x } }_{ 2 }-{ \dot { x } }_{ 1 } \right) -\frac { { k }_{ 1 } }{ { m }_{ 1 } } { x }_{ 1 }-\frac { { d }_{ 1 } }{ { m }_{ 1 } } { \dot { x } }_{ 1 }$

${ \ddot { x } }_{ 2 }=\frac { f(t) }{ { m }_{ 2 } } -\frac { { k }_{ 12 } }{ { m }_{ 2 } } \left( { x }_{ 2 }-{ x }_{ 1 } \right) -\frac { { d }_{ 12 } }{ { m }_{ 2 } } \left( { \dot { x } }_{ 2 }-{ \dot { x } }_{ 1 } \right)$

-
Could you maybe explain a little more what you're trying to do? I understand that you don't want an exact answer, but perhaps giving a simplified problem would be helpful to people trying to answer your question. It also might be a good idea to post what you have for your transfer functions. More explicitly speaking, you have $\frac{U_1(s)}{U_2(s)}$ and $\frac{U_1(s)U_2(s)}{F(s)}$ in equations but don't set them equal to anything or define of your variables. – user758556 Mar 14 '13 at 7:32

Your system has 1 input, the force $F$, and two outputs, the position of $m_1$ and $m_2$. You should end with two equations, one which describes the behavior of $F$ to $m_1$ and one which describes the behavoir from $F$ to $m_2$, i.e. $\frac{U_1(s)}{F(s)}$ and $\frac{U_2(s)}{F(s)}$. These two equations then describe the behavoir of the coupled system.

So your initial laplace transformed equations, I believe, look something like $\alpha U_1(s) + F(s) = \beta U_2(s)$ and $\gamma U_1(s) = \beta U_2(s)$, where $\alpha, \beta, \gamma$ are polynomials of $s$ and the system parameters. Then you just first rewrite the last equations to $U_1(s) = \beta/\gamma U_2(s)$ such that you can subsitute it in the first by which you obtain $U_1(s)/F(s)$. Then rewrite to $U_2(s) = \gamma/\beta U_1(s)$ and subsitute again to obtain $U_2(s)/F(s)$.

You can simulate the behavoir simply with Simulink or Matlab itself e.g. step(tf) or lsim(tf) with tf the transfer function $U_1(s)/F(s)$ or $U_2(s)/F(s)$.

--edit--

Full solution...

We have the following

$$m_1\ddot{x}_1 = k_2(x_2 - x_1) + d_2(\dot{x}_2 - \dot{x}_1) - k_1x_1 - d_1\dot{x}_1 \\ m_2\ddot{x}_2 = F - k_2(x_2 - x_1) - d_2(\dot{x}_2 - \dot{x}_1)$$

Lets transfer this to the laplace domain and group the terms together such $\alpha X_1(s) = \beta X_2(s)$ and $\gamma X_2(s) = F + \zeta X_1(s)$. With $\alpha, \beta, \gamma, \zeta$ polynomials of $s$ and the coefficients $m_i,k_i,d_i$. We get

$$\left(m_1s^2 + (d_1 + d_2)s + (k_1 + k_2)\right)X_1(s) = \left(d_2s + k_2\right)X_2(s) \\ \left(m_2s^2 + d_2s + k_2\right)X_2(s) = F(s) + \left(d_2s + k_2\right)X_1(s) \Leftrightarrow F(s) = \left(m_2s^2 + d_2s + k_2\right)X_2(s) - \left(d_2s + k_2\right)X_1(s)$$

Rewrite the first equation to $\alpha X_1(s) = \beta X_2(s)$ to $X_1(s) = \frac{\beta}{\alpha}X_2(s)$ and $X_2(s) = \frac{\alpha}{\beta}X_1(s)$ we get

$$X_1(s) = \frac{d_2s + k_2}{m_1s^2 + (d_1 + d_2)s + (k_1 + k_2)}X_2(s) \\ X_2(s) = \frac{m_1s^2 + (d_1 + d_2)s + (k_1 + k_2)}{d_2s + k_2}X_1(s)$$

Now lets plugin each one such that we retrieve expressions for $\frac{X_1(s)}{F(s)}$ and $\frac{X_2(s)}{F(s)}$. Lets do $\frac{X_1(s)}{F(s)}$ first, we get

$$F(s) = \left(m_2s^2 + d_2s + k_2\right)\frac{m_1s^2 + (d_1 + d_2)s + (k_1 + k_2)}{d_2s + k_2}X_1(s) - \left(d_2s + k_2\right)X_1(s) \\ F(s) = \frac{\left(m_2s^2 + d_2s + k_2\right)\left(m_1s^2 + (d_1 + d_2)s + (k_1 + k_2)\right) - \left(d_2s + k_2\right)^2}{d_2s + k_2}X_1(s)$$

Now as you can see the numerator term is a quiet nasty one... You have to first multiply everything and then group all terms such that you get a polynomial of $s$ and the coefficients. For this reason it is also better now to start filling in the values of the coefficients which will make this job then much easier. If you want a analytic solution you can 1) compute it yourself by multiplying everything and collecting the terms 2) let Matlab do this for you. I choose to do it with Matlab.

syms m_1 m_2 d_1 d_2 k_1 k_2 s;
simple((m_2*s^2 + d_2*s + k_2)*(m_1*s^2 + (d_1 + d_2)*s + (k_1 + k_2)) - (d_2*s + k_2)^2)


The simple command will call several commands which could simplify the equation. In this case you will always need the answer for the collect function because this will group the coefficients. In this case it will return

$$\kappa = (m_1m_2)s^4 + (d_2m_1 + m_2(d_1 + d_2))s^3 + (m_2(k_1 + k_2) + k_2m_1 - d_2^2 + d_2(d_1 + d_2))s^2 + (k_2(d_1 + d_2) - 2d_2k_2 + d_2(k_1 + k_2))s + k_2(k_1 + k_2) - k_2^2$$

Hence you will get

$$\frac{X_1(s)}{F(s)} = \frac{d_2s + k_2}{\kappa}$$

I hope it is clear now and you can solve $\frac{X_2(s)}{F(s)}$ yourself.

-- edit oh well lets do $X_2(s)$ also... --

$$F(s) = \left(m_2s^2 + d_2s + k_2\right)X_2(s) - \left(d_2s + k_2\right)\frac{d_2s + k_2}{m_1s^2 + (d_1 + d_2)s + (k_1 + k_2)}X_2(s) \\ F(s) = \left(m_2s^2 + d_2s + k_2\right)X_2(s) - \frac{(d_2s + k_2)^2}{m_1s^2 + (d_1 + d_2)s + (k_1 + k_2)}X_2(s) F(s) = \frac{\left(m_1s^2 + (d_1 + d_2)s + (k_1 + k_2)\right)\left(m_2s^2 + d_2s + k_2\right) - \left(d_2s + k_2\right)^2}{m_1s^2 + (d_1 + d_2)s + (k_1 + k_2)}X_2(s) \frac{X_2(s)}{F(s)} = \frac{m_1s^2 + (d_1 + d_2)s + (k_1 + k_2)}{\kappa}$$

-
When I try this, I end up with a transfer function that can't be simplified enough to put into Simulink. – edc1591 Mar 14 '13 at 20:14
Can you please provide the initial equations from which you start? So the ones before you transfer it to the s-domain? – WG- Mar 14 '13 at 20:22
I've editing the original question to add my differential equations. – edc1591 Mar 14 '13 at 21:40
@edc1591: I expect you wound up with a ratio of two polynomials in $s$. Look up the Matlab function residue, this lets you do partial fraction expansion to get a sum of simpler terms. – Retarded Potential Mar 15 '13 at 18:10
@edc1591 I added the solutions. I hope from that it is clear now. The real problem is actually just in the nasty terms you will get when you want to find a analytic solution. – WG- Mar 17 '13 at 1:22

I'm not sure, but from what you've said I don't think what you've got are the transfer functions. You probably wouldn't get a transfer function directly from each equation.

Rather, from two linear second-order differential equations you should get two s-domain equations. A time derivative becomes multiplication by s. So something like:

$$x''(t) = A x'(t) + B x(t) + Cy'(t) + D y(t) + f(t)$$

(where $A,B,C,D$ are constants) becomes:

$$s^2 X(s) = (As+B) X(s) + (Cs+D) Y(s) + F(s)$$

(assuming $x(0)=y(0)=0$ and $x'(0)=y'(0)=0$). Do the same for the other equation, and do some algebra on both equations to get something like:

$$X(s) = \frac{p(s)}{q(s)}F(s)$$

where $p$ and $q$ are polynomials. Then $\frac{p(s)}{q(s)}$ is the transfer function for $x$. From it you can tell how $x(t)$ reacts to the input $f(t)$. You get another transfer function for $y(t)$.

Wikipedia's article has an example without dissipation. Here you have dissipation, you should be able to write each transfer function as a sum of terms $V/(s+\sigma)$, representing exponential decay with time constant $\sigma$, and/or terms $(Us + V)/((s + \sigma)^2 + \omega^2)$, representing a decaying sine wave with frequency $\omega$. Matlab/Simulink do a lot of the work for you.

-