Problem:
 Score : 500 points    Problem Statement There is a building with n rooms, numbered 1 to n .  We can move from any room to any other room in the building.  Let us call the following event a move : a person in some room i goes to another room j~ (i \neq j) .  Initially, there was one person in each room in the building.  After that, we know that there were exactly k moves happened up to now.  We are interested in the number of people in each of the n rooms now. How many combinations of numbers of people in the n rooms are possible?  Find the count modulo (10^9 + 7) .      Constraints  All values in input are integers.  3 \leq n \leq 2 \times 10^5  2 \leq k \leq 10^9         Input Input is given from Standard Input in the following format:  n  k       Output Print the number of possible combinations of numbers of people in the n rooms now, modulo (10^9 + 7) .        Sample Input 1 3 2      Sample Output 1 10  Let c_1 , c_2 , and c_3 be the number of people in Room 1 , 2 , and 3 now, respectively. There are 10 possible combination of (c_1, c_2, c_3) :   (0, 0, 3)  (0, 1, 2)  (0, 2, 1)  (0, 3, 0)  (1, 0, 2)  (1, 1, 1)  (1, 2, 0)  (2, 0, 1)  (2, 1, 0)  (3, 0, 0)   For example, (c_1, c_2, c_3) will be (0, 1, 2) if the person in Room 1 goes to Room 2 and then one of the persons in Room 2 goes to Room 3 .       Sample Input 2 200000 1000000000      Sample Output 2 607923868  Print the count modulo (10^9 + 7) .       Sample Input 3 15 6      Sample Output 3 22583772

Solution (Python):
import sys
import itertools
input = sys.stdin.readline
sys.setrecursionlimit(100000)


def read_values():
    return map(int, input().split())


def read_index():
    return map(lambda x: x - 1, map(int, input().split()))


def read_list():
    return list(read_values())


def read_lists(N):
    return [read_list() for n in range(N)]


def functional(N, mod):
    F = [1] * (N + 1)
    for i in range(N):
        F[i + 1] = (i + 1) * F[i] % mod
    return F


INV = dict()
def inv(a, mod):
    if a in INV:
        return INV[a]
    i = pow(a, mod - 2, mod)
    INV[a] = i
    return i 


def C(F, a, b, mod):
    return F[a] * inv(F[a - b], mod) * inv(F[b], mod) % mod 


def main():
    N, K = read_values()
    mod = 10 ** 9 + 7
    F = functional(5 * 10 ** 5, mod)
    res = 0
    if N <= K:
        res = C(F, 2 * N - 1, N - 1, mod)
    else:
        c = C(F, N - 1, 0, mod)
        for k in range(K + 1):
            res += pow(c, 2, mod) * N * inv(N - k, mod) % mod
            c *= ((N - k - 1) * inv(k + 1, mod)) % mod
            c %= mod
    
    print(res % mod)

if __name__ == "__main__":
    main()

Solution (Python):
import sys
import itertools
input = sys.stdin.readline
sys.setrecursionlimit(100000)


def read_values():
    return map(int, input().split())


def read_index():
    return map(lambda x: x - 1, map(int, input().split()))


def read_list():
    return list(read_values())


def read_lists(N):
    return [read_list() for n in range(N)]


def functional(N, mod):
    F = [1] * (N + 1)
    for i in range(N):
        F[i + 1] = (i + 1) * F[i] % mod
    return F


def inv(a, mod):
    return pow(a, mod - 2, mod)


def C(F, a, b, mod):
    return F[a] * inv(F[a - b], mod) * inv(F[b], mod) % mod 


def main():
    N, K = read_values()
    mod = 10 ** 9 + 7
    F = functional(5 * 10 ** 5, mod)
    res = 0
    if N <= K:
        res = C(F, 2 * N - 1, N - 1, mod)
    else:
        c = C(F, N - 1, 0, mod)
        for k in range(K + 1):
            res += pow(c, 2, mod) * N * inv(N - k, mod) % mod
            c *= ((N - k - 1) * inv(k + 1, mod)) % mod
            c %= mod
    
    print(res % mod)

if __name__ == "__main__":
    main()

Solution (Python):
n,k = map(int,input().split())
P = 10**9+7

ans = 1
nCl = 1
n1Cl = 1
for l in range(1,min(k+1,n)):
    nCl = nCl*(n+1-l)*pow(l,P-2,P)%P
    n1Cl = n1Cl*(n-l)*pow(l,P-2,P)%P
    ans = (ans+nCl*n1Cl)%P

print(ans)

Solution (Python):
def combinations(n,k,mod):
    x,y = 1,1
    for i in range(k):
        x = (x*(n-i))%mod
        y = (y*(i+1))%mod
    return (x*pow(y,mod-2,mod))%mod
    

mod = 10**9+7
n,k = map(int,input().split())

def cmap(n,k):
    x = [1]
    y = [1]
    c = [1]
    for i in range(k+1):
        x.append((x[i]*(n-i))%mod)
        y.append((y[i]*(i+1))%mod)
        c.append((x[i+1]*pow(y[i+1],mod-2,mod))%mod)
    return c
    

nc=cmap(n,n-1)
n1c=cmap(n-1,n-2)

ans = 0
for i in range(min(n-1,k)+1):
    ans = (ans+(nc[i]*n1c[i]))%mod

print(ans)

Solution (Python):
n, k = list(map(int, input().split()))

max_len = 2 * n - 1  
mod = 10**9 + 7


def modinv(x):
    '''
    xの逆元を求める。フェルマーの小定理より、 x の逆元は x ^ (mod - 2) に等しい。計算時間はO(log(mod))程度。
    Python標準のpowは割った余りを出すことも可能。
    '''
    return pow(x, mod-2, mod)





modinv_table = [-1] * (max_len + 1)
modinv_table[0] = None  
factori_table = [1] * (max_len + 1)
factori_inv_table = [1] * (max_len + 1)
for i in range(1, max_len + 1):
    factori_table[i] = factori_table[i-1] * (i) % mod

modinv_table[1] = 1

for i in range(2, max_len + 1):
    modinv_table[i] = (-modinv_table[mod % i] * (mod // i)) % mod
    factori_inv_table[i] = factori_inv_table[i-1] * modinv_table[i] % mod


def binomial_coefficients(n, k):
    '''
    n! / (k! * (n-k)! )
    0 <= k <= nを満たさないときは変な値を返してしまうので、先にNoneを返すことにする。
    場合によっては0のほうが適切かもしれない。
    '''
    if not 0 <= k <= n:
        return None
    return (factori_table[n] * factori_inv_table[k] * factori_inv_table[n-k]) % mod


def binomial_coefficients2(n, k):
    '''
    (n * (n-1) * ... * (n-k+1)) / (1 * 2 * ... * k)
    '''
    ans = 1
    for i in range(k):
        ans *= n-i
        ans *= modinv_table[i + 1]
        ans %= mod
    return ans



if k >= n-1:
    
    print(binomial_coefficients(2 * n - 1, n))
else:
    
    
    
    ans = 0
    for j in range(k+1):
        if j == 0:
            ans += 1
        else:
            ans += binomial_coefficients(n, j) * binomial_coefficients(n-1, j) 
            ans %= mod
    print(ans)

Solution (Python):
n, k = map(int, input().split())
mod = 10**9 + 7

SIZE = 2 * 10 ** 5 + 1
fact = [0] * SIZE
inv = [0] * SIZE
finv = [0] * SIZE
fact[0], fact[1] = 1, 1
inv[1] = 1
finv[0], finv[1] = 1, 1
for i in range(2, SIZE):
    fact[i] = fact[i - 1] * i % mod
    inv[i] = mod - inv[mod % i] * (mod // i) % mod
    finv[i] = finv[i - 1] * inv[i] % mod


def nCr(n, r):
    if n < 0 or r < 0 or n < r:
        return 0
    return fact[n] * (finv[r] * finv[n - r] % mod) % mod


ans = 1
limit = min(k+1, n)
for i in range(1, limit):
    vacant = i
    fill = n-i
    add = i
    pattern = nCr(n, vacant) * nCr(add+fill-1, add)
    ans = (ans + pattern) % mod

print(ans)

Solution (Python):
n, k = map(int, input().split())
mod = 10**9 + 7

SIZE = 2 * 10 ** 5 + 1
fact = [0] * SIZE
inv = [0] * SIZE
finv = [0] * SIZE
fact[0], fact[1] = 1, 1
inv[1] = 1
finv[0], finv[1] = 1, 1
for i in range(2, SIZE):
    fact[i] = fact[i - 1] * i % mod
    inv[i] = mod - inv[mod % i] * (mod // i) % mod
    finv[i] = finv[i - 1] * inv[i] % mod


def nCr(n, r):
    if n < 0 or r < 0 or n < r:
        return 0
    return fact[n] * (finv[r] * finv[n - r] % mod) % mod


ans = 1
limit = min(k+1, n)
for i in range(1, limit):
    vacant = i
    fill = n-i
    add = i
    pattern = nCr(n, vacant) * nCr(add+fill-1, add)
    ans = (ans + pattern) % mod

print(ans)

Solution (Python):



from functools import reduce




def main():

    mod = 1000000007
    n, k = map(int, input().split())

    fact1 = [0] * (n + 1)  
    fact1[0] = 1

    for i in range(1, n + 1):
        fact1[i] = (fact1[i - 1] * i) % mod

    fact2 = [0] * (n + 1)  
    fact2[n] = pow(fact1[n], mod-2, mod)

    for i in range(n, 0, -1):
        fact2[i - 1] = (fact2[i] * i) % mod

    def cmb2(n, r, mod):
        if r < 0 or r > n:
            return 0
        return fact1[n] * fact2[n - r] * fact2[r] % mod

    ans = 0
    for i in range(min(n - 1, k) + 1):
        ans += cmb2(n, i, mod) * cmb2(n - 1, i, mod)
        ans %= mod
    print(ans)


main()

Solution (Python):
n, k = list(map(int, input().split()))

mod = 10**9 + 7


def modinv(x):
    return pow(x, mod-2, mod)





max_len = 2 * n - 1 

factori_table = [1] * (max_len + 1)
factori_inv_table = [1] * (max_len + 1)
for i in range(1, max_len + 1):
    factori_table[i] = factori_table[i-1] * (i) % mod
    factori_inv_table[i] = modinv(factori_table[i])

def binomial_coefficients(n, k):
    
    if k <= 0 or k >= n:
        return 1
    return (factori_table[n] * factori_inv_table[k] * factori_inv_table[n-k]) % mod

if k >= n-1:
    
    print(binomial_coefficients(2 * n - 1, n))
else:
    
    
    
    ans = 0
    for j in range(k+1):
        if j == 0:
            ans += 1
        else:
            ans += binomial_coefficients(n, j) * binomial_coefficients(n-1, j) 
            ans %= mod
    print(ans)

Solution (Python):
class Factorial:
    def __init__(self,n,mod):
        self.f=f=[0]*(n+1)
        f[0]=b=1
        self.mod=mod
        for i in range(1,n+1):f[i]=b=b*i%mod
        self.inv=inv=[0]*(n+1)
        inv[n]=b=pow(self.f[n],mod-2,mod)
        for i in range(n,0,-1):inv[i-1]=b=b*i%mod
    def factorial(self,i):
        return self.f[i]
    def ifactorial(self,i):
        return self.inv[i]
    def comb(self,n,k):
        if n>=k:
            return self.f[n]*self.inv[n-k]%self.mod*self.inv[k]%self.mod
        else:
            return 0
    def homp(self,n,k):
        return self.comb(n+k-1,n-1)
M=10**9+7
n,k=map(int,input().split())
f=Factorial(n*2,M)
c=f.comb
h=f.homp
a=0
for i in range(min(k+1,n)):
    a+=c(n,i)*h(n-i,i)%M
print(a%M)

Solution (Python):
class Factorial:
    def __init__(self,n,mod):
        self.f=f=[0]*(n+1)
        f[0]=b=1
        self.mod=mod
        for i in range(1,n+1):f[i]=b=b*i%mod
        self.inv=inv=[0]*(n+1)
        inv[n]=b=pow(self.f[n],mod-2,mod)
        for i in range(n,0,-1):inv[i-1]=b=b*i%mod
    def factorial(self,i):
        return self.f[i]
    def ifactorial(self,i):
        return self.inv[i]
    def comb(self,n,k):
        if n>=k:
            return self.f[n]*self.inv[n-k]%self.mod*self.inv[k]%self.mod
        else:
            return 0
    def homp(self,n,k):
        return self.comb(n+k-1,n-1)
M=10**9+7
n,k=map(int,input().split())
f=Factorial(n*2,M)
c=f.comb
h=f.homp
a=0
for i in range(min(k+1,n)):
    a+=c(n,i)*h(n-i,i)%M
print(a%M)

Solution (Python):
mod = 10**9 + 7
max_N = 2 * 10**5

def power(b, e):
	if e == 0: return 1

	half = power(b, e // 2)

	if e % 2 == 0:
		return (half * half) % mod
	else:
		return (((half * half) % mod) * b) % mod

def mod_inv(n):
	return power(n, mod - 2)

fac = [1] * (max_N + 1)
for i in range(1, max_N + 1):
	fac[i] = (i * fac[i - 1]) % mod

fac_inv = [1] * (max_N + 1)
fac_inv[max_N] = mod_inv(fac[max_N])
for i in range(max_N - 1, -1, -1):
	fac_inv[i] = (fac_inv[i + 1] * (i + 1)) % mod

def choose(n, k):
	return (((fac[n] * fac_inv[k]) % mod) * fac_inv[n - k]) % mod


if __name__ == "__main__":
	line = input().split(" ")
	N, K = int(line[0]), int(line[1])
	K = min(K, N - 1)

	ans = 0
	for i in range(K + 1):
		ans = (ans + (choose(N, i) * choose(N - i - 1 + i, i)) % mod) % mod

	print(ans)

Solution (Python):
def comb_list(n, r, mod):
    ret = [1]
    for i in range(1, r+1):
        ret.append(ret[-1] * (n+1-i) * pow(i, mod-2, mod) % mod)
    return ret

n,k = map(int,input().split())
mod = 10**9+7

r = min(k,n-1)
ans = 0
nCl = comb_list(n, r, mod)
n1Cl = comb_list(n-1, r, mod)

for i in range(r+1):
    ans = (ans + nCl[i] * n1Cl[i]) % mod

print(ans)

Solution (Python):
def comb_list(n, r, mod):
    ret = [1]
    for i in range(1, r+1):
        ret.append(ret[-1] * (n+1-i) * pow(i, mod-2, mod) % mod)
    return ret

n,k = map(int,input().split())
mod = 10**9+7

r = min(k,n-1)
ans = 0
nCl = comb_list(n, r, mod)
n1Cl = comb_list(n-1, r, mod)

for i in range(r+1):
    ans = (ans + nCl[i] * n1Cl[i]) % mod

print(ans)

Solution (Python):
def modinv(a,m):
    return pow(a,m-2,m)

n,k = map(int,input().split())
P = 10**9+7

ans = 1
nCl = 1
n1Cl = 1
for l in range(1,min(k+1,n)):
    nCl = nCl*(n+1-l)*modinv(l,P)%P
    n1Cl = n1Cl*(n-l)*modinv(l,P)%P
    ans = (ans+nCl*n1Cl)%P

print(ans)

Solution (Python):
def modinv(a,m):
    return pow(a,m-2,m)

n,k = map(int,input().split())
P = 10**9+7

ans = 1
nCl = 1
n1Cl = 1
for l in range(1,min(k+1,n)):
    nCl = nCl*(n+1-l)*modinv(l,P)%P
    n1Cl = n1Cl*(n-l)*modinv(l,P)%P
    ans = (ans+nCl*n1Cl)%P

print(ans)

Solution (Python):
n, k = map(int, input().split())
k = min(k, n - 1)
p = 10 ** 9 + 7
size = n + 100
f_list = [1] * size
if_list = [1] * size
for i in range(1, size):
    f_list[i] = (f_list[i - 1] * i) % p
    if_list[i] = pow(f_list[i], p-2, p)

def cmb(n, r):
    return f_list[n] * if_list[r] * if_list[n-r] % p


ans = 1
for i in range(1, k + 1):
    ans += cmb(n, i) * cmb(n - 1, i)
    ans = ans % p
print(ans)

Solution (Python):
from collections import defaultdict


def combination(a, b):
    if b > a - b:
        return combination(a, a - b)
    return fact[a] * ifact[b] * ifact[a-b]


MOD = 10**9+7
n, k = map(int, input().split())
k = min(k, n-1)

fact = defaultdict(int)
fact[0] = 1
for i in range(1, n+1):
    fact[i] = fact[i-1] * i
    fact[i] %= MOD

ifact = defaultdict(int)
ifact[n] = pow(fact[n], MOD-2, MOD)
for i in reversed(range(1, n + 1)):
    ifact[i-1] = ifact[i] * i
    ifact[i-1] %= MOD

ans = 0
for i in range(k+1):
    ans += combination(n, i) * combination((n-i-1)+i, i)
    ans %= MOD

print(ans % MOD)

Solution (Python):
from collections import defaultdict


def mod_inv(a, m):
    return pow(a, m-2, m)


def combination(a, b):
    if b > a - b:
        return combination(a, a - b)
    ans = fact[a] * ifact[b] * ifact[a-b]
    return ans


MOD = 10**9+7
n, k = map(int, input().split())
k = min(k, n-1)

fact = defaultdict(int)
fact[0] = 1
for i in range(1, n+1):
    fact[i] = fact[i-1] * i
    fact[i] %= MOD

ifact = defaultdict(int)
ifact[n] = mod_inv(fact[n], MOD)
for i in reversed(range(1, n + 1)):
    ifact[i-1] = ifact[i] * i
    ifact[i-1] %= MOD

ans = 0
for i in range(k+1):
    ans += combination(n, i) * combination((n-i-1)+i, i)
    ans %= MOD

print(ans % MOD)

Solution (Python):
from collections import defaultdict


def mod_inv(a, m):
    return pow(a, m-2, m)


def combination(a, b):
    if b > a - b:
        return combination(a, a - b)
    ans = fact[a] * ifact[b] * ifact[a-b]
    return ans


MOD = 10**9+7
n, k = map(int, input().split())
k = min(k, n-1)

fact = defaultdict(int)
fact[0] = 1
ifact = defaultdict(int)
ifact[0] = 1
for i in range(1, n+1):
    fact[i] = fact[i-1] * i
    fact[i] %= MOD
    ifact[i] = mod_inv(fact[i], MOD)
    ifact[i] %= MOD

ans = 0
for i in range(k+1):
    ans += combination(n, i) * combination((n-i-1)+i, i)
    ans %= MOD

print(ans % MOD)

Solution (Python):
n, k = map(int, input().split())
MOD = 1000000007


def pow(x, t):
    if t == 0:
        return 1
    a = pow(x, t >> 1)
    return a * a * (x ** (t & 1)) % MOD


fact = [0] * (n + 1)
fact[0] = 1
for i in range(1, n + 1):
    fact[i] = fact[i - 1] * i
    fact[i] %= MOD
ifact = [0] * (n + 1)
ifact[n] = pow(fact[n], MOD - 2)
for i in range(n, 0, -1):
    ifact[i - 1] = ifact[i] * i
    ifact[i - 1] %= MOD


def nCr(x, y):
    if y < 0 or y > x:
        return 0
    return fact[x] * ifact[y] * ifact[x - y] % MOD


ans = 0
for i in range(min(n - 1, k) + 1):
    ans += nCr(n, i) * nCr(n - 1, i)
    ans %= MOD
print(ans)

Solution (Python):
n,k=map(int,input().split())
def find_power(n,mod):
    
    powlist=[0]*(n+1)
    powlist[0]=1
    powlist[1]=1
    for i in range(2,n+1):
        powlist[i]=powlist[i-1]*i%(mod)
    return powlist

def find_inv_power(n):
    
    powlist=find_power(n,10**9+7)
    check=powlist[-1]
    first=1
    uselist=[0]*(n+1)
    secondlist=[0]*30
    secondlist[0]=check
    secondlist[1]=check**2
    for i in range(28):
        secondlist[i+2]=(secondlist[i+1]**2)%(10**9+7)
    a=format(10**9+5,"b")
    for j in range(30):
        if a[29-j]=="1":
            first=(first*secondlist[j])%(10**9+7)
    uselist[n]=first
    for i in range(n,0,-1):
        uselist[i-1]=(uselist[i]*i)%(10**9+7)
    return uselist
mod=10**9+7
a=find_power(4*10**5+100,mod)
b=find_inv_power(4*10**5+100)
def combi(n,r,mod):
    if n<r:
        return 0
    elif n>=r:
        return (a[n]*b[r]*b[n-r])%(mod)
if n<=k:
    K=combi(2*n-1,n-1,mod)
    print(K)
    
else:
    K=combi(2*n-1,n-1,mod)
    for i in range(k+1,n):
        K-=combi(n-1,n-i-1,mod)*combi(n,i,mod)
    print(K%mod)

Solution (Python):
import sys
import math
import itertools
from collections import deque
import copy
import bisect

MOD = 10 ** 9 + 7
INF = 10 ** 18 + 7
input = lambda: sys.stdin.readline().strip()

n, k = map(int, input().split())

inv_list = [1]  
for i in range(n):
    inv_list.append(int(MOD - inv_list[(MOD % (i + 2)) - 1] * (MOD // (i + 2)) % MOD))
if k >= n:
    ans = 1
    for i in range(2 * n - 1):
        ans *= (i + 1)
        ans %= MOD
    for i in range(n):
        ans *= inv_list[i]
        ans *= inv_list[i]
        ans %= MOD
    ans *= n
    ans %= MOD
    print(ans)
else:
    ans = 1
    nCm = 1
    nmHm = 1
    for i in range(k):
        nCm *= n - i
        nCm *= inv_list[i]
        nCm %= MOD
        nmHm *= (n - 1) - i
        nmHm *= inv_list[i]
        nmHm %= MOD
        ans += nCm * nmHm
        ans %= MOD
    print(ans)

Solution (Python):
import math

def modcomb(n, k, mod, fac, ifac):
    x = fac[n]
    y = ifac[n-k] * ifac[k] % mod
    
    return x * y % mod

def makeTables(n, mod):
    fac = [1, 1]
    ifac = [1, 1]
    inv = [0, 1]

    for i in range(2, n+1):
        fac.append(i * fac[i-1] % mod)
        inv.append(-inv[mod % i] * (mod // i) % mod)
        ifac.append(inv[i] * ifac[i-1] % mod)
    return fac, ifac

def answer(n, k, mod, fac, ifac):
    if k >= n:
        k = n - 1
    ans = 0
    for i in range(k+1):
        c = modcomb(n, i, mod, fac, ifac) * modcomb(n-1, i, mod, fac, ifac) % mod
        ans = (ans + c) % mod

    return ans

n, k = map(int, input().split())
mod = 1000000007

fac, ifac = makeTables(n, mod)
print(answer(n, k, mod, fac, ifac))

Solution (Python):
import sys
import numpy as np
sys.setrecursionlimit(10 ** 7)

N, K = map(int, input().split())
MOD = 10 ** 9 + 7




def cumprod(arr, MOD):
    L = len(arr)
    Lsq = int(L**.5+1)
    arr = np.resize(arr, Lsq**2).reshape(Lsq, Lsq)
    for n in range(1, Lsq):
        arr[:, n] *= arr[:, n-1]
        arr[:, n] %= MOD
    for n in range(1, Lsq):
        arr[n] *= arr[n-1, -1]
        arr[n] %= MOD
    return arr.ravel()[:L]


def make_fact(U, MOD):
    x = np.arange(U, dtype=np.int64)
    x[0] = 1
    fact = cumprod(x, MOD)
    x = np.arange(U, 0, -1, dtype=np.int64)
    x[0] = pow(int(fact[-1]), MOD-2, MOD)
    fact_inv = cumprod(x, MOD)[::-1]
    return fact, fact_inv


U = 10**6
fact, fact_inv = make_fact(N * 2 + 10, MOD)
fact, fact_inv = fact.tolist(), fact_inv.tolist()


def mod_comb_k(n, k, mod):
    return fact[n] * fact_inv[k] % mod * fact_inv[n - k] % mod


ans = 0
for i in range(N):
    if K < i:
        continue

    if N - 1 <= K:
        ans = mod_comb_k(N + N - 1, N - 1, MOD)
        break

    if i == 0:
        ans += 1
        continue

    a = int(mod_comb_k(N - 1, i, MOD)) * int(mod_comb_k(N, i, MOD))
    a %= MOD
    ans += a
    ans %= MOD

    '''
    a = int(fact[N]) * int(fact_inv[i]) % MOD * int(fact_inv[N - 1])
    a = a * int(fact[N-1]) % MOD * int(fact_inv[i]) % MOD * \
        int(fact_inv[N-i-1]) % MOD
    ans = (a + ans) % MOD
    '''

print(ans)

Solution (Python):
import sys
read = sys.stdin.read
readline = sys.stdin.readline
readlines = sys.stdin.readlines

MOD = 10**9+7

fac = [1, 1]
f_inv = [1, 1]
inv = [0, 1]

def prepare(n, mod):
    for i in range(2, n+1):
        fac.append((fac[-1] * i) % mod)
        inv.append((-inv[mod % i] * (mod//i)) % mod)
        f_inv.append((f_inv[-1] * inv[-1]) % mod)


def modcmb(n, r, mod):
    if n < 0 or r < 0 or r > n:
        return 0

    return fac[n] * f_inv[r] * f_inv[n-r] % mod


def main():
    N,K = map(int, readline().split())

    prepare(N, MOD)

    ans = 0
    for i in range(min(N-1, K)+1):
        ans += modcmb(N, i, MOD) * modcmb(N-1, i, MOD)
        ans %= MOD

    print(ans)


if __name__ == "__main__":
    main()

Solution (Python):
n,k=map(int,input().split())
mod=10**9+7
n_max=4*(10**5+1)
F,FI=[0]*(n_max+1),[0]*(n_max+1)
F[0],FI[0]=1,1
for i in range(n_max):
  F[i+1]=(F[i]*(i+1))%mod
FI[n_max-1]=pow(F[n_max-1],mod-2,mod)
for i in reversed(range(n_max-1)):
  FI[i]=(FI[i+1]*(i+1))%mod
def comb(x,y):
  return (F[x]*FI[x-y]*FI[y])%mod

if k>=n-1:
  print(comb(2*n-1,n))
else:
  ans=1
  for i in range(k):
    ans=(ans+comb(n,i+1)*comb(n-1,i+1))%mod

  print(ans)

Solution (Python):
import sys
import numpy as np
sys.setrecursionlimit(10 ** 7)





N, K = map(int, input().split())
MOD = 10 ** 9 + 7




def cumprod(arr, MOD):
    L = len(arr)
    Lsq = int(L**.5+1)
    arr = np.resize(arr, Lsq**2).reshape(Lsq, Lsq)
    for n in range(1, Lsq):
        arr[:, n] *= arr[:, n-1]
        arr[:, n] %= MOD
    for n in range(1, Lsq):
        arr[n] *= arr[n-1, -1]
        arr[n] %= MOD
    return arr.ravel()[:L]


def make_fact(U, MOD):
    x = np.arange(U, dtype=np.int64)
    x[0] = 1
    fact = cumprod(x, MOD)
    x = np.arange(U, 0, -1, dtype=np.int64)
    x[0] = pow(int(fact[-1]), MOD-2, MOD)
    fact_inv = cumprod(x, MOD)[::-1]
    return fact, fact_inv


U = 10**6
fact, fact_inv = make_fact(N*2+10, MOD)


def mod_comb_k(n, k, mod):
    return fact[n] * fact_inv[k] % mod * fact_inv[n - k] % mod


ans = 0
for i in range(N):
    if K < i:
        continue

    if N - 1 <= K:
        ans = mod_comb_k(N + N - 1, N - 1, MOD)
        break

    if i == 0:
        ans += 1
        continue

    a = int(mod_comb_k(N - 1, i, MOD)) * int(mod_comb_k(N, i, MOD))
    a %= MOD
    ans += a
    ans %= MOD

    '''
    a = int(fact[N]) * int(fact_inv[i]) % MOD * int(fact_inv[N - 1])
    a = a * int(fact[N-1]) % MOD * int(fact_inv[i]) % MOD * \
        int(fact_inv[N-i-1]) % MOD
    ans = (a + ans) % MOD
    '''

print(ans)