化学数学 安藤耕司のページ chap03

第 ₃ 章 常微分方程式

化学における1 階常微分方程式と言えば、まずは反応速度式である。本 章では、できるだけ反応速度論から例を引きながら、常微分方程式につ いて議論する。ただし、反応速度論そのものについては詳しくは説明し ない。時刻t における化学種 A の濃度を [A]tのように書く場合が多いが、 ここでは簡単のためA(t) または単に A と書く。

2 階常微分方程式については、強制振動子を例にとる。こちらは、分子 運動のモデルとして用いられることがある。例えば、分子振動や回転の 運動方程式に、媒質からの摩擦や、光電場による振動外力を現象論的に 加えたモデルである。

3.1 ^{変数分離型}

常微分方程式の初等解法で最も基本的なのは、変数分離型である。典 型的な例について計算した後、3.1.5 節で一般論をまとめる。

3.1.1 ^単分子 1 ^次反応

単純な1 次反応 A → P の反応速度は、 dA

dt ^{= −kA (3.1)}

と表される。k は速度定数である。導関数が自身に比例することから、こ の解は指数関数

A(t) = A(0)e^−kt (3.2) であることが直ちに分かるが、丁寧に書くならば次のようになる。

(3.1) ⇒

∫ _dA A ^{= −k}

∫

dt ⇒ ln A(t) = −kt + C

t = 0 とおけば C = ln A(0) なので、

⇒ ln(A(t)/A(0)) = −kt ⇒ (3.2)

補足 上では不定積分の後、積分定数_Cを初期条件_A(0)から決めたが、次のよ うに定積分として考えてもよい。

∫ _A(t)

A(0)

dA^′ A^{′ = −k}

∫ t 0

dt^′ ⇒ [ln A^′]A(t)_A(0)= −k[t^′]t₀

3.1.2 ^単分子 2 ^次反応

２次反応の速度式は、

dA

dt ^{= −kA}

2 _(3.3)

と表される。この場合は、

−

∫ _dA A^{2 = k}

∫

dt ⇒ ¹ A(t) ⁻

1

A(0) ^{= kt}

⇒ A(t) = ^A(0)

1 + ktA(0) ^(3.4) となる。

3.1.3 2 ^分子 2 ^次反応

反応 A + B → P が速度則

−^dA

dt ^{= kAB}

に従うとする。A(t) = A(0) − x(t) とおく。すなわち、x(t) は反応によ る減少分である。この反応では１分子の A と１分子の B が反応するの で、B(t) = B(0) − x(t) と書ける。(これとは異なった分子数比の場合に は、x(t) に適当な係数を付ければ良い。) さらに A(0) = a, B(0) = b と おくと、

dx

dt = k(a − x)(b − x)

⇒

∫ 1

(a − x)(b − x)^{dx = k}

∫ dt

左辺の被積分関数を部分分数に分解する。 1

a − b

∫ ( 1 b − x ⁻

1 a − x

)

dx = k

∫ dt

ただしa ̸= b とする。a = b の場合は後で別に考える。t = 0 で x = 0 を 考慮して積分を実行すると、

1 a − b

( ln ^b

b − x ^{− ln} a a − x

)

= kt x について解けば、

x(t) = ^b(e

κt_{− 1)}

e^κt− (b/a)^, κ = k(a − b) となる。

長時間における振る舞いを調べよう。a < b のときは κ < 0 なので、 t → ∞ で e^κt→ 0 だから、x → a となる。

練習問題 a > b のときはどうなるか。

以上のように、t → ∞ で x は a と b の小さい方に漸近する。これは、初 期濃度の小さい方が無くなるまで反応が進むことを意味する。

濃度A(t), B(t) で表すと、 A(t) = a − x(t) = ^{(a − b)a}

a − be^−κt, B(t) = b − x(t) = ^{(b − a)b} b − ae^κt

対称性から、B(t) の式は A(t) で a と b を交換したものに等しい。また、 a と b の大小により κ の符号が決まり、分母の指数関数が t → ∞ で発散 したり消えたりする。

a = b の場合には、

∫ 1

(a − x)^{2dx = k}

∫ dt 1

a − x ⁻ 1

a ^{= kt ⇒ x(t) = a} akt 1 + akt

⇒ A(t) = B(t) = ^a 1 + akt よって、t → ∞ では x → a すなわち A(t) → 0 となる。

練習問題 (自己触媒反応) 反応 A → P が速度則 dA

dt ^{= −kAP}

を満たすとする。右辺に生成物P の濃度が現れている点が特異的である。 これは、生成物が反応を促進することを意味する。このような反応は自己触 媒反応と呼ばれる。A(0) = a, P (0) = p, A(t) = a − x(t), P (t) = p + x(t) とおいて、x に関する微分方程式を解け。

3.1.4 ^対向反応

次式のような対向反応

A⇌B

の正・逆反応ともに1 次反応とし、速度定数は kf, kb^{とする。このとき、}

速度式は

dA

dt ^{= −kfA + kbB} dB

dt ^{= kfA − kbB} となる。両辺の和をとると、

d

dt(A + B) = 0

なので、和A(t) + B(t) は一定である。これは、物質総量の保存を表す。 初期条件をA(0) = a ̸= 0, B(0) = 0 とする。A + B は不変なので、

A(t) + B(t) = A(0) + B(0) = a これよりB(t) を消去すると、

dA(t)

dt ^{= −kfA(t) + kb(a − A(t))}

= −(kf + kb) (

A(t) − ^kba kf + kb

)

ここで、_{kba/(kf + kb}) は定数だからその時間微分はゼロであることを利 用して、上式を

d dt

(

A(t) − ^kba kf + kb

)

= −(kf + kb) (

A(t) − ^kba kf + kb

)

と書き直す。括弧内を一つの関数と見れば、これは_{1 次反応と同じ形な} ので、解は指数関数

A(t) − ^kba kf + kb

= (

a − ^kba kf + kb

)

e^−(kf+kb)t

となる。整理すると

A(t) = ^A(0) kf + kb

[kfe^−(kf+kb)t+ kb

] (3.5)

ここで、A(0) = a を戻した。上式の t → ∞ の極限を A(∞) = ^kbA(0)

kf + kb

(3.6)

とすれば、

A(t) = (A(0) − A(∞))e^−(kf+kb)t+ A(∞) (3.7) となり、A(0) → A(∞) のつながりが見やすい形になる。

B(t) は B(t) = A(0) − A(t) から計算できる。

3.1.5 ^{変数分離型の一般論}

変数分離型

✓ ✏

dA

dt = f (A) g(t)

のように、右辺がA だけの関数 f (A) と変数 t だけの関数 g(t) の積で 表される場合を変数分離型と呼ぶ。この解は

∫ dA f (A) ⁼

∫

g(t)dt

の両辺を積分することで得られる。

✒ ✑

前節の反応速度論の例では、g(t) は定数だったので、上式の右辺は容 易に積分された。また、f (A) は A の多項式だったので、上式の左辺は 1/f (A) を直接または部分分数分解して積分した。

適当な変数変換により、変数分離型に帰着される場合がある。例えば、 dA

dt ^{= f (A/t)}

のように、右辺がA/t の関数で表される場合は、B = A/t とおけば、 dB

dt ⁼

f (B) − B t に帰着する。

練習問題 p, q, r を定数として、 dA

dt = f (pt + qA + r) の形の場合は、B = pt + qA + r とおけば、

dB

dt = p + qf (B) に帰着することを示せ。

3.2 ^{線形方程式}

3.2.1 ^逐次反応

次式のような逐次反応

A−→ B^ka −→ C^kb の各段階が1 次反応であるとすると、速度式は

dA

dt ^{= −kaA} dB

dt ^{= kaA − kbB} dC

dt ^{= kbB}

(3.8)

となる。両辺の和をとれば、 d

dt(A + B + C) = 0

これは、総物質量の保存である。初期条件をA(0) ̸= 0, B(0) = C(0) = 0 とする。

まず、A の速度式は単分子 1 次反応の場合と同じで、解は A(t) = A(0) e^−kat

となる。これをB の式に代入すれば、 dB

dt ^{+ kbB = kaA(0) e}

−kat

(3.9)

となる。左辺は、B(t) およびその微分の 1 次式になっている。次節で改 めて解説するように、このようなものをB(t) に関する線形微分方程式と 呼ぶ。上式の右辺をゼロとしたものを斉次、上式のように右辺が変数_{t の} 関数になっているものを非斉次の微分方程式と呼ぶ。

斉次の場合は_{B の代りに ˜}B と書くことにすると、 d ˜B

dt ^{+ kbB = 0˜} これは1 次反応と同じで、解は

B(t) = ˜˜ B(0) e^−kbt

となる。この斉次解を参考に、非斉次方程式の解_{B(t) を} B(t) = u(t) e^−kbt

とおいてみる。すると、式(3.9) は u(t) の微分方程式として、 du

dt ^{= kaA(0) e}

(kb−ka)t

を与える。これは直ちに積分できて u(t) − u(0) = ^kaA(0)

kb− ka

(e^(kb−ka)t− 1⁾

となる。B(0) = 0 より u(0) = 0 であるから、結局 B(t) = ^kaA(0)

kb− ka

(e^−kat− e^−kbt)

が得られる。次節で示すように、上記は非斉次の線形微分方程式を解く 一般的な手続きに従っている。

補足 _C(t)については、上の_B(t)を積分してもよいし、あるいは保存則 A(t) + B(t) + C(t) = A(0) + B(0) + C(0) = A(0)

から求めることも出来る。いずれにしても、

C(t) = A(0) (

1 + ¹ ka− kb

(kbe^−kat− kae^−kbt) )

となる。

3.2.2 1 階線形常微分方程式の一般論

f (t), g(t) が与えられたとき、A(t) に関する微分方程式で、A とその微 分に関して_{1 次であるような}

線形常微分方程式 _{(1 階, 非斉次)}

✓ ✏

dA

dt + f (t)A = g(t) (3.10)

✒ ✑

の形のものを線形微分方程式と呼ぶ。上式で、右辺をゼロとした d ˜A

dt ^{+ f (t) ˜A = 0 (3.11)} を斉次(homogeneous) 方程式、上式のように右辺に g(t) が残っているも のを非斉次(inhomogeneous) 方程式と呼ぶ。

斉次方程式の解A(t) が得られたとして、それを用いて˜

A(t) = u(t) ˜A(t) (3.12) とおくと、非斉次方程式(3.10) の左辺は

du dt^A(t)˜ と簡略化される。

実際、斉次方程式(3.11) は変数分離型なので次のように解かれる。

∫ d ˜A A˜ ^{= −}

∫

f (t)dt ⇒ A(t) = ˜^˜ A(0) e^{−∫0tf(τ )dτ} (3.13)

よって、非斉次方程式(3.10) は u(t) の微分方程式として du

dt ⁼ 1 A(0)˜ ^e

∫t

0^f(τ )dτ_g(t)

となる。これを積分すれば、 u(t) = u(0) + ¹

A(0)˜

∫ t 0

e^{∫0sf(τ )dτ}g(s)ds

両辺に式_{(3.13) の ˜}A(t) を掛けることにより、式 (3.10) の解として、 1 階線形常微分方程式の一般解

✓ ✏

A(t) = e^{−∫0tf(τ )dτ} (

A(0) +

∫ t 0

e^{∫0sf(τ )dτ}g(s)ds )

(3.14)

✒ ✑

を得る。前節の逐次反応の場合は_{f (t) = kb} (定数) だったので、右辺に相 当する積分は容易だった。

上式は、1 階線形常微分方程式 (3.10) の解を与える一般的な公式になっ ている。個々の問題において上の公式をそのままあてはめても良いが、前 節で逐次反応について実行したように、上記の手続きに沿って計算を進 めることも容易である。その要点は、(1) まず斉次方程式を解き、(2) そ の解に新たな関数(ここでは u(t)) を掛けたものを非斉次方程式へ代入す ると直接積分できる形に簡略化されるというものである。これは定数変 化法と呼ばれる。

3.2.3 ^定数係数 2 ^{階線形常微分方程式}

2 階の微分方程式の代表例は Newton 方程式であろう。

外力F (t) の下にある質量 m の質点の座標を y とすると、Newton の運 動方程式は

d²y dt^{2 =}

F (t) m

例えば、y を平衡位置からの変位として、ばね定数 k の調和振動子ポテン シャル_{V (y) = ky}²/2 の下にあるとき、

F = −ky

となる。質点がまさつ係数ζ の媒質中にあるときは、速度に比例する摩 擦抵抗

F = −ζ^dy dt

を受ける。これら以外の力をg(t) とおく。例えば、質点が電荷 q を持って おり、外部から電場_{E(t) = E0}cos ωt が掛かる場合には、

g(t) = qE(t) = qE0cos ωt (3.15) 等となる。以上を全て考慮すると、運動方程式は

d²y dt^{2 +}

ζ m

dy dt ⁺

k m^{y =}

g(t)

m ^(3.16)

となる。これは、定数係数で非斉次の2 階線形常微分方程式である。 定数係数2 階線形常微分方程式 (非斉次)

✓ ✏

d²y dt^{2 + a}

dy

dt + by = f (t) (3.17)

✒ ✑

以後、簡単のため上式を

y^′′+ ay^′+ by = f (t) (3.18) と書く。対応する斉次方程式は

y^′′+ ay^′+ by = 0 (3.19) となる。以下に見るように、斉次方程式の一般解は比較的容易に得られ る。2 階の微分方程式なので、一般解は任意定数を 2 つ含むことになる。 補足 証明は省略するが、一般に_n階の斉次線形常微分方程式には、一次独立 な解を_n個作ることができ、それより多く作ることはできない。一般解はそれ らの一次結合で表され、その係数が_n個の任意定数となる。例えば、₂階の微分 方程式を解くことは₂回積分することなので、積分定数が₂つ出る。_Newton方 程式の場合であれば、₂定数を定めることは、初期座標と初期速度の₂つを定め ることに相当する。

いま、式(3.17) の一つの特解 y0(t) が分ったとする。これは y^′′₀ + ay₀^′ + by0 = f (t)

を満たす。これを式(3.17) から引くと、z(t) ≡ y(t) − y0(t) で定義される z(t) は斉次方程式 (3.19) を満たすことが分かる。

z^′′+ az^′+ bz = 0

よって、この一般解が任意定数_{c1, c2}を含むz(t) として求まったとすると、 y(t) = y₀(t) + z(t) (3.20) が非斉次方程式(3.17) の一般解となる。このように、

✓ ✏

非斉次方程式(3.17) の一般解は、斉次方程式 (3.19) の一般解と、非斉 次方程式の特解の和で表される。

✒ ✑

補足 上の議論からも直ちに分かるように、上記はより一般的な_n階の線形微 分方程式₍定数係数でなくてもよい₎についても言える。

定数係数斉次方程式の解法は、以下が標準的である。微分演算子を記 号D = d/dt で表すことにすると、式 (3.19) は

(D²+ aD + b)y = 0 (3.21) と書かれる。ここで、変数λ に関する二次方程式

λ²+ aλ + b = 0 (3.22) を考える。これを、この微分方程式の特性方程式と呼ぶ。この解を_{λ1, λ2} とすると、式_{(3.21) は}

(D − λ1)(D − λ2)y = 0 と変形される。そこで、

(D − λ1)y1 = 0, (D − λ2)y2 = 0

を満たす_{y1, y2}を求めれば、これらの任意の線形結合は式_{(3.21) を満た} す。上式は

(D − λi)y = 0 ⇒ y^′− λiy = 0 ⇒ yi(t) = cie^λit (3.23)

のように、定数係数の1 階斉次線形微分方程式として容易に解ける。 λ1 ̸= λ2^{のときは、式}(3.23) が独立な二つの解を与えるので、

z(t) = c1e^λ1t+ c2e^λ2t (3.24) が斉次方程式の一般解として十分である。

λ1 = λ2 (重解) のときは、式 (3.21) は

(D − λ1)²y = 0 (3.25) となり、上のような考察では独立な解は一つしか得られない。しかし、式 (3.25) に左から e^{−λ1tを掛け、恒等式∗}

e^−λ1t(D − λ1)y = D(e^−λ1ty) (3.26) を2 回用いることにより

e^−λ1t(D − λ1)²y = D(e^−λ1t(D − λ1)y) = D²(e^−λ1ty) = 0 が得られる。2 回微分したものが 0 となるのは 1 次式だから、

⇒ e^−λ1ty = c1+ c2t ⇒ z(t) = (c1+ c2t)e^λ1t (3.27) が一般解となる。

非斉次方程式(3.17) の特解を求めるにも、式 (3.26) を利用できる。 (D − λ1)(D − λ2)y = f

の両辺に左から_e^−λ1tを掛け、式(3.26) を用いれば、 D(e^−λ1t(D − λ2)y) = e^−λ1tf 積分して

e^−λ1t(D − λ2)y =

∫ t 0

e^−λ1sf (s)ds + C

特解を求めればよいので、C = 0 とする。両辺に e^{λ1t と}e^{−λ2t を掛けて再} び式(3.26) を用いれば、

D(e^−λ2ty₀) = e^(λ1−λ2)t

∫ t 0

e^−λ1sf (s)ds 積分して両辺に_e^λ2t を掛ければ、特解_{y0(t) は}

y0(t) = e^λ2t

∫ t 0

e^{(λ1−λ2)τ}

∫ τ 0

e^−λ1sf (s)dsdτ (3.28) と求まる。

∗これは_f^′_{g + f g}^′ _{= (f g)}^′に他ならない。

まとめ 式(3.17) のような非斉次微分方程式を解くには、その特解を一 つ見つけ、対応する斉次方程式の一般解との和をとれば、それが非斉次 方程式の一般解となっている。斉次方程式を解くには、対応する特性方 程式(3.22) を考える。それが重解を持たないならば、指数関数の 1 次結 合(3.24) が一般解となる。重解を持つときには、式 (3.27) が一般解とな る。非斉次方程式の特解は、式(3.28) とすればよい。

3.3 ^{連立線形方程式}

3.3.1 ^逐次反応 ( ^再考 )

3.2.1 節で扱った逐次反応

A−→ B^ka −→ C^kb

を再考する。C(t) は保存則から求まるので、A(t) と B(t) のみを扱えばよ い。速度方程式_{(3.8) は}

d dt

[A B ]

=

[−ka 0 ka −kb

] [A B ]

≡ K [A

B ]

(3.29)

のように行列形式で書くことができる。最右辺で行列_{K を定義した。} 上式は、1 次反応の微分方程式

d

dt^{A = −kA ⇒} A(t) = A(0)e^−kt

を行列形式に拡張した形になっている。これと同様に、式_{(3.29) の解を} [A(t)

B(t) ]

= e^Kt

[A(0) B(0) ]

(3.30)

と書けるかどうか検討する。行列の指数関数は e^Kt = I + Kt +^(Kt)

2

2 ⁺ (Kt)³

3! ^{+ · · · =}

∞

∑

n=0

(Kt)ⁿ n! で定義される。(1.3.3 節、定義 1.7) これを t で微分すれば

d dt^e

Kt _{= K + K}2_{t +} ^K3t2

2! ⁺ K⁴t³

3! ^{+ · · · =}

∞

∑

n=1

Kⁿtⁿ⁻¹

(n − 1)! ^{= Ke}

Kt

であり、確かに式(3.30) は式 (3.29) の解となっている。

よって、あとは_e^Ktを計算すればよいことになる。それには、_{1.3.3 節で} 見たように、K を対角化すればよい。今の例では、それは容易に可能で、

T⁻¹KT =

[−ka 0 0 −kb

]

T =

[kb− ka 0 ka 1 ]

, T⁻¹ = ¹ kb − ka

[ 1 0

−ka kb− ka

]

となる。よって、 e^Kt = T

[e^−kat 0 0 e^−kbt

] T⁻¹

= ¹ kb− ka

[ (kb− ka)e^−kat 0 ka(e^−kat− e^−kbt) (kb − ka)e^−kbt

]

3.2.1 節と同様、初期条件として A(0) ̸= 0, B(0) = 0 とすれば、 [A(t)

B(t) ]

= e^KtA(0) [1

0 ]

= ^A(0) kb− ka

[ (kb− ka)e^−kat ka(e^−kat− e^−kbt)

]

となる。これは3.2.1 節の結果と一致している。

3.3.2 ^{高階微分方程式}

3.2.3 節で考察した 2 階の定数係数線形斉次方程式 (3.19) y^′′+ ay^′+ by = 0

を再考する。_{z ≡ y}^′により新たな関数z(t) を定義すれば、上式は 1 階の 連立微分方程式

d dt

[y z ]

=

[ 0 1

−b −a ] [y

z ]

に帰着する。右辺の2 × 2 行列を K とおくと、微分方程式の解は [y(t)

z(t) ]

= e^Kt [y(0)

z(0) ]

となる。上式のK の固有値 λ は次式から求まる。 

−λ 1

−b −a − λ     

= λ²+ aλ + b = 0

これは、与えられた微分方程式の係数から決まる特性方程式_{(3.22) に他} ならない。

3 階の微分方程式

y^′′′+ ay^′′+ by^′+ cy = 0 であれば、_{z = y}^′ および_{w = z}^′ _{= y}^′′ を導入すれば、

d dt





 y z w





⁼







0 1 0 0 0 1

−c −b −a











 y z w







となり、やはり連立1 階微分方程式に帰着される。

補足 _n階微分方程式であれば、新たな変数を_{n − 1}個導入することで、_n変数 に関する₁階連立微分方程式に帰着される。例えば、_n階微分を_y⁽ⁿ⁾と書くこ とにして

y⁽ⁿ⁾+ an−1y⁽ⁿ⁻¹⁾+ · · · + a1y^′+ a0y = 0 が与えられたとする。

y₀= y, y₁ = y^′, y₂ = y^′′= y₁^′, · · · , yn= y⁽ⁿ⁾ のように_y0, · · · , yn^{を定義すれば、}

d dt









 y₀ y₁ ... yn











=











0 1 0 · · · 0

0 0 1 · · · 0

... ^...

−a₀ −a₁ −a₂ · · · −a_n−1



















 y₀ y₁ ... yn











となる。右辺の_{n × n}行列を_Aとすると、その固有値_λを決める特性方程式は λⁿ+ a_n−1λⁿ⁻¹+ · · · + a₁λ + a₀ = 0 (3.31) である。

n = 3の場合に確認してみる。固有値を決める方程式は|A − λI| = 0^すなわち

−λ 1 0

0 −λ 1

−a₀ −a₁ −a₂− λ      

= 0

である。左辺の₃列目を_λ倍して₂列目に加えて真中の_−λを消し、続いて₂列 目を_λ倍して₁列目に加えて左上の_−λを消せば、

0 1 0

0 0 1

−a0− a1λ − a2λ²− λ³ −a1− a2λ − λ² −a2− λ      

= 0

となる。₁列目について余因子展開すれば、右上の_{2 × 2}部分は単位行列なので λ³+ a₂λ²+ a₁λ + a₀ = 0

が得られる。同様の方法により、一般の_nについて式_(3.31)が導かれる。

3.4 Laplace ^変換

線形常微分方程式の解法の一つにLaplace 変換による方法がある。

✓ ✏

定義 3.1 t ≥ 0 で定義された関数 f (t) の Laplace 変換を L{f (t)} = F (s) =

∫ ∞ 0

e^−stf (t)dt (3.32)

で定義する。s は上の積分が収束する範囲で考える。

✒ ✑

補足 上式の積分がある_s0で収束すれば、_{s > s0}なるすべての_sで収束するこ とが示される。積分が収束する最小の_s0を収束座標という。

よく使う関数の変換をいくつか挙げると、

✓ ✏

例 3.1 a を定数として、 L{a} = ^a

s^{, L{e}

at_{} =} ¹

s − a^{, L{t}

n_{} =} ^n!

s^{n+1 (3.33)} L{sin at} = ^a

s²+ a^2, L{cos at} = ^s

s²+ a^{2 (3.34)} L{sinh at} = ^a

s²− a², L{cosh at} = ^s

s²− a^{2 (3.35)}

✒ ✑

練習問題 式(3.33)–(3.35) を確かめよ^†。

(^ヒント: (3.34)^は(3.33)^第2^式でa → ia^{とし、実部と虚部。}) Laplace 変換の特徴は、次式のような微分の変換にある。

✓ ✏

定理 3.1 微分の Laplace 変換

L{f^′(t)} = sF (s) − f (0) (3.36) L{f^′′(t)} = s²F (s) − sf (0) − f^′(0) (3.37)

✒ ✑

ただし、_f^′_{, f}^′′は1 階および 2 階の微分を表す。上式のように、Laplace 変 換により微分演算は変数s の掛け算に置き換わる。

練習問題 上の_{2 式を確かめよ。}

微分方程式への応用

以上を利用すると、微分方程式が簡単に解けることがある。 例 3.2 定数係数 1 階線形斉次型

f^′(t) + af (t) = 0 Laplace 変換すると、

sF (s) − f (0) + aF (s) = 0 F (s) について解くと、

F (s) = f (0) ¹ s + a

右辺は式(3.33) の 2 番目の式の形なので、直ちに逆変換できて、 f (t) = f (0)e^−at

が得られる。

このように、式(3.36) や (3.37) を用いて微分方程式を Laplace 変換した 式は、F (s) と s の単純な式になる。それを F (s) について解いたものを、 必要ならば部分分数分解し、式(3.33)–(3.35) と見比べて逆変換が見つか れば、微分方程式は直ちに解かれる。

†双曲線関数 sinh x = (e^x− e^−x)/2, cosh x = (e^x+ e^−x)/2

練習問題 式(3.37) と (3.34) を利用して、 f^′′(t) = −a²f (t)

を解け。(調和振動子の Newton 方程式はこの形である。)

次の定理は定義から容易に導かれ、応用において有用となる。

✓ ✏

定理 3.2 F (s) を f (t) の Laplace 変換とすると、F (s + a) は e^−atf (t) のLaplace 変換である。

✒ ✑

証明 式_(3.32)から F (s + a) =

∫ ∞ 0

e^−(s+a)tf (t)dt =

∫ ∞ 0

e^−st{e^−atf (t)}dt = L{e^−atf (t)}

例 _3.3

F (s) = ¹

(s + a)² ⇒ f (t) = te^−at F (s) = ^{s − b}

(s − b)²+ a² ⇒ f (t) = e^btcos at

練習問題 次の微分方程式を解け。

f^′′(t) + f^′(t) − 2f (t) = −3e^−2t, (f (0) = 0, f^′(0) = 1)

3.4.1 ^畳み込み (convolution)

次式の性質は特に重要である。

✓ ✏

定理 3.3 Laplace 変換の積は、畳み込み積分の Laplace 変換である。 L{f (t)} · L{g(t)} = L

{∫ t 0

f (t − τ )g(τ )dτ }

(3.38)

✒ ✑

証明 右辺は定義により、

∫ ∞ 0

dte^−st

∫ t 0

f (t − τ )g(τ )dτ

変数を_tと_τ の組から_τ と_{u ≡ t − τ} に変換すれば、積分は分離されて

∫ ∞ 0

dτ

∫ ∞ 0

due^{−s(τ +u)}f (u)g(τ ) = L{f }L{g}

注意 上の変数変換は (u

τ )

= ^{(1 −1} 0 1

) ( t τ

)

と書けるので、変換の_Jacobian は_{J =}

1 −1 0 1

= 1 である。また、積分範囲は0 ≤ t, 0 ≤ τ ≤ t^なので、(τ, t) 平面の第₁象限で直線_{t = τ}の上側となる。この領域を覆うのに、直線_{t = τ + u} を、_t軸の切片である_uをパラメータとして動かすと考えれば、0 ≤ τ, 0 ≤ u^が 得られる。

例 3.4 定数係数 1 階線形方程式

f^′+ kf = g(t)

を考える。f (t), g(t) の Laplace 変換を F (s), G(s) とおく。両辺を Laplace 変換すると

sF (s) − f (0) + kF (s) = G(s) F (s) について解けば

F (s) = ¹

s + k^{f (0) +} 1

s + k^G(s)

右辺第_{1 項は e}^−ktのLaplace 変換に比例している。第 2 項は e^−ktのLaplace 変換とg(t) の Laplace 変換の積であるから、式 (3.38) が適用される。よっ て、上式の逆変換は

f (t) = e^−ktf (0) +

∫ t 0

e^{−k(t−τ )}g(τ )dτ

となる。これは、一般公式(3.14) と一致している。

補足 定数係数₂階線形方程式の一般形

f^′′+ af^′+ bf = g(t)

に応用してみる。_Laplace変換は、

(s²F (s) − sf (0) − f^′(0)) + a(sF (s) − f (0)) + bF (s) = G(s)

⇒ (s²+ as + b)F (s) = (s + a)f (0) + f^′(0) + G(s)

F (s)に掛かる因子は、以前に考察した特性方程式_(3.22)の左辺に他ならない。

その解_{λ1, λ2}により、 F (s) = ^{s + a}

(s − λ1)(s − λ2)^{f (0) +}

1

(s − λ1)(s − λ2)^(f

′_{(0) + G(s)}⁾

λ₁+ λ₂ = −a^{を用いると、} s + a

(s − λ1)(s − λ2) ⁼ 1 λ1− λ2

( −λ₂ s − λ1

+ ^λ1 s − λ2

)

1

(s − λ₁)(s − λ₂) ⁼ 1 λ₁− λ₂

( 1

s − λ₁ ⁻ 1 s − λ₂

)

と部分分数分解される。これらの逆変換は指数関数となるので、

f (t) = ^{f (0)} λ₁− λ₂

(

−λ₂e^λ1t+ λ₁e^λ2t)+ ^f

′₍₀₎

λ₁− λ₂ (

e^λ1t− e^λ2t)

+ ¹

λ1− λ2

∫ t 0

(e^{λ1(t−τ )}− e^{λ2(t−τ ))}g(τ )dτ

1^行目がf (0)^とf^′(0)を任意定数とする斉次方程式の一般解、₂行目が非斉次方 程式の特解になっている。

応用例_{: 結合振動子}

x(t), y(t) が次の連立運動方程式を満たすとする^‡。

¨

x(t) = −^{∂U (x)}

∂x ^{− cy (3.39)}

¨

y(t) = −ω²y − cx (3.40)

‡簡単のため、自由度x, y ともに質量は 1 とした。あるいは、質量加重座標に変換し たと考えてもよい。

これは、ポテンシャルエネルギー

V (x, y) = U (x) +¹ 2^ω

2_y2_{+ cxy}

から導かれる。第1 項は x に関するポテンシャル、第 2 項は y に関する調 和振動子ポテンシャル、第3 項は x, y 間の相互作用を表し、c は結合の強 さを表す定数とする。y の運動方程式を Laplace 変換すると

Y (s) = ^s

s²+ ω^{2y(0) +} 1

s²+ ω^{2 ˙y(0) −} c

s²+ ω^2X(s) となるから、逆変換して

y(t) = y(0) cos ωt + ^˙y(0)

ω ^{sin ωt −} c ω

∫ t 0

sin ω(t − τ )x(τ )dτ

これを式(3.39) に代入すれば、y(t) は表から消えて x(t) だけを含む運動 方程式となる。上式の最後の積分を部分積分して

∫ t 0

sin ω(t−τ )x(τ )dτ = ¹

ω [cos ω(t − τ )x(τ )]^t_τ=0−¹ ω

∫ t 0

cos ω(t−τ ) ˙x(τ )dτ

と書き換えると、x(t) の運動方程式は、

¨

x(t) = −^{∂U (x)}

∂x ⁻ c² ω²

∫ t 0

cos ω(t − τ ) ˙x(τ )dτ + ^c

2

ω²x(t) + R(t) となる^§。右辺第2 項の積分は速度 ˙x に掛かる摩擦の形をしており、時刻 0 から t までの履歴あるいは記憶に依存している^¶。

3.5 ^{速度場と安定性}

3.5.1 ^{ロジスティ ック方程式}

dy

dt ^{= ay ⇒} y(t) = y(0)e^at

§x(t) を含まない項を R(t) = −cy(0) cos ωt − _ω^c ˙y(0) sin ωt −_ω^c22x(0) cos ωt にまとめ た。

¶Brown 運動を記述する Langevin 方程式に、摩擦の記憶効果を取り入れた「一般化 Langevin 方程式」は、上記と同様にして導かれる。

a > 0 ならば y は指数関数的に増大する。例えば、細胞分裂による増殖は この形の方程式に従う。このとき_{a は増殖率である。}

例えば、細胞数が増加するにつれ環境が悪化し、増殖率が減少すると して、これをa = b − gy と表すことにする。ただし、b > 0, g > 0 とす る。このときy(t) の変化を記述する微分方程式は

dy

dt = (b − gy)y (3.41) となる。これは、ロジスティック方程式と呼ばれる。

練習問題 式(3.41) を解け。

上の解より、t → ∞ で y → b/g に漸近することが分る。

速度場

dy/dt の符号が y の増減の方向を表すことに着目すれば、微分方程式を 解かなくとも解の大体の様子を調べることができる。

まず、式(3.41) の右辺は y = 0 と y = b/g で 0 となる。このとき y の変 化の速度が0 となるので、これらは停留点である。

次に、式(3.41) の右辺は 0 < y < b/g で正、y < 0 と y > b/g で負であ る。よって、前者の領域でy は増加し、後者の領域では減少する。これ らより、y = 0 の両側ではその点から離れていく向きにあり、y = b/g の 両側では引き込まれる向きにあることが分る。すなわち、_{y = 0 は不安定} 点、y = b/g は安定点である。

3.5.2 Lotka-Volterra ^モデル

A と B の変化が次式で与えられるとする。 dA

dt ^{= (p1− q1B)A} dB

dt ^{= (−p2+ q2A)B}

(3.42)

ただし、_{p1, q1, p2, q2}は全て正の定数とする。また、A > 0, B > 0 とす る。これはLotka-Volterra モデルと呼ばれるものの一種で、以下に説明 するように、餌食A と捕食者 B が共存する系のモデルと見なせる。

まず、第一式において_−q1B の項が無いとすると、A は e^{p1tのように増} 大する。しかし、_−q1B の項があるためにこの増大は抑制される。すなわ ち、B の増加は A の増殖を抑制する。第二式について同様に考えると、B は_{A が無いと e}^−p2tのように減衰するが、A の存在がこの減衰を抑制する。

この方程式系の安定点は、_{(p1− q1}B)A = 0 かつ (−p2+ q2A)B = 0 よ り、(A, B) = (0, 0) および (p2/q2, p1/q1) の二点である。また、

0 < B < p1/q1 ^で

dA

dt ^{> 0, B > p1/q1 で} dA

dt ^{< 0} 0 < A < p2/q2 ^で

dB

dt ^{< 0, A > p2/q2 で} dB

dt ^{> 0}

であるから、速度場の向きは図XX のようになる。この図から、軌道は 安定点_{(p2/q2, p1/q1}) の周りを周回すると予想される。すなわち、餌食 A と捕食者B は、一方が増えたら他方が減るという具合に入れ違いつつ振 動する。

補足 式_(3.42)の第一式より、 1 A

dA dt ⁼

d(ln A)

dt ^{= p1− q1B}

周回運動の周期を_Tとすると、ある時刻_tから一周期積分して 0 = ln A(t + T ) − ln A(t) = p₁T − q₁

∫ _t+T

t

B(τ )dτ

これより、_tから_{t + T} の間の_Bの平均個体数が次式のように求まる。 1

T

∫ t+T t

B(τ )dτ = ^p1 q1

同様に、_(3.42)の第二式より、_Aの平均個体数は式

1 T

∫ t+T t

A(τ )dτ = ^p2 q₂

となる。いま、何らかの環境の変化が生じて、式_(3.42)が dA

dt ^{= (p1− q1B)A + αA} dB

dt ^{= (−p2+ q2A)B + βB}

のように変ったとする。α > 0, β > 0^{とすれば、これは}A, B^{両者にとって好都} 合な変化である。式の上では、_p1が_{p1+ α}に、_p2が_{p2− β}に変ったことに相 当する。よって、平均個体数については

p₁+ α q1

> ^p1 q1

, ^{p2− β} q2

< ^p2 q1

となり、_Bについては増え、_Aについては減る。すなわち、両者にとって好都合な変 化であっても、餌食_Aの平均個体数はむしろ減少してしまう。逆に、α < 0, β < 0 とすれば、両者にとって不利な環境変化であり、この場合には捕食者の平均個 体数が減少する。これは、例えば害虫を駆除するために農薬を散布したところ、 害虫の天敵の方に打撃を与えてしまい、害虫がむしろ増加してしまうような現 象を説明する。