単振動

水平な台の上に置かれた一端を壁に固定されたバネを考えます。壁につながれていないばねの端には質量ｍの物体が付いています。ここでは摩擦力は無視します。また、物体は方向を変えずに直線状を伸び縮みすると考えます。

物体の静止位置をつり合いの位置と呼び、そのときのバネの長さを自然長といいます。バネの伸びる向きをｘの正の方向とします。

物体をつり合いの位置からxだけ動かすと、バネはそれをもとの位置(つり合いの位置)に戻そうとして
F=-kx　　　---(1)
の力を及ぼします。これは物体をもとの位置に戻そうとする復元力です。バネがこのような力を及ぼすことをフックの法則といいます。kはばね定数と呼ばれ、復元力を表す定数です。

運動(ニュートン)の第２法則により力Ｆが働くと質点は加速度aをもつので
F=ma ---(2)
という比例関係が成り立ちます。物体の動きと向きを組み合わせたベクトルを速度といい、時刻ｔの瞬間の速度は
$\frac{\text{d}x}{\text{d}t}=v$ ---(3)
で表されます。物体の加速度はこの速度が時間と共に変化していく量を表します。加速度を速度で表すと
$a=\frac{\text{d}v}{\text{d}t}$
となります。この式に(3)を代入すると
$a=\frac{\text{d}^2v}{\text{d}t^2}$ ---(4)
となります。
(4)と(2)を用いて(1)を書き直すと運動方程式が得られます。
$m\frac{\text{d}^2x}{\text{d}t^2}=-kx$ ---(5)
物体はこの運動方程式にしたがい、最大のxと最小のxの間を行ったり来たりします。この最大値と最小値を＋Rと－Rで表し、Rを振幅と呼びます。

つぎにこの運動と円運動との関係を調べます。

質点Pは半径R上の円周上を動いています。そのx座標は
x=Rcos$\theta$
であらわされます。$\theta$はPOAの角度を表します。
t=$\theta$のとき、Pが$\theta$=Φの位置にあり、その後Θが単位時間当たりωという一定の割合で変化するなら
$\theta$=ωt+Φ
なので、
x=Rcos(ωt+Φ)---(6)
となります。このωは角速度と呼ばれます。これを時間で２階微分してみましょう。
$\frac{\text{d}^2x}{\text{d}t^2}=-ω^2Rcos(ωt+Φ)=-ω^2x$ ---(7)
となります。これで(5)式の一般解が(6)式であることが分かります。

(5)と(7)を比べると
$ω^2=\frac{k}{m}$
となります。ωをバネの固有振動数と呼びます。それはωがこのバネ固有の振動数だからです。これで単振動の性質が得られました。

連成振動

つぎにこの固有の意味をさらに詳しく見るために連成振動を見ます。

質量mの物体が２つあり、その左右にバネがつながれています。左の物体の左のバネの左端は壁につながれていて、右の物体の右のバネの右端も壁につながれています。物体は方向を変えずに直線状を動き、バネは方向を変えずに直線状を伸び縮みすると考えます。摩擦力は無視します。左の物体のつり合いの位置からのずれをxであらわし、右の物体のつり合いの位置からのずれをyとします。真ん中のバネの伸びはy-xとなります。左右のバネのばね定数はk、真ん中のバネのばね定数はk`とします。
この連成振動の運動方程式は

$m\frac{\text{d}^2x}{\text{d}t^2}=-kx+k'(y-x)$
$\ \ \ \ \ \ =(-k-k')x+k'y$ ---(8)
$m\frac{\text{d}^2y}{\text{d}t^2}=-k'(y-x)-ky$
$ \ \ \ \ \ \ =k'x+(-k-k')y$ ---(8`)

これらの式に一般解の(6)式を振幅を$A,B$,
x=Acos(ωt+Φ)　--- (9),
y=Bcos(ωt+Φ) --- (9')
と仮定して代入します。(8),(8')も(5)の形をしていないためにあくまでも仮定です。
$-mω^2Acos(ωt+Φ)=(-k-k')Acos(ωt+Φ)+k'Bcos(ωt+Φ)$
$-mω^2Bcos(ωt+Φ)=k'Acos(ωt+Φ)+(-k-k')Bcos(ωt+Φ)$
両辺をcos(ωt+Φ)で割り、
$-mω^2A=(-k-k')A+k'B$ ---(10)
$-mω^2B=k'A+(-k-k')B$ ---(10')
左辺を右辺に移動してまとめると
$(-k-k'+mω^2)A+k'B=0$
$k'A+(-k-k'+mω^2)B=0$
この連立方程式は変数が$\omega$,A,Bと３つあるのに対して、式は２つしかありません。
この連立方程式の解は
$(-k-k'+mω^2)^2-k'^2)B=0$
この連立方程式が自明でない解をもつ、つまりB＝0でない解をもつ条件は
$(-k-k'+mω^2)^2-k'^2)=0$です。
これを解くと
$\omega=\sqrt{\frac{k}{m}}$　---(11)
または
$\omega=\sqrt{\frac{k+2k'}{m}}$ ---(11')

それぞれのωを連立方程式(10)(10')に入れると
$\omega=\sqrt{\frac{k}{m}}$のとき、
-kA=-(k+k')A+K'B は -k'A+k'B=0
-kB=-k'A+(k+K')Bは-k'A+k'B=0
となりどちらも同じ式となります。
これは２つの式が同じものであることを示していて、残る未知数がA,Bと２つあるのに対して、式は１つしかないことと同じことです。したがって、解を比率で表します。
したがってB/A=1

$\omega=\sqrt{\frac{k+2k'}{m}}$のとき、
-(k+2k')A=-(k+k')A+K'B は k'A+k'B=0
-(k+2k')B=-k'A+(k+K')Bはk'A+k'B=0
となりどちらも同じ式となります。したがってB/A=-1
ωは固有値、固有振動数です。A、Bは固有ベクトルですが、比率でしか求まりません。

ここで２つのベクトル(A,B)=(A,A)と(A,B)=(A,-A)の内積はゼロになり２つのベクトルが直交しているのが分かります。固有ベクトルは常に直交しています。

行列式と固有値、固有ベクトル

(10),(10')を行列表現に直して解くと

-m\omega^2
  \left(
    \begin{array}{cc}
      A \\
      B \\
    \end{array}
  \right)=
  \left(
    \begin{array}{cc}
      -k-k' &k' \\
      k'    & -k-k' \\
    \end{array}
  \right)
  \left(
    \begin{array}{cc}
      A \\
      B \\
    \end{array}
  \right)

となります。これは$\omega^2$を定数とすると一般の固有値と固有ベクトルの形をしています。左辺を右辺に移動すると

  \left(
    \begin{array}{cc}
      m\omega^2-k-k' &k' \\
      k'    &m\omega^2-k-k' \\
    \end{array}
  \right)
  \left(
    \begin{array}{cc}
      A \\
      B \\
    \end{array}
  \right)=0

これが解をもつためには左辺の左の行列の行列式がゼロであればよいのです。あとは上述の例と同じです。

したがって、この場合の固有値問題は３つの変数をもち２つの式からなる連立方程式が自明でない解をもつための条件を固有値を求めることからはじめて、固有振動数をもとめ、出来上がった連立方程式が自由解をもつことからAとBの比率として解を得たことになります。

基準座標の利用

つぎに(8)と(8)'を変形して(5)の形にすることを考えてみます。
まず(8)+(8)'を作ってみます。
$m\frac{\text{d}^2(x+y)}{\text{d}t^2}=-k(x+y)$
となります。X=x+yとおくと
$m\frac{\text{d}^2X}{\text{d}t^2}=-kX$
となります。これは(5)と同じ形をしています。この場合の一般解は
$X=A\cos( \theta+\phi)$
となり上式に代入すると
$m\frac{\text{d}^2X}{\text{d}t^2}=-mω^2Acos(ωt+Φ)=-mω^2X$
となり、
$\omega=\sqrt{\frac{k}{m}}$
となります。

また、(8)-(8)'を作ってみます。
$m\frac{\text{d}^2(x-y)}{\text{d}t^2}=-(k+2k')(x-y)$
$Y=x-y$とおくと
$m\frac{\text{d}^2Y}{\text{d}t^2}=-(k+2k')Y$
となり、
この場合の一般解は
$Y=A\cos( \theta+\phi)$
上式に代入すると
$m\frac{\text{d}^2Y}{\text{d}t^2}=-mω^2Acos(ωt+Φ)=-mω^2Y$
よって、$\omega=\sqrt{\frac{k+2k'}{m}}$
となります。
(11),(11')と同様のωが得られました。

例として、k=k'=10N/mとし、物体の質量を100gとするとそれぞれの固有振動数ω_1,ω2と振幅の比率 B1/A1, B2/A2は

from numpy import linalg as LA
import numpy as np
A=[[-200,100],
  [100,-200]]
w,v=LA.eig(A)
np.sqrt(-w),v[0][0]/v[1][0],v[0][1]/v[1][1],v

(array([10.        , 17.32050808]),
 1.0,
 -1.0,
 array([[ 0.70710678, -0.70710678],
        [ 0.70710678,  0.70710678]]))

となります。vを行列として示していることは固有ベクトルを特解として表示していることになります。

バネの動きから固有値・固有ベクトルを理解する。

単振動

単振動をもう一度検討します。(5)と(6)式から
$-kx=-\omega^2x$
を得ます。
一般には両辺に-1を掛けて、マイナス符号を落とし、$x$で割って、両辺から$x$を落とします。
しかし、ここでは
$[-k]x=\lambda x$
と書き換えます。[]は行列を表し、$x$は縦ベクトルを表すと想像します。
$[-k]x-\lambda x=0$
と変形し、$x$が自明でない解をもつためには$-k-\lambda$の行列式が零である必要があります。すると$\lambda=-k$となります。すると
$(-k+k)x=x$
となるので、$x$は自由変数になります。そうするとバネの位置はどこでも可能ということになってしまいます。しかし、(6)式($x=R\cos(ωt+Φ))$より
$-R\le x \le R$
となり最大値と最小値があります。

連成振動

連成振動の方が一般的な例です。(8),(8')を再度書くと
$m\frac{\text{d}^2x}{\text{d}t^2}=-kx+k'(y-x)=(-k-k')x+k'y$
$m\frac{\text{d}^2y}{\text{d}t^2}=-k'(y-x)-ky=k'x+(-k-k')y$
となります。
x=Acos(ωt+Φ)　--- (9),
y=Bcos(ωt+Φ) --- (9')
を２回時間で微分すると

$-m\omega^2 x=-kx+k'(y-x)=(-k-k')x+k'y$
$-m\omega^2y=-k'(y-x)-ky=k'x+(-k-k')y$
$-m\omega^2=\lambda$, $a=-k-k'$,$b=k'$として行列表現に直して解くと

\lambda
  \left(
    \begin{array}{cc}
      x \\
      y \\
    \end{array}
  \right)=
  \left(
    \begin{array}{cc}
      a &b \\
      b   & a \\
    \end{array}
  \right)
  \left(
    \begin{array}{cc}
      x \\
      y \\
    \end{array}
  \right)

書き換えると

  \left(
    \begin{array}{cc}
      a -\lambda &b \\
      b   & a -\lambda \\
    \end{array}
  \right)
  \left(
    \begin{array}{cc}
      x \\
      y \\
    \end{array}
  \right)=0

これが自明でない解をもつためには

  \left|
    \begin{array}{cc}
      a -\lambda &b \\
      b   & a -\lambda \\
    \end{array}
  \right|=0

$(a-\lambda)^2-b^2=(a+b-\lambda)(a-b-\lambda)=0$
よって
$\lambda=a+b$または$\lambda=a-b$
$\lambda=a+b$のとき

    \begin{array}{cc}
      (a -a-b)x +by=-bx+by=0 \\
      bx+ (a -a-b)y=bx-by=0 \\
    \end{array}

よって、$x=y$

$\lambda=a-b$のとき

    \begin{array}{cc}
      (a -a+b)x +by=bx+by=0 \\
      bx+ (a -a+b)y=bx+by=0 \\
    \end{array}

よって、$x=-y$
$y$は自由変数であるので、固有ベクトル$x,y$はどのような値を取ることができる。しかし、
x=Acos(ωt+Φ)　
y=Bcos(ωt+Φ)
より$-A\le x\le A$,$-B\le y\le B$となり、上限と下限がある。固有値$\lambda$は２つのバネの動きの向きを表している。$\lambda=a+b$のときには２つのバネの動きは同じ方向を向き、$\lambda=a-b$のときには逆の方向を向いている。

その他のいろいろな振動

天井につるしたバネ

運動方程式は
$m\frac{\text{d}^2x}{\text{d}t^2}=-kx$

天井につるしたおもり

運動方程式は
$ m\frac{\text{d}^2\phi}{\text{d}t^2}=-\frac{g}{l}\sin\phi=-\frac{g}{l}\phi $

天井につりしたバネでつながれた２つの振子

運動方程式は
$ m\frac{\text{d}^2x}{\text{d}t^2}=-m\frac{g}{l}x+ k(y-x)$
$ m\frac{\text{d}^2y}{\text{d}t^2}=-m\frac{g}{l}-k(y-x)$
となります。行列式で表すと

-\omega^2
  \left(
    \begin{array}{cc}
      A \\
      B \\
    \end{array}
  \right)=
  \left(
    \begin{array}{cc}
      -\frac{g}{l}-\frac{k}{m} &\frac{k}{m} \\
      \frac{k}{m}    & -\frac{g}{l}-\frac{k}{m} \\
    \end{array}
  \right)
  \left(
    \begin{array}{cc}
      A \\
      B \\
    \end{array}
  \right)

天井につるした２重振り子

$x_1=l_1\sin \theta_1$ ---(1)
$y_1=l_1\cos \theta_1$ ---(2)
$x_2=l_1\sin \theta_1+l_1\sin\theta_2$ ---(3)
$y_2=l_1\cos \theta_1+l_1\cos\theta_2$ ---(4)
となります。T_1,T_2を針金の張力とすると、運動方程式は
$ m_1\frac{\text{d}^2x_1}{\text{d}t^2}=-T_1\sin\theta_1+T_2\sin\theta_2$ ---(5)
$ m_1\frac{\text{d}^2y_1}{\text{d}t^2}=m_1g-T_1\cos\theta_1+T_2\cos\theta_2$ ---(6)
$ m_2\frac{\text{d}^2x_2}{\text{d}t^2}=-T_2\sin\theta_2$ ---(7)
$ m_2\frac{\text{d}^2y_2}{\text{d}t^2}=m_2g-T_2\cos\theta_2$ ---(8)

(1)-(4)は$\theta_1$, $\theta_2$が微小であると
$x_1\simeq l_1\theta_1$
$y_1\simeq l_1$
$x_2\simeq l_1\theta_1+l_1\theta_2$
$y_2\simeq -l_1-l_2$
これを２階微分すると
$ m_1\frac{\text{d}^2x_1}{\text{d}t^2}=l_1\frac{\text{d}^2\theta_1}{\text{d}t^2}$
$ m_1\frac{\text{d}^2y_1}{\text{d}t^2}=0$
$ m_2\frac{\text{d}^2x_2}{\text{d}t^2}=l_1\frac{\text{d}^2\theta_1}{\text{d}t^2}+l_2\frac{\text{d}^2\theta_2}{\text{d}t^2}$
$ m_2\frac{\text{d}^2y_2}{\text{d}t^2}=0$
これらを(5)-(8)に代入すると
$ m_1l_1\frac{\text{d}^2\theta_1}{\text{d}t^2}=-T_1\sin\theta_1+T_2\sin\theta_2$
$ 0=m_1g-T_1\cos\theta_1+T_2\cos\theta_2$
$ l_1\frac{\text{d}^2\theta_1}{\text{d}t^2}+l_2\frac{\text{d}^2\theta_2}{\text{d}t^2}=-T_2\sin\theta_2$
$ 0=m_2g-T_2\cos\theta_2$
これらをさらに$\sin \theta \approx \theta$, $\cos \theta\approx 1$とすると
$ m_1l_1\frac{\text{d}^2\theta_1}{\text{d}t^2}=-T_1\theta_1+T_2\theta_2$
$ 0=m_1g-T_1+T_2$
$ l_1\frac{\text{d}^2\theta_1}{\text{d}t^2}+l_2\frac{\text{d}^2\theta_2}{\text{d}t^2}=-T_2\theta_2$
$ 0=m_2g-T_2$
これらを整理すると
$ \frac{\text{d}^2\theta_1}{\text{d}t^2}=-\frac{(m_1+m_2)}{m_1l_1}g\theta_1+ \frac{m_2}{m_1l_1}g\theta_2$
$ \frac{\text{d}^2\theta_2}{\text{d}t^2}=\frac{(m_1+m_2)}{m_1l_2}g\theta_1-\frac{(m_1+m_2)}{m_1l_2}g\theta_2$

行列式で表すと

-\omega^2
  \left(
    \begin{array}{cc}
      \theta_1 \\
      \theta_2 \\
    \end{array}
  \right)=
  \left(
    \begin{array}{cc}
      -\frac{m_1+m_2}{m_1l_1}g &\frac{m_2}{m_1l_1}g \\
      \frac{m_1+m_2}{m_1l_2}g    & -\frac{m_1+m_2}{m_1l_2}g \\
    \end{array}
  \right)
  \left(
    \begin{array}{cc}
      \theta_1 \\
      \theta_2 \\
    \end{array}
  \right)

"天井につるしたバネ","天井につるしたおもり","天井につりしたバネでつながれた２つの振子","天井につるした２重振り子"は同じ形をしているのが分かります。

参考

小形正男　振動・波動　裳華房
【力学】連成振動(n=2)
【力学】連成振動(n=3)

ーーーーーーーーーーーーーーーーーーーーーーーー
Python3ではじめるシステムトレード【第2版】環境構築と売買戦略

「画像をクリックしていただくとpanrollingのホームページから書籍を購入していただけます。

Python3ではじめるシステムトレード改訂版　アマゾンから購入いただけます。