2025年1月26日日曜日

250126

Ruby


Stirling number(2)

|Stirling1(n,k)| のk列目を求める方とわかっている場合、以下の方が高速。(Stirling2(n,k) も同様のことが言える)

def f(n)
  return 1 if n == 0
  (1..n).inject(:*)
end

# m次以下を取り出す
def mul(f_ary, b_ary, m)
  s1, s2 = f_ary.size, b_ary.size
  ary = Array.new(s1 + s2 - 1, 0)
  (0..s1 - 1).each{|i|
    (0..s2 - 1).each{|j|
      ary[i + j] += f_ary[i] * b_ary[j]
    }
  }
  ary[0..m]
end

# m次以下を取り出す
def power(ary, n, m)
  return [1] if n == 0
  k = power(ary, n >> 1, m)
  k = mul(k, k, m)
  return k if n & 1 == 0
  return mul(k, ary, m)
end

# 符号は無視
def stirling(n, k, m = 1)
  return [1] + [0] * n if k == 0
  if m == 1
    ary = [0] + (1..n).map{|i| 1r / i}
  else
    ary = [0] + (1..n).map{|i| 1r / f(i)}
  end
  ary = power(ary, k, n)
  fk = f(k)
  (0..n).map{|i| (ary[i] * f(i) / fk).to_i}
end

n = 10
(0..10).each{|i| p [i, stirling(n, i)]}
p ""
(0..10).each{|i| p [i, stirling(n, i, 2)]}

出力結果
[0, [1, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0]]
[1, [0, 1, 1, 2, 6, 24, 120, 720, 5040, 40320, 362880]]
[2, [0, 0, 1, 3, 11, 50, 274, 1764, 13068, 109584, 1026576]]
[3, [0, 0, 0, 1, 6, 35, 225, 1624, 13132, 118124, 1172700]]
[4, [0, 0, 0, 0, 1, 10, 85, 735, 6769, 67284, 723680]]
[5, [0, 0, 0, 0, 0, 1, 15, 175, 1960, 22449, 269325]]
[6, [0, 0, 0, 0, 0, 0, 1, 21, 322, 4536, 63273]]
[7, [0, 0, 0, 0, 0, 0, 0, 1, 28, 546, 9450]]
[8, [0, 0, 0, 0, 0, 0, 0, 0, 1, 36, 870]]
[9, [0, 0, 0, 0, 0, 0, 0, 0, 0, 1, 45]]
[10, [0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 1]]
""
[0, [1, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0]]
[1, [0, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1]]
[2, [0, 0, 1, 3, 7, 15, 31, 63, 127, 255, 511]]
[3, [0, 0, 0, 1, 6, 25, 90, 301, 966, 3025, 9330]]
[4, [0, 0, 0, 0, 1, 10, 65, 350, 1701, 7770, 34105]]
[5, [0, 0, 0, 0, 0, 1, 15, 140, 1050, 6951, 42525]]
[6, [0, 0, 0, 0, 0, 0, 1, 21, 266, 2646, 22827]]
[7, [0, 0, 0, 0, 0, 0, 0, 1, 28, 462, 5880]]
[8, [0, 0, 0, 0, 0, 0, 0, 0, 1, 36, 750]]
[9, [0, 0, 0, 0, 0, 0, 0, 0, 0, 1, 45]]
[10, [0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 1]]

0 件のコメント:

コメントを投稿

注: コメントを投稿できるのは、このブログのメンバーだけです。