Cylinder rolling in a groove on a spring-mounted cart

A system with two degrees of freedom: a cart on two springs and a cylinder rolling without slipping inside a circular groove on the cart. The coupled nonlinear equations of motion come from Lagrange's equations, then from a Simulink model, and both can be run below.

Problem statement

Derive the equations of motion of the system shown in the figure. A circular cylinder of mass $m$ and radius $r$ rolls without slipping inside a circular groove of radius $R$ cut into a cart of mass $M$. The cart is attached to fixed walls by two springs of equal stiffness $k$. Use Lagrange's equations.

Figure 1. Cart M on two springs k; cylinder (m, r) rolls inside a groove of radius R with centre O′. Generalized coordinates: cart displacement x and angle Θ of the line O′C from the vertical.

NoteThe circles in the figure represent cylinders seen end-on (plane motion), not spheres. This matters for the moment of inertia: $I_C=\tfrac12 mr^2$ for a solid cylinder (a sphere would have $\tfrac25 mr^2$).

Generalized coordinates and kinematics

The system has two degrees of freedom: the cart can move horizontally, and the cylinder can move along the groove. The paper chooses:

$x$ – displacement of the cart from the equilibrium position (both springs undeformed); this is also the spring elongation,
$\Theta$ – angle between the vertical through the groove centre $O'$ and the line $O'C$ to the cylinder centre. With the sign convention used in the paper's velocity equations, positive $\Theta$ moves the cylinder toward negative $x$.

The cylinder centre always lies on a circle of radius $R-r$ around $O'$. Its position is

$$x_C = x - (R-r)\sin\Theta,\qquad y_C = -(R-r)\cos\Theta \;(+\text{const})$$

Rolling without slipping

Rolling without slipping means that the arc travelled on the groove equals the arc travelled on the cylinder's circumference ($\overset{\frown}{CB} = \overset{\frown}{CA}$ in the paper). If $\rho$ is the spin angle of the cylinder measured from the vertical:

$$R\,\Theta = r\,(\rho+\Theta)\quad\Longrightarrow\quad \rho = \frac{R-r}{r}\,\Theta,\qquad \omega = \dot\rho = \frac{R-r}{r}\,\dot\Theta$$

Because the cart only translates (it does not rotate), this angular velocity relative to the cart is also the absolute angular velocity of the cylinder – exactly what the kinetic energy needs.

Velocity of the cylinder centre

The absolute velocity of the centre $C$ is the vector sum of the cart velocity $\dot x$ and the velocity of $C$ relative to the cart, which has magnitude $(R-r)\dot\Theta$ and is tangent to the circle:

$$v_{Cx} = \dot x - (R-r)\dot\Theta\cos\Theta,\qquad v_{Cy} = (R-r)\dot\Theta\sin\Theta$$

$$v_C^2 = \dot x^2 + (R-r)^2\dot\Theta^2 - 2\dot x\dot\Theta(R-r)\cos\Theta$$

The term $-2\dot x\dot\Theta(R-r)\cos\Theta$ is the coupling between the cart and the cylinder: it is what makes the two motions influence each other.

Energies of the system

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

3.1 Attractive (potential) energy

The two springs have the same stiffness and elongations $+x$ and $-x$; after squaring they store the same energy. The gravitational part depends on how high the cylinder centre rises above its lowest position:

$$z = (R-r) - (R-r)\cos\Theta = (R-r)(1-\cos\Theta)$$

$$\mathcal{A} = \frac{kx^2}{2}+\frac{kx^2}{2} + mgz = kx^2 + mg(R-r)(1-\cos\Theta)$$

3.2 Repulsive (kinetic) coenergy

The cart translates with $\dot x$; the cylinder translates with $v_C$ and spins with $\omega$. The moment of inertia of a solid cylinder about its centroidal axis is $I_C = \tfrac12 m r^2$:

$$\mathcal{R}^* = \tfrac12 m v_C^2 + \tfrac12 M\dot x^2 + \tfrac12 I_C\,\omega^2$$

$$\tfrac12 I_C\,\omega^2 = \frac12\cdot\frac12 m r^2\,\frac{(R-r)^2}{r^2}\dot\Theta^2 = \frac14 m (R-r)^2\dot\Theta^2$$

$$\mathcal{R}^* = \frac12 (M+m)\dot x^2 + \frac34 m(R-r)^2\dot\Theta^2 - m\dot x\dot\Theta(R-r)\cos\Theta$$

3.3 Lagrangian

$$\mathcal{L} = \frac12 M\dot x^2 + \frac14 m(R-r)^2\dot\Theta^2 + \frac m2\dot x^2 + \frac m2 (R-r)^2\dot\Theta^2 - m\dot x\dot\Theta(R-r)\cos\Theta - kx^2 - mg(R-r)(1-\cos\Theta)$$

Lagrange equations

There is no dissipation and no external input, so the right-hand side of both Lagrange equations is zero:

$$\frac{d}{dt}\frac{\partial\mathcal{L}}{\partial\dot q_i} - \frac{\partial\mathcal{L}}{\partial q_i} = 0,\qquad q_1=\Theta,\; q_2=x$$

4.1 Equation for $\Theta$

$$\frac{\partial\mathcal{L}}{\partial\dot\Theta} = \frac32 m(R-r)^2\dot\Theta - m\dot x(R-r)\cos\Theta$$

$$\frac{d}{dt}\frac{\partial\mathcal{L}}{\partial\dot\Theta} = \frac32 m(R-r)^2\ddot\Theta - m(R-r)\left(\ddot x\cos\Theta - \dot x\dot\Theta\sin\Theta\right)$$

$$\frac{\partial\mathcal{L}}{\partial\Theta} = m\dot x\dot\Theta(R-r)\sin\Theta - mg(R-r)\sin\Theta$$

The terms $\pm m(R-r)\dot x\dot\Theta\sin\Theta$ cancel. Dividing by $m(R-r)$:

(Eq. 1)$$\frac32 (R-r)\ddot\Theta - \ddot x\cos\Theta + g\sin\Theta = 0$$

4.2 Equation for $x$

$$\frac{\partial\mathcal{L}}{\partial\dot x} = (M+m)\dot x - m\dot\Theta(R-r)\cos\Theta$$

$$\frac{d}{dt}\frac{\partial\mathcal{L}}{\partial\dot x} = (M+m)\ddot x - m(R-r)\left(\ddot\Theta\cos\Theta - \dot\Theta^2\sin\Theta\right),\qquad \frac{\partial\mathcal{L}}{\partial x} = -2kx$$

(Eq. 2)$$(M+m)\ddot x - m(R-r)\left(\ddot\Theta\cos\Theta - \dot\Theta^2\sin\Theta\right) + 2kx = 0$$

These are two coupled nonlinear second-order differential equations. Each one contains the highest derivative of the other coordinate ($\ddot x$ appears in Eq. 1, $\ddot\Theta$ in Eq. 2).

4.3 Form used for Simulink (as in the paper)

(Eq. 1)$$\ddot\Theta = \frac23\frac{\ddot x\cos\Theta}{R-r} - \frac23\frac{g\sin\Theta}{R-r}$$

(Eq. 2)$$\ddot x = \frac{m(R-r)\ddot\Theta\cos\Theta}{M+m} - \frac{m(R-r)\dot\Theta^2\sin\Theta}{M+m} - \frac{2kx}{M+m}$$

Simulink implementation

Each equation feeds a chain of two integrators ($\ddot\Theta\to\dot\Theta\to\Theta$ and $\ddot x\to\dot x\to x$). In the paper's form, the block computing $\ddot\Theta$ needs $\ddot x$, and the block computing $\ddot x$ needs $\ddot\Theta$. This creates an algebraic loop (a signal path with no integrator in it), shown dashed below.

Figure 2. Structure of the Simulink model for the paper's form of the equations. Simulink can solve the algebraic loop iteratively, but this is slower and can be less accurate.

Explicit form (no algebraic loop)

A cleaner approach is to treat Eq. 1 and Eq. 2 as a linear system in the unknowns $\ddot x$ and $\ddot\Theta$ and solve it once:

$$\begin{bmatrix} M+m & -m(R-r)\cos\Theta \\ -\cos\Theta & \tfrac32(R-r)\end{bmatrix}\begin{bmatrix}\ddot x\\ \ddot\Theta\end{bmatrix} = \begin{bmatrix} -m(R-r)\dot\Theta^2\sin\Theta - 2kx \\ -g\sin\Theta\end{bmatrix}$$

$$\ddot x = \frac{\tfrac32\left(-m(R-r)\dot\Theta^2\sin\Theta - 2kx\right) - m g\sin\Theta\cos\Theta}{\tfrac32(M+m) - m\cos^2\Theta}$$

$$\ddot\Theta = \frac{-(M+m)g\sin\Theta + \cos\Theta\left(-m(R-r)\dot\Theta^2\sin\Theta - 2kx\right)}{(R-r)\left[\tfrac32(M+m) - m\cos^2\Theta\right]}$$

The denominator $\tfrac32(M+m) - m\cos^2\Theta \ge \tfrac32 M + \tfrac12 m > 0$ never vanishes, so the explicit form is always well defined. This is the form the simulator below uses.

Simulation results from the paper

The paper shows four simulations, each with the cylinder angle (in rad) and the cart displacement $x$ plotted together. The initial conditions are not stated in the paper. Matching the plots numerically shows that they were $x(0)=0.2$ m (first case) or $0.3$ m (others), with $\Theta(0)=0$ and both velocities zero. The buttons load each case into the simulator.

Case	M	m	r	R	k	x(0)	What the paper observes
A	2	1	0.05	1	100	0.2	Cart oscillates almost sinusoidally (±0.2 m); Θ stays within ±0.29 rad with irregular-looking peaks.
B	2	1	0.05	1	20	0.3	Softer springs: lower frequency, stronger exchange of energy between cart and cylinder (Θ up to 0.5 rad).
C	2	1	0.05	0.5	100	0.3	Smaller groove: the cylinder swings much more (Θ up to ≈0.85 rad).
D	6	3	0.1	1.5	800	1 *	Θ drifts to tens of radians (the cylinder "goes around"); the cart keeps oscillating.

* Case D is described in the paper as "large initial conditions (> 1 m)". This motion is chaotic, so the exact curve depends strongly on the initial condition and the solver. The simulator reproduces the behaviour qualitatively, not point for point.

Paper's conclusions

The system oscillates without damping, which is correct because there are no dissipative forces or external inputs.
For smaller $k$ the oscillation frequency is lower, and the cart's motion depends more on the cylinder. For large $k$ the cart oscillates almost as a pure sinusoid. For small $k$ a periodically repeating change of amplitude appears (beats).
For large masses, large initial conditions and very large $k$, the cylinder angle grows without bound while the cart keeps oscillating. The paper interprets this as the cylinder jumping out of the groove (see section 7).

Review and extensions

The derivation is correctEvery step was checked: the rolling condition, $v_C^2$, the energies, the Lagrangian, the partial derivatives and both final equations. A numerical simulation of these equations reproduces the paper's plots A–C closely, and total energy is conserved to machine precision (see the energy plot in the simulator).

7.1 Why the curves look irregular: two normal modes

Linearizing for small $\Theta$ ($\sin\Theta\approx\Theta$, $\cos\Theta\approx1$, $\dot\Theta^2\approx0$) gives a linear 2-DOF system:

$$\underbrace{\begin{bmatrix} M+m & -m\ell \\ -m\ell & \tfrac32 m\ell^2\end{bmatrix}}_{\mathbf{M}}\begin{bmatrix}\ddot x\\ \ddot\Theta\end{bmatrix} + \underbrace{\begin{bmatrix} 2k & 0 \\ 0 & mg\ell\end{bmatrix}}_{\mathbf{K}}\begin{bmatrix}x\\ \Theta\end{bmatrix} = \mathbf{0},\qquad \ell = R-r$$

The natural frequencies follow from $\det(\mathbf{K}-\omega^2\mathbf{M})=0$. The response is a sum of two harmonics with different frequencies. When they are close, the amplitude rises and falls periodically (beats); this is the "periodically repeating change of amplitude" noted in the paper. For case B ($k=20$) the frequencies are 0.39 Hz and 0.71 Hz; for case A ($k=100$) they are 0.41 Hz and 1.49 Hz. The cart then follows mainly the fast spring mode, which is why it looks almost sinusoidal for large $k$.

7.2 "The cylinder jumps out of the groove" – what the model really says

In the model the groove is a complete circle of radius $R$, and the cylinder is assumed to stay in contact with it at all times. A $\Theta$ that grows without bound therefore means that the cylinder runs around the whole circle (one revolution = $2\pi$ rad), not that it jumps out. The model cannot describe leaving the groove at all. Physically, two things limit validity:

the groove on a real cart is only a short arc, so the cylinder leaves it as soon as $|\Theta|$ exceeds the arc's half-angle;
contact is maintained only while the normal force is positive. From the radial equation of motion of the cylinder:

$$N = m\left(g\cos\Theta + \ddot x\sin\Theta + (R-r)\dot\Theta^2\right) \;>\; 0$$

When $N$ becomes negative, the cylinder separates from the groove and the equations of motion are no longer valid. The simulator computes $N(t)$ and reports when contact is lost. For case D this happens after only ≈0.13 s. Even in case C, $N$ briefly drops below zero (to about −9 N) when the cart decelerates hard while the cylinder is high in the groove. Cases A and B stay in contact the whole time (in case A, $N \ge 4.7$ N).

7.3 Energy check

Without friction, $\mathcal{R}^*+\mathcal{A}$ must stay constant. This is a simple and powerful way to test both the model and the numerical solver. The ragged $\Theta$ curve in the paper's case D may come from a coarse solver step or the iterative solution of the algebraic loop. The explicit form with a fine step conserves energy accurately.

Interactive simulator

The equations are solved numerically in explicit form (4th-order Runge–Kutta, 1 ms step). Angle inputs also accept expressions such as pi/6. An optional viscous friction coefficient $b$ on the cart adds a term $-b\dot x$ to Eq. 2, so you can see how damping changes the response.

Cases from the paper

Also show the linearized model

System parameters

M (cart) [kg]

m (cylinder) [kg]

r [m]

R (groove) [m]

k (each spring) [N/m]

g [m/s²]

b (cart friction) [N·s/m]

Initial conditions

x(0) [m]

ẋ(0) [m/s]

Θ(0) [rad]

Θ̇(0) [rad/s]

Duration [s]

Cart displacement x(t)

Cylinder angle Θ(t)

Normal force N(t)

Contact requires N > 0

Phase portrait Θ̇(Θ)

Energy balance

$\mathcal{R}^*+\mathcal{A}+W_b$ must stay constant ($W_b$ = energy dissipated by friction).

Animation

Results table (every 0.1 s)

Review questions

Why is the cylinder's angular velocity $\frac{R-r}{r}\dot\Theta$ and not $\frac{R}{r}\dot\Theta$?

The arc on the groove ($R\Theta$) equals the arc on the cylinder, $r(\rho+\Theta)$. Part of the turning ($\Theta$) comes from the cylinder's centre going around $O'$, so the spin is $\rho = \frac{R-r}{r}\Theta$. Check: if $r\to R$ the cylinder fills the groove and cannot spin.

Why can the spin $\omega$ be used directly in the kinetic energy, even though the cart moves?

The cart only translates. A translating frame does not rotate, so angular velocities relative to the cart are equal to absolute angular velocities.

What is the physical meaning of the term $-m\dot x\dot\Theta(R-r)\cos\Theta$ in $\mathcal{R}^*$?

It is the dynamic coupling: the cylinder's absolute velocity depends on both the cart's velocity and its own motion in the groove. Without it, the equations would separate into an independent spring–mass system and an independent "pendulum".

What happens if the cart is held fixed ($x\equiv 0$)?

Eq. 1 becomes $\frac32(R-r)\ddot\Theta + g\sin\Theta = 0$: a pendulum with small-oscillation frequency $\omega=\sqrt{\frac{2g}{3(R-r)}}$. The factor $\frac32$ comes from the cylinder's rotational inertia.

Why do the oscillations never decay, and how would you add damping?

There is no dissipation, so energy is conserved. Viscous friction on the cart adds a generalized force $Q_x=-b\dot x$ to the right-hand side of the $x$ equation (Rayleigh dissipation function $\tfrac12 b\dot x^2$). Try the "Damped" preset.

Does $\Theta\to\infty$ in case D mean the cylinder flies out of the groove?

Not in the model: it means the cylinder goes all the way around a full circular groove. In reality the model stops being valid earlier, when the normal force $N$ becomes negative or when the cylinder reaches the edge of the groove.