【清华集训2016】石家庄的工人阶级队伍比较坚强

我作为沙比OI选手果然沙比。

明明都推出来了我都没看到（把减号写个异或？）

直接使用二次剩余版本在这 /kk

---

首先我们手动打个比较表，可以发现，比较是一个减法，我们构造 (0) 为平， (1) 为胜， (2) 为负。

对于两个选手 (x)，(y)，对 (x) 贡献就是 (b_{mathrm{cnt}left(x - y, 1 ight),mathrm{cnt}left(x - y, 2 ight)})。其中 (-) 是三进制不退位减法，(cnt(x,y)) 是 (x) 在三进制下值为 (y) 的数位的个数。

我们记 (B_x = b_{mathrm{cnt}left(x, 1 ight),mathrm{cnt}left(x, 2 ight)})

然后式子就变成了

[egin{align*} f_{i,x} & = sum_{y} f_{i - 1, y} B_{x-y} \ f_{i,z} & = sum_{x + y = z} f_{i - 1, x} B_{y} end{align*} ]

其中 (+) 就是三进制不进位加法，可以当做异或。

那么显然我们构造出了卷积，可以使用高维 DFT 来解决（也就是三进制 FWT）。其中做 (T) 次也就是点值快速幂。

那么我们需要三次单位根。三次单位根有很好的性质 (omega^2 + omega + 1 = 0)，解得一个根 (omega = frac{-1 + sqrt{3} mathrm{i}}{2})。

然后这个 (-3) 在模 (P) 意义下的二次剩余是有的，但是我不会证，也不会excipolla，要学习这个姿势的dalao们请看 yhx-12243 的 blog

所以我们决定扩域，使用 (a + by) 解决，其中 (y^2 = -3)，显然可以用一个二元组实现乘法和加法。

那么我们有了单位根自然可以做 DFT 了。直接套用DFT矩阵 写个高维 DFT 就行了 （FFT什么的不存在的，这个就是直接DFT就行了）。

根据单位根性质，注意最后要除以序列总长度。

复杂度 (Oleft(3^n left(n + log T ight) ight))

Hack数据里有一个假的数据，不满足题目限制，然后 yhx-12243 跑 Excipolla 挂了？ 反正扩域没事情。只是常数大了 (4) 倍，垫底了。

#include <bits/stdc++.h>

const int MAXN = 531441;
int m, T, mod, pw3[13], N;
typedef long long LL;
void reduce(int & x) { x += x >> 31 & mod; }
int mul(int a, int b) { return (LL) a * b % mod; }
int fastpow(int a, int b) {
	int res = 1;
	for (; b; b >>= 1, a = mul(a, a)) if (b & 1) res = mul(res, a);
	return res;
}
int os;
struct comp {
	int r, i;
	comp() { r = i = 0; }
	comp(int a, int b) { r = a, i = b; }
	comp operator * (comp b) {
		return comp(((LL) r * b.r + (LL) i * b.i % mod * os) % mod, ((LL) r * b.i + (LL) i * b.r) % mod);
	}
	comp operator + (comp b) {
		comp t(r + b.r - mod, i + b.i - mod);
		reduce(t.r); reduce(t.i);
		return t;
	}
} A[MAXN], B[MAXN];
comp fastpow(comp a, int b) {
	comp res(1, 0);
	for (; b; b >>= 1, a = a * a) if (b & 1) res = res * a;
	return res;
}

comp W, Ws;
void FWT(comp * A, int typ) {
	comp Wn = typ == 1 ? W : Ws;
	comp Wns = typ == 1 ? Ws : W;
	for (int mid = 1; mid != N; mid *= 3)
		for (int k = 0; k != N; k += mid * 3)
			for (int l = 0; l != mid; ++l) {
				comp & x = A[l + k], & y = A[l + k + mid], & z = A[l + k + mid * 2];
				comp ta = x + y + z;
				comp tb = x + y * Wn + z * Wns;
				comp tc = x + y * Wns + z * Wn;
				x = ta, y = tb, z = tc;
			}
}

int tr[20][20];
int get(int x, int y) {
	int res = 0;
	while (x) res += x % 3 == y, x /= 3;
	return res;
}
int main() {
	std::ios_base::sync_with_stdio(false), std::cin.tie(0);
	std::cin >> m >> T >> mod;
	os = mod - 3;
	pw3[0] = 1;
	for (int i = 1; i != 13; ++i) pw3[i] = pw3[i - 1] * 3;
	N = pw3[m];
	const int iv2 = mod + 1 >> 1;
	const int iv3 = mul(mod % 3 == 1 ? 1 : iv2, mod - mod / 3);
	const int liminv = fastpow(iv3, m);
	W = comp(mod - iv2, iv2);
	Ws = W * W;
	for (int i = 0; i < N; ++i) std::cin >> A[i].r;
	for (int i = 0; i <= m; ++i)
		for (int j = 0; j + i <= m; ++j)
			std::cin >> tr[i][j];
	for (int i = 0; i < N; ++i)
		B[i].r = tr[get(i, 1)][get(i, 2)];
	FWT(A, 1); FWT(B, 1);
	for (int i = 0; i < N; ++i)
		A[i] = A[i] * fastpow(B[i], T);
	FWT(A, -1);
	for (int i = 0; i < N; ++i) A[i].r = mul(A[i].r, liminv);
	for (int i = 0; i < N; ++i) std::cout << A[i].r << '
';
	return 0;
}