7 応用

行列乗は、定数係数線形微分方程式への応用が有名である。 その基礎から少し紹介する。 定数係数の線形微分方程式には、例えば
  $\displaystyle
y''+4y=\cos t\hspace{1zw}(\mbox{$y=y(t)$: 未知関数}) $ (16)
のように単独の方程式 (2 階単独) や、
  $\displaystyle
\left\{\begin{array}{ll}
y_1'' &= 2y_2'+y_1+2y_2\\
y_2'' &= 3y_1'+y_1-y_2+e^{2t}
\end{array}\right. \hspace{1zw}(\mbox{$y_1(t),y_2(t)$: 未知}) $ (17)
のように連立 (2 階 2 未知連立) の方程式もあるが、 未知関数の数を増やせば、すべて 1 階の連立微分方程式
  $\displaystyle
\begin{array}{l}
\left\{\begin{array}{ll}
y_1' &= a_{11} y_1+a...
... [.5zh]
\hspace{1zw}(\mbox{$y_j=y_j(t)$: 未知、$a_{ij},f_j(t)$: 既知})
\end{array}$ (18)
に変形できるし、 逆に未知関数を減らすことで 1 未知関数 1 本の $n$ 階 微分方程式
  $\displaystyle
\begin{array}{l}
y^{(n)}=b_{n-1}y^{(n-1)}+b_{n-2}y^{(n-2)}+\cdo...
...1y'+b_0y+g(t)\\
\hspace{1zw}(\mbox{$y=y(t)$: 未知、$b_j,g(t)$: 既知})
\end{array}$ (19)
の形に帰着させることもできる。

例えば、(16) は、$y'=z$ ($z=z(t)$: 未知) とすれば、 $z'=y''$ より

  $\displaystyle
\left\{\begin{array}{ll}
y' &=z\\
z' &=-4y+\cos t
\end{array}\right.$ (20)
の 1 階連立になるし、(17) は $y_1'=y_3$, $y_2'=y_4$ と すれば
  $\displaystyle
\left\{\begin{array}{ll}
y_1' & =y_3\\
y_2' & =y_4\\
y_3' & =y_1+2y_2 +2y_4\\
y_4' & =y_1-y_2+3y_3 + e^{2t}
\end{array}\right.$ (21)
の 1 階連立となる。 逆に、(17) から $y_2$ (または $y_1$) を 消去して $y_1$ の 4 階微分方程式を作ることもできる。 それには、演算子法を利用すると便利だが、その詳しい説明は 微分方程式の本を参照してもらうことにして、 その計算だけを紹介する。 (17) は微分演算子 $D=d/dt$ を用いて
$\displaystyle \left\{\begin{array}{ll}
(D^2-1)y_1-2(D+1)y_2 &=0\\
-(3D+1)y_1+(D^2+1)y_2 &=e^{2t}\end{array}\right.$
と書けるので、$y_2$ を消去すれば、
$\displaystyle -2(D+1)(3D+1)y_1+(D^2+1)(D^2-1)y_1=2(D+1)e^{2t}=6e^{2t}
$
となるので、
\begin{eqnarray*}\lefteqn{-2(D+1)(3D+1)+(D^2+1)(D^2-1)
\ =\ -2(3D^2+4D+1)+D^4-1}
\\ &=& D^4-6D^2-8D-3\end{eqnarray*}
より、
  $\displaystyle
y_1^{(4)}-6y_1''-8y_1'-3y_1=6e^{2t}$ (22)
が得られる。

よって、1 階連立方程式 (18) を 解くことと $n$ 階単独方程式 (19) を解くことは 実質的に同じで、 微分方程式の本では通常どちらかの解法のみを紹介し、 他方の方程式はそちらに変換する方法を紹介する、というのが普通である。

一般的には、より解法を重視する工学系の教科書では 単独方程式 (19) の方で解法を紹介し、 理論を重視する理学系の教科書では 連立方程式 (18) の方で解法を紹介することが 多いように思う。 それは、連立方程式 (18) の解法には行列の理論、 特に $e^{At}$ の話が必要になるからだと思われる。 単独方程式 (19) の場合は、 行列を使用せずに特性方程式 (連立方程式 (18) の 場合の固有方程式に対応する) を解くことで、代入法や定数変化法、演算子法、 ラプラス変換などにより行列を用いずに解を求めることができる。

なお、(18) や (19) の 係数 $a_{ij},b_j$ が定数でない場合、すなわち $t$ の既知関数である場合は、 方程式の解法は特殊な場合にしか存在せず、問題はかなり難しい。 本節では定数係数の場合のみ扱う。

方程式に現れる既知関数 $f_j(t)$ がすべて 0 である場合の方程式を斉次形、 一つでも 0 でないものがある場合は非斉次形と呼ぶ。

$\displaystyle Y(t) = \left[\begin{array}{c}{y_1(t)}\\ {y_2(t)}\\ \vdots\\ {y_n(...
...{f_2(t)}\\ \vdots\\ {f_n(t)}
\end{array}\right],
\hspace{1zw}A=[a_{ij}]_{i,j}
$
とすると、(18) は
  $\displaystyle
Y'(t)=AY(t)+F(t)$ (23)
と書けるが、その初期値問題、すなわち初期条件
$\displaystyle Y(t_0)=Y_0=\left[\begin{array}{c}{\alpha_1}\\ {\alpha_2}\\ \vdots\\ {\alpha_n}
\end{array}\right]
$
(すなわち $y_j(t_0)=\alpha_j$) を満たす解を考える。

$F(t)=O$ の斉次方程式の場合、$Y'(t)=AY(t)$ の初期値問題の解は、 一意的に

  $\displaystyle
Y(t)=e^{A(t-t_0)}Y_0$ (24)
と書ける。それは、斉次方程式 $Y'=AY$ の両辺に $e^{-At}$ を左からかけると、 定理 5.3 により、
$\displaystyle e^{-At}Y'=e^{-At}AY=-(e^{-At})'Y
$
となり、 $(e^{-At}Y)'=(e^{-At})'Y+e^{-At}Y'=O$ となるので $e^{-At}Y$ は 定数行列であることになる。よって $e^{-At}Y=e^{-At_0}Y_0$ となるので、 両辺に $e^{At}$ をかければ (24) が得られる。

より一般の非斉次方程式の場合は、

$\displaystyle (e^{-At}Y)'=-e^{-At}AY+e^{-At}Y'=e^{-At}(-Ay+Y')=e^{-At}F(t)
$
となるので、$t$ に関して $t_0$ から $t$ まで両辺を積分すると
$\displaystyle e^{-At}Y-e^{-At_0}Y_0=\int_{t_0}^t e^{-As}F(s)ds
$
となり、両辺 $e^{At}$ 倍して移項すれば、
  $\displaystyle
Y
=e^{At}e^{-At_0}Y_0+\int_{t_0}^t e^{At}e^{-As}F(s)ds
=e^{A(t-t_0)}Y_0+\int_{t_0}^t e^{A(t-s)}F(s)ds$ (25)
が得られる。 いずれも、$e^{At}$ の計算によって解がシンプルな形で表現されることになる。

さて、これらを用いて、実際に (20) の 初期値問題

  $\displaystyle
y(0)=y_0,
\hspace{1zw}z(0)=y'(0)=z_0$ (26)
および (21) の初期値問題
  $\displaystyle
y_1(0)=\alpha_1,
\hspace{1zw}y_2(0)=\alpha_2,
\hspace{1zw}y_3(0)=y_1'(0)=\beta_1,
\hspace{1zw}y_4(0)=y_2'(0)=\beta_2$ (27)
の解を求めてみる。

(20) は、

$\displaystyle \left[\begin{array}{c}y\\ z\end{array}\right]'
= \left[\begin{ar...
...\cos t\end{array}\right]
= A\left[\begin{array}{c}y\\ z\end{array}\right]+F(t)
$
と書けるので、初期値問題の解は、
  $\displaystyle
\left[\begin{array}{c}y\\ z\end{array}\right] = e^{At}\left[\beg...
...right]
+\int_0^te^{A(t-s)}\left[\begin{array}{c}0\\ \cos s\end{array}\right]ds$ (28)
となる。$A$ の固有値は、
$\displaystyle \left\vert\begin{array}{cc}\lambda & -1\\ 4&\lambda\end{array}\right\vert
=\lambda^2+4 = 0
$
より $\lambda_1=-2i$, $\lambda_2=2i$、固有ベクトルは
$\displaystyle \mbox{\boldmath$x$}_1=\left[\begin{array}{c}1\\ \lambda_1\end{arr...
...1\\ \lambda_2\end{array}\right]=\left[\begin{array}{c}1\\ 2i\end{array}\right]
$
となる。よって $A$ は対角化可能で、
$\displaystyle Q=[\mbox{\boldmath$x$}_1\ \mbox{\boldmath$x$}_2] = \left[\begin{a...
...ay}\right]
=\frac{1}{4}\left[\begin{array}{cc}2 & i\\ 2 & -i\end{array}\right]
$
により
$\displaystyle Q^{-1}AQ
= \left[\begin{array}{cc}\lambda_1&0\\ 0&\lambda_2\end{...
...Q^{-1}e^{At}Q =\left[\begin{array}{cc}e^{-2it}&0\\ 0&e^{2it}\end{array}\right]
$
となる。よって、
$\displaystyle e^{At}$ $\textstyle =$ $\displaystyle Q\left[\begin{array}{cc}e^{-2it}&0\\  0&e^{2it}\end{array}\right]...
...&e^{2it}\end{array}\right]
\left[\begin{array}{cc}2&i\\  2&-i\end{array}\right]$  
  $\textstyle =$ $\displaystyle \frac{1}{4}\left[\begin{array}{cc}e^{-2it}&e^{2it}\\  -2ie^{-2it}&2ie^{2it}\end{array}\right]
\left[\begin{array}{cc}2&i\\  2&-i\end{array}\right]$  
  $\textstyle =$ $\displaystyle \frac{1}{4}\left[\begin{array}{cc}2e^{-2it}+2e^{2it}& ie^{-2it}-i...
...ft[\begin{array}{cc}4\cos 2t& 2\sin 2t\\  -8\sin 2t& 4\cos 2t\end{array}\right]$  
  $\textstyle =$ $\displaystyle \left[\begin{array}{cc}\cos 2t& \sin 2t/2\\  -2\sin 2t& \cos 2t\end{array}\right]$ (29)
となる。これにより
$\displaystyle e^{At}\left[\begin{array}{c}y_0\\ z_0\end{array}\right]
=
\left...
...array}{c}y_0\cos 2t+(z_0/2)\sin 2t\\ -2y_0\sin 2t+z_0\cos 2t\end{array}\right]
$
であり、また、積分の項は、
\begin{eqnarray*}\lefteqn{\int_0^te^{A(t-s)}\left[\begin{array}{c}0\\ \cos s\end...
...n{array}{c}\cos 3t+3\cos t-4\\ 2\sin 3t+6\sin t\end{array}\right]\end{eqnarray*}
となるが、この最後の項を ${}^T\!{[\alpha\ \beta]}/12$ とし、 $\cos mt=C_m$, $\sin mt=S_m$ と略して書くと、
\begin{eqnarray*}\alpha
&=&
C_2(C_3+3C_1-4)+S_2(S_3+3S_1)
\\ &=&
(C_3C_2+S_...
...-6(S_2C_1-S_1C_2)+8S_2
\ =\
2S_1-6S_1+8S_2
\\ &=&
-4S_1+8S_2\end{eqnarray*}
となるので、結局、
$\displaystyle \int_0^te^{A(t-s)}\left[\begin{array}{c}0\\ \cos s\end{array}\rig...
...1}{3}\left[\begin{array}{c}\cos t-\cos 2t\\ -\sin t+2\sin 2t\end{array}\right]
$
となる。よって、初期値問題 (20), (26) の解は、(28) より
$\displaystyle \left[\begin{array}{c}y\\ z\end{array}\right]
=\left[\begin{array...
...1}{3}\left[\begin{array}{c}\cos t-\cos 2t\\ -\sin t+2\sin 2t\end{array}\right]
$
となることがわかる。この最後の式から $y'=z$ となること、$t=0$ のとき $y=y_0$, $z=z_0$ となることが容易に確認できる。 しかし、ここまでの手順はかなり長く、(16) を 単独方程式のまま特性方程式と代入法などで解く方がはるかに 短く容易である。

次は、(21) を考える。 こちらは最後までは計算しないが、 ジョルダン標準形位までは求めてみる。 この場合は、

$\displaystyle Y=\left[\begin{array}{c}y_1\\ y_2\\ y_3\\ y_4\end{array}\right],
...
...],
\hspace{1zw}
F = \left[\begin{array}{c}0\\ 0\\ 0\\ e^{2t}\end{array}\right]
$
となる。まずは $A$ の固有値から。
\begin{eqnarray*}\lefteqn{\left\vert\begin{array}{cccc}\lambda &0&-1&0\\ 0&\lamb...
...(\lambda^2-2\lambda-3)
% \\ &=&
\ =\
(\lambda+1)^3(\lambda-3)\end{eqnarray*}
より、固有値は $\lambda_1=3$$\lambda_2=-1$ (重複度 3) となる。

$\lambda_1=3$ の固有空間は 1 次元で、固有ベクトル 1 つを求めればよい。 固有ベクトルの方程式は、

$\displaystyle (A-3E)\mbox{\boldmath$x$}
=\left[\begin{array}{cccc}-3&0&1&0\\ 0&...
...t[\begin{array}{c}x_1\\ x_2\\ x_3\\ x_4\end{array}\right]
=\mbox{\boldmath$0$}
$
なので、行を入れかえてから逆順に消去法を行うと、
\begin{eqnarray*}\lefteqn{A-3E
\ =\
\left[\begin{array}{cccc}-3&0&1&0\\ 0&-3...
...}{cccc}0&0&0&0\\ -1&1&0&0\\ -3&0&1&0\\ -3&0&0&1\end{array}\right]\end{eqnarray*}
となるので、$x_1=s$ とすれば $x_2=s$, $x_3=3s$, $x_4=3s$ より 固有ベクトル
$\displaystyle \mbox{\boldmath$v$}_1=\left[\begin{array}{c}1\\ 1\\ 3\\ 3\end{array}\right]s
$
が求まる。$s$$s=1$ と取ればよい。

次は $\lambda=-1$ の固有ベクトルを考える。この場合は、 固有空間の次元は 3 以下だが、

\begin{eqnarray*}\lefteqn{A+E
\ =\
\left[\begin{array}{cccc}1&0&1&0\\ 0&1&0&...
...ay}{cccc}0&0&0&0\\ 1&1&0&0\\ 1&0&1&0\\ -1&0&0&1\end{array}\right]\end{eqnarray*}
なので、$-1$ に対する固有ベクトルは、
$\displaystyle \mbox{\boldmath$v$}_2=\left[\begin{array}{c}1\\ -1\\ -1\\ 1\end{array}\right]s
$
であり、固有空間の次元は 1 となる。これも $s=1$ とすればよい。 固有値 $-1$ の重複度は 3 なので、この場合、
$\displaystyle (A+E)\mbox{\boldmath$v$}_3=\mbox{\boldmath$v$}_2,
\hspace{1zw}
(A+E)\mbox{\boldmath$v$}_4=\mbox{\boldmath$v$}_3
$
となる広義固有ベクトル $\mbox{\boldmath$v$}_3,\mbox{\boldmath$v$}_4$ を取る必要がある。 なお、当然これらも一意に決定するわけではなく、例えば $\mbox{\boldmath$v$}_3$ は、
$\displaystyle \mbox{\boldmath$v$}_3 = \mbox{\boldmath$v$}^{(0)}_3 + s\mbox{\boldmath$v$}_2
$
$s$ の自由度を持つ解が求まり、 $\mbox{\boldmath$v$}_4$ はその $s$ に対してさらに
$\displaystyle \mbox{\boldmath$v$}_4 = \mbox{\boldmath$v$}^{(0)}_4 +s\mbox{\boldmath$v$}^{(0)}_3+t\mbox{\boldmath$v$}_2
$
$t$ の自由度を持つ解が求まるはずである。 これらもとりあえずは $s=t=0$ のように適当に固定してよい。

まずは $\mbox{\boldmath$v$}_3$ を求める。 方程式 $(A+E)\mbox{\boldmath$x$}=\mbox{\boldmath$v$}_2$ の解を求めるには、 拡大係数行列 $[A+E\ \mbox{\boldmath$v$}_2]$ の 消去法を行う。実際には $A+E$ の部分の変形なので、その手順は $\mbox{\boldmath$v$}_2$ の 計算と同じで、違うのは一番右側の列のみとなる。

\begin{eqnarray*}\lefteqn{[A+E\ \vert\ \mbox{\boldmath$v$}_2]
\ =\
\left[\beg...
...0&0&0\\ 1&1&0&0&1/2\\ 1&0&1&0&1\\ -1&0&0&1&-3/2\end{array}\right]\end{eqnarray*}
となるので、$x_1=s$ に対して $x_2=-s+1/2$, $x_3=-s+1$, $x_4=s-3/2$ となり、
$\displaystyle \mbox{\boldmath$v$}_3=s\left[\begin{array}{c}1\\ -1\\ -1\\ 1\end{array}\right]+\frac{1}{2}\left[\begin{array}{c}0\\ 1\\ 2\\ -3\end{array}\right]
$
となる。ここで $s=0$ としたものを $\mbox{\boldmath$v$}_3$ とする。 最後は $\mbox{\boldmath$v$}_4$ を求める。分数を消すために、方程式を 2 倍して、
$\displaystyle (A+E)(2\mbox{\boldmath$v$}_4)=2\mbox{\boldmath$v$}_3
$
と考え、 $2\mbox{\boldmath$v$}_4$ を求める計算を行う。
\begin{eqnarray*}\lefteqn{[A+E\ \vert\ 2\mbox{\boldmath$v$}_3]
\ =\
\left[\be...
...0&0&0&0&0\\ 1&1&0&0&2\\ 1&0&1&0&0\\ -1&0&0&1&-1\end{array}\right]\end{eqnarray*}
なので、 $\mbox{\boldmath$v$}_4$
$\displaystyle \mbox{\boldmath$v$}_4
=\frac{t}{2}\left[\begin{array}{c}1\\ -1\\ ...
...ray}\right]
+\frac{1}{2}\left[\begin{array}{c}0\\ 2\\ 0\\ -1\end{array}\right]
$
となる。これも $t=0$ とすれば $Q=[\mbox{\boldmath$v$}_1\ \mbox{\boldmath$v$}_2\ \mbox{\boldmath$v$}_3\ \mbox{\boldmath$v$}_4]$
$\displaystyle Q
=\left[\begin{array}{cccc}1&1&0&0\\ 1&-1&1/2&1\\ 3&-1&1&0\\ 3&1&-3/2&-1/2\end{array}\right]
$
となり、ジョルダン標準形は、
$\displaystyle J
=\left[\begin{array}{cc}J(3,1)&\raisebox{-.5ex}{\Large$0$}\\ [....
...[\begin{array}{cccc}3&0&0&0\\ 0&-1&1&0\\ 0&0&-1&1\\ 0&0&0&-1\end{array}\right]
$
となる。次は $Q^{-1}$ を求める。これも消去法で計算するが、 $QX=E$ ($X=Q^{-1}$) の両辺を 2 倍して、$(2Q)X=2E$ に対する拡大係数行列の 消去法を行い、なるべく分数計算を避けて計算する。
\begin{eqnarray*}\lefteqn{[Q\ \vert\ 2E]
=
\left[\begin{array}{cccc\vert cccc...
...5&-4\\
0&0&8&0 &5&-2&3&-4\\ 0&0&0&2 &1&2&-1&0\end{array}\right]\end{eqnarray*}
よって、
$\displaystyle Q^{-1}=\frac{1}{32}
\left[\begin{array}{cccc}3&2&5&4\\ 29&-2&-5&-...
... 29/32&-1/16&-5/32&-1/8\\
5/8&-1/4&3/8&-1/2\\ 1/2&1&-1/2&0\end{array}\right]
$
となる。定理 6.1 より、
\begin{eqnarray*}e^{J(3,1)t}
&=&
e^{3t},
\\
e^{J(-1,3)t}
&=&
e^{-t}\left...
...eft[\begin{array}{ccc}1&t&t^2/2\\ 0&1&t\\ 0&0&1\end{array}\right]\end{eqnarray*}
となり、よって
$\displaystyle e^{Jt}
=\left[\begin{array}{cc}e^{J(3,1)t}&\raisebox{-.5ex}{\Larg...
...t}&te^{-t}&t^2e^{-t}/2\\
0&0&e^{-t}&te^{-t}\\ 0&0&0&e^{-t}\end{array}\right]
$
となる。 あとは前と同様にして $e^{At}=Qe^{Jt}Q^{-1}$$e^{At}$ が求まり、 (25) を使えば $Y(t)$ が求まる、 ということになる。 しかし、この先もかなり大変な計算が待っていることが想像できる。

公式 (24), (25) は 一見シンプルな形であり、理論展開には便利で重要だが、 具体的な計算に向くかといえばそうでもなく、 特に大きな $n$ ではあまり実用的ではないことが これらの例からもわかる。

具体的な計算目的なら、むしろ単独方程式の方を特性方程式と代入法や 定数変化法 (やラプラス変換) などで解く方が易しい場合が多いだろう。

竹野茂治@新潟工科大学
2022-05-02