Problem:


Score : 900 points

Problem Statement

In "Takahashi-ya", a ramen restaurant, basically they have one menu: "ramen", but N kinds of toppings are also offered. When a customer orders a bowl of ramen, for each kind of topping, he/she can choose whether to put it on top of his/her ramen or not. There is no limit on the number of toppings, and it is allowed to have all kinds of toppings or no topping at all. That is, considering the combination of the toppings, 2^N types of ramen can be ordered.

Akaki entered Takahashi-ya. She is thinking of ordering some bowls of ramen that satisfy both of the following two conditions:

Do not order multiple bowls of ramen with the exactly same set of toppings.

Each of the N kinds of toppings is on two or more bowls of ramen ordered.

You are given N and a prime number M. Find the number of the sets of bowls of ramen that satisfy these conditions, disregarding order, modulo M. Since she is in extreme hunger, ordering any number of bowls of ramen is fine.

Constraints

2 \leq N \leq 3000

10^8 \leq M \leq 10^9 + 9

N is an integer.

M is a prime number.

Subscores

600 points will be awarded for passing the test set satisfying N ≤ 50.

Input

Input is given from Standard Input in the following format:

N M

Output

Print the number of the sets of bowls of ramen that satisfy the conditions, disregarding order, modulo M.

Sample Input 1

2 1000000007

Sample Output 1

2

Let the two kinds of toppings be A and B. Four types of ramen can be ordered: "no toppings", "with A", "with B" and "with A, B". There are two sets of ramen that satisfy the conditions:

The following three ramen: "with A", "with B", "with A, B".

Four ramen, one for each type.

Sample Input 2

3 1000000009

Sample Output 2

118

Let the three kinds of toppings be A, B and C. In addition to the four types of ramen above, four more types of ramen can be ordered, where C is added to the above four. There are 118 sets of ramen that satisfy the conditions, and here are some of them:

The following three ramen: "with A, B", "with A, C", "with B, C".

The following five ramen: "no toppings", "with A", "with A, B", "with B, C", "with A, B, C".

Eight ramen, one for each type.

Note that the set of the following three does not satisfy the condition: "'with A', 'with B', 'with A, B'", because C is not on any of them.

Sample Input 3

50 111111113

Sample Output 3

1456748

Remember to print the number of the sets modulo M. Note that these three sample inputs above are included in the test set for the partial score.

Sample Input 4

3000 123456791

Sample Output 4

16369789

This sample input is not included in the test set for the partial score.



Sample input:
2 1000000007


Sample output:
2

Submission status: Wrong Answer

Solution code:
import numpy as np

def prepare(n, m):
    f = 1
    for i in range(1, n + 1):
        f = f * i % m
    fn = f
    inv = [1] * (n + 1)
    f = pow(f, m - 2, m)
    inv[n] = f
    for i in range(n, 0, -1):
        f = f * i % m
        inv[i - 1] = f
    return fn, inv

def a_x(a, x, m):
    ret = 1
    yield ret
    for _ in range(x):
        ret = ret * a % m
        yield ret


N, M = map(int, input().split())

fn, inv = prepare(N, M)


stir2 = np.zeros(N + 2, dtype=np.int64)
stir2[0] = 1
upd = np.arange(2, N + 3, dtype=np.int64)

ex2 = [2]
for i in range(N):
    ex2.append(ex2[-1] ** 2 % M)

ans = 0
si = 1

for i in range(N+1):
    nCi = fn * inv[i] * inv[N-i] % M
    i_with = np.fromiter(a_x(pow(2, N-i, M), i, M), dtype=np.int64)
    i_on = (stir2[ :i+1] * i_with % M).sum() % M
    ans = (ans + nCi * (i_on % M) * (ex2[N-i] % M) * si) % M
    stir2[1 : i+2] = (stir2[1 : i+2] * upd[ :i+1] + stir2[ :i+1]) % M
    si *= -1


print(ans)

Submission status: Wrong Answer

Solution code:
from collections import defaultdict, deque, Counter
from heapq import heappush, heappop, heapify
import math
import bisect
import random
from itertools import permutations, accumulate, combinations, product
import sys
import string
from bisect import bisect_left, bisect_right
from math import factorial, ceil, floor
from operator import mul
from functools import reduce
from pprint import pprint


INF = 10 ** 20
def LI(): return list(map(int, sys.stdin.readline().split()))
def I(): return int(sys.stdin.readline())
def LS(): return sys.stdin.readline().split()
def S(): return sys.stdin.readline().strip()
def IR(n): return [I() for i in range(n)]
def LIR(n): return [LI() for i in range(n)]
def SR(n): return [S() for i in range(n)]
def LSR(n): return [LS() for i in range(n)]
def SRL(n): return [list(S()) for i in range(n)]
def MSRL(n): return [[int(j) for j in list(S())] for i in range(n)]


n, mod = LI()


fac = [1] * (n + 1)
inv = [1] * (n + 1)
for j in range(1, n + 1):
    fac[j] = fac[j-1] * j % mod


inv[n] = pow(fac[n], mod-2, mod)
for j in range(n-1, -1, -1):
    inv[j] = inv[j+1] * (j+1) % mod


def comb(n, r):
    if r > n or n < 0 or r < 0:
        return 0
    return fac[n] * inv[n - r] * inv[r] % mod



sterling = [[0] * (n + 1) for _ in range(n + 1)]
sterling[0][0] = 1
# n個の区別できるものをk個の(n個以下)の区別不可能なグループに分けることに(i, j)が対応
# n >= k
for i in range(1, n + 1):
    for j in range(i + 1):
        sterling[i][j] = sterling[i - 1][j] * (j + 1)
        if j:
            sterling[i][j] += sterling[i - 1][j - 1]
        sterling[i][j] %= mod



ans = pow(2, pow(2, n, mod), mod)
for l in range(1, n + 1):
    ret = 0
    for k in range(l + 1):
        cumsum = (sterling[l - 1][k - 1] + (k + 1) * sterling[l - 1][k]) % mod
        cumsum = cumsum * pow(pow(2, n - l, mod), k, mod) % mod
        cumsum = cumsum * pow(2, pow(2, n - l, mod), mod) % mod
        cumsum = cumsum * comb(n, l) % mod
        ret += cumsum
    if l % 2:
        ans -= ret
    else:
        ans += ret
    ans %= mod

print(ans)

Submission status: Wrong Answer

Solution code:
from collections import defaultdict, deque, Counter
from heapq import heappush, heappop, heapify
import math
import bisect
import random
from itertools import permutations, accumulate, combinations, product
import sys
import string
from bisect import bisect_left, bisect_right
from math import factorial, ceil, floor
from operator import mul
from functools import reduce
from pprint import pprint


INF = 10 ** 20
def LI(): return list(map(int, sys.stdin.readline().split()))
def I(): return int(sys.stdin.readline())
def LS(): return sys.stdin.readline().split()
def S(): return sys.stdin.readline().strip()
def IR(n): return [I() for i in range(n)]
def LIR(n): return [LI() for i in range(n)]
def SR(n): return [S() for i in range(n)]
def LSR(n): return [LS() for i in range(n)]
def SRL(n): return [list(S()) for i in range(n)]
def MSRL(n): return [[int(j) for j in list(S())] for i in range(n)]


n, mod = LI()


fac = [1] * (n + 1)
inv = [1] * (n + 1)
for j in range(1, n + 1):
    fac[j] = fac[j-1] * j % mod


inv[n] = pow(fac[n], mod-2, mod)
for j in range(n-1, -1, -1):
    inv[j] = inv[j+1] * (j+1) % mod


def comb(n, r):
    if r > n or n < 0 or r < 0:
        return 0
    return fac[n] * inv[n - r] * inv[r] % mod



sterling = [[0] * (n + 2) for _ in range(n + 2)]
sterling[1][1] = 1
# n個の区別できるものをk個の(n個以下)の区別不可能なグループに分けることに(i, j)が対応
# n >= k
for i in range(2, n + 2):
    for j in range(1, i + 1):
        sterling[i][j] = (j * sterling[i - 1][j] % mod + sterling[i - 1][j - 1]) % mod




ans = pow(2, pow(2, n, mod), mod)
for l in range(1, n + 1):
    ret = 0
    for k in range(0, l + 1):
        cumsum = (sterling[l][k] + (k + 1) * sterling[l][k + 1]) % mod
        cumsum = cumsum * pow(pow(2, n - l, mod), k, mod) % mod
        cumsum = cumsum * pow(2, pow(2, n - l, mod), mod) % mod
        cumsum = cumsum * comb(n, l) % mod
        ret += cumsum
    if l % 2:
        ans -= ret
    else:
        ans += ret
    ans %= mod

print(ans)

Submission status: Wrong Answer

Solution code:
n, m = map(int, input().split())
MOD = m
list_size = 10000

f_list = [1] * list_size
f_r_list = [1] * list_size

for i in range(list_size - 1):
	f_list[i + 1] = int((f_list[i] * (i + 2)) % MOD)

f_r_list[-1] = pow(f_list[-1], MOD - 2, MOD)

for i in range(2, list_size + 1):
	f_r_list[-i] = int((f_r_list[-i + 1] * (list_size + 2 - i)) % MOD)

def comb(N, r):
	if N < r:
		return 0
	elif N == 0 or r == 0 or N == r:
		return 1
	else:
		return (((f_list[N - 1] * f_r_list[N - r - 1]) % MOD) * f_r_list[r - 1]) % MOD 

st = [[0 for _ in range(n+1)] for _ in range(n+1)]
st[0][0] = 1
for i in range(1, n+1):
	st[i][0] = 1
	for j in range(1, i+1):
		st[i][j] = (st[i-1][j-1] + (j+1) * st[i-1][j]) % MOD

ans = 0
for i in range(n+1):
	res = (comb(n, i) * pow(2, pow(2, n-i, MOD), MOD)) % MOD
	#print(res)
	tmp = 0
	for j in range(i+1):
		tmp += st[i][j] * pow(2, (n-i) * j, MOD)
		tmp %= MOD
	res *= tmp
	res %= MOD

	if i%2 == 0:
		ans += res
	else:
		ans -= res
	ans %= MOD

print(ans)

Submission status: Wrong Answer

Solution code:
from collections import defaultdict, deque, Counter
from heapq import heappush, heappop, heapify
import math
import bisect
import random
from itertools import permutations, accumulate, combinations, product
import sys
import string
from bisect import bisect_left, bisect_right
from math import factorial, ceil, floor
from operator import mul
from functools import reduce
from pprint import pprint
sys.setrecursionlimit(2147483647)


INF = float('inf')
def LI(): return list(map(int, sys.stdin.readline().split()))
def I(): return int(sys.stdin.readline())
def LS(): return sys.stdin.readline().split()
def S(): return sys.stdin.readline().strip()
def IR(n): return [I() for i in range(n)]
def LIR(n): return [LI() for i in range(n)]
def SR(n): return [S() for i in range(n)]
def LSR(n): return [LS() for i in range(n)]
def SRL(n): return [list(S()) for i in range(n)]
def MSRL(n): return [[int(j) for j in list(S())] for i in range(n)]
mod = 1000000007



from collections import defaultdict, deque, Counter
from heapq import heappush, heappop, heapify
import math
import bisect
import random
from itertools import permutations, accumulate, combinations, product
import sys
import string
from bisect import bisect_left, bisect_right
from math import factorial, ceil, floor
from operator import mul
from functools import reduce


INF = float('inf')
def LI(): return list(map(int, sys.stdin.readline().split()))
def I(): return int(sys.stdin.readline())
def LS(): return sys.stdin.readline().split()
def S(): return sys.stdin.readline().strip()
def IR(n): return [I() for i in range(n)]
def LIR(n): return [LI() for i in range(n)]
def SR(n): return [S() for i in range(n)]
def LSR(n): return [LS() for i in range(n)]
def SRL(n): return [list(S()) for i in range(n)]
def MSRL(n): return [[int(j) for j in list(S())] for i in range(n)]
mod = 1000000007


n, mod = LI()

def mod_combination(n, r, mod=mod):
    if r > n:
        return 0
    fact = 1
    for i in range(n - r + 1, n + 1):
        fact = fact * i % mod
    divisor = 1
    for j in range(1, r + 1):
        divisor = divisor * j % mod
    return fact * pow(divisor, mod - 2, mod) % mod



sterling_number = [[0] * (n + 2) for _ in range(n + 2)]
for i in range(1, n + 2):
    sterling_number[i][1] = 1
# n個の区別できるものをk個の(n個以下)の区別不可能なグループに分けることに(i, j)が対応
# n >= k
for i in range(1, n + 2):
    for j in range(2, i + 1):
        sterling_number[i][j] = (j * sterling_number[i - 1][j] % mod + sterling_number[i - 1][j - 1]) % mod




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




n_combination_list = [1] * (n + 1)
for i in range(1, n + 1):
    n_combination_list[i] = factorial_list[n] * pow(factorial_list[n - i], mod - 2, mod) % mod * pow(factorial_list[i], mod - 2, mod) % mod



ans = pow(2, pow(2, n, mod), mod)


for l in range(1, n + 1):
    ret = 0
    for k in range(0, l + 1):
        cumsum = (sterling_number[l][k] + (k + 1) * sterling_number[l][k + 1] % mod) % mod
        cumsum = cumsum * pow(pow(2, n - l, mod), k, mod) % mod
        cumsum = cumsum * pow(2, pow(2, n - l, mod), mod) % mod
        ret = (ret + cumsum) % mod
    ret = ret * n_combination_list[l] % mod
    if l % 2:
        ret *= -1
    ans = (ans + ret) % mod



print(ans)

Submission status: Wrong Answer

Solution code:
import sys
input = sys.stdin.readline
import numpy as np

N,MOD = map(int,input().split())

"""
余事象を調べる。包除の原理を使う。
A[n] = （1,2,...,n）が1杯以下、他は何でも良い
B[n,l] : (1,2,...,n) をlグループに分ける方法の個数
A[n]
・0杯のグループあり
・なし
"""

B = np.zeros((N+1,N+1), dtype=np.int64)
B[0,0] = 1
for n in range(1,N+1):
    # 1番を単独で使う
    B[n,1:] = B[n-1,:-1]
    # 1番をどこかに混ぜてもらう
    B[n,1:] += B[n-1,1:] * np.arange(1,N+1) % MOD
    B[n] %= MOD

# 2^{kl}
pow_2 = np.ones((N+1,N+1), dtype=np.int64)
for n in range(1,N+1):
    pow_2[1,n] = 2 * pow_2[1,n-1] % MOD
for n in range(2,N+1):
    pow_2[n] = pow_2[n-1] * pow_2[1] % MOD

A = np.zeros(N+1, dtype=np.int64)
for n in range(N+1):
    A[n] = (pow(2,pow(2,N-n,MOD-1),MOD) * B[n,1:] % MOD * (pow_2[N-n,1:] + pow_2[N-n,:-1] * np.arange(1,N+1) % MOD) % MOD).sum() % MOD

comb = np.zeros((N+1,N+1),dtype = np.int64)
comb[:,0] = 1
for n in range(1,N+1):
    comb[n,1:] = (comb[n-1,1:] + comb[n-1,:-1]) % MOD

A[::2] *= (-1)
A *= comb[N]
answer = pow(2,pow(2,N,MOD-1),MOD) - A.sum()
answer %= MOD
print(answer)

Submission status: Wrong Answer

Solution code:
N,M = list(map(int,input().split()))

table = [1,1]
while len(table) <= N:
    temp = 1
    for i in range(len(table)-1):
        table[i+1] += temp
        temp = table[i+1]- temp
        table[i+1] = M
    table.append(1)

S = [1]

rev2 = pow(2, M-2, M)
base = pow(2, N, M)
ans = 0
S = [1]
for K in range(N+1):
    res = table[K] % M
    res = (res * pow(2, pow(2, N - K, M-1), M)) % M
    b = 1
    v = 0
    T = [0]*(K+2)
    for L in range(K):
        T[L+1] = s = (S[L] + (L+1)*S[L+1]) % M
        v += s * b
        b = (b * base) % M
    v += b
    T[K+1] = 1
    S = T
    res = (res * v) % M
    if K % 2:
        ans -= res
    else:
        ans += res
    ans %= M

    base = (base * rev2) % M
print(ans)

Submission status: Wrong Answer

Solution code:
N, M = map(int, input().split())
fact = [1]*(N+1)
rfact = [1]*(N+1)
for i in range(1, N+1):
    fact[i] = r = (i * fact[i-1]) % M
    rfact[i] = pow(r, M-2, M)


memo = {}
def f(N, K):
    if (N, K) in memo:
        return memo[N, K]
    if K == 0:
        return N == 0
    elif K == 1 or N == K:
        return 1
    r = memo[N, K] = (f(N-1, K-1) + K*f(N-1, K)) % M
    return r

S = [1]

rev2 = pow(2, M-2, M)
base = pow(2, N, M)
ans = 0
S = [1]
for K in range(N+1):
    r = (fact[N] * rfact[K] * rfact[N-K]) % M
    r = (r * pow(2, base, M)) % M
    b = 1
    v = 0
    T = [0]*(K+2)
    for L in range(K):
        T[L+1] = s = (S[L] + (L+1)*S[L+1]) % M
        v += s * b
        b = (b * base) % M
    v += b
    T[K+1] = 1
    S = T
    r = (r * v) % M
    if K % 2:
        ans -= r
    else:
        ans += r
    ans %= M

    base = (base * rev2) % M
print(ans)

Submission status: Time Limit Exceeded

Solution code:
n, m = map(int, input().split())
MOD = m
list_size = n+2

f_list = [1] * list_size
f_r_list = [1] * list_size

for i in range(list_size - 1):
	f_list[i + 1] = int((f_list[i] * (i + 1)) % MOD)

f_r_list = [pow(x, MOD-2, MOD) for x in f_list]

def comb(N, r):
	if N < r or r < 0:
		return 0
	else:
		return (((f_list[N] * f_r_list[N - r]) % MOD) * f_r_list[r]) % MOD 

st = [[0 for _ in range(n+1)] for _ in range(n+1)]
st[0][0] = 1
for i in range(1, n+1):
	st[i][0] = 1
	for j in range(1, i+1):
		st[i][j] = (st[i-1][j-1] + (j+1) * st[i-1][j]) % MOD

ans = 0
for i in range(n+1):
	res = (comb(n, i) * pow(2, pow(2, n-i, MOD-1), MOD)) % MOD
	tmp = 0
	for j in range(i+1):
		tmp += st[i][j] * pow(2, (n-i) * j, MOD)
		tmp %= MOD
	res *= tmp
	res %= MOD

	if i%2 == 0:
		ans += res
	else:
		ans -= res
	ans %= MOD

print(ans)

Submission status: Time Limit Exceeded

Solution code:
import sys

sys.setrecursionlimit(10 ** 6)
input = sys.stdin.readline

def main():
    def com(com_n, com_r):
        return fac[com_n] * inv[com_r] * inv[com_n - com_r] % md

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

    # combinationの準備
    n_max = n+3
    fac = [1]
    inv = [1] * (n_max + 1)
    k_fac_inv = 1
    for i in range(1, n_max + 1):
        k_fac_inv = k_fac_inv * i % md
        fac.append(k_fac_inv)
    k_fac_inv = pow(k_fac_inv, md - 2, md)
    for i in range(n_max, 1, -1):
        inv[i] = k_fac_inv
        k_fac_inv = k_fac_inv * i % md

    # スターリング数発展の作成
    stir = [[1] * (n + 1) for _ in range(n + 1)]
    si1 = stir[0]
    for i in range(1, n + 1):
        si = stir[i]
        for j in range(1, i):
            si[j] = si1[j] * (j + 1) + si1[j - 1]
        si1 = si

    # 本体
    ans = 0
    pm = 1
    for m in range(n + 1):
        a = 0
        e2 = pow(2, n - m)
        for k in range(m + 1):
            a = (a + stir[m][k] * pow(e2, k, md)) % md
        ans = (ans + pm * a * pow(2, e2, md) * com(n, m)) % md
        pm *= -1
    print(ans)

main()

Submission status: Time Limit Exceeded

Solution code:
# coding: utf-8
# Your code here!

import sys
sys.setrecursionlimit(10**6)
readline = sys.stdin.readline #文字列入力のときは注意

n,MOD = [int(i) for i in readline().split()]

SIZE=3001; #MOD=10**9+7 #998244353 #ここを変更する

inv = [0]*(SIZE+1)# inv[j] = j^{-1} mod MOD
fac = [0]*(SIZE+1)# fac[j] = j! mod MOD
finv = [0]*(SIZE+1)# finv[j] = (j!)^{-1} mod MOD
inv[1] = 1
fac[0] = fac[1] = 1
finv[0] = finv[1] = 1
for i in range(2,SIZE+1):
    inv[i] = MOD -(MOD//i)*inv[MOD%i]%MOD
    fac[i] = fac[i-1]*i%MOD
    finv[i]= finv[i-1]*inv[i]%MOD

def choose(n,r): # nCk mod MOD の計算
    if 0 <= r <= n:
        return (fac[n]*finv[r]%MOD)*finv[n-r]%MOD
    else:
        return 0


"""
make the table of Sterling numbers of the second kind
Sterling[ball][box]
"""
SIZE = n
Sterling2 = [[0]*(SIZE+1) for _ in range(SIZE+1)]
Sterling2[0][0] = Sterling2[1][1] = 1
for i in range(2,SIZE+1):
    for k in range(1,i+1):
        Sterling2[i][k] = (Sterling2[i-1][k-1] + k*Sterling2[i-1][k])%MOD





SIZE = n
f = [[0]*(SIZE+1) for _ in range(SIZE+1)]
f[0][0] = 1
for i in range(1,SIZE+1):
    for k in range(i+1):
        f[i][k] = (f[i-1][k] + f[i-1][k-1] + k*f[i-1][k])%MOD

#print(f)

def size(s):
    res = 0
    for b in range(n+1):
        res += f[s][b]*pow(2,b*(n-s),MOD)
    return res*pow(2,pow(2,n-s,MOD-1),MOD)%MOD

#for i in range(n+1):
#    print(size(i))

ans = 0
sgn = 1
for j in range(n+1):
    ans += sgn*choose(n,j)*size(j)
    ans %= MOD
    sgn *= -1

print(ans)

Submission status: Time Limit Exceeded

Solution code:
# coding: utf-8
# Your code here!

import sys
sys.setrecursionlimit(10**6)
readline = sys.stdin.readline #文字列入力のときは注意

n,MOD = [int(i) for i in readline().split()]

SIZE=3001; #MOD=10**9+7 #998244353 #ここを変更する

inv = [0]*(SIZE+1)# inv[j] = j^{-1} mod MOD
fac = [0]*(SIZE+1)# fac[j] = j! mod MOD
finv = [0]*(SIZE+1)# finv[j] = (j!)^{-1} mod MOD
inv[1] = 1
fac[0] = fac[1] = 1
finv[0] = finv[1] = 1
for i in range(2,SIZE+1):
    inv[i] = MOD -(MOD//i)*inv[MOD%i]%MOD
    fac[i] = fac[i-1]*i%MOD
    finv[i]= finv[i-1]*inv[i]%MOD

def choose(n,r): # nCk mod MOD の計算
    if 0 <= r <= n:
        return (fac[n]*finv[r]%MOD)*finv[n-r]%MOD
    else:
        return 0


"""
make the table of Sterling numbers of the second kind
Sterling[ball][box]
"""
SIZE = n
Sterling2 = [[0]*(SIZE+1) for _ in range(SIZE+1)]
Sterling2[0][0] = Sterling2[1][1] = 1
for i in range(2,SIZE+1):
    for k in range(1,i+1):
        Sterling2[i][k] = (Sterling2[i-1][k-1] + k*Sterling2[i-1][k])%MOD

"""
def M(s,u):
    res = 0
    for b in range(s+1):
        res += Sterling2[u][b]*pow(2,b*(n-s))
    return res%MOD
"""    

def size(s):
    res = 0
    for j in range(s+1):
        res2 = 0
        for b in range(s+1):
            res2 += Sterling2[j][b]*pow(2,b*(n-s),MOD)%MOD
        res2 %= MOD
        res += choose(s,j)*res2%MOD
    return res*pow(2,pow(2,n-s,MOD-1),MOD)%MOD

#for i in range(n+1):
#    print(size(i))

ans = 0
sgn = 1
for j in range(n+1):
    ans += sgn*choose(n,j)*size(j)
    ans %= MOD
    sgn *= -1

print(ans)

Submission status: Time Limit Exceeded

Solution code:
# coding: utf-8
# Your code here!

import sys
sys.setrecursionlimit(10**6)
readline = sys.stdin.readline #文字列入力のときは注意

n,MOD = [int(i) for i in readline().split()]

SIZE=3001; #MOD=10**9+7 #998244353 #ここを変更する

inv = [0]*(SIZE+1)# inv[j] = j^{-1} mod MOD
fac = [0]*(SIZE+1)# fac[j] = j! mod MOD
finv = [0]*(SIZE+1)# finv[j] = (j!)^{-1} mod MOD
inv[1] = 1
fac[0] = fac[1] = 1
finv[0] = finv[1] = 1
for i in range(2,SIZE+1):
    inv[i] = MOD -(MOD//i)*inv[MOD%i]%MOD
    fac[i] = fac[i-1]*i%MOD
    finv[i]= finv[i-1]*inv[i]%MOD

def choose(n,r): # nCk mod MOD の計算
    if 0 <= r <= n:
        return (fac[n]*finv[r]%MOD)*finv[n-r]%MOD
    else:
        return 0


"""
make the table of Sterling numbers of the second kind
Sterling[ball][box]
"""
SIZE = n
Sterling2 = [[0]*(SIZE+1) for _ in range(SIZE+1)]
Sterling2[0][0] = 1
for i in range(1,SIZE+1):
    for k in range(i+1):
        Sterling2[i][k] = (Sterling2[i-1][k-1] + k*Sterling2[i-1][k])%MOD

#sumS =[sum(Sterling2[u]) for u in range(n+1)]
"""
def M(s,u):
    res = 0
    for b in range(s+1):
        res += Sterling2[u][b]*pow(2,b*(n-s))
    return res%MOD
"""    

def size(s):
    res = 0
    for j in range(s+1):
        res2 = 0
        for b in range(s+1):
            res2 += Sterling2[j][b]*pow(2,b*(n-s))
        res2 %= MOD
        res += choose(s,j)*res2%MOD
    return res*pow(2,pow(2,n-s,MOD-1),MOD)%MOD

#for i in range(n+1):
#    print(size(i))

ans = 0
sgn = 1
for j in range(n+1):
    ans += sgn*choose(n,j)*size(j)
    ans %= MOD
    sgn *= -1

print(ans)

Submission status: Time Limit Exceeded

Solution code:
# coding: utf-8
# Your code here!

import sys
sys.setrecursionlimit(10**6)
readline = sys.stdin.readline #文字列入力のときは注意

n,MOD = [int(i) for i in readline().split()]

SIZE=3001; #MOD=10**9+7 #998244353 #ここを変更する

inv = [0]*(SIZE+1)# inv[j] = j^{-1} mod MOD
fac = [0]*(SIZE+1)# fac[j] = j! mod MOD
finv = [0]*(SIZE+1)# finv[j] = (j!)^{-1} mod MOD
inv[1] = 1
fac[0] = fac[1] = 1
finv[0] = finv[1] = 1
for i in range(2,SIZE+1):
    inv[i] = MOD -(MOD//i)*inv[MOD%i]%MOD
    fac[i] = fac[i-1]*i%MOD
    finv[i]= finv[i-1]*inv[i]%MOD

def choose(n,r): # nCk mod MOD の計算
    if 0 <= r <= n:
        return (fac[n]*finv[r]%MOD)*finv[n-r]%MOD
    else:
        return 0


"""
make Sterling_table
Sterling[ball][box]
"""
SIZE = n+5
Sterling2 = [[0]*(SIZE+1) for _ in range(SIZE+1)]
Sterling2[0][0] = 1
for i in range(1,SIZE+1):
    for k in range(i+1):
        Sterling2[i][k] = (Sterling2[i-1][k-1] + k*Sterling2[i-1][k])%MOD

def M(s,u):
    res = 0
    p = pow(2,n-s,MOD-1)
    for b in range(n+1):
        res += Sterling2[u][b]*pow(2,b*(n-s)+p,MOD)%MOD
    return res%MOD
    

def size(s):
    res = 0
    for j in range(s+1):
        res += choose(s,j)*M(s,j)%MOD
    return res%MOD

#for i in range(n+1):
#    print(size(i))

ans = 0
sgn = 1
for j in range(n+1):
    ans += sgn*choose(n,j)*size(j)
    ans %= MOD
    sgn *= -1

print(ans)

Submission status: Time Limit Exceeded

Solution code:

N, M = map(int, input().split())
fact = {i: None for i in range(N+1)}
fact_inverse = {i: None for i in range(N+1)}

fact[0] = fact[1] = 1
fact_inverse[0] = fact_inverse[1] = 1
for i in range(2,N+1):
    fact[i] = (fact[i-1]*i) % M
    fact_inverse[i] = pow(fact[i], M-2, M)
    
def binom(n, r):
    '''
    Calculate nCr = n! /(r! (n-r)!) % mod
    '''
    if n < r or n < 0 or r < 0:
        return 0
    else:
        return fact[n] * (fact_inverse[r] * fact_inverse[n-r] % M) % M

ans = 0

# Preprocess
# calculate 2**(2**n)
pow_pow2 = [0] * (N+1)
pow_pow2[0] = 2
for i in range(1, N+1):
    pow_pow2[i] = pow(pow_pow2[i-1], 2, M)

# (#ways to order i ramens, with each of n toppings chosen at most once)
ways2 = {n: {i: 0 for i in range(N+1)} for n in range(N+1)}
for n in range(N+1):
    ways = 0 # (#ways to choose n toppings, with each topping chosen at most once)
    temp = 0
    if n >= 1:
        for i in range(n+1):
            ways2[n][i] = (ways2[n][i] +  ways2[n-1][i]) % M # n-th topping is not used
            if i >= 1: ways2[n][i] = (ways2[n][i] + ways2[n-1][i-1]) % M # only n-th topping is used in an order
            ways2[n][i] = (ways2[n][i] + (i*ways2[n-1][i]) % M) % M # n-th topping is chosen with another topping
            temp = (temp + (ways2[n][i] * pow(2, (N-n)*i, M)) % M) % M # degree of freedom on each order is 2**(N-n)
    else: # n = 0
        ways2[n][0] = 1
        temp += ways2[n][0]
    ways = (temp * pow_pow2[N-n]) % M
    if n % 2 == 0:
        ans = (ans + (binom(N, n) * ways) % M) % M
    else:
        ans = (ans - (binom(N, n) * ways) % M) % M
print(ans)

Submission status: Time Limit Exceeded

Solution code:
class Combinatorics:
    def __init__(self, N, mod):
        '''
        Preprocess for calculating binomial coefficients nCr (0 <= r <= n, 0 <= n <= N)
        over the finite field Z/(mod)Z.
        Input:
            N (int): maximum n
            mod (int): a prime number. The order of the field Z/(mod)Z over which nCr is calculated.
        '''
        self.mod = mod
        self.fact = {i: None for i in range(N+1)}         # n!
        self.inverse = {i: None for i in range(1, N+1)}   # inverse of n in the field Z/(MOD)Z
        self.fact_inverse = {i: None for i in range(N+1)} # inverse of n! in the field Z/(MOD)Z
        
        # preprocess
        self.fact[0] = self.fact[1] = 1
        self.fact_inverse[0] = self.fact_inverse[1] = 1
        self.inverse[1] = 1
        for i in range(2, N+1):
            self.fact[i] = i * self.fact[i-1] % self.mod
            q, r = divmod(self.mod, i)
            self.inverse[i] = (- (q % self.mod) * self.inverse[r]) % self.mod
            self.fact_inverse[i] = self.inverse[i] * self.fact_inverse[i-1] % self.mod
    
    def perm(self, n, r):
        '''
        Calculate nPr = n! / (n-r)! % mod
        '''
        if n < r or n < 0 or r < 0:
            return 0
        else:
            return (self.fact[n] * self.fact_inverse[n-r]) % self.mod
    
    def binom(self, n, r):
        '''
        Calculate nCr = n! /(r! (n-r)!) % mod
        '''
        if n < r or n < 0 or r < 0:
            return 0
        else:
            return self.fact[n] * (self.fact_inverse[r] * self.fact_inverse[n-r] % self.mod) % self.mod
        
    def hom(self, n, r):
        '''
        Calculate nHr = {n+r-1}Cr % mod.
        Assign r objects to one of n classes.
        Arrangement of r circles and n-1 partitions:
            o o o | o o | | | o | | | o o | | o
        '''
        if n == 0 and r > 0:
            return 0
        if n >= 0 and r == 0:
            return 1
        return self.binom(n + r - 1, r)

N, M = map(int, input().split())
com = Combinatorics(N, M)

# Preprocess
# calculate 2**n and 2**(2**n)
pow2 = [0] * (N*N + 1)
pow_pow2 = [0] * (N+1)
pow2[0] = 1; pow_pow2[0] = 2
for i in range(1, N*N + 1):
    pow2[i] = (pow2[i-1] * 2) % M
for i in range(1, N+1):
    pow_pow2[i] = pow(pow_pow2[i-1], 2, M)

# (#ways to choose n toppings, with each topping chosen at most once)
ways = [0]*(N+1)
# (#ways to order i ramens, with each of n toppings chosen at most once)
ways2 = {n: {i: 0 for i in range(N+1)} for n in range(N+1)}
for n in range(N+1):
    temp = 0
    if n >= 1:
        for i in range(n+1):
            ways2[n][i] = (ways2[n][i] +  ways2[n-1][i]) % M # n-th topping is not used
            if i >= 1: ways2[n][i] = (ways2[n][i] + ways2[n-1][i-1]) % M # only n-th topping is used in an order
            ways2[n][i] = (ways2[n][i] + (i*ways2[n-1][i]) % M) % M # n-th topping is chosen with another topping
            temp = (temp + (ways2[n][i] * pow2[(N-n)*i]) % M) % M # degree of freedom on each order is 2**(N-n)
    else: # n = 0
        ways2[n][0] = 1
        temp += ways2[n][0]
    ways[n] = (temp * pow_pow2[N-n]) % M

ans = 0
for n in range(N+1):
    if n % 2 == 0:
        ans = (ans + (com.binom(N, n) * ways[n]) % M) % M
    else:
        ans = (ans - (com.binom(N, n) * ways[n]) % M) % M
print(ans)

Submission status: Time Limit Exceeded

Solution code:
class Combinatorics:
    def __init__(self, N, mod):
        '''
        Preprocess for calculating binomial coefficients nCr (0 <= r <= n, 0 <= n <= N)
        over the finite field Z/(mod)Z.
        Input:
            N (int): maximum n
            mod (int): a prime number. The order of the field Z/(mod)Z over which nCr is calculated.
        '''
        self.mod = mod
        self.fact = {i: None for i in range(N+1)}         # n!
        self.inverse = {i: None for i in range(1, N+1)}   # inverse of n in the field Z/(MOD)Z
        self.fact_inverse = {i: None for i in range(N+1)} # inverse of n! in the field Z/(MOD)Z
        
        # preprocess
        self.fact[0] = self.fact[1] = 1
        self.fact_inverse[0] = self.fact_inverse[1] = 1
        self.inverse[1] = 1
        for i in range(2, N+1):
            self.fact[i] = i * self.fact[i-1] % self.mod
            q, r = divmod(self.mod, i)
            self.inverse[i] = (- (q % self.mod) * self.inverse[r]) % self.mod
            self.fact_inverse[i] = self.inverse[i] * self.fact_inverse[i-1] % self.mod
    
    def perm(self, n, r):
        '''
        Calculate nPr = n! / (n-r)! % mod
        '''
        if n < r or n < 0 or r < 0:
            return 0
        else:
            return (self.fact[n] * self.fact_inverse[n-r]) % self.mod
    
    def binom(self, n, r):
        '''
        Calculate nCr = n! /(r! (n-r)!) % mod
        '''
        if n < r or n < 0 or r < 0:
            return 0
        else:
            return self.fact[n] * (self.fact_inverse[r] * self.fact_inverse[n-r] % self.mod) % self.mod
        
    def hom(self, n, r):
        '''
        Calculate nHr = {n+r-1}Cr % mod.
        Assign r objects to one of n classes.
        Arrangement of r circles and n-1 partitions:
            o o o | o o | | | o | | | o o | | o
        '''
        if n == 0 and r > 0:
            return 0
        if n >= 0 and r == 0:
            return 1
        return self.binom(n + r - 1, r)

N, M = map(int, input().split())
com = Combinatorics(N, M)

# Preprocess
ways = [0]*(N+1)
# (#ways to order i ramens, with each of n toppings chosen at most once)
ways2 = {n: {i: 0 for i in range(N+1)} for n in range(N+1)}
for n in range(N+1):
    temp = 0
    if n >= 1:
        for i in range(n+1):
            ways2[n][i] = (ways2[n][i] +  ways2[n-1][i]) % M # n-th topping is not used
            if i >= 1: ways2[n][i] = (ways2[n][i] + ways2[n-1][i-1]) % M # only n-th topping is used in an order
            ways2[n][i] = (ways2[n][i] + (i*ways2[n-1][i]) % M) % M # n-th topping is chosen with another topping
            temp = (temp + (ways2[n][i] * pow(2, (N-n) * i, M)) % M) % M # degree of freedom on each order is 2**(N-n)
    else: # n = 0
        ways2[n][0] = 1
        temp += ways2[n][0]
    ways[n] = temp * pow(2, 2**(N-n), M)

ans = 0
for n in range(N+1):
    if n % 2 == 0:
        ans = (ans + (com.binom(N, n) * ways[n]) % M) % M
    else:
        ans = (ans - (com.binom(N, n) * ways[n]) % M) % M
print(ans)

Submission status: Time Limit Exceeded

Solution code:
N,M = list(map(int,input().split()))

table = [1,1]
while len(table) <= N:
    temp = table[0]
    for i in range(len(table)-1):
        table[i+1] += temp
        temp = table[i+1]- temp
        table[i+1] %= M
    table.append(1)


S = [1]

rev2 = pow(2, M-2, M)
base = pow(2, N, M)
ans = 0
S = [1]
for K in range(N+1):
    res = table[K] % M
    res = (res * pow(2, pow(2, N - K, M-1), M)) % M
    b = 1
    v = 0
    T = [0]*(K+2)
    for L in range(K):
        T[L+1] = s = (S[L] + (L+1)*S[L+1]) % M
        v += s * b
        b = (b * base) % M
    v += b
    T[K+1] = 1
    S = T
    res = (res * v) % M
    if K % 2:
        ans -= res
    else:
        ans += res
    ans %= M

    base = (base * rev2) % M
print(ans)