Xmas Contest 2021 选做

Xmas Contest 2021

Xmas Contest 其实已经在 AtCoder 举办了几年了，不过不能直接通过首页进入，而且题目没有英文。$$

C - Count Me

题目大意

给定一个 0、1、? 构成的字符串。将每个 ? 枚举其填 0 或 1 后得到一个 01 串 $S$，对所有 $S$ 的 $f(S)$ 求和，$f(S)$ 为如下序列的数量：

• $S_0,S_1,\dots ,S_N=S$ 皆为 01 序列，其中第 $k$ 项的长度为 $k$。
• $S_k$ 是 $S_{k+1}$ 的子序列。

对 $998244353$ 取模。

• $N\le 2.5\times {10}^5$

经过一些特殊情况的打表我们猜想，答案等价于如下 $N+1$ 阶排列计数：将 0 替换为 <、将 1 替换为 > 之后满足相邻不等关系的排列。

证明是很简单的，首先我们发现操作序列等价于每次把某个全 0 段或者全 1 段的长度减去 $1$，考虑每次去掉最大的那个数，这能和删除一个 0 或 1 的操作序列进行一一映射。

那么我们直接把 LOJ575 不等关系 的代码蒯过来就好了。

复杂度 $O(N{\log }^{2}N)$ 或 $O\left(\frac{N{\log }^{2}N}{\log \log N}\right)$ 或更低。

点击查看代码

#include <bits/stdc++.h>

using namespace std;

typedef long long ll;
typedef unsigned long long ull;

const int L = 20, N = 1 << L, P = 998244353, R = 3;

int n;
int inv[N], fac[N], ifac[N], f[N];
bool key[N];
char s[N], t[N];

int norm(int x) { return x >= P ? (x - P) : x; }

int reduce(int x) { return x < 0 ? x + P : x; }

void add(int& x, int y) {
	if ((x += y) >= P)
		x -= P;
}

int mpow(int x, int k) {
	int ret = 1;
	while (k) {
		if (k & 1)
			ret = ret * (ll)x % P;
		k >>= 1;
		x = x * (ll)x % P;
	}
	return ret;
}

void prepare(int n) {
	inv[1] = 1;
	for (int i = 2; i <= n; ++i) inv[i] = -(P / i) * (ll)inv[P % i] % P + P;
	fac[0] = 1;
	for (int i = 1; i <= n; ++i) fac[i] = fac[i - 1] * (ll)i % P;
	ifac[0] = 1;
	for (int i = 1; i <= n; ++i) ifac[i] = ifac[i - 1] * (ll)inv[i] % P;
}

struct NTT {
	int L;
	vector<int> root;

	NTT() : L(-1) {}

	void prepRoot(int l) {
		L = l;
		root.resize((1 << L) + 1);
		int i, n = 1 << L;
		int *w2 = root.data();
		*w2 = 1;
		w2[1 << l] = mpow(31, 1 << (21 - l));

		for (i = l; i; --i)
			w2[1 << (i - 1)] = (ull) w2[1 << i] * w2[1 << i] % P;

		for (i = 1; i < n; ++i)
			w2[i] = (ull) w2[i & (i - 1)] * w2[i & -i] % P;
	}

	void DIF(int *a, int l) {
		int *j, *k, n = 1 << l, len = n >> 1, r, *o;

		for (; len; len >>= 1)
			for (j = a, o = root.data(); j != a + n; j += len << 1, ++o)
				for (k = j; k != j + len; ++k) {
					r = (ull) *o * k[len] % P;
					k[len] = reduce(*k - r);
					add(*k, r);
				}
	}

	void DIT(int *a, int l) {
		int *j, *k, n = 1 << l, len = 1, r, *o;

		for (; len != n; len <<= 1)
			for (j = a, o = root.data(); j != a + n; j += len << 1, ++o)
				for (k = j; k != j + len; ++k) {
					r = norm(*k + k[len]);
					k[len] = ull(*k - k[len] + P) * *o % P;
					*k = r;
				}
	}

	void fft(int *a, int lgn, int d = 1) {
		if (L < lgn) prepRoot(lgn);
		int n = 1 << lgn;
		if (d == 1) DIF(a, lgn);
		else {
			DIT(a, lgn);
			reverse(a + 1, a + n);
			ull nv = P - (P - 1) / n;
			for (int i = 0; i < n; ++i) a[i] = a[i] * nv % P;
		}
	}
} ntt;

const int BRUTE_N2_LIMIT = 50;

int pool[N << 2];

void dc(int l, int r, int* top) {
	if (l == r) {
		if (l == 0) return;
		if (key[l])
			f[l] = norm(P - f[l]);
		else
			f[l] = 0;
		return;
	}
	if (r - l <= BRUTE_N2_LIMIT) {
		for (int i = l; i <= r; ++i) {
			for (int j = l; j < i; ++j)
				f[i] = (f[i] + f[j] * (ll)ifac[i - j]) % P;
			dc(i, i, top);
		}
		return;
	}
	int lg = 31 - __builtin_clz(r - l);
	int d = (r - l) / lg + 1;
	int lgd = 0;
	while ((1 << lgd) < d) ++lgd;
	++lgd;

	for (int i = 0; i < lg; ++i) {
		copy(ifac + i * d, ifac + (i + 2) * d, top + (i << lgd));
		fill(top + (i << lgd) + 2 * d, top + ((i + 1) << lgd), 0);
		ntt.fft(top + (i << lgd), lgd, 1);
	}
	fill(top + (lg << lgd), top + (lg << (lgd + 1)), 0);
	for (int i = 0; i < lg; ++i) {
		if (i) ntt.fft(top + ((lg + i) << lgd), lgd, -1);
		for (int j = 0; j < min(d, r - l - i * d + 1); ++j)
			f[l + i * d + j] = norm(f[l + i * d + j] + top[((lg + i) << lgd) + d + j]);
		dc(l + i * d, min(l + (i + 1) * d - 1, r), top + (lg << (lgd + 1)));
		if (i + 1 < lg) {
			copy(f + l + i * d, f + l + (i + 1) * d, top + ((lg + i) << lgd));
			fill(top + ((lg + i) << lgd) + d, top + ((lg + i + 1) << lgd), 0);
			ntt.fft(top + ((lg + i) << lgd), lgd, 1);
		}
		for (int j = i + 1; j < lg; ++j) {
			for (int k = 0; k < (1 << lgd); ++k)
				top[((lg + j) << lgd) + k] = (top[((lg + j) << lgd) + k] + top[((j - i - 1) << lgd) + k] * (ll)top[((lg + i) << lgd) + k]) % P;
		}
	}
}

int solve(int n) {
	for (int i = 1; i <= n; ++i)
		key[i] = s[i] == '>';
	key[0] = true;
	key[n + 1] = true;
	fill(f, f + n + 2, 0);
	f[0] = 1;
	dc(0, n + 1, pool);
	int ans = f[n + 1];
	bool par = false;
	for (int i = 1; i <= n + 1; ++i)
		if (key[i])
			par ^= true;
	if (par)
		ans = norm(P - ans);
	return ans;
}

int main() {
	prepare((1 << L) - 1);
	int n; scanf("%d%s", &n, t + 1);
	int lst = 0;
	int ans = fac[n + 1];
	t[n + 1] = '?';
	for (int i = 1; i <= n + 1; ++i) {
		if (t[i] == '?') {
			if (i - lst > 1) {
				for (int j = lst + 1; j < i; ++j) s[j - lst] = t[j] == '0' ? '>' : '<';
				ans = ans * (ll)solve(i - lst - 1) % P;
			}
			lst = i;
		}
	}
	printf("%d\n", ans);
	return 0;
}

D - Determinant?

题目大意

给 $N$ 个 $K\times K$ 的矩阵 $A_i$，计算 $\sum _{\sigma \in {\mathfrak{S}}_N}\mathrm{sgn}(\sigma )A_{\sigma (1)}\cdots A_{\sigma (N)}$，答案对 $998244353$ 取模。

• $1\le N\le 32$
• $1\le K\le 8$

我们先从简单的问题想起。假设 $K=1$，$A$ 退化为标量，这时，当 $N\ge 2$ 时，奇偶排列配对得知答案必为 $0$。

考虑计算最终矩阵的 $M_{i_0{i}_N}$ 项，将其贡献拆解为 $A_{\sigma (1)i_0{i}_1}\cdots A_{\sigma (N)i_{N-1}{i}_N}$。此时，对于一条具体的路径 $i_0\to i_1\to \cdots \to i_N$，设 $B_{jk}=A_{j{i}_{k-1}{i}_k}$，那么这部分贡献的求和就是

$$\sum _{\sigma \in {\mathfrak{S}}_N}\mathrm{sgn}(\sigma )B_{\sigma (1)1}\cdots B_{\sigma (N)N}=detB$$

于答案 $M_{i_0{i}_N}$，我们应该可以将其写作

$$\sum _{S\subseteq E}det{A}_{[S]}\cdot (\sum _{\sigma \in {\mathfrak{S}}_N}\mathrm{sgn}(\sigma )[S_{\sigma (1)}\cdots S_{\sigma (N)}\text{ is }{i}_0⇝i_N])$$

那么问题来了，对于一个给定的边集，后面这个值有什么样的性质呢？我们注意到，如果一个点有两次从它出发又回到自身，那么这两个回路的边顺序交换之后，$\sigma$>那么问题来了，对于一个给定的边集，后面这个值有什么样的性质呢？我们注意到，如果一个点有两次从它出发又回到自身，那么这两个回路的边顺序交换之后，$\sigma$ 的奇偶性改变。因此，起点的入度 $\ge 2$ 的话就会贡献抵消，其它点的入度 $\ge 3$ 的话就会贡献抵消，因此边数不能超过 $2K-1$。

上面这个证明好像漏了环长是偶数的情况, 我有时间再修修, 但是好像有些代数证明比如说 这个.

所以，当 $N\ge 2K$ 的时候答案一定是 $O$，否则我们直接状压 DP 就好了，复杂度为 $O(2^{N}N{K}^3)=O(4^K{K}^4)$。

C T T 2 0 1 9 情 景 再 现

点击查看代码

#include <bits/stdc++.h>

using namespace std;

using ull = unsigned long long;
using ll = long long;

const int P = 998244353;

int norm(int x) { return x >= P ? x - P : x; }
int reduce(int x) { return x < 0 ? x + P : x; }
int neg(int x) { return x ? P - x : 0; }
void add(int& x, int y) { if ((x += y - P) < 0) x += P; }
void sub(int& x, int y) { if ((x -= y) < 0) x += P; }
void fam(int& x, int y, int z) { x = (x + y * (ull)z) % P; }
int mpow(int x, unsigned k) {
	if (k == 0) return 1;
	int ret = mpow(x * (ull)x % P, k >> 1);
	if (k & 1) ret = ret * (ull)x % P;
	return ret;
}
int quo2(int x) { return ((x & 1) ? x + P : x) >> 1; }

const int _ = 8;

int dp[1 << 15][_][_];
int A[15][_][_];

int main() {
	ios::sync_with_stdio(false); cin.tie(NULL);

	int N, K; cin >> N >> K;
	if (N >= K * 2) {
		for (int i = 0; i != K; ++i)
			for (int j = 0; j != K; ++j)
				cout << 0 << " \n"[j == K - 1];
		return 0;
	}
	for (int i = 0; i != N; ++i)
		for (int j = 0; j != K; ++j)
			for (int k = 0; k != K; ++k)
				cin >> A[i][k][j];
	for (int i = 0; i != K; ++i) dp[0][i][i] = 1;
	for (int s = 1; s != 1 << N; ++s)
		for (int i = 0; i != N; ++i) if ((s >> i) & 1) {
			if (__builtin_parity(s >> (i + 1))) {
				for (int j = 0; j != K; ++j)
					for (int k = 0; k != K; ++k)
						for (int l = 0; l != K; ++l)
							fam(dp[s][j][k], dp[s ^ (1 << i)][j][l], P - A[i][k][l]);
			} else {
				for (int j = 0; j != K; ++j)
					for (int k = 0; k != K; ++k)
						for (int l = 0; l != K; ++l)
							fam(dp[s][j][k], dp[s ^ (1 << i)][j][l], A[i][k][l]);
			}
		}
	for (int i = 0; i != K; ++i)
		for (int j = 0; j != K; ++j)
			cout << dp[(1 << N) - 1][i][j] << " \n"[j == K - 1];

	return 0;
}

E - E and PI

题目大意

对于正整数 $j$，数轴上点 $P_j$ 的坐标为 $j\mathrm{e}mod1$，其中 $\mathrm{e}$ 是自然对数的底。对于正整数 $n$，考虑点 $P_1,\dots ,P_{n-1}$ 将线段 $[0,1]$ 分割成 $n$ 段，记 $f(n)$ 为这 $n$ 段中的最大值。

记 $n_K$ 是第 $K$ 个满足 $f(n)>f(n+1)$ 的正整数，给定 $K$，求出 $n_K$ 和 $f(n_K)$。显然 $f(n_K)$ 可以表示成 $a+b\mathrm{e}$ 的形式，只需要输出 $n_K,a,b$ 同余 $998244353$ 的结果。

• $1\le K\le {10}^{11}$

由于本人不太懂具体道理，所以这里主要讲讲大概找到什么规律。

经过打表和 OEIS 搜索之后我们可以发现，$n_K=\mathtt{\text{A006259}}(K+3)-1$，这个数列由 $\mathrm{e}$ 的最佳有理逼近给出。也就是，将 $\mathrm{e}$ 在 Stern-Brocot 树上走的字符串 $\mathrm{e}={\mathrm{RL}}^0{\mathrm{RLR}}^2{\mathrm{LRL}}^4{\mathrm{RLR}}^6\dots$ 截断到前 $K+2$ 个字符，得到的数 $-1$ 就是答案：

$$\frac{2}{1},\frac{3}{1},\frac{5}{\mathbf{2}},\frac{8}{\mathbf{3}},\frac{11}{\mathbf{4}},\frac{19}{\mathbf{7}},\frac{30}{\mathbf{11}},\frac{49}{\mathbf{18}},\frac{68}{\mathbf{25}},\frac{87}{\mathbf{32}},\frac{106}{\mathbf{39}},\dots$$

显然我们可以在 $O(1)$ 的时间内完成 ${\mathrm{L}}^n$ 或 ${\mathrm{R}}^n$ 的计算，因此我们就可以在 $O(\sqrt{K})$ 的时间内暴力计算出 $n_K$ 了。

然后经过打表，我们又发现对应的 $f(n_K)$ 是这样的：

$$1+0\mathrm{e},-2+\mathrm{e},-5+2\mathrm{e},3-\mathrm{e},-8+3\mathrm{e},11-4\mathrm{e},30-11\mathrm{e},49-18\mathrm{e},68-25\mathrm{e},-19+7\mathrm{e},\dots$$

我们发现除了第一个，对于剩下的 $a+b\mathrm{e}$，${\displaystyle \frac{-a}{b}}$ 的值依然都是出现在有理逼近里的！接着我们更具体地看，容易想到 $a,b$ 何者为负只取决于 ${\displaystyle \frac{a}{b}}$ 和 $\mathrm{e}$ 的大小关系，也就是说，如果 ${\displaystyle \frac{a}{b}}$ 的下一个字符是 $\mathrm{L}$，那么说明 ${\displaystyle \frac{a}{b}}>\mathrm{e}$，出现在答案里的应该是 $a-b\mathrm{e}$，反之则是 $-a+b\mathrm{e}$。再仔细看，我们发现有些数出现的位置相比于原序列是变得靠后了：

$$1+0\mathrm{e},-2+\mathrm{e},-5+2\mathrm{e},\boxed{{\displaystyle 3-\mathrm{e}}},-8+3\mathrm{e},11-4\mathrm{e},30-11\mathrm{e},49-18\mathrm{e},68-25\mathrm{e},\boxed{{\displaystyle -19+7\mathrm{e}}},\dots$$

看起来规律还是不够显然，我们再把 $\mathrm{LR}$ 标出来，并且用红色表示位置不对的数，灰色表示这些数在逼近序列里的位置。

$$\begin{gathered} \frac{1}{0}\stackrel{\mathrm{R}}{\to }\frac{2}{1}\stackrel{\mathrm{L}}{\to }\boxed{{\displaystyle \frac{3}{1}\stackrel{\mathrm{R}}{\to }}}\frac{5}{2}\stackrel{\mathrm{L}}{\to }\boxed{{\displaystyle \frac{3}{1}\stackrel{\mathrm{R}}{\to }}}\frac{8}{3}\stackrel{\mathrm{L}}{\to }\frac{11}{4}\stackrel{\mathrm{R}}{\to }\boxed{{\displaystyle \frac{19}{7}\stackrel{\mathrm{L}}{\to }}}\frac{30}{11}\stackrel{\mathrm{R}}{\to }\frac{49}{18}\stackrel{\mathrm{R}}{\to }\frac{68}{25}\stackrel{\mathrm{R}}{\to }\boxed{{\displaystyle \frac{19}{7}\stackrel{\mathrm{L}}{\to }}}\frac{87}{32}\stackrel{\mathrm{R}}{\to }\frac{106}{39}\stackrel{\mathrm{L}}{\to } \\ \boxed{{\displaystyle \frac{193}{71}\stackrel{\mathrm{R}}{\to }}}\frac{299}{110}\stackrel{\mathrm{L}}{\to }\frac{492}{181}\stackrel{\mathrm{L}}{\to }\frac{685}{252}\stackrel{\mathrm{L}}{\to }\frac{878}{323}\stackrel{\mathrm{L}}{\to }\frac{1071}{394}\stackrel{\mathrm{L}}{\to }\boxed{{\displaystyle \frac{193}{71}\stackrel{\mathrm{R}}{\to }}}\frac{1264}{465}\stackrel{\mathrm{L}}{\to }\dots \end{gathered}$$

如果用每个数接下来的字符代表它，那么我们发现这个情况就是它恰好把前面提到的 ${\mathrm{LR}}^{2n}/{\mathrm{RL}}^{2n}$，替换成了 ${\mathrm{L}\mathrm{R}}^{2n-1}\mathrm{L}\mathrm{R}/{\mathrm{R}\mathrm{L}}^{2n-1}\mathrm{R}\mathrm{L}$。

这也是好处理的，于是我们就可以 $O(\sqrt{K})$ 解决这道题了。

（那么如果我们再用类似整式递推的方法加个速，可以做到 $O(K^{1/4}\log K)$，这题就能出到 ${10}^{18}$ 了）

点击查看代码

#include <bits/stdc++.h>

using namespace std;

using ull = unsigned long long;

const int P = 998244353;

int norm(int x) { return x >= P ? x - P : x; }
int reduce(int x) { return x < 0 ? x + P : x; }
int neg(int x) { return x ? P - x : 0; }
void add(int& x, int y) { if ((x += y - P) < 0) x += P; }
void sub(int& x, int y) { if ((x -= y) < 0) x += P; }
void fam(int& x, int y, int z) { x = (x + y * (ull)z) % P; }
int mpow(int x, unsigned k) {
	if (k == 0) return 1;
	int ret = mpow(x * (ull)x % P, k >> 1);
	if (k & 1) ret = ret * (ull)x % P;
	return ret;
}
int quo2(int x) { return ((x & 1) ? x + P : x) >> 1; }

void cal1(ull K) {
	int a = 1, b = 0;
	auto apply = [&](bool f, int N) {
		N = min((ull)N, K);
		K -= N;
		if (!f) fam(a, b, N);
		else fam(b, a, N);
		return K == 0;
	};
	for (int i = 0; ; ++i) {
		if (apply(i & 1, 1)) break;
		if (apply(~i & 1, i * 2)) break;
		if (apply(i & 1, 1)) break;
	}
	cout << norm(a + b - 1);
}

int main() {
	ios::sync_with_stdio(false);
	cin.tie(NULL); cout.tie(NULL);

	ull K; cin >> K;
	cal1(K + 2);
	if (K == 1) { cout << " 1 0\n"; return 0; }
	if (K == 2) { cout << " 998244351 1\n"; return 0; }
	--K;
	int x = 0, y = 1, a = 1, b = 0;
	auto apply = [&](bool f, int N) {
		N = min((ull)N, K);
		K -= N;
		if (!f) { fam(a, b, N); fam(x, y, N); }
		else { fam(b, a, N); fam(y, x, N); }
		return K == 0;
	};
	for (int i = 0; ; ++i) {
		if (K >= (i + 1) * 2) {
			apply(i & 1, 1); apply(~i & 1, i * 2), apply(i & 1, 1);
			continue;
		}
		if (K < i * 2 - 1) ++K;
		else if (K == i * 2 - 1) K = 0;
		bool flag = (K == 0 || K == i * 2 + 1) ^ (i & 1);
		int a1, a2;
		if (!apply(i & 1, 1))
			if (!apply(~i & 1, i * 2))
				apply(i & 1, 1);
		a1 = norm(x + y); a2 = norm(a + b);
		if (flag) a1 = reduce(-a1);
		else a2 = reduce(-a2);
		cout << ' ' << a1 << ' ' << a2 << '\n';
		break;
	}

	return 0;
}

H - Homework from Zhejiang

题目大意

给定大小分别为 $MN$ 的有向带权图 $G_1,G_2$，权值矩阵分别为 $A,B$。令 $G=G_1◻G_2$，也即点集由 $V_1\times V_2$ 构成，其中 $(i,j)$ 到 $(k,j)$ 有边，权值为 $A_{ik}$，另外 $(i,j)$ 到 $(i,l)$ 有边，权值为 $B_{jl}$。

对于每个 $(i,j)$，求以 $(i,j)\in G_1◻G_2$ 为根的外向树权值和，权值定义为边权之积。答案对 $998244353$ 取模。

• $2\le M,N\le 500$

这就是经典的图的 Cartesian 积。令 $L^{(1)},L^{(2)}$ 分别为 $G_1,G_2$ 的 Laplacian 矩阵，易得 $G_1◻G_2$ 的 Laplacian 矩阵为 $L=L^{(1)}\otimes I_N+I_M\otimes L^{(2)}$。我们知道，计算以一个点 $u$ 为根的生成树数量就是在 $L_{uu}$ 处的主余子式，只需计算 Laplacian 矩阵的伴随矩阵。由于 Laplacian 矩阵不满秩，伴随矩阵能写作 $x{y}^{\mathsf{T}}$ 的形式，其中 $Lx=0$，$y^{\mathsf{T}}L=0$，根据 Laplacian 矩阵的性质，显然 $y$ 只需取全 $1$ 向量，而我们可以在 $O(N^3)$ 时间内消元得到 $x$。得到了 $x^{(1)},x^{(2)}$ 之后，易见 $L(x^{(1)}\otimes x^{(2)})=0$，可知 $x^{(1)}\otimes x^{(2)}$ 就是答案了。但其实还有点问题，求出来的向量其实整体差一个常数倍，我们还需要定出这个常数倍。

设 $L^{(1)},L^{(2)}$ 的特征值分别为 ${\lambda }_1,\dots ,{\lambda }_M$ 和 ${\mu }_1,\dots ,{\mu }_N$，我们有 $L$ 的特征值为 ${\lambda }_i+{\mu }_j$（证明见 [1]），由于一定有 $0$ 这个特征值，不妨设 ${\lambda }_1={\mu }_1=0$，那么 $G$ 以每个点为根的生成树数量之和，也就是特征多项式的 $1$ 次项系数，也即求出除 ${\lambda }_1+{\mu }_1$ 以外的特征值乘积。这可以用求出 $L^{(1)}$ 和 $L^{(2)}$ 的特征多项式之后利用结式计算，答案记为 $T$。这样，我们乘以常数 $c$ 使得 $c(x^{(1)}\otimes x^{(2)}{)}^{\mathsf{T}}\cdot 1=T$ 之后，答案就对了。

由于答案之和有可能刚好是 $998244353$ 的倍数，上述做法会寄，我们需要稍微仔细一点处理。

设 $x^{(1)},x^{(2)}$ 分别为 $G_1,G_2$ 的答案向量，这是容易算的，首先得到他们的常数倍，然后对一个非零位置求一次行列式就能得到。那么注意到他们的和分别是 ${\lambda }_2\cdots {\lambda }_M$ 和 ${\mu }_2\cdots {\mu }_N$，所以此时 $G$ 的答案向量一定是

$$(\prod _{i=2}^M\prod _{j=2}^N({\lambda }_i+{\mu }_j))\cdot x^{(1)}\otimes x^{(2)}$$

而此时的乘积就是 $\mathrm{res}(\frac{|tI-A|}{t},\frac{|(-t)I-B|}{-t})$。

复杂度 $O(M^3+N^3)$。

点击查看代码

#include <bits/stdc++.h>

using namespace std;

using ull = unsigned long long;
using ll = long long;

const int P = 998244353;

int norm(int x) { return x >= P ? x - P : x; }
int reduce(int x) { return x < 0 ? x + P : x; }
int neg(int x) { return x ? P - x : 0; }
void add(int& x, int y) { if ((x += y - P) < 0) x += P; }
void sub(int& x, int y) { if ((x -= y) < 0) x += P; }
void fam(int& x, int y, int z) { x = (x + y * (ull)z) % P; }
int mpow(int x, unsigned k) {
	if (k == 0) return 1;
	int ret = mpow(x * (ull)x % P, k >> 1);
	if (k & 1) ret = ret * (ull)x % P;
	return ret;
}
int quo2(int x) { return ((x & 1) ? x + P : x) >> 1; }

namespace PinkRabbit {

	const int MN = 510, Mod = P;
	typedef long long LL;

	int N, A[MN][MN], H[MN][MN], Ans[MN];

	void Hessenberg() {
		for (int i = 1; i <= N; ++i)
			for (int j = 1; j <= N; ++j) {
				H[i][j] = A[i][j];
			}
		for (int i = 1; i <= N - 2; ++i) {
			if (!H[i + 1][i]) for (int j = i + 2; j <= N; ++j) if (H[j][i]) {
				for (int k = i; k <= N; ++k) std::swap(H[i + 1][k], H[j][k]);
				for (int k = 1; k <= N; ++k) std::swap(H[k][i + 1], H[k][j]);
				break;
			}
			if (!H[i + 1][i]) continue;
			int val = mpow(H[i + 1][i], P - 2);
			for (int j = i + 2; j <= N; ++j) {
				int coef = (LL)val * H[j][i] % Mod;
				for (int k = i; k <= N; ++k) H[j][k] = (H[j][k] + (LL)H[i + 1][k] * (Mod - coef)) % Mod;
				for (int k = 1; k <= N; ++k) H[k][i + 1] = (H[k][i + 1] + (LL)H[k][j] * coef) % Mod;
			}
		}
	}

	void Characteristic() {
		Hessenberg();
		for (int i = 1; i <= N; ++i)
			for (int j = 1; j <= N; ++j)
				H[i][j] = Mod - H[i][j];
		static int P[MN][MN];
		P[0][0] = 1;
		for (int i = 1; i <= N; ++i) {
			P[i][0] = 0;
			for (int j = 1; j <= i; ++j) P[i][j] = P[i - 1][j - 1];
			int val = 1;
			for (int j = i - 1; j >= 0; --j) {
				int coef = (LL)val * H[j + 1][i] % Mod;
				for (int k = 0; k <= j; ++k) P[i][k] = (P[i][k] + (LL)P[j][k] * coef) % Mod;
				if (j) val = (LL)val * (Mod - H[j + 1][j]) % Mod;
			}
		}
		for (int i = 0; i <= N; ++i) Ans[i] = P[N][i];
	}

}

const int _ = 510;

namespace Gauss {
	int bas[_][_], mat[_][_];
	void solve(int N, int adj[]) {
		memset(bas, 0, sizeof(bas));
		for (int i = 1; i <= N; ++i) bas[i][i] = 1;
		for (int i = 1; i <= N; ++i) {
			int piv = 1;
			while (piv <= N && !mat[i][piv]) ++piv;
			if (piv > N) {
				memcpy(adj, bas[i], sizeof(bas[i]));
				break;
			}
			int c = mpow(mat[i][piv], P - 2);
			for (int j = piv; j <= N; ++j)
				mat[i][j] = mat[i][j] * (ull)c % P;
			for (int j = 1; j <= i; ++j)
				bas[i][j] = bas[i][j] * (ull)c % P;
			for (int j = i + 1; j <= N; ++j) {
				c = P - mat[j][piv];
				for (int k = piv; k <= N; ++k)
					fam(mat[j][k], mat[i][k], c);
				for (int k = 1; k <= i; ++k)
					fam(bas[j][k], bas[i][k], c);
			}
		}
	}
}

struct G {
	int N;
	int mat[_][_];
	int adj[_];
	vector<int> charp;

	int det() {
		int ret = 1;
		for (int i = 1; i <= N; ++i) {
			int pivot = i;
			while (!mat[pivot][i] && pivot <= N) ++pivot;
			if (pivot > N) return 0;
			if (pivot != i) {
				ret = P - ret;
				for (int j = i; j <= N; ++j) std::swap(mat[i][j], mat[pivot][j]);
			}
			ret = ret * (ull)mat[i][i] % P;
			int nv = mpow(mat[i][i], P - 2);
			for (int j = i; j <= N; ++j) mat[i][j] = mat[i][j] * (ull)nv % P;
			for (int j = i + 1; j <= N; ++j) for (int k = N; k >= i; --k)
				fam(mat[j][k], P - mat[j][i], mat[i][k]);
		}
		return ret;
	}

	void solve() {
		for (int i = 1; i <= N; ++i)
			for (int j = 1; j <= N; ++j) {
				int x; cin >> x;
				add(mat[j][j], x);
				sub(mat[i][j], x);
			}
		PinkRabbit::N = N;
		memcpy(PinkRabbit::A, mat, sizeof(mat));
		PinkRabbit::Characteristic();
		charp.resize(N + 1);
		memcpy(charp.data(), PinkRabbit::Ans, sizeof(int) * (N + 1));
		for (int i = 1; i <= N; ++i)
			for (int j = 1; j != i; ++j)
				swap(mat[i][j], mat[j][i]);
		memcpy(Gauss::mat, mat, sizeof(mat));
		Gauss::solve(N, adj);
		int piv = 1; while (!adj[piv]) ++piv;
		for (int i = 1; i <= N; ++i) mat[piv][i] = mat[i][piv] = 0;
		mat[piv][piv] = 1;
		int r = det() * (ull)mpow(adj[piv], P - 2) % P;
		for (int i = 1; i <= N; ++i) adj[i] = adj[i] * (ull)r % P;
	}
} A, B;

int res(vector<int> p, vector<int> q) {
	int ret = 1, N = p.size() - 1, M = q.size() - 1;
	if (N < M) {
		swap(N, M); swap(p, q);
		if ((N * M) % 2) ret = reduce(-ret);
	}
	while (M) {
		for (int i = N - M; i >= 0; --i) {
			int c = reduce(-p[M + i]);
			for (int j = 0; j <= M; ++j)
				fam(p[j + i], q[j], c);
		}
		while (N && p[N] == 0) --N; p.resize(N + 1);
		if (p[N] == 0) return 0;
		ret = ret * (ull)mpow(p[N], M) % P;
		int c = mpow(p[N], P - 2);
		for (int i = 0; i <= N; ++i) p[i] = p[i] * (ull)c % P;
		swap(N, M); swap(p, q);
		if ((N * M) % 2) ret = reduce(-ret);
	}
	return ret;
}

int main() {
	ios::sync_with_stdio(false); cin.tie(NULL);

	cin >> A.N >> B.N;
	A.solve(); B.solve();

	vector<int> p = A.charp, q = B.charp;
	for (int i = p.size() - 2; i >= 0; i -= 2)
		p[i] = reduce(-p[i]);
	p.erase(p.begin()); q.erase(q.begin());
	int r = res(p, q);
	for (int i = 1; i <= A.N; ++i)
		for (int j = 1; j <= B.N; ++j)
			cout << (A.adj[i] * (ull)B.adj[j] % P * r % P) << " \n"[j == B.N];

	return 0;
}

[1] 潘佳奇，浅谈线性代数与图论的关系，IOI 2021 中国国家集训队论文