A pendulum with a horizontally oscillating fulcrum can be described using the Lagrange formalism (related to a setting, as it is shown e.g. here) by the following system of nonlinear ordinary differential equations (ODE):

$$\begin{cases} m\ddot x+kx+ml\left(\ddot \theta \cos \left(\theta\right)-\dot \theta^2 \sin \left(\theta\right)\right) & \text{= 0} \\ \cos \left(\theta\right) \ddot x+g\sin \left(\theta\right)+l \ddot \theta & \text{= 0} \end{cases},$$

with $m$ being the pendulum's point mass, $l$ being the length of the pendulum's massless, inextensible and always taut cord, $k$ being the spring constant, $g \approx 9.81 \frac{m}{s^2}$ being Earth's average gravitational acceleration, $x$ being the horizontal deflection of the spring, and $\theta$ being the angle of rotation of the pendulum.

After I had comprehended the derivation, I asked myself how the equations of motion would look like in concrete terms. Of course, it's hardly possible to find analytical solutions for these equations.

But what is the easiest way to generate numerical approximations? Is it possible online, with appropriate programs (e.g. Scilab; if so, what would be the code?) or at all? Does the system have to be somehow manipulated, linearized, simplified ... before (if so, how?)? Which and how many boundary conditions are necessary?

Thank you and all the best!

  • $\begingroup$ The simple pendulum also has no closed-form solution in general. The usual closed-form solution assumes a very small angle of deflection. A similar assumption here would be that $\theta$ is small and the range of motion of $x$ is also small compared to $l$, and then we have the simplifying assumptions $\cos\theta \approx 1$ and $\sin\theta\approx\theta.$ $\endgroup$
    – David K
    Feb 24, 2020 at 18:15
  • $\begingroup$ A similar question is Approximating spring-cart-pendulum system, with python and javascript codes for RK4 in the answers. $\endgroup$ Feb 25, 2020 at 9:53
  • $\begingroup$ This is quite the strange situation with the mass $M=0$ of the fulcrum. Imagine an initial situation where you pull the spring but let the pendulum hang straight down. Then you let loose. In the initial moment the equations say that the system is not solvable. Or more physically, that the acceleration of the fulcrum is infinite, which translates to an infinite acceleration of the pendulum. $\endgroup$ Feb 25, 2020 at 12:39

2 Answers 2


Assume that the fulcrum has a positive mass $M$, as small as it may be. Then you get Lagrangian terms $\frac{M}2\dot x^2-\frac{k}2x^2$ for the spring.

The pendulum moves along the curve $x_P(t)=x(t)+l\sinθ(t)$, $y_P(t)=l(1-\cosθ(t))$. Its velocity and kinetic energy are thus \begin{align} \dot x_P&=\dot x+l\dotθ\cosθ\\ \dot y_P&=l\dotθ\sinθ\\ \frac{m}2\left(\dot x_P^2+\dot y_P^2\right)&=\frac{m}2\left(\dot x^2+2l\cosθ\,\dotθ\dot x+l^2\dotθ^2\right) \end{align} so that the full Lagrangian is $$ L=\frac{M+m}2\dot x^2+\frac{m}2\left(2l\cosθ\,\dotθ\dot x+l^2\dotθ^2\right) -\frac{k}2x^2-lg(1-\cosθ) $$ Now define the impulses $\newcommand{\pd}[2]{\frac{\partial#1}{\partial#2}}$ \begin{align} p&=\pd{L}{\dot x}=(M+m)\dot x+ml\cosθ\,\dotθ\\ \pi&=\pd{L}{\dot θ}=ml^2\dotθ+ml\cosθ\,\dot x \end{align} so that the Euler-Lagrange equations now read \begin{align} \dot p&=\pd{L}{x}=-kx\\ \dot\pi&=\pd{L}{θ}=-ml\sinθ\,\dotθ\dot x-lg\sinθ \end{align}

The first system is a linear system for $\dot x,\dotθ$ with determinant $$ml^2(M+m\sin^2θ).$$ One can see that this becomes problematic if $M$ is small or zero and the pendulum goes through the lower part of its arc, $θ\approx 0$.

In detail, the solution of this system is \begin{align} \dot x &= \frac{lp-\cos θ\pi}{l(M+m\sin^2θ)}\\[.5em] \dot θ &= \frac{(M+m)\pi-ml\cosθ\,p}{ml^2(M+m\sin^2θ)} \end{align}

This first order system can now be given to any numerical solver function, set $M=m/1000$ to get the effect of a very lightweight fulcrum.


Considering a fulcrum mass $M$, the Lagrangian modeling is straightforward

Kinetic energy

$$ T = \frac{1}{2} \left(m \left(l^2 \theta'^2 \sin ^2(\theta)+\left(l \theta' \cos (\theta)+x_1'\right)^2\right)+M x_1'^2\right) $$

Potential energy

$$ V = g l m (1-\cos (\theta))+\frac{1}{2} k x_1^2 $$

with the lagrangian

$$ L = T-V $$

giving the movement equations

$$ \left\{ \begin{array}{rcl} x_1''&=&\frac{2 m \sin (\theta) \left(l \theta'^2+g \cos (\theta)\right)-2 k x_1}{m+2 M-\cos (2 \theta) m} \\ \theta ''&=&\frac{\sec (\theta) \left(k x_1-\tan (\theta) \left(l m \cos (\theta) \theta'^2+g (m+M)\right)\right)}{l \left((m+M) \tan^2(\theta)+M\right)} \\ \end{array} \right. $$

as can be observed, if $M=0$ the denominators for $x_1''$ and $\theta''$ become $0$ as $\theta\to 0$. Considering $M<<m$ we can simulate the beautiful patterns appearing for $x_p = x_1+l\sin(\theta), y_p = -l\cos(\theta)$ executing the following MATHEMATICA script.

xp[t] = x1[t] + l Sin[theta[t]];
yp[t] = -l Cos[theta[t]];
T = 1/2 (M D[x1[t], t]^2 + m (D[xp[t], t]^2 + D[yp[t], t]^2));
V = 1/2 k x1[t]^2 + l (1 - Cos[theta[t]]) m g;
L = T - V;
vars = {x1[t], theta[t]};
equs = Grad[L, vars] - D[Grad[L, D[vars, t]], t];
sols = Solve[equs == 0, D[vars, {t, 2}]][[1]] // FullSimplify;
equs = Thread[D[vars, {t, 2}] == (D[vars, {t, 2}] /. sols)];
parms = {l -> 1, m -> 1, M -> 1/10, k -> 100, g -> 10};
tmax = 10;
equs0 = equs /. parms;
cinits = {x1[0] == 0, x1'[0] == -4, theta[0] == 0, theta'[0] == 0};
solx = NDSolve[Join[equs0, cinits], {x1, theta}, {t, 0, tmax}];
ParametricPlot[Evaluate[{xp[t], yp[t]} /. parms /. solx], {t, 0, tmax},  PlotRange -> All]

enter image description here


You must log in to answer this question.

Not the answer you're looking for? Browse other questions tagged .