📔

Botter自主ゼミノート 2.9 ウィーナ過程(2)

2022/12/11に公開

やること

https://www.amazon.co.jp/dp/4254209444
を読んで、確率微分方程式による最適化問題を解けるようになることです。

これまでのBotter自主ゼミノート

Botter自主ゼミノート 1.2 確定システムの制御の回顧
Botter自主ゼミノート 2.1 確率過程とは?, 2.2 確率過程の数学的表現
Botter自主ゼミノート 1.2 数式導出
Botter自主ゼミノート 2.3 確率モーメント
Botter自主ゼミノート 2.4 確率過程の分類
Botter自主ゼミノート 2.5 エルゴード性
Botter自主ゼミノート 2.6 確率過程の周波数表現
Botter自主ゼミノート 2.7 マルコフ過程
Botter自主ゼミノート 2.8 正規型確率過程
Botter自主ゼミノート 2.9 ウィーナ過程(1)

2.9 ウィーナ過程

ブラウン運動とウィーナ過程の関係性について

ウィーナ過程は物理現象としてよく知られたブラウン運動を数学的に記述したものです。教科書では、ブラウン運動過程が式(2.58)で与えられるような正規型過程であることを、液体中の粒子の運動方程式を考えることで示しています。

p(t,w) = \frac{1}{\sqrt{2\pi t}\sigma}\exp\left(-\frac{w^2}{2\sigma^2 t}\right) \tag{2.58}

まず、液体中の粒子は質量mの質点とし、その位置をw(t)(ただし一次元)とします。この粒子の速度をu(t)=\frac{dw(t)}{dt}とします。

すると、ニュートンの第二法則(Wikipedia)より、式(2.70)が成り立ちます。

m\frac{du(t)}{dt} = -\xi u(t) + \gamma_0(t)\quad (t \ge 0) \tag{2.70}

式(2.70)の右辺第一項は液体による粘性抵抗力で、\xiはその係数です。また、\gamma_0(t)は粒子に働く外力です。粒子が不規則に動く原因として、外力\gamma_0(t)が確率的に不規則に動くものと考えます。外力\gamma_0(t)は粒子に衝突する液体分子による力であると考え、白色雑音と仮定します。(\gamma_0(t)\gamma_0(\tau) (t \ne \tau)は互いに相関を持たないものとします)

式(2.70)のような運動方程式をランジュヴァン方程式といい、ブラウン運動の解析モデルとして用いられています。この式を書き直すと式(2.71)のようになります。

\frac{du(t)}{dt} + \beta u(t) = \gamma(t) \tag{2.71}

ただし、\beta = \frac{\xi}{m}, \gamma(t) = \frac{\gamma_0(t)}{m}とします。また、\gamma(t)は平均値と相関関数が式(2.72)のように与えられる確率過程であると仮定します。ここで、Dは正の定数、\delta(t)はディラックのデルタ関数です。これは、液体の分子が粒子に不規則に衝突し、その相関は\delta(t-s)と仮定するのが適切という考え方から来ています。

\mathcal{E}\{ \gamma(t) \} = 0, \quad \mathcal{E}\{ \gamma(t)\gamma(s)\} = 2D\delta(t-s) \tag{2.72}

式(2.71)は一階線形常微分方程式なので、これを解くと式(2.73)となります。

u(t) = u_0 e^{-\beta t} + \int_{0}^{t} e^{-\beta(t-s)} \gamma(s) ds \tag{2.73}

ここで、u(t) = \frac{dw(t)}{dt}, w(0) = 0に留意して式(2.73)を積分します。

\begin{aligned} w(t) &= \int_{0}^{t} u(s) ds\\ &= \int_{0}^{t} \left[ u_0 e^{-\beta s} + \int_{0}^{s} e^{-\beta (s-\tau)} \gamma(\tau) d\tau \right] ds \tag{2.74} \end{aligned}

式(2.74)を変形して以下の式(2.75)を得ます。

\begin{aligned} w(t) - u_0 \int_{0}^{t} e^{-\beta s} ds &= \int_{0}^{t} e^{-\beta s} \left[\int_{0}^{s} e^{\beta \tau} \gamma(\tau) d\tau \right] ds \tag{2.75} \end{aligned}

左辺の積分を実行して…

\begin{aligned} w(t) - \frac{1}{\beta} u_0 (1 - e^{-\beta t}) &= \int_{0}^{t} e^{-\beta s} \left[\int_{0}^{s} e^{\beta \tau} \gamma(\tau) d\tau \right] ds \end{aligned}

ここで教科書では右辺の部分積分を実行して以下の式(2.76)を得るとされています。

\begin{aligned} w(t) - \frac{1}{\beta} u_0 (1 - e^{-\beta t}) &= \frac{1}{\beta} \int_{0}^{t} (1 - e^{-\beta (t-s)}) \gamma(s) ds \tag{2.76} \end{aligned}

式(2.76)が得られたので、W(t) := w(t) - \frac{1}{\beta}u_0(1-e^{-\beta t})の平均値と分散を計算すると…

\begin{aligned} \mathcal{E}\{W(t)\} &= \frac{1}{\beta} \int_{0}^{t} \left( 1-e^{-\beta(t-s)}\right) \mathcal{E}\{\gamma(s)\} ds \\ &= \frac{1}{\beta} \int_{0}^{t} \left( 1-e^{-\beta(t-s)}\right) \cdot 0 \: ds \\ &= 0 \end{aligned}
\begin{aligned} \mathcal{E}\{W^2(t)\} &= \frac{1}{\beta^2} \int_{0}^{t} \int_{0}^{t} \left( 1-e^{-\beta(t-s_1)}\right) \left( 1-e^{-\beta(t-s_2)}\right) \mathcal{E}\{\gamma(s_1)\gamma(s_2)\} ds_1 ds_2 \\ &= \frac{1}{\beta^2} \int_{0}^{t} \int_{0}^{t} \left( 1-e^{-\beta(t-s_1)}\right) \left( 1-e^{-\beta(t-s_2)}\right) 2D \delta(s_1 - s_2) ds_1 ds_2 \\ \end{aligned}

\delta(s_1 - s_2)s_1 = s_2の時のみ値を持ち、その積分した値は1であるので…

\begin{aligned} &= \frac{1}{\beta^2} \int_{0}^{t} \left( 1-e^{-\beta(t-s)}\right)^2 2D \: ds \\ &= \frac{2D}{\beta^2} \int_{0}^{t} 1 - 2e^{-\beta t} e^{\beta s} + e^{-2 \beta t} e^{2 \beta s} ds \\ &= \frac{2D}{\beta^2} \left[s - \frac{2}{\beta} e^{-\beta t} e^{\beta s} + \frac{1}{2 \beta} e^{-2 \beta t} e^{2 \beta s} \right]_{0}^{t}\\ &= \frac{2D}{\beta^2} \left(t - \frac{2}{\beta} + \frac{1}{2 \beta} \right) - \left(0 - \frac{2}{\beta} e^{-\beta t} + \frac{1}{2 \beta} e^{-2 \beta t} \right) \\ &= \frac{2D}{\beta^2} \left(t - \frac{3}{2 \beta} + \frac{2}{\beta} e^{-\beta t} - \frac{1}{2 \beta} e^{-2 \beta t} \right) \\ \mathcal{E}\{W^2(t)\} &= \frac{D}{\beta^3} \left(2 \beta t - 3 + 4 e^{-\beta t} - e^{-2 \beta t} \right) \\ \end{aligned}

となります。

w(t) = \int_{0}^{t} u(s) dsであり、u(t)は式(2.73)のように白色雑音の積分和で与えられるので、w(t)は中心極限定理によって正規分布します。そのことから、W(t)過程の分布は以下のようになります。

\begin{aligned} p(t, w) &= \frac{1}{\sqrt{2 \pi \sigma^2}} \exp{-\frac{(w-\mu)^2}{2 \sigma^2}} \\ &= \sqrt{ \frac{\beta^3}{2 \pi D (2 \beta t - 3 + 4 e^{-\beta t} - e^{-2 \beta t} ) } } \exp \left( -\frac{\beta^3 |w - \frac{1}{\beta} u_0(1-e^{-\beta t})|^2}{2D(2\beta t - 3 + 4e^{-\beta t}-e^{-2\beta t}} \right) \tag{2.77}\\ \end{aligned}

ここで、t \gg \frac{1}{\beta}となる十分長い時刻tにおいては、指数関数e^{-\beta t}, e^{-2 \beta t}は微小で、定数項は2\beta tに対して無視できます。これに基づいて上式を変形すると…。

\begin{aligned} p(t, w) &= \sqrt{ \frac{\beta^2}{4 \pi D t} } \exp \left( -\frac{\beta^2 |w - \frac{1}{\beta} u_0|^2}{4 D t} \right) \tag{2.77}\\ \end{aligned}

上式で得られたW(t)の分布から、W(t) \sim \frac{1}{\sqrt{t}}であることが分かります。であればw(t) \sim \frac{1}{\sqrt{t}}であるので、\frac{1}{\beta}u_0(1-e^{\beta t})w(t)と比較すると無視できることが分かります。

これを利用して式(2.77)は式(2.78)に近似できます。

\begin{aligned} p(t, w) &\simeq \sqrt{ \frac{\beta^2}{4 \pi D t} } \exp \left( -\frac{\beta^2 w^2}{4 D t} \right) \tag{2.78}\\ \end{aligned}

ここで、D_0 = \frac{D}{\beta^2}とすると、以下の式(2.79)となります。

\begin{aligned} p(t, w) &\simeq \frac{1}{\sqrt{4 \pi D_0 t}} \exp \left( -\frac{w^2}{4 D_0 t} \right) \tag{2.79}\\ \end{aligned}

ここで、\sigma^2 = 2D_0 \: (= \frac{2D}{\beta^2} = \frac{2Dm^2}{\xi^2})とすると、式(2.79)は前回のノートで出てきた、ウィーナ過程の増分は正規分布するものとする、という定義から求めた確率密度分数式(2.58)と同じになります。

p(t,w) = \frac{1}{\sqrt{2\pi t}\sigma}\exp\left(-\frac{w^2}{2\sigma^2 t}\right) \tag{2.58}

ブラウン運動の外力Dの大きさについて

つぎに、外力の大きさD (あるいはD_0)を求めてみます。式(2.73)から、U(t) := u(t) - u_0 e^{-\beta t}の平均と分散を求めます。

教科書によると計算はW(t)の平均と分散の計算と同様だそうです。

\begin{aligned} \mathcal{E}\{U(t)\} &= \int_{0}^{t} e^{-\beta(t-s)} \mathcal{E}\{\gamma(s)\} ds\\ \mathcal{E}\{U(t)\} &= \int_{0}^{t} e^{-\beta(t-s)} \cdot 0 ds\\ \mathcal{E}\{U(t)\} &= 0 \end{aligned}
\begin{aligned} \mathcal{E}\{U^2(t)\} &= \int_{0}^{t} \int_{0}^{t} e^{-\beta(t-s_1)} e^{-\beta(t-s_2)} \mathcal{E}\{\gamma(s_1) \gamma(s_2)\} ds_1 ds_2\\ &= \int_{0}^{t} \int_{0}^{t} e^{-\beta(t-s_1)} e^{-\beta(t-s_2)} 2D \delta(s_1 - s_2) ds_1 ds_2\\ &= 2D \int_{0}^{t} e^{-2\beta(t-s)} ds\\ &= 2D e^{-2 \beta t} \int_{0}^{t} e^{2\beta s} ds\\ &= 2D e^{-2 \beta t} \frac{1}{2 \beta} \left[ e^{2\beta s} \right]_{0}^{t}\\ &= \frac{D}{\beta} e^{-2 \beta t} (e^{2\beta t} - 1)\\ &= \frac{D}{\beta} (1 - e^{-2 \beta t}) \end{aligned}

中心極限定理により、U(t)の分布も正規型となるので

\begin{aligned} p(t, u) &= \frac{1}{\sqrt{2 \pi \sigma^2}} \exp{-\frac{(u-\mu)^2}{2 \sigma^2}} \\ &= \sqrt{ \frac{\beta}{2 \pi D(1-e^{-2\beta t})} } \exp\left( -\frac{\beta |u-u_0 e^{-\beta t}|^2}{2D(1-e^{-2 \beta t})}\right) \tag{2.80} \end{aligned}

W(t)の時と同様にt \to \inftyの安定状態での分布を考えると、

lim_{t \to \infty} p(t,u) = \sqrt{\frac{\beta}{2 \pi D}} \exp\left(-\frac{\beta u^2}{2D}\right) \tag{2.81}

となる。

教科書によると、統計熱力学で知られていることとして、熱平衡状態における粒子の速度分布は式(2.82)のマックスウェル分布で与えられるそうです。ここで、kはボルツマン定数, Tは絶対温度だそうです。

p(t) \propto e^{-\frac{mu^2}{2kT}} \tag{2.82}

式(2.81)と式(2.82)は一致しなければいけないので、両者を見比べて、式(2.83)を得ます。

\begin{aligned} D &= \frac{\beta k T}{m} = \frac{\xi k T}{m^2}\\ D_0 &= \frac{kT}{\beta m} = \frac{kT}{\xi} \tag{2.83} \end{aligned}

もともと\xiは粘性抵抗力として導入されていましたが、ストークスの法則により粘性抵抗力は-6\pi a \eta u (\eta: 液体の粘性係数)となることから、\xi = 6\pi a \eta uで与えられます。

ここまでの計算により、白色雑音\sigma(t)の強さを表すものとして導入された係数Dは、式(2.83)のような物理的意味をもつことが分かりました。また、Dは実験によって求めることができ、式(2.83)はアインシュタインの関係式と呼ばれています。

ここまでの計算では、外力\gamma(t)は時間的に相関のない白色雑音という仮定をしただけで、その分布が正規型である過程はしていないことに注意するようにと教科書にはありました。

Discussion