Four masses, four springs and two dampers

A four-degree-of-freedom translational system with damping. The coupled equations of motion come from Lagrange's equations, with spring elongations as generalized coordinates and generalized forces for the dampers. The symmetry of the mass matrix and a Rayleigh dissipation function check the derivation. Normal modes, modal damping, and the simulator complete the picture.

Problem statement

Determine the equations of motion of the system shown in the figure.

Figure 1. m₁ is tied to the left wall by spring k₁ and damper B₁. m₂ and m₃ slide on a fixed shelf, chained to m₁ through k₂ and k₃. m₄ is linked to m₁ by k₄ and to the right wall by damper B₂.

The masses slide without friction. Energy is lost only in the two viscous dampers $B_1$ and $B_2$. These are non-conservative forces, so they cannot be included in the Lagrangian. They enter Lagrange's equations as generalized forces.

Generalized coordinates

The system has four degrees of freedom. The paper chooses the elongations of the four springs $x_1, x_2, x_3, x_4$ as generalized coordinates; $x_i$ belongs to spring $k_i$. They form a complete and independent set (geometric admissibility). With them, the absolute displacements of the masses are

$$y_1 = x_1,\qquad y_2 = x_1 + x_2,\qquad y_3 = x_1+x_2+x_3,\qquad y_4 = x_1 + x_4$$

$$v_1 = \dot x_1,\qquad v_2 = \dot x_1+\dot x_2,\qquad v_3 = \dot x_1+\dot x_2+\dot x_3,\qquad v_4 = \dot x_1+\dot x_4$$

This choice makes the potential energy very simple (each spring depends on one coordinate). The price is that the kinetic energy couples all coordinates. The alternative choice, absolute displacements $y_i$, does the opposite: the kinetic energy becomes diagonal and the springs couple the coordinates (see section 7).

Energies and generalized forces

TerminologyAs in the course: attractive energy $\mathcal{A}$ = potential energy, repulsive coenergy $\mathcal{R}^*$ = kinetic coenergy, $\mathcal{L}=\mathcal{R}^*-\mathcal{A}$.

$$\mathcal{R}^* = \tfrac12 m_1v_1^2 + \tfrac12 m_2v_2^2 + \tfrac12 m_3v_3^2 + \tfrac12 m_4v_4^2,\qquad \mathcal{A} = \tfrac12k_1x_1^2 + \tfrac12k_2x_2^2 + \tfrac12k_3x_3^2 + \tfrac12k_4x_4^2$$

$$\mathcal{L} = \tfrac12 m_1\dot x_1^2 + \tfrac12 m_2(\dot x_1+\dot x_2)^2 + \tfrac12 m_3(\dot x_1+\dot x_2+\dot x_3)^2 + \tfrac12 m_4(\dot x_1+\dot x_4)^2 - \tfrac12\sum_{i=1}^4 k_ix_i^2$$

In the paper's written Lagrangian the spring terms for $k_2, k_3, k_4$ appear with a plus sign, and the velocities are not squared. These are typos; the derivatives that follow use the correct form.

Damper forces

The Lagrange equations with non-conservative forces are

$$\frac{d}{dt}\left(\frac{\partial\mathcal{L}}{\partial\dot x_i}\right) - \frac{\partial\mathcal{L}}{\partial x_i} = \Xi_i,\qquad i=1,\dots,4$$

The paper takes

$$\Xi_1 = -B_1\dot x_1,\qquad \Xi_2 = 0,\qquad \Xi_3 = 0,\qquad \Xi_4 = -B_2(\dot x_1+\dot x_4)$$

IncompleteThe damper $B_2$ acts on mass $m_4$, whose displacement $y_4 = x_1+x_4$ depends on both $x_1$ and $x_4$. Its force therefore does virtual work for a variation of $x_1$ too, so $\Xi_1$ must also contain $-B_2(\dot x_1+\dot x_4)$. See section 7.

Lagrange equations (paper)

Coordinate $x_1$

$$\frac{d}{dt}\Big[m_1\dot x_1 + m_2(\dot x_1+\dot x_2) + m_3(\dot x_1+\dot x_2+\dot x_3) + m_4(\dot x_1+\dot x_4)\Big] + k_1x_1 = -B_1\dot x_1$$

(1)$$(m_1+m_2+m_3+m_4)\ddot x_1 + (\color{#c27a00}{m_1}+m_3)\ddot x_2 + m_3\ddot x_3 + m_4\ddot x_4 + k_1x_1 = -B_1\dot x_1$$

Coordinates $x_2, x_3, x_4$

(2)$$(m_2+m_3)\ddot x_1 + (m_2+m_3)\ddot x_2 + m_3\ddot x_3 + k_2x_2 = 0$$

(3)$$m_3\ddot x_1 + m_3\ddot x_2 + m_3\ddot x_3 + k_3x_3 = 0$$

(4)$$m_4\ddot x_1 + m_4\ddot x_4 + k_4x_4 = -B_2\dot x_1 - B_2\dot x_4$$

The result is a system of four coupled second-order linear differential equations. Equations (2)–(4) are correct. Equation (1) contains two errors (highlighted, and discussed in section 7).

Simulink model

The paper's block diagram has four integrator chains, one per coordinate. Each summing junction collects the spring term, the damper terms, and the second derivatives of the other coordinates multiplied by mass ratios. Every equation contains the accelerations of the others (the mass matrix is full), so the diagram is full of algebraic loops.

Writing the system in matrix form makes the structure clear and gives a loop-free implementation:

$$\mathbf{M}\ddot{\mathbf{x}} + \mathbf{C}\dot{\mathbf{x}} + \mathbf{K}\mathbf{x} = \mathbf{0}\quad\Longrightarrow\quad \ddot{\mathbf{x}} = -\mathbf{M}^{-1}\left(\mathbf{C}\dot{\mathbf{x}} + \mathbf{K}\mathbf{x}\right)$$

$$\mathbf{M} = \begin{bmatrix} m_1+m_2+m_3+m_4 & m_2+m_3 & m_3 & m_4\\ m_2+m_3 & m_2+m_3 & m_3 & 0\\ m_3 & m_3 & m_3 & 0\\ m_4 & 0 & 0 & m_4\end{bmatrix},\quad \mathbf{C} = \begin{bmatrix} B_1+B_2 & 0 & 0 & B_2\\ 0&0&0&0\\ 0&0&0&0\\ B_2 & 0 & 0 & B_2\end{bmatrix},\quad \mathbf{K} = \operatorname{diag}(k_1,k_2,k_3,k_4)$$

Figure 2. Vectorized Simulink model: one gain for M⁻¹, two for C and K, and two vector integrators. There are no algebraic loops.

Simulation results from the paper

The paper shows one simulation, about 80 s long, of all four elongations. They are damped oscillations: large initially (x₁ up to about ±0.4) and decaying within roughly 30–40 s. A small, fast oscillation then remains visible until the end. The parameter values and initial conditions are not given, so the plot cannot be reproduced exactly. The simulator uses illustrative values that produce the same picture: a dominant 4 s oscillation that dies out and a lightly damped 1.3 s oscillation that persists.

The paper concludes that the oscillations of $m_1$ and $m_4$ decay faster because these masses are attached to the walls through dampers, while those of $m_2$ and $m_3$ decay more slowly because they are connected only through springs.

Critical review

Wrong coefficient in equation (1): the mass matrix must be symmetric

Expanding $\frac{d}{dt}\partial\mathcal{L}/\partial\dot x_1$ gives $\ddot x_2$ the coefficient $m_2+m_3$, not $m_1+m_3$. There is a quick check for this kind of slip. The kinetic energy is a quadratic form $\tfrac12\dot{\mathbf{x}}^{\mathsf T}\mathbf{M}\dot{\mathbf{x}}$, so the mass matrix from Lagrange's equations is always symmetric. The coefficient of $\ddot x_2$ in equation (1) must equal the coefficient of $\ddot x_1$ in equation (2), which is $m_2+m_3$. With $m_1+m_3$ the matrix is not symmetric, and the model is no longer physical. Its energy can rise even with dampers present (see the energy plot in the simulator with the paper model).

The damper B₂ is missing from equation (1)

The cleanest way to obtain the damper terms is the Rayleigh dissipation function, half of the dissipated power written in the generalized velocities:

$$\mathcal{D} = \tfrac12 B_1 v_1^2 + \tfrac12 B_2 v_4^2 = \tfrac12 B_1\dot x_1^2 + \tfrac12 B_2(\dot x_1+\dot x_4)^2,\qquad \Xi_i = -\frac{\partial\mathcal{D}}{\partial\dot x_i}$$

$$\Xi_1 = -B_1\dot x_1 - B_2(\dot x_1+\dot x_4),\qquad \Xi_4 = -B_2(\dot x_1+\dot x_4)$$

Like the mass matrix, the damping matrix $\mathbf{C}$ must be symmetric. The paper's version has a $B_2$ term in row 4, column 1, but none in row 1, column 4.

Corrected equation (1)

$$(m_1+m_2+m_3+m_4)\ddot x_1 + (m_2+m_3)\ddot x_2 + m_3\ddot x_3 + m_4\ddot x_4 + (B_1+B_2)\dot x_1 + B_2\dot x_4 + k_1x_1 = 0$$

Equations (2)–(4) stay as in the paper.

Cross-check with absolute coordinates

With absolute displacements $y_i$, Newton's law for each mass gives the same model in a more familiar form. The mass matrix is diagonal and the springs couple the coordinates:

$$\begin{aligned} m_1\ddot y_1 &= -k_1y_1 - B_1\dot y_1 + k_2(y_2-y_1) + k_4(y_4-y_1)\\ m_2\ddot y_2 &= -k_2(y_2-y_1) + k_3(y_3-y_2)\\ m_3\ddot y_3 &= -k_3(y_3-y_2)\\ m_4\ddot y_4 &= -k_4(y_4-y_1) - B_2\dot y_4\end{aligned}$$

Substituting $y = x_1, x_1+x_2, \dots$ and combining equations reproduces the corrected Lagrange equations. Adding all four equations gives the $x_1$ equation, which shows directly that $B_2$ must appear in it.

Which masses really decay faster?

The paper's explanation, that the masses attached to dampers lose energy faster, is intuitive but not quite how a coupled system works. Damping acts on modes, not on individual masses. Every coordinate is a mix of all four modes. A mode decays quickly if the damped masses ($m_1$, $m_4$) move strongly in it, and slowly if they hardly move. The persistent fast oscillation in the plot is a lightly damped mode in which $m_2$ and $m_3$ oscillate against each other while $m_1$ and $m_4$ stay almost at rest. It is therefore most visible in $x_2$ and $x_3$, but it is present in all coordinates. The simulator shows the frequencies, approximate damping ratios and decay times of the modes.

Interactive simulator

The equations are written as $\ddot{\mathbf x} = -\mathbf M^{-1}(\mathbf C\dot{\mathbf x} + \mathbf K\mathbf x)$ and integrated with 4th-order Runge–Kutta (5 ms step). Initial velocities are zero. The default values are illustrative, because the paper does not give its own.

Cases

Model

Also show the other model

Show

Masses [kg]

m₁

m₂

m₃

m₄

Springs [N/m] and dampers [N·s/m]

k₁

k₂

k₃

k₄

B₁

B₂

Initial elongations [m]

x₁(0)

x₂(0)

x₃(0)

x₄(0)

Duration [s]

Spring elongations x₁…x₄

Normal modes (corrected model, undamped frequencies + approximate modal damping)

The mode shapes show how much each mass moves (absolute displacement, normalized). Small entries for m₁ and m₄ mean the dampers barely act, so the mode decays slowly.

Energy

Physical energy $\mathcal{R}^*+\mathcal{A}$ of the selected model's motion. With dampers it must never increase; the stat above reports any increase.

Animation (selected model, displacements exaggerated)

Results table (every 0.5 s)

Review questions

Why can the dampers not be included in the Lagrangian?

Their forces depend on velocity and dissipate energy; they cannot be derived from a potential. They enter Lagrange's equations as generalized forces $\Xi_i$, obtained from virtual work or from the Rayleigh dissipation function $\mathcal{D}$.

How do you obtain a generalized force correctly?

Write the virtual work of the physical force in terms of the variations of the generalized coordinates. For example, $\delta W = -B_2 v_4\,\delta y_4 = -B_2 v_4(\delta x_1 + \delta x_4)$, so the force appears in both $\Xi_1$ and $\Xi_4$.

Why must the mass matrix be symmetric?

It comes from the kinetic energy $\tfrac12\dot{\mathbf x}^{\mathsf T}\mathbf M\dot{\mathbf x}$. Only the symmetric part of a matrix contributes to a quadratic form, and the Lagrange equations produce exactly that symmetric part. A non-symmetric coefficient pattern signals an algebra error.

What changes if we use absolute displacements instead of spring elongations?

The mass matrix becomes diagonal, while the stiffness matrix becomes non-diagonal (springs couple neighbours). The physics is identical; only where the coupling appears changes. For Simulink, a diagonal mass matrix avoids algebraic loops.

Can a single mass have its "own" damping in a coupled system?

No. Damping belongs to the modes. Each mass's motion is a combination of all modes, each decaying at its own rate. A mass can look weakly damped because it takes part strongly in a mode that the dampers hardly affect.

Does the motion always die out completely?

Only if every mode makes at least one damper move. Here an undamped mode would need $m_1$ and $m_4$ to stay at rest, which in turn forces $m_2$ and $m_3$ to stay at rest. So all modes decay, though some can decay very slowly.

Explanation of Method