It is agreed that in an actual experiment the initial motion out of the vertical equilibrium state is likely to be dominated by a rolling motion. Presentations of the corresponding stability analyses for chains of two and three balls can be found in the supplementary information below.
The critical value of alpha = 3 for N = 2 balls is confirmed. For N = 3 balls, the critical value is alpha = 18.55 being only about 6% below the sliding case, for which the corresponding critical value is 19.68. This indicates that for larger values of N the stability thresholds for the sliding and rolling cases are likely to be very close to one another. Therefore, the idealized analysis in the original paper still provides a reasonable approximation and the results can be compared to experimental measurements for chains of sufficient length.
We use equation (2.4) of the paper for the dimensionless potential energy
\[
E=\frac{\alpha}{6}\sum_{i< j}\frac{(\mathbf{n}_i\cdot\mathbf{n}_j)x_{ij}^2-3(\mathbf{n}_i\cdot\mathbf{x}_{ij})(\mathbf{n}_j\cdot\mathbf{x}_{ij})}{x_{ij}^5}
+\mathbf{e}_2\cdot\sum_i\mathbf{w}_i\,.
\]
With the angle \(t\) defined as in Figure 1, and
\[\mathbf{w}_1=(0,0)\,,\]
\[\mathbf{n}_1=\mathbf{e}_2=(0,1)\,,\]
\[\mathbf{w}_2=\mathbf{x}_{12}=(\sin t,\cos t)\,,\]
\[\mathbf{n}_2=(\sin(2t),\cos(2t))\,,\]
\[ x_{12}=1\,,\]
we then obtain
\[
E=-\frac{\alpha}{12}\Big[\cos(2t)+3\Big]+\cos t\,,
\]
which leads to the critical \(\alpha\) for stability
\[\frac{d^2 E}{d t^2}=\frac{\alpha}{3}\cos(2t)-\cos t\,,\]
\[\frac{d^2 E}{d t^2}\bigg|_{t=0}=0\quad\Leftrightarrow\quad\alpha_{\sf c}=3\,.\]
With the same energy as used above and with the angles \(t\) and \(s\) defined as in Figure 3 we have
\[
E=-\frac{\alpha}{12}\big[M_{12}+M_{23}+M_{13}\big]+G_2+G_3
\]
with
\[
M_{12}=\cos(2t)+3\,,
\]
\[
M_{23}=\cos(2s)+3\,,
\]
\[
M_{13}=\frac{3\cos\big(\frac{t+3s}{2}\big)\cos\big(\frac{3t+s}{2}\big)-\cos(2t+2s)}{4\cos^3\big(\frac{t+s}{2}\big)}\,,
\]
\[
G_2=\cos t\,,
\]
\[
G_3=\cos t+\cos(2t+s)\,.
\]
The determinant of the Hessian matrix \(\mathbf{H}=\begin{pmatrix}E_{tt}&E_{ts}\\E_{st}&E_{ss}\end{pmatrix}\) evaluated at \(t=s=0\) is
\[
\det \mathbf{H}\big|_{t=s=0}=2-65\alpha/24+323\alpha^2/2304\,,
\]
the two roots of this polynomial are given by
\[
\alpha_{1,2}=(65\pm\sqrt{3579})\,48/323\,,
\]
the larger root being \(\alpha_1\approx 18.55\). After the vertical solution (\(t=s=0\)) destabilizes below \(\alpha_1\), there are two symmetric (left/right) bend solutions which then destabilize around \(\alpha\approx 17.58\) and the tower crashes down.