古典版：Vlasov 方程式と Dobrushin の安定性

近似の数理　第6章 §6.2　★★★ 平均場近似が N→∞ で厳密になることを、省略なしに証明できる最も簡単な場合

この節では、粒子の位置と速度 $z=(x,v)\in\mathbb{R}^d\times\mathbb{R}^d$ を持つ古典 $N$ 粒子系を扱う。運動方程式は

$$\dot x_i=v_i,\qquad \dot v_i=-\frac1N\sum_{j=1}^N\nabla W(x_i-x_j)\qquad(i=1,\dots,N)\tag{6.2}$$

とする。$W$ は偶関数（$W(-x)=W(x)$）で、$\nabla W$ はリプシッツ連続、つまりある定数 $L$ について$\lvert \nabla W(x)-\nabla W(y)\rvert \le L\lvert x-y\rvert$ を満たすとする。$W$ が偶関数なので $\nabla W(0)=0$ であり、和に $j=i$ の項を含めても含めなくても同じになる。

経験測度：N 個の点を1つの分布として見る

$N$ 体の状態 $(z_1,\dots,z_N)$ と一体の分布関数 $f(x,v)$ を比べるには、両者を同じ種類の対象にそろえる必要がある。そのための道具が測度である。

経験測度にすると、粒子に番号を付けた情報が消え、「どこに何割の粒子がいるか」だけが残る。同種粒子の系ではこれで十分であり、しかも粒子数 $N$ が違う系どうしや、連続な分布 $f$ と同じ土俵で比べられる。相空間の測度 $\mu$ から位置の分布 $\rho_\mu$（$v$ について積分したもの）を作り、平均の力を

$$(\nabla W*\rho_\mu)(x)=\int\nabla W(x-y)\,\mathrm{d}\rho_\mu(y)$$

と書く（$*$ は畳み込み）。経験測度なら $(\nabla W*\rho_{\mu_N})(x_i)=\frac1N\sum_j\nabla W(x_i-x_j)$ で、まさに (6.2) の右辺の力である。

Vlasov 方程式と、経験測度がその解であること

一体の分布関数 $f_t(x,v)$ に対するVlasov 方程式は

$$\frac{\partial f}{\partial t}+v\cdot\nabla_xf-(\nabla W*\rho_{f})\cdot\nabla_vf=0,\qquad \rho_f(x)=\int f(x,v)\,\mathrm{d} v\tag{6.3}$$

である。これは「各粒子は、全体がつくる平均の力 $-\nabla W*\rho_f$ の中を独立に運動する」という一体の運動方程式を、分布関数の連続の式に書き直したものである。力が $f$ 自身に依存するので非線形の偏微分方程式になる。

経験測度は滑らかでないので (6.3) に直接代入できない。そこで滑らかな試験関数 $g(x,v)$ を掛けて部分積分した形

$$\frac{\mathrm{d} }{\mathrm{d} t}\int g\,\mathrm{d}\mu_t=\int\bigl[v\cdot\nabla_xg-(\nabla W*\rho_{\mu_t})\cdot\nabla_vg\bigr]\mathrm{d}\mu_t\tag{6.4}$$

がすべての $g$ で成り立つことを、測度 $\mu_t$ が Vlasov 方程式を（弱い意味で）満たすことの定義とする。$\mu_t$ が滑らかな密度 $f_t$ を持てば、部分積分を戻して (6.3) と同値になる。

証明

$\int g\,\mathrm{d}\mu_N(t)=\frac1N\sum_ig(x_i(t),v_i(t))$ を時間微分すると、連鎖律と (6.2) により

$$\frac{\mathrm{d} }{\mathrm{d} t}\int g\,\mathrm{d}\mu_N=\frac1N\sum_i\Bigl[\nabla_xg(z_i)\cdot v_i-\nabla_vg(z_i)\cdot\frac1N\sum_j\nabla W(x_i-x_j)\Bigr].$$

内側の和は $(\nabla W*\rho_{\mu_N})(x_i)$ に等しいので、右辺は $\int[v\cdot\nabla_xg-(\nabla W*\rho_{\mu_N})\cdot\nabla_vg]\,\mathrm{d}\mu_N$ である。

この命題の意味は大きい。$N$ 体の問題と一体の問題は同じ方程式の異なる初期値に対する解であり、「平均場近似が正しいか」という問いは「Vlasov 方程式の解は初期値に連続に依存するか」という安定性の問いに置き換わる。初期値が $N$ 個の点の経験測度 $\mu_N(0)$ か、滑らかな $f_0$ か、の違いだけである。

Wasserstein 距離：点の集合と滑らかな分布を比べる

安定性を述べるには測度どうしの距離が要る。ここで素朴な距離は役に立たない。たとえば密度の差の $L^1$ ノルム$\int\lvert f-g\rvert \,\mathrm{d} z$ を測度へ拡張した全変動距離では、経験測度 $\mu_N$ と滑らかな $f$ の距離は、$N$ をどれだけ大きくしても最大値の 2 のままである（$\mu_N$ は $N$ 個の点に集中し、$f$ はその点に確率 0 しか置かないため）。ここで要る距離は、「点の位置が少しずれたら距離も少しだけ変わる」もの、つまり確率の質量を運ぶ距離を測るものである。

カップリング $\pi$ は「$\mu$ の質量のどれだけを $\nu$ のどこへ運ぶか」という輸送計画であり、$W_1$ はその最小の輸送コストである（最適輸送の Monge–Kantorovich 問題）。点の数が同じ経験測度どうしなら、$W_1$ は点の対応付け（置換 $\sigma$）の最小値$\min_\sigma\frac1N\sum_i\lvert z_i-z'_{\sigma(i)}\rvert$ に一致する。点 $z_0$ と $z_0'$ の Dirac 測度なら $W_1=\lvert z_0-z_0'\rvert$ で、点が少し動けば距離も少しだけ変わる。

Dobrushin の安定性定理

一般の初期測度 $\mu_0$ に対して、Vlasov 方程式の解は次の非線形の特性曲線で作る。$\mu_t$ が分かっているとして、

$$\dot X=V,\qquad \dot V=-(\nabla W*\rho_{\mu_t})(X),\qquad (X,V)(0)=z\tag{6.5}$$

の解を $Z_t^\mu(z)=(X_t^\mu(z),V_t^\mu(z))$ と書き、$\mu_t$ を「$\mu_0$ の各点を $Z_t^\mu$ で運んだ測度」（押し出し $\mu_t=(Z_t^\mu)_\#\mu_0$、すなわち $\int g\,\mathrm{d}\mu_t=\int g(Z_t^\mu(z))\,\mathrm{d}\mu_0(z)$）とする。力が $\mu_t$ に依存するので自己無撞着な条件になるが、$\nabla W$ がリプシッツなら Picard の逐次近似（第2章の存在定理と同じ論法）で大域的に一意に解けることが知られている。$\mu_0$ が経験測度なら、$Z_t^\mu$ は (6.2) の各粒子の軌道そのものである。

証明

証明の考え方は「二つの測度の粒子をカップリングで組にし、組ごとに軌道のずれを追う」ことである。

(i) 初期のカップリング。 $\varepsilon >0$ を任意にとり、$\iint\lvert z-z'\rvert \,\mathrm{d}\pi_0\le W_1(\mu_0,\nu_0)+\varepsilon$ となる $\pi_0\in\Pi(\mu_0,\nu_0)$ を選ぶ。組 $(z,z')$ を $\mu$ 側の特性曲線 $Z_t^\mu(z)$ と $\nu$ 側の特性曲線 $Z_t^\nu(z')$ で運び、

$$D(t)=\iint\lvert Z_t^\mu(z)-Z_t^\nu(z')\rvert \,\mathrm{d}\pi_0(z,z')$$

とおく。組を運んだ測度 $\pi_t$（$\iint G\,\mathrm{d}\pi_t=\iint G(Z_t^\mu(z),Z_t^\nu(z'))\,\mathrm{d}\pi_0$）の第1成分の分布は $(Z_t^\mu)_\#\mu_0=\mu_t$、第2成分は $\nu_t$ なので $\pi_t\in\Pi(\mu_t,\nu_t)$ であり、$W_1$ の定義から

$$W_1(\mu_t,\nu_t)\le\iint\lvert z-z'\rvert \,\mathrm{d}\pi_t=D(t),\qquad D(0)\le W_1(\mu_0,\nu_0)+\varepsilon .\tag{6.6}$$

(ii) 力の差の評価。 組 $(z,z')$ を固定し、$X=X_t^\mu(z)$、$X'=X_t^\nu(z')$ と略記する。$\rho_{\mu_t}$ は $\mu_0$ を運んだ測度の位置成分なので$(\nabla W*\rho_{\mu_t})(X)=\int\nabla W(X-X_t^\mu(w))\,\mathrm{d}\mu_0(w)$ と書け、$\nu$ 側も同様である。$\pi_0$ の周辺分布が $\mu_0,\nu_0$ であることを使って一つの積分にまとめると

$$(\nabla W*\rho_{\mu_t})(X)-(\nabla W*\rho_{\nu_t})(X') =\iint\bigl[\nabla W(X-X_t^\mu(w))-\nabla W(X'-X_t^\nu(w'))\bigr]\mathrm{d}\pi_0(w,w').$$

リプシッツ性から被積分関数の大きさは $L\lvert X-X'\rvert +L\lvert X_t^\mu(w)-X_t^\nu(w')\rvert$ 以下なので、

$$\bigl|(\nabla W*\rho_{\mu_t})(X)-(\nabla W*\rho_{\nu_t})(X')\bigr|\le L\lvert X-X'\rvert +L\,D(t).$$

ここが証明の要である。平均の力のずれは「自分の位置のずれ」と「他の粒子の位置のずれの平均」で抑えられ、後者が再び $D$ になる。

(iii) 軌道のずれ。 (6.5) を積分形に書き、$V=V_t^\mu(z)$、$V'=V_t^\nu(z')$ として

$$\lvert X(t)-X'(t)\rvert \le\lvert x-x'\rvert +\int_0^t\lvert V-V'\rvert \,\mathrm{d} s,\qquad \lvert V(t)-V'(t)\rvert \le\lvert v-v'\rvert +\int_0^t\bigl(L\lvert X-X'\rvert +L\,D(s)\bigr)\mathrm{d} s.$$

足し合わせ、$\lvert V-V'\rvert +L\lvert X-X'\rvert \le(1+L)\lvert Z-Z'\rvert$ を使うと

$$\lvert Z_t^\mu(z)-Z_t^\nu(z')\rvert \le\lvert z-z'\rvert +\int_0^t\bigl[(1+L)\lvert Z_s^\mu(z)-Z_s^\nu(z')\rvert +L\,D(s)\bigr]\mathrm{d} s.$$

両辺を $\pi_0$ で積分すると $D(t)\le D(0)+(1+2L)\int_0^tD(s)\,\mathrm{d} s$ を得る。

(iv) Grönwall の不等式で閉じる。 $G(t)=D(0)+(1+2L)\int_0^tD\,\mathrm{d} s$ とおくと $D\le G$ なので $G'=(1+2L)D\le(1+2L)G$、したがって $(\mathrm{e}^{-(1+2L)t}G)'\le0$ から $D(t)\le G(t)\le\mathrm{e}^{(1+2L)t}D(0)$。(6.6) と合わせ、$\varepsilon \to0$ とすれば主張を得る。

証明

命題 6.2 により $\mu_N(t)$ は Vlasov 方程式の解なので、定理 6.4 を $\mu_0=\mu_N(0)$、$\nu_0=f_0$ に使えばよい。初期の収束 $W_1(\mu_N(0),f_0)\to0$ は経験測度の大数の法則（有界連続関数の標本平均が期待値に収束し、1次モーメントも収束する）による。

これが「平均場近似は $N\to\infty$ で厳密になる」の最初の完全な例である。物理屋の直観と照らし合わせると、直観の「力の揺らぎは平均されて消える」は初期時刻の大数の法則に対応し、直観では暗黙にされている「消えた揺らぎが時間発展で増幅されない」ことを保証するのが Dobrushin の安定性である。後者は自明でない。カオス的な力学系では初期のずれが指数関数的に増え、定理の $\mathrm{e}^{(1+2L)t}$ はまさにその最悪の増え方を表している。したがって近似の保証は有限時間 $T$ の間に限られ、$N$ を固定して $t\to\infty$ とすると保証は失われる（$t\sim\log N$ の時間で破綻しうる）。

数値で見る

1次元（$d=1$）で $W(x)=-\cos x$（$\nabla W=\sin x$、$L=1$ の引力的な相互作用）をとり、初期分布 $f_0$ を$x\sim\mathcal N(0,1)$、$v\sim\mathcal N(0,0.3^2)$ の独立な正規分布とした。Vlasov 方程式は、$f_0$ を Gauss–Hermite 求積の重みつき点（$220^2$ 個）で表して特性曲線 (6.5) を解いた（$150^2$ 個との差は $0\le t\le8$ で $5\times10^{-4}$ 未満）。図 6.1 は、凝集の度合い $m(t)=\int\cos x\,\mathrm{d}\rho_t$ を Vlasov 解と、$f_0$ から抽出した $N=100$、$2000$ 粒子の系とで比べたものである。$N$ 粒子系の $m_N(t)$ は Vlasov 解のまわりに揺らぎ、その揺らぎは $N$ とともに小さくなる。12 回の抽出で平均した $\sup_{0\le t\le8}\lvert m_N(t)-m(t)\rvert$ は $N=10^2,10^3,10^4$ でそれぞれ $0.075$、$0.020$、$0.0052$ であり、両対数の傾きは $-0.58$ で $N^{-1/2}$ に近い。同じ系で二つの $N=150$ 粒子配置の $W_1$ を最適割当で直接計算すると、$W_1(\mu_t,\nu_t)/W_1(\mu_0,\nu_0)$ は $t=3$ で $1.12$ にとどまり、評価 $\mathrm{e}^{3t}\approx8100$ よりはるかに小さい。Dobrushin の評価は最悪の場合に対するもので、典型的な配置では緩い（検算スクリプト）。