Problem:


Score : 600 points

Problem Statement

There are 2N people numbered 1 through 2N.
The height of Person i is h_i.

How many ways are there to make N pairs of people such that the following conditions are satisfied?
Compute the answer modulo 998,244,353.

Each person is contained in exactly one pair.

For each pair, the heights of the two people in the pair are different.

Two ways are considered different if for some p and q, Person p and Person q are paired in one way and not in the other.

Constraints

1 \leq N \leq 50,000

1 \leq h_i \leq 100,000

All values in input are integers.

Input

Input is given from Standard Input in the following format:

N
h_1
:
h_{2N}

Output

Print the answer.

Sample Input 1

2
1
1
2
3

Sample Output 1

2

There are two ways:

Form the pair (Person 1, Person 3) and the pair (Person 2, Person 4).

Form the pair (Person 1, Person 4) and the pair (Person 2, Person 3).

Sample Input 2

5
30
10
20
40
20
10
10
30
50
60

Sample Output 2

516



Sample input:
2
1
1
2
3


Sample output:
2

Submission status: Runtime Error

Solution code:
import os
import sys
from heapq import heappop, heappush

import numpy as np


# https://github.com/atcoder/ac-library/blob/master/atcoder/convolution.hpp

def solve(inp):
    def ceil_pow2(n):
        x = 0
        while 1 << x < n:
            x += 1
        return x

    def bit_scan_forward(n):
        x = 0
        while n & 1 == 0:
            n >>= 1
            x += 1
        return x

    def pow_mod(x, n, m):
        r = 1
        y = x % m
        while n:
            if n & 1:
                r = (r * y) % m
            y = y * y % m
            n >>= 1
        return r

    def get_primitive_root(m):
        if m == 2: return 1
        if m == 167772161: return 3
        if m == 469762049: return 3
        if m == 754974721: return 11
        if m == 998244353: return 3
        divs = [2]
        x = (m - 1) // 2
        while x & 1 == 0:
            x >>= 1
        i = 3
        while i * i <= x:
            if x % i == 0:
                divs.append(i)
                x //= i
                while x % i == 0:
                    x //= i
            i += 2
        if x > 1:
            divs.append(x)
        g = 2
        while True:
            ok = True
            for d in divs:
                if pow_mod(g, (m - 1) // d, m) == 1:
                    ok = False
                    break
            if ok:
                return g

    def butterfly_prepare(mod, primitive_root):
        sum_e = np.zeros(30, np.int64)
        sum_ie = np.zeros(30, np.int64)
        es = np.zeros(30, np.int64)
        ies = np.zeros(30, np.int64)
        cnt2 = bit_scan_forward(mod - 1)
        e = pow_mod(primitive_root, (mod - 1) >> cnt2, mod)
        ie = pow_mod(e, mod - 2, mod)
        for i in range(cnt2, 1, -1):
            es[i - 2] = e
            ies[i - 2] = ie
            e = e * e % mod
            ie = ie * ie % mod
        now_e = 1
        now_ie = 1
        for i in range(cnt2 - 1):
            sum_e[i] = es[i] * now_e % mod
            sum_ie[i] = ies[i] * now_ie % mod
            now_e = now_e * ies[i] % mod
            now_ie = now_ie * es[i] % mod
        return sum_e, sum_ie

    def butterfly(aaa, mod, sum_e):
        n = aaa.size
        h = ceil_pow2(n)

        for ph in range(1, h + 1):
            w = 1 << (ph - 1)
            p = 1 << (h - ph)
            now = 1
            for s in range(w):
                offset = s << (h - ph + 1)
                for i in range(p):
                    l = aaa[i + offset]
                    r = aaa[i + offset + p] * now % mod
                    aaa[i + offset] = (l + r) % mod
                    aaa[i + offset + p] = (l - r) % mod
                now = now * sum_e[bit_scan_forward(~s)] % mod

    def butterfly_inv(aaa, mod, sum_ie):
        n = aaa.size
        h = ceil_pow2(n)

        for ph in range(h, 0, -1):
            w = 1 << (ph - 1)
            p = 1 << (h - ph)
            inow = 1
            for s in range(w):
                offset = s << (h - ph + 1)
                for i in range(p):
                    l = aaa[i + offset]
                    r = aaa[i + offset + p]
                    aaa[i + offset] = (l + r) % mod
                    aaa[i + offset + p] = ((l - r) * inow) % mod
                inow = inow * sum_ie[bit_scan_forward(~s)] % mod

    MOD = 998244353
    primitive_root = get_primitive_root(MOD)
    sum_e, sum_ie = butterfly_prepare(MOD, primitive_root)

    def convolution(aaa, bbb, MOD):
        n = aaa.size
        m = bbb.size
        z = 1 << ceil_pow2(n + m - 1)
        raaa = np.zeros(z, np.int64)
        rbbb = np.zeros(z, np.int64)
        raaa[:n] = aaa
        rbbb[:m] = bbb
        butterfly(raaa, MOD, sum_e)
        butterfly(rbbb, MOD, sum_e)
        ccc = raaa * rbbb % MOD
        butterfly_inv(ccc, MOD, sum_ie)
        iz = pow_mod(z, MOD - 2, MOD)
        result = ccc[:n + m - 1] * iz % MOD
        return result

    def mod_pow(x, a, MOD):
        ret = 1
        cur = x
        while a:
            if a & 1:
                ret = ret * cur % MOD
            cur = cur * cur % MOD
            a >>= 1
        return ret

    def precompute_factorials(n, MOD):
        factorials = np.ones(n + 1, dtype=np.int64)
        for m in range(2, n + 1):
            factorials[m] = factorials[m - 1] * m % MOD
        inversions = np.ones(n + 1, dtype=np.int64)
        inversions[n] = mod_pow(factorials[n], MOD - 2, MOD)
        for m in range(n, 2, -1):
            inversions[m - 1] = inversions[m] * m % MOD
        return factorials, inversions

    def counter(arr):
        cnt = {0: 0}
        for a in arr:
            if a in cnt:
                cnt[a] += 1
            else:
                cnt[a] = 1
        if cnt[0] == 0:
            del cnt[0]
        return cnt

    n = inp[0]
    hhh = inp[1:]
    hhh_cnt = counter(hhh)
    same_heights = []
    max_pair_cnt = 0
    max_same_height = 0
    for c in hhh_cnt.values():
        if c <= 1:
            continue
        max_pair_cnt += c // 2
        max_same_height = max(max_same_height, c)
        same_heights.append(c)

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

    facts, finvs = precompute_factorials(max_same_height, MOD)

    same_heights.sort()
    pc = -1
    cmbs = []
    q = []
    for i, c in enumerate(same_heights):
        if pc == c:
            cmb = cmbs[-1].copy()
        else:
            cmb = np.zeros(c // 2 + 1, np.int64)
            for j in range(c // 2 + 1):
                cmb[j] = facts[c] * finvs[j * 2] % MOD * finvs[c - j * 2] % MOD * pairs[j] % MOD
            pc = c
        cmbs.append(cmb)
        q.append((cmb.size, i))

    while len(q) > 1:
        _, i = heappop(q)
        _, j = heappop(q)
        cmbs[i] = convolution(cmbs[i], cmbs[j], MOD)
        heappush(q, (cmbs[i].size, i))

    result = cmbs[q[0][1]]

    ans = 0
    for i, v in enumerate(result):
        if i & 1:
            ans = (ans - v * pairs[n - i]) % MOD
        else:
            ans = (ans + v * pairs[n - i]) % MOD

    return ans


if sys.argv[-1] == 'ONLINE_JUDGE':
    from numba.pycc import CC

    cc = CC('my_module')
    cc.export('solve', '(i8[:],)')(solve)
    cc.compile()
    exit()

if os.name == 'posix':
    # noinspection PyUnresolvedReferences
    from my_module import solve
else:
    from numba import njit

    solve = njit('(i8[:],)', cache=True)(solve)
    print('compiled', file=sys.stderr)

inp = np.fromstring(sys.stdin.read(), dtype=np.int64, sep=' ')
ans = solve(inp)
print(ans)

Submission status: Runtime Error

Solution code:

# 998244353 = 119*2**23+1

mod = 998244353
primitive_root = 3  # mod の原始根
roots = [pow(primitive_root,(mod-1) >> i,mod) for i in range(24)]
inv_roots = [pow(r,mod-2,mod) for r in roots]
# roots[i] = 1 の 2**i 乗根、inv_roots[i] = 1 の 2**i 乗根の逆元

# 順番は変わる


def ntt(A,n):
    for i in range(n):
        m = 1 << (n-i-1)
        for start in range(1 << i):
            w = 1
            start *= m*2
            for j in range(m):
                A[start+j],A[start+j+m] = (A[start+j]+A[start+j+m]) % mod,(A[start+j]-A[start+j+m])*w % mod
                w *= roots[n-i]
                w %= mod
    return A


def inv_ntt(A,n):
    for i in range(n):
        m = 1 << i
        for start in range(1 << (n-i-1)):
            w = 1
            start *= m*2
            for j in range(m):
                A[start+j],A[start+j+m] = (A[start+j]+A[start+j+m]*w) % mod,(A[start+j]-A[start+j+m]*w) % mod
                w *= inv_roots[i+1]
                w %= mod
    a = pow(2,n*(mod-2),mod)
    for i in range(1 << n):
        A[i] *= a
        A[i] %= mod
    return A


def convolution(A,B):
    a,b = len(A),len(B)
    deg = a+b-2
    n = deg.bit_length()
    N = 1 << n
    A += [0]*(N-a)  # A の次数を 2冪-1 にする
    B += [0]*(N-b)  # B の次数を 2冪-1 にする
    A = ntt(A,n)
    B = ntt(B,n)
    C = [(A[i]*B[i]) % mod for i in range(N)]
    C = inv_ntt(C,n)
    return C[:deg+1]


import sys
sys.setrecursionlimit(10**7)
def I(): return int(sys.stdin.readline().rstrip())


N = I()
A = [0]*100001  # A[i] = 身長 i の人が何人いるか
for _ in range(2*N):
    A[I()] += 1

B = [[1,-A[i]*(A[i]-1)//2] for i in range(1,100001) if A[i] >= 2]
b = len(B)

r = 1
while r < b:
    for i in range(0,b,r*2):
        if i+r < b:
            B[i] = convolution(B[i],B[i+r])
    r *= 2

C = B[0]  # B 内の多項式の積

D = [1]  # D[i] = 1*3*…*(2*i-1)
for i in range(N):
    D.append((D[-1]*(2*i+1)) % mod)

D.reverse()

ans = 0
for i in range(b+1):
    ans += C[i]*D[i]
    ans %= mod

print(ans)

Submission status: Runtime Error

Solution code:

# 998244353 = 119*2**23+1

mod = 998244353
primitive_root = 3  # mod の原始根
roots = [pow(primitive_root,(mod-1) >> i,mod) for i in range(24)]
inv_roots = [pow(r,mod-2,mod) for r in roots]
# roots[i] = 1 の 2**i 乗根、inv_roots[i] = 1 の 2**i 乗根の逆元

# 順番は変わる


def ntt(A,n):
    for i in range(n):
        m = 1 << (n-i-1)
        for start in range(1 << i):
            w = 1
            start *= m*2
            for j in range(m):
                A[start+j],A[start+j+m] = (A[start+j]+A[start+j+m]) % mod,(A[start+j]-A[start+j+m])*w % mod
                w *= roots[n-i]
                w %= mod
    return A


def inv_ntt(A,n):
    for i in range(n):
        m = 1 << i
        for start in range(1 << (n-i-1)):
            w = 1
            start *= m*2
            for j in range(m):
                A[start+j],A[start+j+m] = (A[start+j]+A[start+j+m]*w) % mod,(A[start+j]-A[start+j+m]*w) % mod
                w *= inv_roots[i+1]
                w %= mod
    a = pow(2,n*(mod-2),mod)
    for i in range(1 << n):
        A[i] *= a
        A[i] %= mod
    return A


def convolution(A,B):
    a,b = len(A),len(B)
    deg = a+b-2
    n = deg.bit_length()
    N = 1 << n
    A += [0]*(N-a)  # A の次数を 2冪-1 にする
    B += [0]*(N-b)  # B の次数を 2冪-1 にする
    A = ntt(A,n)
    B = ntt(B,n)
    C = [(A[i]*B[i]) % mod for i in range(N)]
    C = inv_ntt(C,n)
    return C[:deg+1]


import sys
sys.setrecursionlimit(10**7)
def I(): return int(sys.stdin.readline().rstrip())


N = I()
A = [0]*100001  # A[i] = 身長 i の人が何人いるか
for _ in range(2*N):
    A[I()] += 1

B = [[1,-A[i]*(A[i]-1)//2] for i in range(1,100001) if A[i] >= 2]
b = len(B)

r = 1
while r < b:
    for i in range(0,b,r*2):
        if i+r < b:
            B[i] = convolution(B[i],B[i+r])
    r *= 2

C = B[0]  # B 内の多項式の積

D = [1]  # D[i] = 1*3*…*(2*i-1)
for i in range(N):
    D.append((D[-1]*(2*i+1)) % mod)

D.reverse()


ans = 0
for c,d in zip(C,D):
    ans += c*d
    ans %= mod

print(ans)

Submission status: Runtime Error

Solution code:
import numpy as np
from collections import Counter
from heapq import *
import sys
input = sys.stdin.buffer.readline
read = sys.stdin.buffer.read


def convolve(a, b):
    l = len(a) + len(b)
    return np.round(np.real(np.fft.irfft(np.fft.rfft(a, l) * np.fft.rfft(b, l), l))).astype(np.int64)


class Arr:

    def __init__(self, arr):
        self.l = len(arr)
        self.arr = arr

    def __lt__(self, other):
        return self.l < other.l


N = int(input())
hs = list(map(int, read().split()))

ch = Counter(hs)
arrays = []
for k, c in ch.items():
    arr = [1]
    while c >= 2:
        arr.append(c * (c - 1) // 2)
        c -= 2
    arrays.append(Arr(np.array(arr)))

heapify(arrays)
while len(arrays) >= 2:
    a1 = heappop(arrays)
    a2 = heappop(arrays)
    heappush(arrays, Arr(convolve(a1.arr, a2.arr) % MOD))
    
final = arrays[0].arr

MOD = 998244353
table_len = 100010

fac = [1, 1]
for i in range(2, table_len):
    fac.append(fac[-1] * i % MOD)

finv = [0] * table_len
finv[-1] = pow(fac[-1], MOD - 2, MOD)
for i in range(table_len-1, 0, -1):
    finv[i - 1] = finv[i] * i % MOD
    

ans = 0

for i in range(N + 1):
    tmp = pow(2, MOD - 1 - N + i, MOD) * fac[2 * (N - i)] % MOD * finv[N - i] % MOD
    # print("t", tmp)
    if len(final) > i:
        tmp *= int(final[i])
        tmp %= MOD
    else:
        break
    if i % 2 == 1:
        ans -= tmp
    else:
        ans += tmp
    ans %= MOD
    # print("a", ans)

print(ans % MOD)

Submission status: Runtime Error

Solution code:
ROOT = 3
MOD = 998244353
roots  = [pow(ROOT,(MOD-1)>>i,MOD) for i in range(24)] # 1 の 2^i 乗根
iroots = [pow(x,MOD-2,MOD) for x in roots] # 1 の 2^i 乗根の逆元

def untt(a,n):
    for i in range(n):
        m = 1<<(n-i-1)
        for s in range(1<<i):
            w_N = 1
            s *= m*2
            for p in range(m):
                a[s+p], a[s+p+m] = (a[s+p]+a[s+p+m])%MOD, (a[s+p]-a[s+p+m])*w_N%MOD
                w_N = w_N*roots[n-i]%MOD

def iuntt(a,n):
    for i in range(n):
        m = 1<<i
        for s in range(1<<(n-i-1)):
            w_N = 1
            s *= m*2
            for p in range(m):
                a[s+p], a[s+p+m] = (a[s+p]+a[s+p+m]*w_N)%MOD, (a[s+p]-a[s+p+m]*w_N)%MOD
                w_N = w_N*iroots[i+1]%MOD
            
    inv = pow((MOD+1)//2,n,MOD)
    for i in range(1<<n):
        a[i] = a[i]*inv%MOD

def convolution(a,b):
    la = len(a)
    lb = len(b)
    deg = la+lb-2
    n = deg.bit_length()
    if min(la, lb) <= 50:
        if la < lb:
            la,lb = lb,la
            a,b = b,a
        res = [0]*(la+lb-1)
        for i in range(la):
            for j in range(lb):
                res[i+j] += a[i]*b[j]
                res[i+j] %= MOD
        return res

    N = 1<<n
    a += [0]*(N-len(a))
    b += [0]*(N-len(b))
    untt(a,n)
    untt(b,n)
    for i in range(N):
      a[i] = a[i]*b[i]%MOD
    iuntt(a,n)
    return a[:deg+1]

def convolution_all(polys):
    if not polys: return [1]
    polys.sort(key=lambda x:len(x))
    N = len(polys)
    height = (N-1).bit_length()
    N0 = 1<<(height)
    data = [None]*(N0*2)
    data[N0:N0+N] = polys
    data[N0+N:] = [[1] for _ in range(N0-N)]
    
    #print(data)
    #print(polys)
    for k in range(N0-1,0,-1):
        data[k] = convolution(data[2*k], data[2*k+1])        

    return data[1]

SIZE=10**2 + 10 #998244353 #ここを変更する

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
fac[0] = fac[1] = 1
finv[0] = finv[1] = 1
for i in range(2,SIZE):
    fac[i] = fac[i-1]*i%MOD
finv[-1] = pow(fac[-1],MOD-2,MOD)
for i in range(SIZE-1,0,-1):
    finv[i-1] = finv[i]*i%MOD
    inv[i] = finv[i]*fac[i-1]%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

dfac = [0]*SIZE  # fac[j] = j! mod MOD
dfac[0] = dfac[1] = 1
for i in range(2,SIZE):
    dfac[i] = dfac[i-2]*i%MOD

import sys
readline = sys.stdin.readline
read = sys.stdin.read

n,*h = map(int, read().split())
from collections import Counter

def get(k): # k 個の頂点からiペア作れる場合の数
    res = [1]
    for i in range(k//2):
        res.append(res[-1]*choose(k-2*i,2)%MOD)
    for i in range(2,k//2 + 1):
        res[i] = res[i]*finv[i]%MOD
    return res        


d = Counter(h)
res = []
for k,v in d.items():
    if v >= 2: res.append(get(v))

bad = convolution_all(res)

#print(bad)
#bad = [1,5,7,3]
ans = 0
sgn = 1
for i,ai in enumerate(bad):
    ans += sgn*dfac[2*n-2*i-1]%MOD*ai  
    ans %= MOD
    sgn *= -1
    

print(ans%MOD)

Submission status: Runtime Error

Solution code:
mod=998244353
n,q=map(int,input().split())
lrd=[list(map(int,input().split())) for _ in range(q)]
s=0
ss=[]
tmp=1
for i in range(n):
  ss.append(tmp)
  s+=tmp
  s%=mod
  tmp*=10
  tmp%=mod
ss.reverse()
#print(ss,s)
for l,r,d in lrd:
  for i in range(l-1,r):
    s-=ss[i]
    ss[i]=d*pow(10,n-i-1,mod)
    s+=ss[i]
    s%=mod
  print(s)