非斉次微分方程式の初期値問題の公式

非斉次方程式

非斉次（常微分）方程式とは、以下の方程式である： $˙ 𝒙 = 𝐴 𝑡 𝒙 + 𝒋 𝑡 (1)$ これは、前章で扱った斉次方程式 $˙ 𝒙 = 𝐴 𝑡 𝒙 (2)$ に、時刻 $𝑡$ の関数 $𝒋 𝑡$ を加えたものである。 $𝒋 𝑡$ は、斉次方程式の解に擾乱を与えるものなので、ソース項という。例えばバネの振動であれば、バネに外から与える力が $𝒋 𝑡$ に含まれる。

非斉次方程式( $1$ )の初期値問題の解の公式は、斉次方程式( $2$ )の時間発展行列 $𝐺 𝑡$ が分かれば得られる。この章では、これについて述べる。

6.1初期値問題の解の公式

非斉次方程式( $1$ )の初期値問題の解は、斉次方程式( $2$ )の時間発展行列 $𝐺 𝑡$ から求めることができる。この節では、その公式()を示す。

6.1.1デュアメルの原理

斉次方程式( $2$ )の場合、時間発展行列 $𝐺 𝑡$ が分かれば、初期値問題の解は $𝒙 𝑡 = 𝐺 𝑡 𝒙 0 (3)$ となるのであった。

斉次方程式の場合と同様に、 $𝒙 𝛿 𝑡$ を1次近似すると $𝒙 𝛿 𝑡 ≐ 𝒙 0 + ˙ 𝒙 0 𝛿 𝑡 = 𝒙 0 + (𝐴 0 𝒙 0 + 𝒋 0) 𝛿 𝑡 ≐ 𝑒 𝐴 0 𝛿 𝑡 𝒙 0 + 𝒋 0 𝛿 𝑡 𝒙 (1) \equiv 𝑒 (0) 𝒙 (0) + 𝒋 (0) 𝛿 𝑡$ となる。ただし、最後の式では、以下のように定義している：（ $𝑖 = 0, 1, \dots$ ） $𝑒 (𝑖) \equiv 𝑒 𝐴 𝑖 𝛿 𝑡 𝛿 𝑡 𝒙 (𝑖) \equiv 𝒙 𝑖 𝛿 𝑡$

この記法を使って、 $𝒙 2 𝛿 𝑡, 𝒙 3 𝛿 𝑡$ まで計算してみると、以下のようになる： $𝒙 2 𝛿 𝑡 ≐ 𝑒 (1) 𝒙 (1) + 𝒋 (1) 𝛿 𝑡 \approx 𝑒 (1) 𝑒 (0) ⏟ \equiv 𝑒 (1, 0) 𝒙 (0) + (𝑒 (1) 𝒋 (0) + 𝒋 (1)) 𝛿 𝑡 𝒙 3 𝛿 𝑡 \approx 𝑒 (2, 1, 0) 𝒙 (0) + (𝑒 (2, 1) 𝒋 (0) + 𝑒 (1) 𝒋 (1) + 𝒋 (2)) 𝛿 𝑡 = 𝑒 (2, 1, 0) [𝒙 (0) + (𝑒 - (0) 𝒋 (0) + 𝑒 - (1, 0) 𝒋 (1) + 𝑒 - (2, 1, 0) 𝒋 (2)) 𝛿 𝑡]$ ここまで展開すると、規則性が見える。即ち、以下が成り立つと予想できる： $𝒙 𝑡 = 𝐺 𝑡 (𝒙 0 + \int 𝑡 𝑡' = 0 𝐺 - 1 𝑡' 𝒋 𝑡')$ この式をデュアメルの原理といい、実際、正しいことが示せる（以下の【6.1-注1】）。よって、前章の方法を使うなどして斉次方程式の時間発展行列 $𝐺 𝑡$ （グリーン行列）を求めておけば、デュアメルの原理によって、非斉次方程式の初期値問題が解ける。

【6.1-注1】非斉次方程式の解の公式：デュアメルの原理

非斉次方程式( $1$ )において、初期値 $𝒙 0$ が与えられたときの解 $𝒙 𝑡$ は、斉次方程式( $2$ )の時間発展行列 $𝐺 𝑡$ を用いて、以下のように書ける（デュアメルの原理）： $𝒙 𝑡 = 𝐺 𝑡 (𝒙 0 + \int 𝑡 𝑡' = 0 𝐺 - 1 𝑡' 𝒋 𝑡') (4)$ $𝐺 𝑡$ のことを、（非斉次方程式の初期値問題における）グリーン行列と呼ぶ。

係数行列 $𝐴$ が、時間に依存しない定数であれば、 $𝐺 𝑡 = 𝑒 𝐴 𝑡$ と書けるので、以下のようになる：

$𝒙 𝑡 = 𝑒 𝐴 𝑡 𝒙 0 + \int 𝑡 𝑡' = 0 𝑒 𝐴 (𝑡 - 𝑡') 𝒋 𝑡' (5)$

証明

式( $4$ )において、 $𝑡 = 0$ とすれば、積分が消えるので、両辺が等しくなり、初期条件を満たすことが分かる。非斉次方程式( $1$ )を満たすことは、式( $4$ )を辺々微分すれば分かる： $𝑑 𝑑 𝑡 (左辺) = ˙ 𝒙 𝑡 𝑑 𝑑 𝑡 (右辺) = ˙ 𝐺 𝑡 (𝒙 0 + \int 𝑡 𝑡' = 0 𝐺 - 1 𝑡' 𝒋 𝑡') + 𝐺 𝑡 (𝟎 + 𝐺 - 1 𝑡 𝒋 𝑡) ∣ ˙ 𝐺 𝑡 = 𝐴 𝐺 𝑡 = 𝐴 𝑡 𝒙 𝑡 + 𝒋 𝑡 ◼$

後の章で扱うが、この議論は非斉次の偏微分方程式についても成り立つ。デュアメルの原理という用語はそちらを指すことが多い。上で述べた常微分方程式の場合の式( $4$ )は、定数変化法と呼ばれることが多い。定数変化法と呼ぶのは、「斉次方程式の解 $𝒙 = 𝐺 𝑡 𝒙 0$ において、定数 $𝒙 0$ を変数 $𝒚$ に置き換えたもの」 $𝒙 = 𝐺 𝑡 𝒚$ を非斉次方程式( $1$ )に代入して $𝒚$ を求めるという手法による。実際に代入してみると $˙ 𝒙 = 𝐴 𝑡 𝒙 + 𝒋 𝑡 ˙ 𝐺 𝑡 𝒚 + 𝐺 𝑡 ˙ 𝒚 = 𝐴 𝑡 𝐺 𝑡 𝒚 + 𝒋 𝑡 𝐺 𝑡 ˙ 𝒚 = 𝒋 𝑡 ∴ 𝒚 = \int 𝑡 𝑡' = 0 𝐺 - 1 𝑡' 𝒋 𝑡' + 𝑪$ となり、積分定数 $𝑪$ は、初期条件 $𝒚 0 = 𝒙 0$ より $𝑪 = 𝒙 0$ なので、 $𝒚$ が式( $4$ )の $(\dots)$ 部分に確かに一致することが分かる。

6.2例題

例題として、放物運動と強制振動を扱う。

6.2.1【例題】放物運動をデュアメルの原理で解く

放物運動の運動方程式 $¨ 𝒙 = 𝒈$ を解くことを考える（あまり面白くない例だが感覚をつかむため）。これは、非斉次微分方程式である。実際、非斉次方程式( $1$ )の形で書くと $𝑑 𝑑 𝑡 [𝒙 ˙ 𝒙] = [0100] ⏟ 𝐴 [𝒙 ˙ 𝒙] + [0 𝒈] ⏟ 𝒋 𝑡$ となる。 $𝐴$ はすでにジョルダン標準形（固有値が $0$ ）になっているので、 $𝑒 𝐴 𝑡$ は容易に求まる： $𝑒 𝐴 𝑡 = [1 𝑡 01]$ これらをデュアメルの原理( $5$ )に代入すると、よく知られた初期値問題の解が得られる： $[𝒙 ˙ 𝒙] 𝑡 = [1 𝑡 01] [𝒙 0 ˙ 𝒙 0] + \int 𝑡 𝑡' = 0 [1 𝑡 - 𝑡' 01] [0 𝒈] = ⎡ ⎢ ⎢ ⎣ 𝒙 0 + 𝑡 ˙ 𝒙 0 + 12 𝑡 2 𝒈 ˙ 𝒙 0 + 𝑡 𝒈 ⎤ ⎥ ⎥ ⎦$

6.2.2【例題】時間がたった後の強制振動：式( $7$ )

強制振動とは、減衰振動の方程式に、時間に依存する外力 $𝑓 𝑡$ を加えたものである： $¨ 𝑥 = - 𝜔 2 𝑥 - 2 𝛼 ˙ 𝑥 + 𝑚 - 1 𝑓 𝑡$ これを、非斉次方程式( $1$ )の形で書けば $𝑑 𝑑 𝑡 [𝑥 ˙ 𝑥] = [01 - 𝜔 2 - 2 𝛼] ⏟__⏟__⏟ 𝐴 [𝑥 ˙ 𝑥] + [0 𝑚 - 1 𝑓 𝑡] ⏟ 𝒋 𝑡$ となるので、その解は、デュアメルの原理( $5$ )により以下のようになる：（ $𝑒 𝐴 𝑡$ は前章の5.4節で既に求めている） $[𝑥 ˙ 𝑥] 𝑡 = 𝑒 𝐴 𝑡 [𝒙 0 ˙ 𝒙 0] + \int 𝑡 𝑡' = 0 𝑒 𝐴 (𝑡 - 𝑡') [0 𝑚 - 1 𝑓 𝑡'] (6)$

ここでは、 $𝑓 𝑡$ として、周期的なものを考える： $𝑓 𝑡 \equiv 𝑓 0 c o s 𝜔 0 𝑡$ この時、時間が経った後でどのような状態に落ち着くかを知りたい。式( $6$ )の右辺第1項は $𝑒 - 𝛼 𝑡$ が含まれるため、 $𝑡 𝛼 ≫ 1$ を満たす時刻では無視してよい。また、第2項について、 $𝑥 𝑡$ に効いてくるのは、 $𝑒 𝐴 (𝑡 - 𝑡')$ の右上成分 $[𝑒 𝐴 (𝑡 - 𝑡')] 12$ だけである。この成分は、前章の【5.4-注1】で与えた減衰振動の解において、 $˙ 𝑥 0$ の係数である： $[𝑒 𝐴 (𝑡 - 𝑡')] 12 = 𝑒 𝜆 + (𝑡 - 𝑡') - 𝑒 𝜆 - (𝑡 - 𝑡') 𝜆 + - 𝜆 - 𝜆 \pm \equiv - 𝛼 \pm \sqrt 𝛼 2 - 𝜔 2$ （後の計算を簡単にするため、場合分けをまとめて1つの式にしている。）。これを使うと、 $𝑥 𝑡$ は、 $𝑡 ≫ 𝛼 - 1$ の時、以下のようになる： $𝑥 𝑡 \approx 0 + \int 𝑡 𝑡' = 0 [𝑒 𝐴 (𝑡 - 𝑡')] 12 𝑓 0 𝑚 c o s 𝜔 0 𝑡' = 𝑓 0 𝑚 R e \int 𝑡 𝑡' = 0 [𝑒 𝐴 (𝑡 - 𝑡')] 12 𝑒 𝑖 𝜔 0 𝑡' = 𝑓 0 𝑚 R e \int 𝑡 𝑡' = 0 𝑒 𝑖 𝜔 0 𝑡' + 𝜆 + (𝑡 - 𝑡') - 𝑒 𝑖 𝜔 0 𝑡' + 𝜆 - (𝑡 - 𝑡') 𝜆 + - 𝜆 - = 𝑓 0 𝑚 R e [𝑒 𝑖 𝜔 0 𝑡' + 𝜆 + (𝑡 - 𝑡') 𝑖 𝜔 0 - 𝜆 + - 𝑒 𝑖 𝜔 0 𝑡' + 𝜆 - (𝑡 - 𝑡') 𝑖 𝜔 0 - 𝜆 -] 𝑡 𝑡' = 0 𝜆 + - 𝜆 - = 𝑓 0 𝑚 R e 𝑒 𝑖 𝜔 0 𝑡 𝑖 𝜔 0 - 𝜆 + - 𝑒 𝑖 𝜔 0 𝑡 𝑖 𝜔 0 - 𝜆 - 𝜆 + - 𝜆 - - (𝑡' = 0) ⏟ \approx 0 ∣ 分母の実数化 \approx 𝑓 0 𝑚 R e 𝜔 2 - 𝜔 20 - 𝑖 2 𝛼 𝜔 0 (𝜔 2 - 𝜔 20) 2 + 4 𝛼 2 𝜔 20 𝑒 𝑖 𝜔 0 𝑡 ∣ ∣ ∣ ∣ ∣ 分数部分の複素数を指数関数で表す。 t a n 𝜃 = 2 𝛼 𝜔 0 𝜔 20 - 𝜔 2 = 𝑓 0 𝑚 R e 𝑒 𝑖 𝜃 \sqrt (𝜔 2 - 𝜔 20) 2 + 4 𝛼 2 𝜔 20 𝑒 𝑖 𝜔 0 𝑡 = 𝑓 0 𝑚 1 \sqrt (𝜔 2 - 𝜔 20) 2 + 4 𝛼 2 𝜔 20 c o s (𝜔 0 𝑡 + 𝜃) (7)$ よって、角周波数 $𝜔 0$ で振動することが分かる。ただし、外力の振動に対して $𝜃$ だけずれる。