乱数の質その2

モンテカルロ法だとよくわからなかったglibcとMelsenne Twisterの違い。ならば単純にヒストグラムを取ったらどうなるという訳でテストです。
一様乱数だと処理が面倒なので、0から9999までの乱数を発生させて、各数値の発生回数をカウントします。

全体で1010回の試行の結果は次の通り。

分散の値はMTの方が大きくなっているけど、そういうものなのかな??

実際的にはホワイトノイズ(例のシャーってノイズ)にA/Dコンを接続したものや量子力学的現象を利用するとかいろいろな方法が提案されていますが、何が一番なんだろうな。というか、乱数の質ってなんだか哲学的な話になってきたので簡単には結論は出せないんでしょうけど。

乱数の質?

乱数を使用してあれこれという場面は結構ありますが、質がどうなのかのチェックはしていなかったと思い、ちょっとテストした結果です。

比較にはLinux標準のglibcに含まれるrand()と Mersenne Twisterと呼ばれる手法によるgenrand_real1()を利用して一様乱数を発生させてどんな感じなのかを調べます。ちなみにMersenne Twisterは非常に質が良いとされている最近出てきた乱数生成手法。
ただヒストグラムを取っても面白くないので、モンテカルロ法で円周率を計算させてみます。

プログラム自体は単純なので省略。各10億(109)点(20億組)から円周率を計算することを1000回繰り返すので、1兆(1012)組の乱数を発生させることになります。
MT.hの中身を見ていないのですが、けっこう面倒なことをしているようで、rand()が360分かかったのに対して、genrand()の方は415分が所要時間。

で、結果はどうなのかというと、
円周率計算結果
と大差無しの結果が得られました。なんか肩透かしだな。というよりrand()もなかなかやるじゃないかという訳で安心して使用できそうだという結論。それともあと4桁くらい余計に生成させると差が出るんだろうか?  ただ4桁増やすと、結果出るのが4150000分後って288日もかかるんでやる気は無いけど。

追記:
もしかして、109個ってのが多すぎだったかも。3桁くらい減らしてみたら、違いが出るか?

多角形からの円周率

円周率の算出方法として内接する多角形の一辺の長さの合計から割り出す方法があります。

基本から考えると、こんな感じかな。
まず円に内接する多角形を考える。この多角形は中心を頂点とするn個の二等辺三角形に分割でき、底辺の長さの合計が円周に近似される。
さて角数をnとすると、各三角形の頂角θは

であり、各三角形の底角θ1

頂点から底辺の中央に向けて垂線を引くと、底辺の1/2の点で交わるので、円の半径をrとすれば、底辺の1/2の長さr1

となるので、n角形の辺の長さの合計Rは

と表される。
元の円の半径rを0.5とすれば、
となる。

例えば、4角形ならば4cos(90-180/4)=2.828…、5角形なら5cos(90-180/5)=2.938…
120角形ならば、120cos(90-180/120)=3.141233…
といった具合。

つまり、

とまぁ、Texを使えば数式をそれらしく書けるのが嬉しくて下らないことをだらだら書きました。

現代では三角関数が使えるので、それこそ電卓でも持ち出せば1000角形だろうが一瞬で計算OKと。それでは意味はないので、幾何的に処理することを考えないと。
幾何的に処理となると、三平方の定理を使ってやればいいのかな? ちょっと図を書いて考えてみよう。