2016年1月3日日曜日

160103(3)

Ruby


A190945(100)(3)

念のため、オンライン整数列大辞典の
A190945(http://oeis.org/A190945/list)
と比較し、答え合わせしてみる。

def C(n, r)
  r = [r, n - r].min
  return 1 if r == 0
  return n if r == 1
  numerator = (n - r + 1..n).to_a
  denominator = (1..r).to_a
  (2..r).each{|p|
    pivot = denominator[p - 1]
    if pivot > 1
      offset = (n - r) % p
      (p - 1).step(r - 1, p){|k|
        numerator[k - offset] /= pivot
        denominator[k] /= pivot
      }
    end
  }
  result = 1
  (0..r - 1).each{|k|
    result *= numerator[k] if numerator[k] > 1
  }
  return result
end

=begin
隣接をm個除去しておいてd + m個隣接を増やす入れ方を考える。
つまり、
隣接箇所(k)からm箇所選んでに数字を入れ
非隣接箇所(s - k)からmm(= i - d - 2m)箇所数字を入れ
入れた数字(m + mm)にd + m個数字を重複的に追加する方法の数である。

調べる範囲を整理しておく。
上記より、
① 0 ≦ m ≦ k
② 0 ≦ mm ≦ s - k
③ m + mm ≠ 0
④ d + m ≧ 0
となるが、②については、①および④より、
s - k
≧ (i - 1) * i / 2 + 1 - (i - 1) * (i - 2) / 2 = i
≧ i - (d + m) - m = mm
となるので、0 ≦ mm でよい。
=end
def A190945(n)
  a_ary = [1]
  # iまでの数字を使った並びで隣接がj箇所ある並び方がa[j]
  a = [1]
  (2..n).each{|i|
    b = []
    s = (i - 1) * i / 2 + 1
    (0..s - 1).each{|j|
      x = 0
      (0..(i - 1) * (i - 2) / 2).each{|k|
        d = j - k
        (0..k).each{|m|
          if d + m >= 0
            mm = i - d - 2 * m
            # a[k] * C(k, m) * C(s - k, mm) * H(m + mm, d + m)
            x += a[k] * C(k, m) * C(s - k, mm) * C(i - 1, d + m) if 0 <= mm && m + mm != 0
          end
        }
      }
      b[j] = x
    }
    a = b
    a_ary << a[0]
  }
  a_ary
end
ary = A190945(11)

# OEIS A190945のデータ
ary0 =
[1,1,10,1074,1637124,45156692400,
 27230193578558160,420296434943941609215120,
 190200071567439616748736269178720,
 2843464512159537301384360263178682136716160,
 1562137388408002436396705025296003247844758163480828800]
# 一致の確認
p ary == ary0

0 件のコメント:

コメントを投稿

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