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: Time Limit Exceeded

Solution code:
# https://arc096.contest.atcoder.jp/tasks/arc096_c

import itertools

n, m = input().split(" ")
n = int(n)
m = int(m)

p = 2 ** n

pt = []
for i in range(p):
    l = format(i, "b").zfill(n)
    pt.append(list(map(int, list(l))))

# print(pt)

ok = 0

for i in range(1, p + 1):
    # print("i = " + str(i))
    for c in itertools.combinations(pt, i):
        eat = [0] * n

        # print(c)
        li = list(map(list, c))
        # print(li)
        for ramen in li:
            for e_index in range(len(ramen)):
                eat[e_index] += ramen[e_index]
        # print("eat = " + str(eat))
        enough = True
        for v in eat:
            if v < 2:
                enough = False
                break
        if enough:
            ok += 1

print(ok % m)

Submission status: Memory 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 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

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

# Preprocess
# calculate 2**n and 2**(2**n)
pow2 = [0] * (N*N // 4 + 1)
pow_pow2 = [0] * (N+1)
pow2[0] = 1; pow_pow2[0] = 2
for i in range(1, N*N // 4 + 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 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] * 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 = (temp * pow_pow2[N-n]) % M
    if n % 2 == 0:
        ans = (ans + (com.binom(N, n) * ways) % M) % M
    else:
        ans = (ans - (com.binom(N, n) * ways) % M) % M
print(ans)

Submission status: Memory 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)
ans = 0

# 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 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] * 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 = (temp * pow_pow2[N-n]) % M
    if n % 2 == 0:
        ans = (ans + (com.binom(N, n) * ways) % M) % M
    else:
        ans = (ans - (com.binom(N, n) * ways) % M) % M
print(ans)

Submission status: Accepted

Solution code:
import numpy as np


def solve(n, m):
    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

    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

    return ans


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

Submission status: Accepted

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: Accepted

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 #ここを変更する

SIZE += 1
inv = [0]*SIZE # inv[j] = j^{-1} mod MOD
fac = [0]*SIZE # fac[j] = j! mod MOD
finv = [0]*SIZE# 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):
    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+1
f = [[0]*(SIZE) for _ in range(SIZE)]
f[0][0] = 1
for i in range(1,SIZE):
    f[i][0] = 1
    for k in range(1,i+1):
        f[i][k] = (f[i-1][k] + f[i-1][k-1] + k*f[i-1][k])%MOD

#print(f)

#pow(2,n-s,MOD-1)
pp = [1]*SIZE
for i in range(1,SIZE):
    pp[i] = 2*pp[i-1]%(MOD-1)

def size(s):
    """
    for b in range(n+1):
        res += f[s][b]*pow(2,b*(n-s),MOD)
    """
    res = 0
    p = 1
    e = pow(2,(n-s),MOD)
    for b in range(SIZE):
        res += f[s][b]*p
        p *= e
        p %= MOD
    return res*pow(2,pp[n-s],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: Accepted

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

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

A = (B[:,1:] * (pow_2[::-1,1:] + pow_2[::-1,:-1] * np.arange(1,N+1) % MOD) % MOD).sum(axis = 1)
A %= MOD
A *= pow_pow[::-1]
A %= MOD

fact = [1] * (N+1)
fact_inv = [1] * (N+1)
for n in range(1,N+1):
    fact[n] = fact[n-1] * n % MOD
fact_inv[N] = pow(fact[N],MOD-2,MOD)
for n in range(N,0,-1):
    fact_inv[n-1] = fact_inv[n] * n % MOD
fact = np.array(fact, dtype = np.int64)
fact_inv = np.array(fact_inv, dtype = np.int64)

comb = fact_inv * fact_inv[::-1] % MOD * fact[N] % MOD

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

Submission status: Accepted

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]
A %= MOD
answer = pow(2,pow(2,N,MOD-1),MOD) - A.sum()
answer %= MOD
print(answer)

Submission status: Accepted

Solution code:
import numpy as np


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

        return fn, invs

    def nik(n2, k, m):
        ret = 1
        yield ret
        for _ in range(k):
            ret = ret * n2 % m
            yield ret

    fn, invs = prepare(n, m)
    np_range = np.arange(2, n + 3, dtype=np.int64)
    dp = np.zeros(n + 2, dtype=np.int64)
    dp[0] = 1

    nn2 = [2]
    for i in range(n):
        nn2.append(nn2[-1] ** 2 % m)

    ans = 0
    parity = 1

    for i in range(n + 1):
        coef = np.fromiter(nik(pow(2, n - i, m), i, m), dtype=np.int64)
        ncr = fn * invs[i] * invs[n - i] % m
        left = (dp[:i + 1] * coef % m).sum() % m
        ans = (ans + ncr * left % m * nn2[n - i] % m * parity) % m
        dp[1:i + 2] = (dp[1:i + 2] * np_range[:i + 1] + dp[:i + 1]) % m
        parity *= -1

    return ans


n, m = map(int, input().split())
print(solve(n, m))

Submission status: Accepted

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
    f = pow(f, m - 2, m)
    invs = [1] * (n + 1)
    invs[n] = f
    for i in range(n, 0, -1):
        f = f * i % m
        invs[i - 1] = f

    return fn, invs


def nik(n2, k, m):
    ret = 1
    yield ret
    for _ in range(k):
        ret = ret * n2 % m
        yield ret


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

fn, invs = prepare(n, m)
np_range = np.arange(2, n + 3, dtype=np.int64)
dp = np.zeros(n + 2, dtype=np.int64)
dp[0] = 1

nn2 = [2]
for i in range(n):
    nn2.append(nn2[-1] ** 2 % m)

ans = 0
parity = 1

for i in range(n + 1):
    n2 = pow(2, n - i, m)
    coef = np.fromiter(nik(n2, i, m), dtype=np.int64)
    ncr = fn * invs[i] * invs[n - i] % m
    left = (dp[:i + 1] * coef % m).sum() % m
    ans = (ans + ncr * left % m * nn2[n - i] % m * parity) % m
    dp[1:i + 2] = (dp[1:i + 2] * np_range[:i + 1] + dp[:i + 1]) % m
    parity *= -1

print(ans)

Submission status: Accepted

Solution code:
import numpy as np


def prepare(n, m):
    f = 1
    facts = [1]
    for i in range(1, n + 1):
        f = f * i % m
        facts.append(f)

    f = pow(f, m - 2, m)
    invs = [1] * (n + 1)
    invs[n] = f
    for i in range(n, 0, -1):
        f = f * i % m
        invs[i - 1] = f

    return facts, invs


def nik(n2, k, m):
    ret = 1
    yield ret
    for _ in range(k):
        ret = ret * n2 % m
        yield ret


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

facts, invs = prepare(n, m)
fn = facts[n]

np_range = np.arange(2, n + 4, dtype=np.int64)
dp = np.zeros(n + 2, dtype=np.int64)
dp[0] = 1

nn2 = [2]
for i in range(n):
    nn2.append(nn2[-1] ** 2 % m)

ans = 0
parity = 1

for i in range(n + 1):
    n2 = pow(2, n - i, m)
    coef = np.fromiter(nik(n2, i, m), dtype=np.int64)
    ncr = fn * invs[i] * invs[n - i] % m
    left = (dp[:i + 1] * coef % m).sum() % m
    right = nn2[n - i]
    ans = (ans + ncr * left % m * right % m * parity) % m
    dp[1:i + 2] = (dp[1:i + 2] * np_range[:i + 1] + dp[:i + 1]) % m
    parity *= -1

print(ans)

Submission status: Accepted

Solution code:
from collections import defaultdict

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 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

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

# Preprocess
# calculate 2**n and 2**(2**n)
pow2 = [0] * (N*N // 4 + 1)
pow_pow2 = [0] * (N+1)
pow2[0] = 1; pow_pow2[0] = 2
for i in range(1, N*N // 4 + 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 order i ramens, with each of n toppings chosen at most once)
ways2 = {n: defaultdict(int) 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] * 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 = (temp * pow_pow2[N-n]) % M
    if n % 2 == 0:
        ans = (ans + (com.binom(N, n) * ways) % M) % M
    else:
        ans = (ans - (com.binom(N, n) * ways) % M) % M
print(ans)

Submission status: Accepted

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)

S = [1]

rev2 = pow(2, M-2, M)
base = pow(2, N, M)
ans = 0
S = [1]
for K in range(N+1):
    res = (fact[N] * rfact[K] * rfact[N-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: Accepted

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)

S = [1]

rev2 = pow(2, M-2, M)
base = pow(2, N, M)
ans = 0
S = [1]
for K in range(N+1):
    res = (fact[N] * rfact[K] * rfact[N-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: Accepted

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, 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
    r = (r * v) % M
    if K % 2:
        ans -= r
    else:
        ans += r
    ans %= M

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

Submission status: Runtime Error

Solution code:
import numpy as np
from collections import defaultdict

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 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

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

# Preprocess
# calculate 2**n and 2**(2**n)
pow2 = np.power(2, np.array(range(N*N // 4 + 1)))
pow_pow2 = [0] * (N+1)
pow_pow2[0] = 2
for i in range(1, N*N // 4 + 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 order i ramens, with each of n toppings chosen at most once)
ways2 = {n: defaultdict(int) 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] * 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 = (temp * pow_pow2[N-n]) % M
    if n % 2 == 0:
        ans = (ans + (com.binom(N, n) * ways) % M) % M
    else:
        ans = (ans - (com.binom(N, n) * ways) % M) % M
print(ans)

Submission status: Runtime Error

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 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

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

# Preprocess
# calculate 2**n and 2**(2**n)
pow2 = [0] * (N*N // 4 + 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 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] * 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 = (temp * pow_pow2[N-n]) % M
    if n % 2 == 0:
        ans = (ans + (com.binom(N, n) * ways) % M) % M
    else:
        ans = (ans - (com.binom(N, n) * ways) % M) % M
print(ans)

Submission status: Runtime Error

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 + (com.binom(N, n) * ways) % M) % M
    else:
        ans = (ans - (com.binom(N, n) * ways) % M) % M
print(ans)