binomial()の実装方法を調べてみる

$ @which binomial(5, 2)

code:jl

Base.@assume_effects :terminates_locally function binomial(n::T, k::T) where T<:Integer

n0, k0 = n, k

# 定義に合うように処理

k < 0 && return zero(T)

sgn = one(T)

if n < 0

n = -n + k - one(T)

if isodd(k)

sgn = -sgn

end

end

k > n && return zero(T)

(k == 0 || k == n) && return sgn

k == 1 && return sgn*n

# 半分で折り返せる

if k > (n>>1)

k = (n - k)

end

x = nn = n - k + one(T)

nn += one(T)

rr = T(2)

while rr <= k

xt = div(widemul(x, nn), rr)

x = xt % T

x == xt || throw(OverflowError(LazyString("binomial(", n0, ", ", k0, ") overflows")))

rr += one(T)

nn += one(T)

end

copysign(x, sgn)

end

漸化式の形式で表現しておく

$ a_0 = n-k+1

$ a_r = \frac{a_{r-1}\cdot(n-k+r)}{r}

code:jl

function binomial_(n, k)

if k > (n ÷ 2)

k = (n - k)

end

x = n = n - k + 1

n += 1

r = 2

while r <= k

x = div(x * n, r)

r += 1

n += 1

end

return x

end