$ \DeclareMathOperator{\arccosh}{arccosh} \DeclareMathOperator{\arcsinh}{arcsinh} \DeclareMathOperator{\rank}{rank} \DeclareMathOperator{\rot}{rot} \DeclareMathOperator{\grad}{grad} \DeclareMathOperator{\diver}{div} $

二重振り子

通常の振り子のおもりにさらに振り子を付けた装置を二重振り子といいます。
今回は二重振り子の運動方程式を Lagrange の運動方程式を用いて立ててみたいと思います。

1装置

二重振り子の装置
座標やおもりの重さなどは上の図のようにとることとします。
また、
\begin{equation} \left(\begin{matrix} x_1 \\ y_1 \end{matrix}\right) = l_1 \left(\begin{matrix} \sin\theta_1 \\ \cos\theta_1 \end{matrix}\right) \end{equation}
(1)
および、
\begin{equation} \left(\begin{matrix} x_2 \\ y_2 \end{matrix}\right) = l_1 \left(\begin{matrix} \sin\theta_1 \\ \cos\theta_1 \end{matrix}\right) + l_2 \left(\begin{matrix} \sin\theta_2 \\ \cos\theta_2 \end{matrix}\right) \end{equation}
(2)
に注意します。

2Lagrangian の導出

それでは Lagrangian を計算しましょう。
まず、運動エネルギー $T$ は、
\begin{eqnarray*} T & = & \frac{1}{2} m_1 \left((l_1\dot\theta_1\cos\theta_1)^2 + (l_1\dot\theta_1\sin\theta_1)^2\right) + \\ & & \frac{1}{2} m_2 \left((l_1\dot\theta_1\cos\theta_1 + l_2\dot\theta_2\cos\theta_2)^2 + (l_1\dot\theta_1\sin\theta_1 + l_2\dot\theta_2\sin\theta_2)^2\right) \\ & = & \frac{1}{2} m_1 l_1^2 \dot\theta_1^2 + \frac{1}{2} m_2 \left(l_1^2 \dot\theta_1^2 + l_2^2 \dot\theta_2^2 + 2l_1l_2\dot\theta_1\dot\theta_2\cos(\theta_1-\theta_2)\right) \end{eqnarray*}
(3)
一方、ポテンシャルエネルギー $V$ は、
\begin{equation} V = - m_1 g \times l_1 \cos\theta_1 - m_2 g \times (l_1 \cos\theta_1 + l_2 \cos\theta_2) \end{equation}
(4)
です。
従って Lagrangian $L$ は、
\begin{equation} L = \frac{1}{2} m_1 l_1^2 \dot\theta_1^2 + \frac{1}{2} m_2 \left(l_1^2 \dot\theta_1^2 + l_2^2 \dot\theta_2^2 + 2l_1l_2\dot\theta_1\dot\theta_2\cos(\theta_1-\theta_2)\right) + m_1 g l_1 \cos\theta_1 + m_2 g (l_1 \cos\theta_1 + l_2 \cos\theta_2) \end{equation}
(5)
となります。

3運動方程式の導出

それでは運動方程式を求めましょう。 座標変数として $\theta_1, \theta_2$ があるので、運動方程式は二本得ることができます。

3.1$\theta_1$ から

\begin{eqnarray*} \frac{\partial L}{\partial \dot\theta_1} & = & m_1l_1^2\dot\theta_1 + \frac{1}{2}m_2\left(2l_1^2\dot\theta_1+2l_1l_2\dot\theta_2\cos(\theta_1-\theta_2)\right) \\ & = & m_1l_1^2\dot\theta_1 + m_2\left(l_1^2\dot\theta_1+l_1l_2\dot\theta_2\cos(\theta_1-\theta_2)\right) \\ \therefore \frac{d}{dt}\frac{\partial L}{\partial \dot\theta_1} & = & m_1l_1^2\ddot\theta_1 + m_2\left(l_1^2\ddot\theta_1 + l_1l_2\ddot\theta_2\cos(\theta_1-\theta_2) - l_1l_2\dot\theta_2\sin(\theta_1-\theta_2) \times (\dot\theta_1-\dot\theta_2)\right) \end{eqnarray*}
(6)
また、
\begin{eqnarray*} \frac{\partial L}{\partial \theta_1} & = & \frac{1}{2}m_2\left(2l_1l_2\dot\theta_1\dot\theta_2(-\sin(\theta_1-\theta_2))\right) - m_1gl_1\sin\theta_1 - m_2gl_1\sin\theta_1 \end{eqnarray*}
(7)
なので、 運動方程式は、
\begin{equation} \frac{d}{dt}\frac{\partial L}{\partial \dot\theta_1} = \frac{\partial L}{\partial \theta_1} \end{equation}
(8)
\begin{equation} \Leftrightarrow (m_1+m_2)l_1^2\ddot\theta_1 + m_2l_1l_2\cos(\theta_1-\theta_2)\ddot\theta_2 = - m_2l_1l_2\dot\theta_2^2\sin(\theta_1-\theta_2) - (m_1+m_2)gl_1\sin\theta_1 \end{equation}
(9)
となります。

3.2$\theta_2$ から

次にもう一つの運動方程式を求めます。
\begin{eqnarray*} \frac{\partial L}{\partial \dot\theta_2} & = & m_2l_2^2\dot\theta_2 + m_2l_1l_2\dot\theta_1\cos(\theta_1-\theta_2) \\ \therefore \frac{d}{dt}\frac{\partial L}{\partial \dot\theta_2} & = & m_2l_2^2\ddot\theta_2 + m_2l_1l_2\left(\ddot\theta_1\cos(\theta_1-\theta_2) - \dot\theta_1\sin(\theta_1-\theta_2) \times (\dot\theta_1-\dot\theta_2)\right) \end{eqnarray*}
(10)
また、
\begin{equation} \frac{\partial L}{\partial \theta_2} = m_2l_1l_2\dot\theta_1\dot\theta_2\sin(\theta_1-\theta_2) - m_2gl_2\sin\theta_2 \end{equation}
(11)
なので、 運動方程式は、
\begin{equation} \frac{d}{dt}\frac{\partial L}{\partial \dot\theta_2} = \frac{\partial L}{\partial \theta_2} \end{equation}
(12)
\begin{equation} \Leftrightarrow l_1\cos(\theta_1-\theta_2)\ddot\theta_1 + l_2\ddot\theta_2 = l_1\dot\theta_1^2\sin(\theta_1-\theta_2) - g\sin\theta_2 \end{equation}
(13)
となります。

3.3連立方程式を解く

さて、上で運動方程式が一応出たわけですが、この形だと $\ddot\theta_1$$\ddot\theta_2$ がごちゃ混ぜになっていて、とても役に立つ形ではありません。 なので上の(9)式、(13)式を $\ddot\theta_1, \ddot\theta_2$ の連立一次方程式だと思って解き、 $\ddot\theta_1 = ???, \ddot\theta_2 = ???$ という形に変形しましょう。
ただの連立一次方程式です。 やることは簡単だよね。 途中で記号がたくさん出てきてヤバイことになるますがー。
とりあえず頑張って解くとこんな感じになると思います[1]
\begin{eqnarray*} \ddot\theta_1 & = & \frac{-m_1g\sin\theta_1 - m_2\left(g\sin\theta_1+l_2\dot\theta_2^2\sin(\theta_1-\theta_2)+\left(l_1\dot\theta_1^2\sin(\theta_1-\theta_2)-g\sin\theta_2\right)\cos(\theta_1-\theta_2)\right)}{l_1\left(m_1+m_2\left(\sin(\theta_1-\theta_2)\right)^2\right)} \\ \ddot\theta_2 & = & \frac{(m_1+m_2)\left(l_1\dot\theta_1^2\sin(\theta_1-\theta_2)-g\sin\theta_2+g\sin\theta_1\cos(\theta_1-\theta_2)\right) + m_2l_2\dot\theta_2^2\cos(\theta_1-\theta_2)\sin(\theta_1-\theta_2)}{l_2\left(m_1+m_2\left(\sin(\theta_1-\theta_2)\right)^2\right)} \end{eqnarray*}
(14)
ヤバイですね、この式。目眩がします。
もちろん、計算が合っている保証も自信もありません。
ええ、ありませんとも。
…でも下のシミュレーションがこの式で「それっぽく」動いているのでたぶん合っているんでしょうねw[2]
まぁ間違っていたとしても、少なくとも、某物性化学の試験[3]の3時間前とかいう、よく分からないテンションになる時に計算しないと心が折れそうになる、ということが理解できたのでよしとします。

[1]2009/08/05 修正
[2]検算してください、誰か!
[3]化学、苦手なのです

4シミュレーション

4.1アプレット

4.2説明と使い方

s1 は $\theta_1$ を、 s2 は $\theta_2$ をそれぞれ表します。
おもりの重さはともに $m_1 = m_2 = 1\text{kg}$、棒の長さは $l_1 = 1.0\text{m}, l_2 = 0.8\text{m}$ だったと思います。
また開始時のおもりの位置を指定したいときには、画面上を一度クリックすると、クリックした方向におもり $m_1$ が向き、 さらにもう一度クリックすると、その方向におもり $m_2$ が向いた状態でシミュレーションが開始されます。
また、後ろに赤茶色で軌跡を描いていますが、ばね振り子と違ってあまりきれいな軌跡にはならない模様です。 残念。