from __future__ import annotations input class CombLUT: MOD = 998_244_353 _max = 0 F = [1] InvF = [1] _fibo = [0, 1] _lucas = [2, 1] @classmethod def expand(cls, n): if n <= cls._max: return n = n*4//3 cls.F.extend([1]*(n - cls._max)) for i in range(cls._max + 1, n+1): cls.F[i] = cls.F[i-1] * i % cls.MOD cls.InvF.extend([1]*(n - cls._max)) cls.InvF[n] = pow(cls.F[n], -1, cls.MOD) for i in range(n-1, cls._max, -1): cls.InvF[i] = cls.InvF[i+1] * (i+1) % cls.MOD cls._max = n @classmethod def comb(cls, n, r): if not (0 <= r <= n): return 0 if n > cls._max: cls.expand(n) return cls.F[n] * cls.InvF[r] % cls.MOD * cls.InvF[n-r] % cls.MOD @classmethod def hyper_comb(cls, n, r): if (r < 0): return 0 if (n < 0): return 0 if n == r == 0: return 1 if not (0 <= r <= n-1+r): return 0 return cls.comb(n-1+r, r) @classmethod def neg_hyper_comb(cls, n, r): if (r < 0): return 0 if n >= 0: return cls.hyper_comb(n, r) sgn = -1 if r&1 else 1 return sgn * cls.comb(-n, r) @classmethod def factorial(cls, n): assert 0 <= n if n > cls._max: cls.expand(n) return cls.F[n] @classmethod def inv_factorial(cls, n): assert 0 <= n if n > cls._max: cls.expand(n) return cls.InvF[n] def subset_sum(A, *, func=None, initial=None): if (func is None): if initial is None: initial = 0 N = len(A) dp = [initial] * (1< list[int]: if maxlen is None: if len(A) == 0 or len(B) == 0: return [] maxlen = len(A) + len(B) - 1 assert isinstance(maxlen, int) if maxlen <= 0: return [] if MOD is None: res = [0]*(maxlen) for i, a in enumerate(A): for j, b in enumerate(B): if not (i+j < maxlen): break res[i+j] += a*b return res else: res = [0]*(maxlen) for i, a in enumerate(A): for j, b in enumerate(B): if not (i+j < maxlen): break res[i+j] = (res[i+j] + a*b) % MOD return res def dummy_frac(MOD=998_244_353, U=10**6): dp = [0, 1] dp += [0]*(U+1 - len(dp)) for x in range(2, U+1): q, r = divmod(MOD, x) if dp[r] == 0: dp[x] = pow(x, -1, MOD) else: dp[x] = -(q) * dp[r] % MOD memo = {} def F(a=0, b=1): if not 1 <= b <= U: b %= MOD if b == 0: raise ZeroDivisionError if b == 1: return a if 1 <= b <= U: return a * dp[b] % MOD inv = memo.get(b) if inv is None: memo[b] = inv = pow(b, -1, MOD) return a*inv%MOD return F Fraction = dummy_frac() class ExpPoly: def __init__(self, lamda, coef: list[Fraction]): self.lamda = lamda self.coef = coef def __add__(self, other): raise NotImplementedError def __mul__(self, other: "int|Fraction|ExpPoly"): if isinstance(other, int): return ExpPoly(self.lamda, [c*other for c in self.coef]) else: assert isinstance(other, ExpPoly) nshift = self.lamda + other.lamda ncoef = convolution_naive(self.coef, other.coef, MOD=998_244_353) # type: ignore return ExpPoly(nshift, ncoef) def __rmul__(self, other: "int|Fraction"): return self.__mul__(other) def integral(self): res = Fraction(0) for n in range(len(self.coef)): c = self.coef[n] res += c * CombLUT.factorial(n) % MOD * pow(self.lamda, -(n+1), MOD) % MOD res %= MOD return res def bit_count(x): x = x - ((x>>1)&0x5555_5555_5555_5555) x = (x&0x3333_3333_3333_3333) + ((x>>2)&0x3333_3333_3333_3333) x = (x + (x>>4))&0x0f0f_0f0f_0f0f_0f0f x += (x>>8) x += (x>>16) x += (x>>32) return x&0x0000_0000_0000_007f def solve(N, K, A): def resolve(k, a): lamda = Fraction(1, a) coef = [Fraction(0)] * k for n in range(k): coef[n] = pow(lamda, n, MOD) * CombLUT.inv_factorial(n) return ExpPoly(lamda, coef) F = [resolve(k, a) for k, a in zip(K, A)] dp = subset_sum(F, func=lambda x, y: x*y, initial=1) I = [f.integral() if isinstance(f, ExpPoly) else 0 for f in dp] popcnt = [Fraction(0)] * (N+1) for S in range(1<