random は Mersenne Twister による疑似乱数発生関数です。書式

で用い、0 以上「正の数」未満の疑似乱数を発生させます。「正の数」が整数の場合は整数に、小数の場合は小数になります。

(%i1) random(10);

(%o1)                                  1


(%i2) random(10);

(%o2)                                  9


(%i3) random(10);

(%o3)                                  5


(%i4) random(10);

(%o4)                                  8


(%i5) random(10);

(%o5)                                  3


(%i6) random(10.0);

(%o6)                          4.155627227214829


(%i7) random(10.0);

(%o7)                          6.125738892538246


(%i8) random(10.0);

(%o8)                          8.418126553359214


(%i9) random(10.0);

(%o9)                          5.055894446438119


(%i10) random(10.0);

(%o10)                         6.784619950185144

［利用例 1］　モンテカルロ法を用いて円周率πの近似値を求めてみましょう。関数 random を用いて、正方形領域 [0, 1) × [0, 1) 内の点 (x, y) をランダムに n 個発生させ（x = random(1.0)、y = random(1.0)）、四分円

に入るものの個数が cであれば、

によりπの近似値が求まります。

(%i11) n: 10000$


(%i12) c: 0$


(%i13) for i: 1 thru n do if (random(1.0)^2 + random(1.0)^2 < 1) then c: c + 1;

(%o13)                               done


(%i14) float(4 * c / n);

(%o14)                              3.1348

［利用例 2］　Viswanath 定数を求めてみます。Viswanath 定数とは、漸化式

で定義される Random Fibonacci 数列 {f_n} に対して、絶対値 |f_n| の n 乗根 の極限値です。まず手はじめに、Random Fibonacci 数の最初の数項を求めてみます。

(%i15) f[1]: 1; f[2]: 1; for i: 3 thru 10 do 
(f[i]: (2 * random(2) - 1) * f[i - 1] + f[i - 2], print(f[i]));

(%o15)                                 1

(%o16)                                 1

0

1

- 1

0

- 1

1

0

1

(%o17)                               done


(%i18) f[1]: 1; f[2]: 1; for i: 3 thru 10 do 
(f[i]: (2 * random(2) - 1) * f[i - 1] + f[i - 2], print(f[i]));

(%o18)                                 1

(%o19)                                 1

0

1

- 1

2

1

3

- 2

5

(%o20)                               done

つまり、+1 または -1 をとる疑似乱数は 2 * random(2) - 1 によって発生させることが出来ます。なお、当然ながら上記実験は実行の度ごとに出力が異なります。大きな n に対して、同じようにして f[n] を計算するのはメモリーの無駄遣いですから、次のように改良します。

(%i21) a: 1; b: 1; for i: 3 thru 10 do 
(c: (2 * random(2) - 1) * b + a, a: b, b: c, print(c));

(%o21)                                 1

(%o22)                                 1

0

- 1

- 1

- 2

- 1

- 3

- 4

- 1

(%o23)                               done

いよいよ本題です。第 100000 項を求め、絶対値の 100000 乗根を計算してみましょう。

(%i24) a: 1; b: 1; for i: 3 thru 10^5 do 
(c: (2 * random(2) - 1) * b + a, a: b, b: c)$

(%o24)                                 1

(%o25)                                 1


(%i27) bfloat(abs(c)^(1/10^5));

(%o27)                        1.132021884715663b0

ちなみに実際の値は 1.13198824..... だそうです。

［リンク］：

1. Mersenne Twister Home Page
2. Divakar Viswanath, Random Fibonacci sequences and the number 1.13198824.....,
   Math. Comp. 69 (2000), no. 231, 1131-1155