読者です 読者をやめる 読者になる 読者になる

倭算数理研究所

科学・数学・学習関連の記事を、「倭マン日記」とは別に書いていくのだ!

2次元の調和振動子

古典力学 1粒子系 2次元 解析力学 調和振動子

古典力学のいろいろな系で運動方程式を解いていくシリーズ(目次)。 以前の記事で1次元の調和振動子運動方程式を解きましたが、今回は2次元の調和振動子を解きます。 直交座標と極座標でそれぞれ解いてみます。

この記事の内容

2次元調和振動子のポテンシャル

ばね定数が  { k } の(等方的な)2次元の調和振動子は、ポテンシャル  { V(x,\,y) }

  { \displaystyle\begin{align*}
  V(x,\,y) = \frac{1}{2}k(x^2 + y^2) = \frac{1}{2}kr^2
    \qquad \left(r = \sqrt{x^2 + y^2}\right)
\end{align*}}

で与えられます。 振動する質点の位置ベクトルを  { \textbf{r} }、質量を  { m } としましょう。

直交座標で解く

 { \nabla r^2 = 2r \nabla r = 2 \textbf{r} } より*1、直交座標での運動方程式

  { \displaystyle\begin{align*}
  m\frac{d^2 \textbf{r}}{dt^2} &= -k \textbf{r}
\end{align*}}

もしくは、成分ごとに書くと

  { \displaystyle\begin{align*}
  m\frac{d^2 x}{dt^2} &= -k x \\
  m\frac{d^2 y}{dt^2} &= -k y
\end{align*}}

となります。 それぞれの方程式は1次元の調和振動子運動方程式と等しいので、解は1次元の結果を使って簡単に書き下せます。 初期位置、初速度のベクトルをそれぞれ  { \textbf{r}_0,\,\textbf{v}_0 } とすると、 { \omega = \sqrt{\frac{k}{m}} } として

  { \displaystyle\begin{align*}
  \textbf{r}(t) &= \textbf{r}_0\cos \omega t + \frac{\textbf{v}_0}{\omega}\sin\omega t \\
  \textbf{v}(t) &= \textbf{v}_0\cos \omega t - \textbf{r}_0\omega\sin\omega t 
\end{align*}}

極座標で解く

調和振動子のポテンシャルは動径  { r } のみに依存するので『中心力ポテンシャル中での質点の運動方程式 ~2次元~』での解法が使えます。 なので、この解法の具体例として、2次元の調和振動子を解いてみましょう。 ただ積分するだけではありますが、初期条件を課すのも含めて結構面倒くさい計算が必要です。

動径成分
 { r,\,t } の関係は以下で与えられます:

  { \displaystyle\begin{align*}
  t - t_0 = \sqrt{\frac{m}{2}}\int \frac{dr}{ \sqrt{E - \frac{1}{2}kr^2 - \frac{\ell^2}{2mr^2}}}
\end{align*}}

ここで  { E,\,t_0 }積分定数です。 右辺の積分 { R } とおいて積分を実行すると

  { \displaystyle\begin{align*}
  R
    &=  \int \frac{dr}{ \sqrt{E - \frac{1}{2}kr^2 - \frac{\ell^2}{2mr^2}}} \\
    &= \int \frac{rdr}{ \sqrt{Er^2 - \frac{1}{2}kr^4 - \frac{\ell^2}{2m}}} \\
    &= \sqrt{\frac{2}{k}} \int \frac{rdr}{ \sqrt{\frac{2E}{k}r^2 - r^4 - \frac{\ell^2}{mk}}} \\
    &= \sqrt{\frac{2}{k}} \int \frac{rdr}{ \sqrt{\frac{E^2}{k^2} - \frac{\ell^2}{mk} - \left(r^2 - \frac{E}{k}\right)^2}}
\end{align*}}

平方根の中が  { 1-\sin^2\varphi } に比例するように、以下で定義される新しい積分変数  { \varphi } を導入しましょう:

  { \displaystyle\begin{align*}
  r^2 - \frac{E}{k} = \sqrt{\frac{E^2}{k^2} - \frac{\ell^2}{mk}}\;\sin\varphi \qquad
  \left(2rdr = \sqrt{\frac{E^2}{k^2} - \frac{\ell^2}{mk}}\;\cos\varphi d\varphi\right)
\end{align*}}

このとき積分  { R }

  { \displaystyle\begin{align*}
  R &= \sqrt{\frac{1}{2k}} \int \frac{\cos\varphi d\varphi}{\sqrt{1-\sin^2\varphi}} \\
     &= \sqrt{\frac{1}{2k}} \int d\varphi \\
     &= \sqrt{\frac{1}{2k}} \varphi \\
     &= \sqrt{\frac{1}{2k}} \sin^{-1} \left(\frac{r^2 - \frac{E}{k}}{\sqrt{\frac{E^2}{k^2} - \frac{\ell^2}{mk}}}\right)
\end{align*}}

 { t - t_0 = \sqrt{\frac{m}{2}}R } だったので、結局

  { \displaystyle\begin{align*}
  t - t_0 = \frac{1}{2}\sqrt{\frac{m}{k}} \sin^{-1} \left(\frac{r^2 - \frac{E}{k}}{\sqrt{\frac{E^2}{k^2} - \frac{\ell^2}{mk}}}\right)
\end{align*}}

となります。 これを  { r^2 } について解くと

  { \displaystyle\begin{align*}
  r^2 &= \sqrt{\frac{E^2}{k^2} - \frac{\ell^2}{mk}}\;\sin\left(2\sqrt{\frac{k}{m}}t +\phi\right) + \frac{E}{k} \qquad
    \left(\phi = -2\sqrt{\frac{k}{m}}\;t_0\right)
\end{align*}}

式を少し簡単にするために、 { m } の代わりに1次元の調和振動子にも出てきた角振動数  { \omega = \sqrt{\frac{k}{m}} } を使って書き換えましょう:

  { \displaystyle\begin{align*}
  r^2 &= \frac{\sqrt{E^2 - \ell^2\omega^2}}{k}\;\sin\left(2\omega t +\phi\right) + \frac{E}{k}
\end{align*}}

両辺の平方根をとると  { r } { t } の関数として求まりますが、その前にこの時点で初期条件を課しておきましょう。

初期条件
中心力ポテンシャル中での質点の運動方程式 ~2次元~』に従って、中心に最も近づく回帰点(の1つ)に質点が位置する時刻が  { t=0 } となるようにしましょう。 これは

  { \displaystyle\begin{align*}
  \left.\frac{dr}{dt}\right|_{t=0} = 0
\end{align*}}

から概ね決められますが、今の場合は  { r }三角関数で書かれているので「質点が中心に最も近い位置にいる時刻を  { t = 0 } とする」という条件から決めた方が簡単です。 このときの動径を  { r_0 } として、初期条件として与えられるとしましょう:

  { \displaystyle\begin{align*}
  r(0) = r_0
\end{align*}}

質点が回帰点にいるとき動径方向の速度は0ですが、偏角方向の速度は持ち、これは角運動量  { \ell } と関係しています。 初速度の偏角成分(角速度ではなく速度の動径方向に直交する成分)  { v_0 } が初期条件として与えられるとすると、 { \ell = r_0 \cdot mv_0 } より

  { \displaystyle\begin{align*}
  v_0 = \frac{\ell}{mr_0}
\end{align*}}

となります。

 { E }積分定数だったので、これも定数で書き下しておきましょう。  { t = 0 } { E } を考えて

  { \displaystyle\begin{align*}
  E 
    &= \left[\frac{1}{2}m\left(\frac{dr}{dt}\right)^2 +  \frac{\ell^2}{2mr^2} + \frac{1}{2}kr^2\right]_{t=0} \\
    &= \frac{1}{2}mv_0^2 + \frac{1}{2}kr_0^2 \\
    &= \frac{kv_0^2}{2\omega^2} + \frac{kr_0^2}{2} \qquad
      \left(\because \omega = \sqrt{\frac{k}{m}}\right) \\
    &= \frac{k}{2}\left(\frac{v_0^2}{\omega^2} + r_0^2\right)
\end{align*}}

ついでに、上記で得られた  { r^2 } の表式に含まれている因子も書き換えておきましょう。

  { \displaystyle\begin{align*}
  \frac{E}{k} &= \frac{1}{2}\left(\frac{v_0^2}{\omega^2} + r_0^2\right)
\end{align*}}

また

  { \displaystyle\begin{align*}
  \frac{E^2 - \ell^2\omega^2}{k^2}
    &= \frac{1}{4}\left(\frac{v_0^2}{\omega^2} + r_0^2\right)^2 - \frac{(mr_0v_0)^2\omega^2}{k^2} \qquad
      \left(\ell = mr_0v_0\right)\\
    &= \frac{1}{4}\left(\frac{v_0^2}{\omega^2} + r_0^2\right)^2 - \frac{r_0^2v_0^2}{\omega^2} \\
    &= \frac{1}{4}\left(\frac{v_0^2}{\omega^2} - r_0^2\right)^2 \\[2mm]
  \frac{\sqrt{E^2 - \ell^2\omega^2}}{k}
    &= \frac{1}{2}\left(\frac{v_0^2}{\omega^2} - r_0^2\right)
\end{align*}}

最後の行では  { \frac{v_0^2}{\omega^2} - r_0^2 > 0 } を使っていますが、これは中心に最も近い回帰点では

  { \displaystyle\begin{align*}
  \frac{\ell^2}{2mr_0^2} > \frac{1}{2}kr_0^2
\end{align*}}

となっていることから導かれます( { \ell = mr_0v_0 })。 初期条件として与える  { r_0,\,v_0 } の値によってはこの条件が満たされない場合がありますが、このとき( { r_0^2 > \frac{v_0^2}{\omega^2} } となっているとき)は、この回帰点が中心から最も遠いものになります。 そして、 { r_0^2 = \frac{v_0^2}{\omega^2} } と等号が成り立つときは、言うまでもないかも知れませんが円軌道になります。

ちなみに、ここで与えた初期条件を直交座標での初期条件として書くと

  { \displaystyle\begin{align*}
  \textbf{r}_0 &= \begin{pmatrix} r_0 \\ 0 \end{pmatrix} \qquad
  \textbf{v}_0 = \begin{pmatrix} 0 \\ v_0 \end{pmatrix}
\end{align*}}

です。 

動径成分 再び
さて、前節の初期条件を  { r^2 } の式に課すと  { \phi = \pm\frac{\pi}{2} } となります。 よって( { t = 0 } { r } が最小になる方を選んで)

  { \displaystyle\begin{align*}
  r^2
    &= -\frac{\sqrt{E^2 - \ell^2\omega^2}}{k}\;\cos 2\omega t + \frac{E}{k} \\
    &= -\frac{1}{2}\left(\frac{v_0^2}{\omega^2} - r_0^2\right) \cos 2\omega t
      + \frac{1}{2}\left(\frac{v_0^2}{\omega^2} + r_0^2\right) \\
    &= \frac{r_0^2}{2} (1 + \cos 2\omega t) + \frac{v_0^2}{2\omega^2}(1 - \cos 2\omega t) \\
    &= r_0^2 \cos^2 \omega t + \frac{v_0^2}{\omega^2} \sin^2\omega t \\[2mm]
  r(t) &= \sqrt{r_0^2 \cos^2 \omega t + \frac{v_0^2}{\omega^2} \sin^2\omega t}
\end{align*}}

を得ます。 これは直交座標の解に今の場合の初期条件を課して  { |\textbf{r}(t)| } を計算したものと一致します。

偏角成分
偏角成分は  { \theta } に対しての正準方程式  { \frac{d\theta}{dt} = \frac{\ell}{mr^2} } から求まります:

  { \displaystyle\begin{align*}
  \frac{d\theta}{dt}
    &= \frac{r_0v_0}{kr^2} \qquad\left(\because \omega = \sqrt{\frac{k}{m}},\,\ell = mr_0v_0\right) \\
  \frac{d\theta}{dt}
    &= \frac{r_0v_0}{r_0^2 \cos^2 \omega t + \frac{v_0^2}{\omega^2} \sin^2\omega t}
\end{align*}}

これを積分するのはちょっと大変なので、直交座標の解(に今の初期条件を課したもの)から解の形に見当をつけて、以下の関数を考えてみます( { A } は定数):

  { \displaystyle\begin{align*}
  f(t) = \tan^{-1}\left(A\tan \omega t\right)
\end{align*}}

両辺を  { t }微分すると

  { \displaystyle\begin{align*}
  \frac{df(t)}{dt}
    &= \frac{1}{1+A^2\tan^2 \omega t}\cdot A \cdot \frac{1}{\cos^2 \omega t} \cdot \omega \\
    &= \frac{A\omega}{\cos ^2 \omega t + A^2\sin^2 \omega t}
\end{align*}}

となるので、 { A = \frac{v_0}{r_0\omega} } とすれば今の場合にちょうど当てはまって

  { \displaystyle\begin{align*}
  \frac{d\theta}{dt}
    &= \frac{\frac{v_0}{r_0\omega}\cdot \omega}{\cos^2 \omega t + \left(\frac{v_0}{r_0\omega}\right)^2 \sin^2\omega t} \\
  \theta(t)
    &= \tan^{-1}\left(\frac{v_0}{r_0\omega}\tan \omega t\right) \\
  \therefore \tan \theta(t)
    &= \frac{v_0}{r_0\omega}\tan \omega t
\end{align*}}

を得ます。 やはりこれも直交座標の解に今の場合の初期条件を課して  { \frac{y(t)}{x(t)} } を計算したものと一致します。

2次元の調和振動子極座標で素直に解こうとすると思いの外大変なのが分かりました。 予定では同じことを  { \frac{1}{r} } に比例するポテンシャルについて行って惑星の運動を導くつもりですが、果たして解けるのかどうか。

惑星の運動などを導く際には、ここで行った方法ではなく、角運動量以外の保存量を使って積分する方法が一般的に使われているので、そちらをやるかもしれません。 ただ、いきなりその方法をやると有難味が分からない気がするので今回は頑張って泥臭い積分を計算してみました。

【追記】

  • 軌道の方程式を追記しました。
  • 別の初期条件の場合を追記しました。
  • 初期条件を通常の中心力ポテンシャルの場合に合わせました。

古典力学〈上〉 (物理学叢書)

古典力学〈上〉 (物理学叢書)

  • 作者: ゴールドスタイン,サーフコ,ポール,Herbert Goldstein,John Safko,Charles Poole,矢野忠,渕崎員弘,江沢康生
  • 出版社/メーカー: 吉岡書店
  • 発売日: 2006/07
  • メディア: 単行本
  • クリック: 7回
  • この商品を含むブログ (8件) を見る
古典力学 下 (物理学叢書 105)

古典力学 下 (物理学叢書 105)

*1:成分ごとに計算した方が簡単ですが。