線形代数II/連立線形微分方程式 の履歴(No.20)

基本†

通常の１階線形微分方程式

$$\frac{dx(t)}{dt}=ax(t)$$

を解くには、まず両辺を $x(t)$ で割り、

$$\frac{\ \frac{dx(t)}{dt}\ }{x(t)}=\frac{d}{dt}\log x(t)=a$$

中辺、右辺を積分して（$C'$: 積分定数）、指数関数を取ることにより

$$\log x(t)=at + C'$$

$$x(t)=Ce^{at}$$

とすればよい。ここで $C=e^{C'}$ であるが、$t=0$ を代入すれば

$$x(0)=C$$

となることに注意すると、

$$x(t)=e^{at}\,x(0)$$

と表すことができる。

連立方程式†

未知の関数 $x_1(t),x_2(t)$ に対して 以下の連立１階線形微分方程式が与えられているとする。

$$\left\{ \begin{matrix} \displaystyle\frac{dx_1(t)}{dt}=\dot x_1(t)=ax_1(t)+bx_2(t) \\ \displaystyle\frac{dx_2(t)}{dt}=\dot x_2(t)=cx_1(t)+dx_2(t) \end{matrix} \right.$$

この方程式は

$$\bm x(t)=\begin{pmatrix}x_1(t)\\x_2(t)\end{pmatrix}, \ \ A=\begin{pmatrix}a&b\\c&d\end{pmatrix}$$

を用いれば

$$\frac{d\bm x(t)}{dt}=\dot{\bm x}(t)=A\bm x(t)$$

と表せるから、行列の関数を用いることで方程式の解を

$$ \bm x(t)=e^{At}\bm x(0) $$

として求められる！！！（ドヤッ

で終わってしまっても良いのだけれど念のため、$A$ が対角化可能な場合について行列の関数を用いない方法との対比を確認しておこう。

$A$ を $P^{-1}AP=D=\begin{pmatrix}\lambda_1&\\&\lambda_2\end{pmatrix}$ の形で対角化できるなら、１年生の二次形式のところでやったのと同様に、

$$ \bm x(t)=P\bm y(t)=P\begin{pmatrix}y_1(t)\\y_2(t)\end{pmatrix} $$

を満たす $\bm y(t)$ へと変数変換することで、

$$ \dot{\bm x}(t)=\frac{d}{dt}\{P\bm y(t)\}=P\dot{\bm y}(t)=A\bm x(t)=AP\bm y(t) $$

すなわち、

$$\dot{\bm y}(t)=P^{-1}AP\bm y(t)=D\bm y(t)$$

のように非対角項がゼロとなる方程式が得られる。これは

$$\left\{\begin{matrix} \dot y_1(t)=\lambda_1y_1(t)\\ \dot y_2(t)=\lambda_2y_2(t) \end{matrix} \right.$$

を表わすから、それぞれ $y_1(t),y_2(t)$ に対して独立に解くことができて、

$$\left\{\begin{matrix} y_1(t)=y_1(0)e^{\lambda_1t}\\ y_2(t)=y_2(0)e^{\lambda_2t}\end{matrix} \right.$$

あるいは $$\bm y(t)=\begin{pmatrix}e^{\lambda_1t}\\&e^{\lambda_2t}\end{pmatrix}\bm y(0)$$

のように解が得られる。

このように、方程式を単純化する基底を取って解を表わした座標を「基準座標」と呼ぶ。

最終的な問題は基準座標に対する解 $\bm y(t)$ を得ることではなく、 $\bm x(t)$ を求めることなので、 $\bm y(0)=P^{-1}\bm x(0)$ を使って変形すると、

$$\bm x(t)=P\bm y(t) =P\begin{pmatrix}e^{\lambda_1t}\\&e^{\lambda_2t}\end{pmatrix}\bm y(0) =P\begin{pmatrix}e^{\lambda_1t}\\&e^{\lambda_2t}\end{pmatrix}P^{-1}\bm x(0) $$

が得られ、これが初期条件 $\bm x(0)$ が与えられた際に $\bm x(t)$ を求めるための式となる。

先にやったように $e^{At}=P\begin{pmatrix}e^{\lambda_1t}\\&e^{\lambda_2t}\end{pmatrix}P^{-1}$ であることを用いてこの式を、

$$\bm x(t)=e^{At}\bm x(0)$$

と書き直すと、通常の線形微分方程式の解とほとんど同じに解ける、というわけだ。

$$ \begin{array}{lll} \dot x(t)=ax(t) &\to& x(t)=x(0)e^{at} \\ \dot{\bm x}(t)=A\bm x(t) &\to& \bm x(t)=e^{At}\,\bm x(0) \\ \end{array} $$

またこの $\bm x(t)$ が与えられた微分方程式を満たすことも容易に確認できる。

$$\dot{\bm x}(t)=\frac{d}{dt}\left\{e^{At}\bm x(0)\right\} =\left\{\frac{d}{dt}e^{At}\right\}\bm x(0) =Ae^{At}\bm x(0) =A\bm x(t) $$

$A$ が対角化可能でない場合にも $\bm x(t)=e^{At}\bm x(0)$ と書いて問題ないが、 $e^{At}$ を固有値を用いてどのように書いたらよいるかを理解するには 発展：ジョルダン標準形 について学ぶ必要がある。（加えて、固有値にゼロを含むときも特別扱いが必要になる）

例題†

$\left\{\begin{matrix} \frac{dx_1(t)}{dt}&=3x_1(t)+\phantom 3x_2(t)\\[3mm] \frac{dx_2(t)}{dt}&=\phantom 3x_1(t)+3x_2(t) \end{matrix}\right.$ 　を初期条件 $x_1(0)=6,x_2(0)=2$ の下で解け。

$\bm x(t)=\begin{pmatrix}x_1(t)\\x_2(t)\end{pmatrix}$、 $A=\begin{pmatrix}3&1\\1&3\end{pmatrix}$ と置けば、 $\dot{\bm x}(t)=A\bm x(t)$ と書けるから、

$\bm x(t)=e^{tA}\bm x(0)$ と表せる。

$P^{-1}AP=\begin{pmatrix}\lambda_1\\&\lambda_2\end{pmatrix}$ と対角化できるならさらに変形して、

$$\bm x(t)=P\begin{pmatrix}e^{\lambda_1t}\\&e^{\lambda_2t}\end{pmatrix}P^{-1}\bm x(0)$$

である（$\lambda_1\ne 0,\lambda_2\ne 0$ の時）。

$A$ の固有値を求め、対角化しよう。

$$|A-\lambda E|=\begin{vmatrix}3-\lambda&1\\1&3-\lambda\end{vmatrix}=(3-\lambda)^2-1=(4-\lambda)(2-\lambda)$$

より、$\lambda_1=4,\lambda_2=2$

それぞれに対応する固有ベクトルは、

$(A-4 E)\bm x=\begin{pmatrix}-1&1\\1&-1\end{pmatrix}\bm x=\bm 0$ より $\bm x=s\begin{pmatrix}1\\1\end{pmatrix}$

$(A-2 E)\bm x=\begin{pmatrix}1&1\\1&1\end{pmatrix}\bm x=\bm 0$ より $\bm x=s\begin{pmatrix}1\\-1\end{pmatrix}$

したがって、$P=\begin{pmatrix}1&1\\1&-1\end{pmatrix}$、$P^{-1}=\frac{1}{2}\begin{pmatrix}1&1\\1&-1\end{pmatrix}$、$P^{-1}AP=\begin{pmatrix}4\\&2\end{pmatrix}$ となって、

$$\begin{aligned} \bm x(t)&=\begin{pmatrix}1&1\\1&-1\end{pmatrix}\begin{pmatrix}e^{4t}\\&e^{2t}\end{pmatrix}\frac{1}{2}\begin{pmatrix}1&1\\1&-1\end{pmatrix}\begin{pmatrix}6\\2\end{pmatrix}\\ &=\begin{pmatrix}e^{4t}&e^{2t}\\e^{4t}&-e^{2t}\end{pmatrix}\begin{pmatrix}1&1\\1&-1\end{pmatrix}\begin{pmatrix}3\\1\end{pmatrix}\\ &=\begin{pmatrix}e^{4t}&e^{2t}\\e^{4t}&-e^{2t}\end{pmatrix}\begin{pmatrix}4\\2\end{pmatrix}\\ &=\begin{pmatrix}4e^{4t}+2e^{2t}\\4e^{4t}-2e^{2t}\end{pmatrix}\\ \end{aligned}$$

が求める解である。

演習問題†

以下の微分方程式を境界条件 $x_1(0)=1,x_2(0)=2$ の下で解け。

$$ \begin{cases} \frac{d}{dt}x_1(t)=-x_1(t)+x_2(t)\\ \frac{d}{dt}x_2(t)=\phantom{-}x_1(t)-x_2(t)\\ \end{cases} $$

基本†

通常の２階線形微分方程式が次のように与えられたとする。

$$\frac{d^2x(t)}{dt^2}=-\omega^2x(t)\hspace{10mm}(\omega>0)$$

$$x''=-\omega^2x$$

これを解くには*1参考：http://w3e.kanazawa-it.ac.jp/m...、 まず全体に $x'$ をかけてから積分して、

$$x'x''+\omega^2xx'=0$$

$$\frac12x'^2+\frac12\omega^2x^2=C_0\hspace{1cm}C_0:\text{ 積分定数}$$

を得る(この関係は物理学ではエネルギー保存則と関係が深い)。

$x$ が実数関数であれば明らかに $C_0\ge 0$ であるから $2C_0=\omega^2A^2\ge 0$ と書くことにすると、

$$x'^2=\omega^2A^2-\omega^2x^2$$

すなわち、

$$x'=\pm\omega\sqrt{A^2-x^2}$$

を得る（三角関数の形が見えてきた）。 変数分離して両辺を積分すると、

$$\pm\frac{x'}{\sqrt{A^2-x^2}}=\omega$$

$$\pm\int \frac{x'}{\sqrt{A^2-x^2}}dt=\pm\int \frac{dx}{\sqrt{A^2-x^2}}=\pm\cos^{-1}\frac xA=\omega t+\theta\hspace{1cm}\theta:\text{ 積分定数}$$

両辺の $\cos$ を取ると $\cos x=\cos(-x)$ より $\pm$ は消えて、

$$ x(t)=A\cos(\omega t+\theta)=C\cos\omega t+C'\sin\omega t $$

ただし、$C=A\cos\theta,C'=-A\sin\theta$ と表せる。

$$\dot x(t)=-\omega C\sin\omega t+\omega C'\cos\omega t$$

より、

$$x(0)=C, \dot x(0)=\omega C'$$

であるから $t=0$ における境界条件を用いれば

$$x(t)=x(0)\cos\omega t+\frac{\dot x(0)}{\omega}\sin\omega t$$

と表せる。

連立方程式†

未知の関数 $x_1(t),x_2(t)$ が次の２階線形連立微分方程式を満たす時、

$$\left\{ \begin{matrix} \displaystyle\frac{d^2x_1(t)}{dt^2}=\ddot x_1(t)=-ax_1(t)-bx_2(t) \\ \displaystyle\frac{d^2x_2(t)}{dt^2}=\ddot x_2(t)=-cx_1(t)-dx_2(t) \end{matrix} \right.$$

この方程式はベクトル関数 $\bm x(t)=\begin{pmatrix}x_1(t)\\x_2(t)\end{pmatrix}$、 行列 $A=\begin{pmatrix}a&b\\c&d\end{pmatrix}$ を用いて、

$$\frac{d^2\bm x(t)}{dt^2}=\ddot{\bm x}(t)=-A\bm x(t)$$

と表せる。

もし $A$ を $P^{-1}AP=D=\begin{pmatrix}\lambda_1&\\&\lambda_2\end{pmatrix}$ の形で対角化できるなら、

$$ \bm x(t)=P\bm y(t)=P\begin{pmatrix}y_1(t)\\y_2(t)\end{pmatrix} $$

の変数変換により一次の時と同様に、

$$\ddot{\bm y}(t)=-P^{-1}AP\bm y(t)=-D\bm y(t)$$

が得られる。

$\lambda_1>0,\ \lambda_2>0$ であれば $\omega_1=\sqrt{\lambda_1},\ \omega_2=\sqrt{\lambda_2}$ と置けて、これは

$$\left\{\begin{matrix} \ddot y_1(t)=-\lambda_1y_1(t)=-\omega_1^2y_1(t)\\ \ddot y_2(t)=-\lambda_2y_2(t)=-\omega_2^2y_2(t) \end{matrix} \right.$$

という方程式を表わすから、その解を

$$\left\{\begin{matrix} y_1(t)=y_1(0)\cos\omega_1t+\frac{\dot y_1(0)}{\omega_1}\sin\omega_1t\\ y_2(t)=y_2(0)\cos\omega_2t+\frac{\dot y_2(0)}{\omega_2}\sin\omega_2t \end{matrix}\right.$$

あるいは

$$\bm y(t)=\begin{pmatrix}\cos\omega_1t\\&\cos\omega_2t\end{pmatrix}\bm y(0)+ \begin{pmatrix}\frac{\sin\omega_1t}{\omega_1}\\&\frac{\sin\omega_1t}{\omega_1}\end{pmatrix}\dot{\bm y}(0)$$

と表せる。$\bm y(0)=P^{-1}\bm x(0), \dot{\bm y}(0)=P^{-1}\dot{\bm x}(0)$ を用いて書き直すと、

$$\bm x(t)=P\bm y(t)=P\begin{pmatrix}\cos\omega_1t\\&\cos\omega_2t\end{pmatrix}P^{-1}\bm x(0)+ P\begin{pmatrix}\frac{\sin\omega_1t}{\omega_1}\\&\frac{\sin\omega_1t}{\omega_1}\end{pmatrix}P^{-1}\dot{\bm x}(0)$$

これが、初期条件 $\bm x(0),\dot{\bm x}(0)$ が与えられた際に $\bm x(t)$ を求める式になる。

$\omega=\sqrt\lambda$ であったことを念頭に置きつつ「行列の関数」をイメージすれば、

$$\bm x(t)=\Bigg[\cos\left(\sqrt{A}\,t\right)\Bigg]\bm x(0)+\Bigg[\frac{\sin\left(\sqrt{A}\,t\right)}{\sqrt{A}}\Bigg]\dot{\bm x}(0)$$

の形であり、通常の調和振動の式

$$x(t)=x(0)\cos\omega t+\frac{\dot x(0)}{\omega}\sin\omega t$$

との類似性が際立つ。

ここで表れた $\sqrt A$ を含む「行列の関数」は先の授業でやったような「固有値の整数次のみの多項式」の形では表せないが、多項式が分数次の項を含んでも良いとして定義を拡張すれば、ほぼ同じ議論でこの書き方を正当化できる。

例題†

次の連立微分方程式を考える。

$$\left\{\begin{matrix} m\ddot x_1(t)=-kx_1(t)-k\{x_1(t)-x_2(t)\}\\ m\ddot x_2(t)=-kx_2(t)-k\{x_2(t)-x_1(t)\}\\ \end{matrix}\right.$$

これは右図のように、３つのバネで繋がれた２つの質点の運動を記述する方程式を考えているのであるが、 物理モデルから上記方程式を組み立てる部分は数学の範疇ではないので、線形代数学においては 上記方程式をどのように解くかに焦点を絞る。

$\bm x(t)=\begin{pmatrix}x_1(t)\\x_2(t)\end{pmatrix}$、 $\displaystyle A=\frac{k}{m}\begin{pmatrix}2&-1\\-1&2\end{pmatrix}$ と置けば、上式は

$$\ddot{\bm x}(t)=-A\bm x(t)$$

と書けるから、$A$ の固有値がすべて正であるとき、その解を形式的に

$$\bm x(t)=\Bigg[\cos\sqrt{A}t\Bigg]\bm x(0)+\Bigg[\frac{\sin\sqrt{A}t}{\sqrt A}\Bigg]\dot{\bm x}(t)$$

と書ける。

今の場合、$|(m/k)A-\lambda' E|=\begin{vmatrix}2-\lambda'&-1\\-1&2-\lambda'\end{vmatrix}=(2-\lambda')^2-1=(1-\lambda')(3-\lambda')=0$ より、

$\lambda'=1,3$ すなわち $\lambda_1=\frac{k}{m}>0,\lambda_2=\frac{3k}{m}>0$ であり、

対応する固有ベクトルはそれぞれ $\begin{pmatrix}1\\1\end{pmatrix},\begin{pmatrix}1\\-1\end{pmatrix}$ であるから、

$P=\begin{pmatrix}1&1\\1&-1\end{pmatrix}$ と置けば、

$$P^{-1}AP=\begin{pmatrix}\lambda_1&0\\0&\lambda_2\end{pmatrix}$$

として対角化できる。

このとき、

$$\bm x(t)= P\begin{pmatrix}\cos\sqrt{\lambda_1}t&0\\0&\cos\sqrt{\lambda_2}t\end{pmatrix}P^{-1}\bm x(0)+ P\begin{pmatrix}\frac{\sin\sqrt{\lambda_1}t}{\sqrt{\lambda_1}}&0\\0&\frac{\sin\sqrt{\lambda_2}t}{\sqrt{\lambda_2}}\end{pmatrix}P^{-1}\dot{\bm x}(0)$$

と表せ、

$$ \begin{aligned} P\begin{pmatrix}a\\&b\end{pmatrix}P^{-1}&=\begin{pmatrix}1&1\\1&-1\end{pmatrix}\begin{pmatrix}a\\ &b\end{pmatrix}\frac{1}{2}\begin{pmatrix}1&1\\1&-1\end{pmatrix}\\ &=\frac{1}{2}\begin{pmatrix}a&b\\a&-b\end{pmatrix}\begin{pmatrix}1&1\\1&-1\end{pmatrix}\\ &=\frac{1}{2}\begin{pmatrix}a+b&a-b\\a-b&a+b\end{pmatrix} \end{aligned} $$

を用いて計算を進めると、

$$ \bm x(t)= \frac{1}{2}\begin{pmatrix} \cos{\omega_1}t+\cos{\omega_2}t&\cos{\omega_1}t-\cos{\omega_2}t\\ \cos{\omega_1}t-\cos{\omega_2}t&\cos{\omega_1}t+\cos{\omega_2}t \end{pmatrix}\bm x(0)+ \frac{1}{2\omega_1}\begin{pmatrix} \sin{\omega_1}t+\tfrac1{\sqrt3}\sin{\omega_2}t& \sin{\omega_1}t-\tfrac1{\sqrt3}\sin{\omega_2}t\\ \sin{\omega_1}t-\tfrac1{\sqrt3}\sin{\omega_2}t& \sin{\omega_1}t+\tfrac1{\sqrt3}\sin{\omega_2}t\\ \end{pmatrix}\dot{\bm x}(0) $$

が $\bm x(0),\dot{\bm x}(0)$ で表わした $\bm x(t)$ の表式となる。 ここで、$\omega_1=\sqrt{\lambda_1}=\sqrt{\frac{k}{m}},\omega_2=\sqrt{\lambda_2}=\sqrt3\omega_1$ である。

演習問題†

$x$ 軸上のみを動ける３つの質点が、ばね定数 $k$ を持つ２つのばねで連結された図のような系を考える。ただし、真ん中の質点の質量を $m$、両端の質点の質量をそれぞれ $m/2$ とする。このとき運動方程式は次のようになる。以下の問いに従って、この方程式の一般解を基準運動に分離せよ。

$$ \begin{cases} (m/2)\ddot x_1(t)=\phantom{-}k\{x_2(t)-x_1(t)\}\\ \hspace{6.2mm}m\ddot x_2(t)=-k\{x_2(t)-x_1(t)\}+k\{x_3(t)-x_2(t)\}\\ (m/2)\ddot x_3(t)=\phantom{-k\{x_2(t)-x_1(t)\}}-k\{x_3(t)-x_2(t)\}\\ \end{cases} $$

(1) $\bm x(t)=\begin{pmatrix}x_1(t)\\x_2(t)\\x_3(t)\end{pmatrix}$ に対して、与えられた方程式を $\ddot \bm x=-A\bm x$ と表した際の係数行列 $A$ を求めよ。

(2) $A$ の固有値 $\lambda_1,\lambda_2,\lambda_3$ を求めよ。また、対応する固有ベクトルを１つずつ求め $\bm p_1,\bm p_2,\bm p_3$ とせよ。ただし、$\lambda_1<\lambda_2<\lambda_3$ とする。

(3) $P=\Big(\bm p_1\ \bm p_2\ \bm p_3\Big)$ に対して $P^{-1}$ を求めよ。

(4) $\bm y=\begin{pmatrix}y_1(t)\\y_2(t)\\y_3(t)\end{pmatrix}$ が $\bm x=P\bm y$ を満たすとき、 $y_1,y_2,y_3$ に対する運動方程式を求め、それぞれの一般解 $y_1^*,y_2^*,y_3^*$ を適当なパラメータを含む形で書き下せ。

(5) $\bm x$ の一般解は $\bm x(t)=\bm p_1y_1^*+\bm p_3y_3^*+\bm p_3y_3^*$ で与えられることを念頭に置き、 $y_1,y_2,y_3$ で表される基準運動がどのような運動を表すか説明せよ。

発展：$n$ 粒子の系†

両端を壁で挟まれた質量 $m$ の２つの質点が３つのバネ（ばね定数 $k$）で繋がれた系を上で解いたが、
これを一般化して 両端を壁で挟まれた $n$ 個の質点が $n+1$ 個のバネで繋がれた系について考える。

例†

電荷密度分布 $\rho(\bm x)$ の作る電位 $\phi(\bm x)$ を求めたい。

基本方程式はポアソン方程式

$$-\nabla^2\phi(\bm x)=\frac1\epsilon_0\rho(\bm x) $$

で与えられるが、この方程式は

$$ \underbrace{-\epsilon_0\nabla^2}_T\phi(\bm x)=\rho(\bm x) $$

としたときの線形演算子 $T$ の逆 $T^{-1}$ を求めれば $\phi(x)=T^{-1}\rho(x)$ として解ける。

$\bm x'$ に置かれた単位点電荷の作る電位を無限遠で $\phi(\bm x)=0$ となる境界条件の下で解けば、

$$ \phi_1(\bm x;\bm x')=-\frac{1}{4\pi\epsilon_0}\frac{1}{|\bm x-\bm x'|} $$

となるが、これは

$$ T\phi_1(\bm x;\bm x')=\delta^2(\bm x-\bm x') $$

を表しており、$T$ に対する同境界条件を満たすグリーン関数が $G(\bm x,\bm x')=\phi_1(\bm x;\bm x')$ であることを表している。

したがって、任意の電荷分布 $\rho(\bm x)$ の作る電位を

$$ \begin{aligned} \phi(\bm x)&=T^{-1}\rho(\bm x)=\iiint G(\bm x,\bm x')\rho(\bm x')d^3\bm x'\\ &=-\frac{1}{4\pi\epsilon_0}\iiint \frac{\rho(\bm x')}{|\bm x-\bm x'|}d^3\bm x'\\ \end{aligned} $$

のようにして求められるのである。

この式の物理的な意味は、各点における電荷密度が作る電位を重ね合わせることで全体としての電位が求まる、と読めるが、重ね合わせることが可能な理由は $T$ の線形性にあることを理解せよ。

非線形な方程式ではこのような「重ね合わせ」は通用しない。

前の単元 ＜＜＜ 　　　　　　　　　　　　　　　線形代数II 　　　　　　　　　　　　　　　＞＞＞ 全体のまとめ

無題†

( 2018-08-23 (木) 15:02:03)

むずかっしー

• 勉強になります -- 2019-06-08 (土) 18:11:21
• 特性根が重解　と　複素数解の場合も言ってほしいけど -- 2019-07-15 (月) 00:39:43
• 参考になります -- 雑魚 2020-07-26 (日) 23:25:58