【清华集训】石家庄的工人阶级队伍比较坚强
给定 (m),然后令 (n=3^m),给定数组 (f_{0,0},f_{0,1}...f_{0,3^m-1})
给定数组 (b),然后定义卷积:
其中 (u) 定义为 (W(x,y)),(v) 定义为 (L(x,y))
其中 (W(x,y)) 被定义为将三进制的每一位进行比较,出现一组 10/21/02
就累加一次贡献,(L(x,y)) 定义为 (W(y,x))
你需要计算 (f_{t}),(t) 给定。
答案对 (p) 取模,保证不存在 (x,y) 满足 (frac{1}{x}+frac{1}{y}=frac{3}{p})
(mle 12,tle 10^9,ple 10^9)
Solution
又来水题解了。
考虑三进制不退位减法:
我们发现结果为 (1) 等价于 win,结果为 (2) 等价于 lose
设 (bit_k(x)) 表示 (x) 中的 (k) 的数量。
那么我们有转移形如:
不难发现贡献实际上只和 (xoplus y) 相关。
设 (odot) 表示三进制不进位加法,那么不难发现:
于是相当于计算 3 进制下异或卷积的 (k) 次幂。
我们搞出 FWT 的点值即可处理了。
注意到 (k) 进制 FWT 相当于在做循环卷积,这样我们只需要将给 (k) 进制单位根代入进去即可。
具体式子:
这样递归即可,由于只需要乘以单项式,所以复杂度是 (mathcal O(k^{m+3})) 的。
由于要次幂一下,会掉精度,所以只能扩域,我们找一个代数单位来表示 (omega_3)
根据单位根的性质,我们注意到 (omega_3^2+omega_3+1=0),这样的话就有 (omega_3^2=-omega_3-1)
这样我们先相当于拿着 (ax^2+bx+c) 类型的多项式在模(?)上做运算,最后我们再把 (a) 表示成 (b,c) 相关消一消。
- 坑:这个模数保证了 (3) 存在逆元,但是不能直接套质数情况的求法...
然后发现这个题卡这个做法,然后你会发现没必要 FWT 的时候取模,所以可以最后取模。
然后手动循环展开一下,就可以过了。最后放一下比较丑的代码
(Code:)
#include<bits/stdc++.h>
using namespace std ;
#define Next( i, x ) for( register int i = head[x]; i; i = e[i].next )
#define rep( i, s, t ) for( register int i = (s); i <= (t); ++ i )
#define drep( i, s, t ) for( register int i = (t); i >= (s); -- i )
#define re register
#define int long long
int gi() {
char cc = getchar() ; int cn = 0, flus = 1 ;
while( cc < '0' || cc > '9' ) { if( cc == '-' ) flus = - flus ; cc = getchar() ; }
while( cc >= '0' && cc <= '9' ) cn = cn * 10 + cc - '0', cc = getchar() ;
return cn * flus ;
}
const int N = 1e6 + 5 ;
const int M = 20 + 5 ;
int lim, m, t, P, b[M][M] ;
struct gf {
int a[3] ;
void init() { a[0] = a[1] = a[2] = 0 ; }
} f[N], g[N], Nx[N] ;
int fpow(int x, int k) {
int ans = 1, base = x ;
while(k) {
if(k & 1) ans = 1ll * ans * base % P ;
base = 1ll * base * base % P, k >>= 1 ;
} return ans ;
}
void inc(int &x, int y) { (x += y) ; }
int bit1(int x) { int ans = 0 ; while( x ) ans += (x % 3 == 1), x /= 3 ; return ans ; }
int bit2(int x) { int ans = 0 ; while( x ) ans += (x % 3 == 2), x /= 3 ; return ans ; }
gf Move(gf x, int l) {
gf ans ; ans.init() ;
ans.a[(0 + l) % 3] = x.a[0], ans.a[(1 + l) % 3] = x.a[1], ans.a[(2 + l) % 3] = x.a[2] ;
return ans ;
}
gf operator * (gf x, gf y) {
gf ans ; ans.init() ;
ans.a[0] = ( x.a[0] * y.a[0] + x.a[1] * y.a[2] + x.a[2] * y.a[1] ) % P ;
ans.a[1] = ( x.a[0] * y.a[1] + x.a[1] * y.a[0] + x.a[2] * y.a[2] ) % P ;
ans.a[2] = ( x.a[0] * y.a[2] + x.a[1] * y.a[1] + x.a[2] * y.a[0] ) % P ;
return ans ;
}
gf operator * (gf x, int y) { gf ans ; rep(i, 0, 2) ans.a[i] = x.a[i] * y % P ; return ans ; }
gf operator + (gf x, gf y) { gf ans ; rep(i, 0, 2) ans.a[i] = (x.a[i] + y.a[i]) ; return ans ; }
gf fpow(gf x, int k) {
gf ans, base = x ; ans.init(), ans.a[0] = 1 ;
while(k) {
if(k & 1) ans = ans * base ;
base = base * base, k >>= 1 ;
} return ans ;
}
gf operator % (gf x, int p) { rep( i, 0, 2 ) x.a[i] %= p ; return x ; }
void FWT(gf *a) {
for(re int k = 1; k < lim; k *= 3 ) for(re int i = 0; i < lim; i += (k * 3) )
for(re int j = i; j < i + k; ++ j) {
rep( l, 0, 2 ) Nx[l].init() ;
rep( l, 0, 2 ) rep( o, 0, 2 )
Nx[l] = (Nx[l] + Move(a[j + o * k], o * l % 3)) ;
rep( l, 0, 2 ) a[j + l * k] = Nx[l] ;
}
}
void IFWT(gf *a) {
for(re int k = 1; k < lim; k *= 3 ) for(re int i = 0; i < lim; i += (k * 3) )
for(re int j = i; j < i + k; ++ j) {
rep( l, 0, 2 ) Nx[l].init() ;
rep( l, 0, 2 ) rep( o, 0, 2 )
Nx[l] = (Nx[l] + Move(a[j + o * k], 3 - o * l % 3)) ;
rep( l, 0, 2 ) a[j + l * k] = Nx[l] ;
}
}
int phi(int x) {
int ans = x ;
for(re int i = 2; i <= sqrt(x); ++ i) {
if(x % i) continue ;
ans = (ans - ans / i) ;
while(x % i == 0) x /= i ;
} if( x != 1 ) ans = ans - ans / x ;
return ans ;
}
signed main()
{
m = gi(), t = gi(), P = gi(), lim = 1 ;
rep( i, 1, m ) lim = lim * 3 ;
rep( i, 0, lim - 1 ) f[i].a[0] = gi() ;
rep( i, 0, m ) rep( j, i, m ) b[i][j - i] = gi() ;
rep( i, 0, lim - 1 ) g[i].a[0] = b[bit1(i)][bit2(i)] ;
FWT(f), FWT(g) ;
rep( i, 0, lim - 1 ) f[i] = f[i] % P, g[i] = g[i] % P ;
rep( i, 0, lim - 1 ) f[i] = f[i] * fpow(g[i], t) ;
int iv = fpow(lim, phi(P) - 1) ; IFWT(f) ;
rep( i, 0, lim - 1 ) f[i] = f[i] % P ;
rep( i, 0, lim - 1 ) f[i] = f[i] * iv ;
rep( i, 0, lim - 1 ) {
int d = (f[i].a[0] - f[i].a[2] + P) % P ;
printf("%lld
", d ) ;
}
return 0 ;
}