メカトロニックなカメ

メカトロニクス技術者になりたいカメです

1自由度振動系のお話 その3

1自由度振動系のお話 その1 - メカトロニックなカメ(無次元化について)

1自由度振動系のお話 その2 - メカトロニックなカメ(周波数応答・過渡応答について)

の続き

復習

質量を\(m\)[kg]、ばね定数を\(k\)[N/m]、減衰係数を\(c\)[Ns/m]、位置を\(x(t)\)[m]、外力を\(f(t)\)[N]とすると、一自由度マスばねダンパ系の運動方程式は下記のようになります。

\begin{align} m\cfrac{d^2x}{dt^2}+c\cfrac{dx}{dt}+kx(t)=f(t) \end{align}

 上式を無次元化した運動方程式は次式のようになります。

\begin{align} \cfrac{d^2\tilde{x}}{d\tilde{t}^2}+2\zeta\cfrac{d\tilde{x}}{d\tilde{t}}+\tilde{x}(\tilde{t})=\tilde{f}(\tilde{t}) \end{align} 

この系の過渡応答には次式のようになります。

 \begin{align}\tilde{x}(\tilde{t})=C_1e^{-(\zeta+\sqrt{\zeta^2-1})\tilde{t}}+C_2e^{-(\zeta-\sqrt{\zeta^2-1})\tilde{t}}\end{align}

\(\zeta<1\)の時、この系の過渡応答には次式のように振動的な応答が発生します。

 \begin{align}\tilde{x}(\tilde{t})=X'e^{-\zeta\tilde{t}}\sin\left(\sqrt{1-\zeta^2}\tilde{t}+\Theta'\right)\end{align}

この過渡応答における振動の周波数比\(\alpha_{damp}\)は無次元化減衰固有周波数と呼ばれ次式で表される。

  \begin{align}\alpha_{damp}=\sqrt{1-\zeta^2}\end{align}

状態方程式(State equation)

制御工学において、高次の微分方程式を1階の微分方程式の形で表した式を状態方程式と呼ばれ、非常に基礎的な考え方である。また、非線形な挙動を含む複雑なシステムは一般に厳密な解を求めることができないため、ソフトウェアを用いたシミュレーションの分野では、複雑なシステムを1階の微分方程式に変換する。そしてルンゲ=クッタ法等の数値積分手法を用いて、シミュレーションを行う。というように1階の微分方程式に変換することは非常に重要であるので、その手順を簡単に示す。

まずは、一自由度マスばねダンパ系の2階の運動方程式を1階の微分方程式に変換するために、次のような新しい変数\(\tilde{v}(\tilde{t})\)を導入します。

 \begin{align}\tilde{v}(\tilde{t})&=\cfrac{d\tilde{x}}{d\tilde{t}}\end{align}

上式を見てわかる通り、\(\tilde{v}(\tilde{t})\)は単なる無次元化速度である。この無次元化速度を用いると、運動方程式は次式のようになる。

\begin{align} \cfrac{d\tilde{v}}{d\tilde{t}}+2\zeta\tilde{v}(\tilde{t})+\tilde{x}(\tilde{t})=\tilde{f}(\tilde{t}) \end{align} 

これで1階の微分方程式になりました。しかしこの方程式には2つの変数\(\tilde{v}, \tilde{x}\)があり、多変数な微分方程式となってしまい扱いが困難になります。そこで、また新しい変数のベクトル\(\tilde{\boldsymbol{w}}\)を次式のように用意します。

 \begin{align}\tilde{\boldsymbol{w}}(\tilde{t})&=\left(\begin{array}{c} \tilde{x}(\tilde{t}) \\\tilde{v}(\tilde{t})\end{array}\right)\end{align}

この新しい変数ベクトル(状態変数ともいう)を用いると、次式のような状態方程式となります。

\begin{align} \cfrac{d\tilde{\boldsymbol{w}}}{d\tilde{t}}=\left[\begin{array}{cc} 0&1\\-1&-2\zeta\end{array}\right]\tilde{\boldsymbol{w}}(\tilde{t})+\left(\begin{array}{c} 0 \\\tilde{f}(\tilde{t})\end{array}\right)\end{align}

 さて、この状態方程式の行列部の固有値を求めてみましょう。固有値の求め方は、次式を満たすような\(\lambda\)を求めることです。

\begin{align} \left[\begin{array}{cc} 0&1\\-1&-2\zeta\end{array}\right]\boldsymbol{z}=\lambda\boldsymbol{z}\end{align}

求められた固有値\(\lambda_1,\lambda_2\)は次式のようになる。

\begin{align} \lambda_{1,2}=-\zeta\pm\sqrt{\zeta^2-1}\end{align}

この値は、まさに過渡応答における指数部に相当します。

対角化(Diagonalization)

固有値には非常に重要な性質を持つことが分かった。そこで対角化して固有値を抜き出す作業を行ってみる。まずは固有値に対応した固有ベクトルを求めてよう。固有ベクトルは次式を満たす\(\boldsymbol{z}_0\)である。

\begin{align} \left[\begin{array}{cc} 0&1\\-1&-2\zeta\end{array}\right]\boldsymbol{z}_0=\lambda_{1,2}\boldsymbol{z}_0\end{align}

\begin{align} \boldsymbol{z}_{1,2}=\left(\begin{array}{c} 1\\\lambda_{1,2}\end{array}\right)=\left(\begin{array}{c} 1\\-\zeta\pm\sqrt{\zeta^2-1}\end{array}\right)\end{align}

 状態方程式の対角化のために、固有ベクトルを横に並べた変換行列\(\boldsymbol{A}\)を導入します。

\begin{align} \boldsymbol{A}&=\left[\begin{array}{cc} 1&1\\\lambda_1&\lambda_2\end{array}\right]=\left[\begin{array}{cc} 1&1\\-\zeta+\sqrt{\zeta^2-1}&-\zeta-\sqrt{\zeta^2-1}\end{array}\right]\\\boldsymbol{A}^{-1}&=-\cfrac{1}{\sqrt{\zeta^2-1}}\left[\begin{array}{cc} \lambda_2&-1\\-\lambda_1&1\end{array}\right]\end{align}

対角化のために状態変数を次のように変換する。

\begin{align} \tilde{\boldsymbol{w}}(\tilde{t})=\boldsymbol{A}\tilde{\boldsymbol{w}}_d(\tilde{t})\end{align}

この状態変数を用いると状態方程式は次のようになります。

\begin{align} \cfrac{d\tilde{\boldsymbol{w}}_d}{d\tilde{t}}&=\boldsymbol{A}^{-1}\left[\begin{array}{cc} 0&1\\-1&-2\zeta\end{array}\right]\boldsymbol{A}\tilde{\boldsymbol{w}}_d(\tilde{t})+\boldsymbol{A}^{-1}\left(\begin{array}{c} 0 \\\tilde{f}(\tilde{t})\end{array}\right)\\&=\left[\begin{array}{cc} \lambda_1&0\\0&\lambda_2\end{array}\right]\tilde{\boldsymbol{w}}_d(\tilde{t})-\cfrac{1}{\sqrt{\zeta^2-1}}\left(\begin{array}{c} -1 \\1\end{array}\right)\tilde{f}(\tilde{t})\end{align}

上式は固有値でまとめることができて、非常にすっきりした形となったといえる。しかし、減衰比\(\zeta\)が1より小さい場合、各係数に虚数複素数が出てきてしまい、扱いが非常に複雑になる。そこで対角化ではなくモード形式への変換を行ってみよう。

モード形式(Mode form)

減衰比が1未満\(\zeta<1\)の時、固有値は次のように表すこともできる。

\begin{align} \lambda_{1,2}=-\zeta\pm i \alpha_{damp}\end{align}

次のような変換行列\(\boldsymbol{B}\)を導入します。 

\begin{align} \boldsymbol{B}&=\left[\begin{array}{cc} \textrm{Re}(\boldsymbol{z}_{1,2})&\textrm{Im}(\boldsymbol{z}_{1,2})\end{array}\right]\\&=\left[\begin{array}{cc} 1&0\\-\zeta&\alpha_{damp}\end{array}\right]\\\boldsymbol{B}^{-1}&=\cfrac{1}{\alpha_{damp}}\left[\begin{array}{cc} \alpha_{damp}&0\\\zeta&1\end{array}\right]\end{align}

先ほどと同様の手順で状態変数を次のように変換する。

\begin{align} \tilde{\boldsymbol{w}}(\tilde{t})=\boldsymbol{B}\tilde{\boldsymbol{w}}_m(\tilde{t})\end{align}

この状態変数を用いると状態方程式は次のようになります。

\begin{align} \cfrac{d\tilde{\boldsymbol{w}}_m}{d\tilde{t}}&=\boldsymbol{B}^{-1}\left[\begin{array}{cc} 0&1\\-1&-2\zeta\end{array}\right]\boldsymbol{B}\tilde{\boldsymbol{w}}_m(\tilde{t})+\boldsymbol{B}^{-1}\left(\begin{array}{c} 0 \\\tilde{f}(\tilde{t})\end{array}\right)\\&=\left[\begin{array}{cc} -\zeta&\alpha_{damp}\\-\alpha_{damp}&-\zeta\end{array}\right]\tilde{\boldsymbol{w}}_m(\tilde{t})+\cfrac{1}{\alpha_{damp}}\left(\begin{array}{c} 0 \\1\end{array}\right)\tilde{f}(\tilde{t})\end{align}

上式のように1自由度振動系で重要な性質を表す減衰比\(\zeta\)と無次元化減衰固有周波数\(\alpha_{damp}\)を抜き出すことに成功しました。

まとめ

以上のように、1自由度振動系を表す2階の微分方程式を1階の微分方程式に置き換える方法にはいくつかある。最後にこれらの特徴をまとめると、

  • 変換なし:元の運動方程式から求めやすい。\(\tilde{\boldsymbol{w}}(\tilde{t})\)に物理的意味(位置、速度)があるため、わかりやすい。1階の微分方程式の形が直感にそぐわない
  • 対角化:固有値を用いた対角行列が得られてすっきりする。\(\tilde{\boldsymbol{w}}_d(\tilde{t})\)に物理的意味がない。減衰比が1未満の時は虚数が現れる。
  • モード形式:減衰比が1未満の時に、減衰比と無次元化減衰固有周波数に分けて表記できて、見通しが良い。\(\tilde{\boldsymbol{w}}_m(\tilde{t})\)に物理的意味がない。