ばねの運動を題材に、減衰振動の微分方程式の解法を解説します。
減衰振動に考える前に、まずは、減衰や強制振動の無いばねの運動を考えます。
このような運動を単振動と呼びます。
このとき、ばねの運動方程式（＝微分方程式）は、以下のように表せます。
\begin{eqnarray}
m\ff{\diff^2 x(t)}{\diff t^2} &=& – kx(t)
\end{eqnarray}
ただし、時刻$t$での質点の位置を$x(t)$、ばね定数を$k$、質点の質量を$m$とします。
この微分方程式を解いていきましょう。特性方程式を考えると、以下のように$\lambda$を求められます。
\begin{eqnarray}
\lambda^2 &=& -\, \ff{k}{m} \EE
\therefore\, \lambda &=& \pm i\sqrt{\ff{k}{m}}
\end{eqnarray}
$\lambda$が虚数となりましたが、オイラーの公式を利用すると、$x(t)$は次のようにできます。
\begin{eqnarray}
x(t) &=& C_1 \cos\sqrt{\ff{k}{m}} t + iC_2 \sin \sqrt{\ff{k}{m}} t
\end{eqnarray}
物理の解として、虚数が含まれるのは具合が悪いため、$C_2=-iC’_2$として、三角関数の合成を行うと、
\begin{eqnarray}
x(t) &=& C’_1 \sin\left( \sqrt{\ff{k}{m}} t + \alpha \right)
\end{eqnarray}
となります。
このとき、$\displaystyle\sqrt{\ff{k}{m}}$を角振動数と見なせます。そこで、$\omega_n = \displaystyle\sqrt{\ff{k}{m}}$と表すことにしましょう。
すると、ばねの変位は
\begin{eqnarray}
x(t) &=& A \sin\left( \omega_n t + \alpha \right)
\end{eqnarray}
と表せます。（$A, \alpha$は定数です）
また、$\omega_n$は、系により決まる固有の定数のため、固有角振動数と呼ばれます。ばねの変位が振り子の角度変化と同じように表せることが分かります。
違う現象でも、解が同じ様に表せる理由は、微分方程式が同じ形であるためです。
ばねの振動に戻りましょう。ばねの振動は実際には空気抵抗などにより、徐々に振動が小さくなっていきます。このように、振動が減衰していく振動のことを減衰振動と呼びます。
減衰振動
ばねの振動がいつまでも続くと仮定して問題を解きましたが、現実には空気抵抗やダンパー等により、エネルギーが失われていきます。
この場合、ばねの振動は徐々に小さくなります。
このようなばねの振動を、減衰振動と呼びます。
減衰振動の例として、速度に比例した抵抗力を受ける場合を考えます。このとき、減衰振動の運動方程式は次のように表せます。
\begin{eqnarray}
m\ff{\diff^2 x}{\diff t^2} &=& -\, kx -\, c \ff{\diff x}{\diff t}
\end{eqnarray}
特性方程式から、$\lambda$を計算できます。
\begin{split}
m\lambda^2 &+ c\lambda + k = 0 \EE
\therefore \, \lambda &= \ff{-c\pm\sqrt{c^2 \,- 4mk}}{2m}
\end{split}
この微分方程式の解は、平方根の正負により三通りに分類できます。
① $c^2 \,- 4mk > 0$
② $c^2 \,- 4mk = 0$
③ $c^2 \,- 4mk < 0$
過減衰：$c^2 \,- 4mk > 0$ の場合
$c^2 \,- 4mk \geq 0$のとき、$\lambda$は実数となります。
このとき、
\begin{split}
\lambda_1 &= \ff{-c-\sqrt{c^2 \,- 4mk}}{2m} \EE
\lambda_2 &= \ff{-c+\sqrt{c^2 \,- 4mk}}{2m} \EE
\end{split}
とすると、$x$は、
\begin{split}
x(t) &= C_1e^{\lambda_1 t} + C_2e^{\lambda_2 t}
\end{split}
と表せます。
$C_1=C_2=1$、$c=2, m=k=1$としてグラフの形を描画すると、以下のようになります。
振動せず、変位がすぐに$0$になることが分かります。このような現象を過減衰と呼びます。
※無周期振動とも呼びます。
臨界減衰：$c^2 \,- 4mk = 0$ の場合
$c^2 \,- 4mk = 0$のとき、すなわち、$c=2\sqrt{mk}$ となる場合を考えます。
このとき、$\lambda = \displaystyle{-\ff{c}{2m}}$となります。したがって$x$は、
\begin{split}
x(t) &= Ce^{-\ff{c}{2m}}
\end{split}
と求められます。初期変位が先ほどの過減衰のグラフと重なるように設定して、グラフを描画すると、以下のような赤線になります。
このときの運動は、後ほど紹介する振動状態との境になるため、臨界減衰と呼びます。このときの$c$を$c_c$で表し、臨界減衰係数と呼びます。$c_c=2\sqrt{mk}$ となります。
また、$\omega_n = \displaystyle{\sqrt{\ff{k}{m}}}$とすると、$c_c=2m\omega_n$とできます。
減衰振動：$c^2 \,- 4mk < 0$ の場合
$c^2 \,- 4mk \leq 0$のとき、$\lambda$は複素数となります。
このとき、
\begin{split}
\lambda_1 &= \ff{-c-i\sqrt{4mk \,- c^2}}{2m} \EE
\lambda_2 &= \ff{-c+i\sqrt{4mk \,- c^2}}{2m} \EE
\end{split}
とできます。
したがって、$x(t)$は、
\begin{eqnarray}
x(t) = e^{-\ff{c}{2m}t}\left( C_1e^{-\ff{i\sqrt{4mk \,- c^2}}{2m}t} + C_2e^{\ff{i\sqrt{4mk \,- c^2}}{2m}t} \right) \tag{1}
\end{eqnarray}
と求められます。
オイラーの公式を利用し、式(1)を変形すると、
\begin{eqnarray}
x(t) = Ae^{-\ff{c}{2m}t}\sin\left( \ff{\sqrt{4mk \,- c^2}}{2m}\, t + \alpha \right) \tag{2}
\end{eqnarray}
となります。
このような振動を減衰振動と呼びます。（$A, \alpha$は初期変位に依存する定数です。）
\begin{eqnarray}
x(t) = Ae^{-\ff{c}{2m}t}\sin\left( \ff{\sqrt{4mk \,- c^2}}{2m}\, t + \alpha \right) \\
\,
\end{eqnarray}
減衰振動のグラフを描画すると、以下のようになります。
三角関数で表される振動が指数関数の影響を受けて、減少していく様子が確認できます。
減衰比・対数減衰率
減衰振動の微分方程式に関する解説は前節で以上になりますが、機械力学や制御工学で題材として良く取り上げられます。
今後の準備のために、いくつかの用語について解説します。
手始めに、定数を次のように定め、式(2)を簡略化しましょう。
さて、$\omega_n = \displaystyle{\sqrt{\ff{k}{m}}}$, $\zeta=\displaystyle{\ff{c}{2\sqrt{mk}}}$として式(2)を変形すると、
\begin{eqnarray}
x(t) &=& Ae^{-\ff{c}{2m}t}\sin\left( \ff{\sqrt{4mk \,- c^2}}{2m}\, t + \alpha \right) \EE
&=& Ae^{-\ff{c}{2m}t}\sin\left( \sqrt{\ff{k}{m} \,- \left(\ff{c}{2m}\right)^2} \, t + \alpha \right) \EE
&=& Ae^{-\ff{c}{2\sqrt{mk}}\cdot\sqrt{\ff{k}{m}}t}\sin\left( \sqrt{\ff{k}{m}}\sqrt{1 \,- \ff{c^2}{4mk}} \, t + \alpha \right) \EE
\therefore \, x(t) &=& Ae^{-\zeta\omega_n t}\sin\left( \sqrt{1 \,- \zeta^2} \,\omega_n\, t + \alpha \right) \tag{3} \EE
\end{eqnarray}
式(3)から、減衰振動の固有角振動数$\omega_d$は、$\sqrt{1 \,- \zeta^2} \omega_n$となることが分かります。
減衰の無い振動と比較して、角振動数が小さくなることが分かります。（$\sqrt{1 \,- \zeta^2}< 1$のため）
減衰比
さて、いきなり$\zeta$（ゼータ）なる定数を導入しましたが、結果として式(3)のように簡単にできました。
そこで、$\zeta$を減衰比と呼びましょう。
\begin{eqnarray}
\zeta = \ff{c}{2\sqrt{mk}} \\
\,
\end{eqnarray}
下の図から分かる通り、振動の振幅が$\zeta$を指数に持つ指数関数に影響を受けながら、減少していくことが分かります。
振動状態における振幅の減衰の程度を表すため、$\zeta$を減衰比と呼ぶのです。
また、$\zeta$の正負により、次のように振動の状態を分類できます。
① $\zeta > 0$ のとき、過減衰
② $\zeta = 0$ のとき、臨界減衰
③ $\zeta < 0$ のとき、減衰振動
対数減衰率
減衰比が分かっていない場合、測定データから減衰比を求める必要があります。
その際に使われるのが、対数減衰率です。対数減衰率$\delta$ は次のように表せます。
\begin{eqnarray}
\delta = \ln \ff{x_p}{x_{p+1}} \\
\,
\end{eqnarray}
$x_p, x_{p+1}$は$p, p+1$番目の振幅の変位の極大値を表します。
このとき、各変位は、以下のようになります。
\begin{eqnarray}
x_p &=& Ae^{-\zeta\omega_n t_p}\sin\left( \sqrt{1 \,- \zeta^2} \,\omega_n\, t_p + \alpha \right) \EE
x_{p+1} &=& Ae^{-\zeta\omega_n t_{p+1}}\sin\left( \sqrt{1 \,- \zeta^2} \,\omega_n\, t_{p+1} + \alpha \right) \EE
\end{eqnarray}
$x_p, x_{p+1}$は極大値のため、$\sin$の角度の中身は$\displaystyle{\ff{\pi}{2}}$と近似できます。２式の商を考えると、以下のようになります。
\begin{eqnarray}
\ff{x_p}{x_{p+1}} &=& e^{\zeta\omega_n(t_{p+1} \,- t_p)}
\end{eqnarray}
さらに、$T’ = t_{p+1} \,- t_p$とすると、$T’ = \displaystyle{\ff{2\pi}{\omega_d}}$となるので、
\begin{eqnarray}
\ff{x_p}{x_{p+1}} &=& e^{\zeta\omega_n T’} = e^{\zeta\omega_n \ff{2\pi}{\omega_d}}
\end{eqnarray}
となります。両辺の対数をとると、
\begin{eqnarray}
\delta &=& \ln\ff{x_p}{x_{p+1}} \EE
&=& \zeta\omega_n \ff{2\pi}{\omega_d} \EE
&=& \ff{2\pi\zeta}{\sqrt{1-\zeta^2}}
\end{eqnarray}
なり、これを$\zeta$に関して整理すると、
\begin{eqnarray}
\zeta &=& \sqrt{\ff{\delta^2}{\delta^2 + (2\pi)^2}}
\end{eqnarray}
と求められます。
\begin{eqnarray}
\zeta &=& \sqrt{\ff{\delta^2}{\delta^2 + (2\pi)^2}} \\
\,
\end{eqnarray}
強制振動
工学的な応用を考える場合、外部から力を受けて振動する強制振動が興味の対象となります。
外力を$f(t)$として強制振動の運動方程式は次のように表せます。
\begin{eqnarray}
m\ff{\diff^2 x}{\diff t^2} + c \ff{\diff x}{\diff t} + kx &=& f(t)
\end{eqnarray}
外力$f(t)$が周期的に与えられる場合を考えしょう。
簡単のため、外力を角振動数$\omega$の正弦関数でモデル化しましょう。つまり、以下のような微分方程式を考えていきます。
\begin{eqnarray}
m\ff{\diff^2 x}{\diff t^2} + c \ff{\diff x}{\diff t} + kx &=& f\sin \omega t \\
\,
\end{eqnarray}
それでは、特性方程式に代入してみましょう。すると、
\begin{split}
m\lambda^2 &+ c\lambda + k = f\sin \omega t
\end{split}
となりますが、この二次方程式を解いても正しい答えを導けそうにありません。このように困ってしまった原因は、解きたい微分方程式の右辺が$0$でないためです。
※正確に言えば、強制振動の微分方程式が非斉次微分方程式であるためです。
微分方程式が非斉次微分方程式である場合、特性方程式を解いても答えを導くことはできません。
非斉次微分方程式は今までの微分方程式と比べて、難易度が上がります。微分方程式の入門は終わり、本格的に微分方程式と格闘していくことになります。
強制振動は共鳴現象や危険速度と関連して、工学的に興味深い対象です。強制振動の微分方程式の解法については、別の機会に考えることにしましょう。
