スポンサーサイト

上記の広告は1ヶ月以上更新のないブログに表示されています。
新しい記事を書く事で広告が消せます。

[Ruby][アルゴリズム入門]乱数

第二章「数値計算」 目次

総合目次

2-1 乱数
2-2 数値積分
2-3 テイラー展開
2-4 非線形方程式の解法
2-5 補間
2-6 多桁計算
2-7 長いπ
2-8 連立方程式の解法
2-9 線形計画法
2-10 最小2乗法

一様乱数(線形合同法)

線形合同法を用いて、一様乱数を発生させる。

線形合同法は、適当なX0から初めて、
Xn=(Ax_(n-1) + C) % M
という式を使って、次々に0~Mの範囲の値を発生させる。

今回、A=109, C=1021, M=32768, X0=13を用いた。
$rndnum = 13 # X0の値

def irand
  $rndnum = ($rndnum * 109 + 1021) % 32768
  return $rndnum
end

100.times do |i|
  puts irand
end


一様性の検定

上記で作った乱数が、どのくらい均一にばらまかれているかをX^2(カイ2乗)検定の手法で計算する。

乱数が1~Mの範囲で発生するとき、iという値の発生回数をfi、iという値の発生期待値をFiとすると、
一様性の検定式
という値を計算する。このX^2の値が小さいほど均一にばらまかれていることを示す。
N = 1000
M = 10
F = N / M
SCALE = 40.0 / F

$rndnum = 13

def irnd
  $rndnum = ($rndnum * 109 + 1021) % 32768
  return $rndnum
end

def rnd
  return irnd / 32767.1
end

e = 0.0
hist = Array.new(M + 1, 0)
N.times do
  rank = (M * rnd).to_i + 1
  hist[rank] += 1
end

1.upto(M) do |i|
  printf "%3d:%3d"%[i, hist[i]]
  ((hist[i] * SCALE).to_i).times do |j|
    print "*"
  end
  puts

  e += (hist[i] - F) * (hist[i] - F) / F.to_f
end

puts " X2 = #{e}"


正規乱数(ボックス・ミュラー法)

正規乱数をボックス・ミュラー法により発生する。

def brnd (sig, m, rn)
  r1 = rand
  r2 = rand
  rn.x = sig * Math.sqrt(-2 * Math.log(r1)) * Math.cos(2 * 3.14159 * r2) + m
  rn.y = sig * Math.sqrt(-2 * Math.log(r1)) * Math.sin(2 * 3.14159 * r2) + m
end

st = Struct.new(:x, :y)
rn = st.new(0, 0)
hist = Array.new(100, 0)
1000.times do |i|
  brnd(2.5, 10.0, rn)
  hist[rn.x.to_i] += 1
  hist[rn.y.to_i] += 1
end

0.upto(20) do |i|
  printf "%3d : I "%[i]
  1.upto(hist[i] / 10) do |j|
    print "*"
  end
  puts
end

テーマ : プログラミング - ジャンル : コンピュータ

コメント
コメントの投稿
管理者にだけ表示を許可する

上記広告は1ヶ月以上更新のないブログに表示されています。新しい記事を書くことで広告を消せます。