5 モンテカルロ積分法

5.1 単純な場合(２次元の面積)

ここでは，モンテカルロ法(Monte Carlo method)呼ばれる乱数を使った一風変わった数値積分の方法を示す．モンテカルロとは有名なモナコのカジノで，博打をするところである．この方法を少し学習すると，名前の由来が分かったような気になるはずである．．

まずは，図4の様な円の面積を求めることにする．もちろん，この円の面積は $\pi$ と分かっているが，これを乱数を用いて数値計算しようと言うのである．円を囲む正方形内部には，ランダムな点がある．正方形内部にある点の個数を，円の内部にある個数をとする．そして，正方形の面積を，円をとする．すると，これらには

$\displaystyle \frac{S_c}{S_s}\simeq\frac{N_c}{N_s}$

(16)

の関係がある．ランダムな点の数増加させると，これらは，ほとんど等式で結ばれるであろう．形状が単純な正方形の面積

は分かっているので，円の面積

は

$\displaystyle S_c\simeq\frac{N_c}{N_s}S_s$

(17)

と計算できる．点の数を増やすと，それがランダムに分布する限り，右辺は円の面積に近づく．正方形の面積と点の数を計算するだけで円の面積が分かるのである．このようにランダムな点(乱数)を用いて，積分する方法をモンテカルロ積分と言う．

円の面積を求めても大したことはないが，例えば図5のような複雑な形状の面積が計算できるとなるとぐっと御利益はある．先ほどの円の面積と同じようにして，計算できるのである．計算したい部分内の点の数 $N_{in}$ は，

$\displaystyle N_{in}$	$\displaystyle =($ 大円 $\displaystyle \Omega_0$ 内部の点数 $\displaystyle ) -($ 大円 $\displaystyle \Omega_0$ 内部かつ小円 $\displaystyle \Omega_1$ 内部の点数 $\displaystyle )$
	$\displaystyle \qquad\qquad -($ 大円 $\displaystyle \Omega_0$ 内部かつ楕円 $\displaystyle \Omega_2$ 内部の点数 $\displaystyle ) -($ 大円 $\displaystyle \Omega_0$ 内部かつ正方形 $\displaystyle \Omega_3$ 内部の点数 $\displaystyle )$	(18)

こうやって，内部の点の数を調べて，後は，式(18)のようにすれば面積を求めることができる．

このように複雑な領域の積分でも適用でき，応用範囲の広い方法である．積分の精度を上げようとすると，サンプル数を増やすことになる．しかし，その精度は，たかだかサンプル数の平方根に比例するのみである．誤差を1/10にしたければ，サンプル数を100倍にする必要がある．計算時間も100倍になる．これが，この方法の最大の弱点である．

精度がサンプル数の平方根に比例するのは，完全な乱数を仮定しているからである．それを計算機で発生させることは不可能なので，その生成も問題となる．できるだけ良い乱数を使わなくてはならないが，この問題は奥が深い．ここではこれ以上立ち入らないが，興味のある者は調べてみると良い．

話は変わるが，ひとつ注意を与えておく．実際のプログラムでは，乱数を使ってを発生させて，それが内部にあるか否かをその場で判断する．配列に座標データを入れて，最後に検査するとメモリーがいくらあっても足りなくなるからである．

**図 4:** モンテカルロ積分．正方形内にランダムな点をばらまく．ランダムな点の数を増やせば，正方形内と円内の点の数の比が，正方形と円の面積の比になる．
$\includegraphics[keepaspectratio, scale=1.0]{figure/MonteCarlo.eps}$

**図 5:** この面積を解析的に計算するのは，大変面倒である．
$\includegraphics[keepaspectratio, scale=1.0]{figure/MonteCarlo_complex.eps}$

5.2 複雑な積分

もう少し複雑な積分が必要な場合を考えてみよう．ある関数 $f(x_1,x_2,x_3,\cdot,x_M)$ の体積積分を考える．この体積積分は

$\displaystyle \int_\Omega f(x_1,x_2,x_3,\cdots,x_M)dx_1dx_2dx_3\cdots dx_M= V\langle f \rangle$

(19)

となる．ここで，

はM次元体の領域 $\Omega$ の体積， $\langle f \rangle$ はその内部の関数の平均値である．これは3次元体の質量を考えると簡単である．その質量は，密度の体積積分となる．これは体積に平均密度を乗じた値に等しい．当たり前の式である．このことから，式(20)の積分の値が欲しければ，体積と平均値が分かればよいことになる．

積分を計算するために，まずは体積である．先ほどは，二次元問題であるが，はかなり複雑な形状の面積をモンテカルロ積分で計算する方法を示した．2次元と同じ考えで，3次元， 4次元， $\cdots$ ，M次元の体積も比較的，容易に求めることができる．分かり切った体積に，領域 $\Omega$ を包み込み，その内部にランダムに配置されたサンプル点の数を数えれば良いのである．体積の計算は簡単である．

残りは，体積内部の平均 $\langle f \rangle$ である．これも簡単で，領域 $\Omega$ 内部にあるサンプル点の平均より求めることができる．即ち，

$\displaystyle \langle f \rangle\simeq\frac{1}{N_{in}}\sum_i f(\boldsymbol{r}_i)$ 領域 $\Omega$ の内部のみ

(20)

となる．ここで， $\boldsymbol{r}_i$ は

番目のサンプル点の $(x_1,x_2,x_3,\cdot,x_N)$ 座標である．また， $N_{in}$ は領域内部のサンプル数である．この計算も簡単で，コンピューターは難なく，平均値に近似値を求めることができる．

以上より，モンテカルロ法を用いると，体積と平均値 $\langle f \rangle$ の近似値が計算できることが分かる．従って，式(20)の近似値を求めることができる．

数値計算で有名な本「NUMERICAL RECIPES in C」 [#!NRC!#]には，図 6の様な形状の重心を求める問題が出ている．このように複雑な積分が必要なときにモンテカルロ法は威力を発揮する．

**図 6:** 「NUMERICAL RECIPES in C」の例題．一部を切り取られたトーラスの重心を求める．
$\includegraphics[keepaspectratio, scale=0.5]{figure/NRC_example.eps}$

ホームページ: Yamamoto's laboratory
著者: 山本昌志

Yamamoto Masashi
2006-01-30