「2018 集训队互测 Day 1」完美的旅行 (容斥)(BM)(多项式)

本文介绍了一种结合多项式求逆和矩阵乘法的算法,用于解决特定类型的路径计数问题。通过使用生成函数和多项式求逆,算法能够在O(mnlog(m))的时间复杂度内找到所有可能路径的方案数。同时,通过矩阵乘法,文章展示了如何在O(n^4+n^2m)的时间复杂度内找到从任意起点到任意终点的路径计数。
摘要由CSDN通过智能技术生成

传送门

假设我们已经求出一条路径与为 x x x 的方案数,发现多条路径放在一起并不好做、

考虑记录超集, f [ a ] [ x ] f[a][x] f[a][x] 表示与是 a a a 的超集,走了 x x x 步的方案数

F a ( x ) F_a(x) Fa(x) 为它的生成函数,设 G a ( x ) G_a(x) Ga(x) 为任意多条路径的生成函数,枚举条数

G a ( x ) = 1 + F a ( x ) + F a ( x ) 2 + . . . = 1 1 − F a ( x ) G_a(x)=1+F_a(x)+F_a(x)^2+...=\frac{1}{1-F_a(x)} Ga(x)=1+Fa(x)+Fa(x)2+...=1Fa(x)1

可以多项式求逆做,复杂度 O ( m n l o g ( m ) ) O(mnlog(m)) O(mnlog(m))

考虑求 f [ t ] [ i ] [ j ] f[t][i][j] f[t][i][j] 表示起点为 i i i 终点为 j j j 的方案数,可以 O ( n 3 m ) O(n^3m) O(n3m)

发现 f [ t ] [ i ] [ j ] f[t][i][j] f[t][i][j] 可以用矩阵乘法转移,那么存在线性递推式

先暴力前 n n n 项,用 B M BM BM 解出线性递推式,然后递推即可

复杂度 O ( n 4 + n 2 m ) O(n^4+n^2m) O(n4+n2m)


#include<bits/stdc++.h>
#define cs const
using namespace std;
cs int Mod = 998244353;
int add(int a, int b){ return a + b >= Mod ? a + b - Mod : a + b; }
int mul(int a, int b){ return 1ll * a * b % Mod; }
int ksm(int a, int b){ int ans = 1; for(;b;b>>=1,a=mul(a,a)) if(b&1) ans = mul(ans, a); return ans; }
int dec(int a, int b){ return a - b < 0 ? a - b + Mod : a - b; }
void Add(int &a, int b){ a = add(a, b); }
cs int N = 64, M = 2e4 + 5;
int n, m, a[N][N], f[N];
#define poly vector<int> 
poly F[N], ans[N];
namespace NTT{
	cs int C = 18;
	poly w[C+1]; int up, bit; poly rev;
	void init(int len){
		up = 1, bit = 0; while(up < len) up <<= 1, ++bit; 
		rev.resize(up); for(int i = 0; i < up; i++) rev[i] = (rev[i>>1]>>1)|((i&1)<<(bit-1));
	} 	
	void prework(){
		for(int i = 1; i <= C; i++) w[i].resize(1 << i - 1);
		int wn = ksm(3, (Mod-1)/(1<<C)); w[C][0] = 1;
		for(int i = 1; i < (1 << C-1); i++) w[C][i] = mul(w[C][i-1], wn);
		for(int i = C-1; i; i--) for(int j = 0; j < (1 << i-1); j++) w[i][j] = w[i+1][j<<1];
	}
	void NTT(poly &a, int typ){
		for(int i = 0; i < up; i++) if(i < rev[i]) swap(a[i], a[rev[i]]);
		for(int i = 1, l = 1; i < up; i <<= 1, ++l)
		for(int j = 0; j < up; j += (i<<1))
		for(int k = 0; k < i; k++){
			int x = a[k+j], y = mul(w[l][k], a[k+j+i]);
			a[k+j] = add(x, y); a[k+j+i] = dec(x, y);
		}
		if(typ == -1){
			reverse(a.begin() + 1, a.end());
			for(int i = 0, iv = ksm(up, Mod - 2); i < up; i++) a[i] = mul(a[i], iv);
		}
	}
	poly inv(poly a){
		int deg = a.size(); poly b(1, ksm(a[0], Mod - 2)), c;
		for(int len = 4; len < (deg << 2); len <<= 1){
			init(len); c = a; c.resize(len >> 1);
			c.resize(up); b.resize(up);
			NTT(b, 1); NTT(c, 1);
			for(int i = 0; i < up; i++) b[i] = mul(b[i], dec(2, mul(b[i],c[i])));
			NTT(b, -1); b.resize(len >> 1);
		} b.resize(deg); return b;
	}
}
poly operator + (poly a, poly b){
	int deg = max(a.size(), b.size()); a.resize(deg); b.resize(deg);
	for(int i = 0; i < deg; i++) Add(a[i], b[i]); return a; 
} 
namespace BM{
	poly best, cur; int fail = -1, val = 0;
	void extend(int n){
		int t = f[n];
		for(int i = 1; i < cur.size(); i++) Add(t, Mod-mul(cur[i], f[n - i]));
		if(!t) return;
		if(cur.empty()){ cur.resize(n + 1); fail = n; val = t; return; }
		int coef = 	mul(t, ksm(val, Mod - 2));
		poly nxt(n - fail);
		nxt.push_back(coef);
		for(int i = 1; i < best.size(); i++) nxt.push_back(mul(coef, dec(0, best[i])));
		nxt = nxt + cur; 
		if((int)cur.size() - n <= (int)best.size() - fail){ best = cur; fail = n; val = t; }
		cur = nxt;
	}
}
int main(){
	NTT::prework();
	scanf("%d%d", &n, &m);
	for(int i = 0; i < n; i++)
		for(int j = 0; j < n; j++)
			scanf("%d", &a[i][j]);
	int len = (n << 1) + 5; 
	static int way[N << 2][N][N];
	static int ret[N][M]; 
	for(int i = 0; i < n; i++) way[0][i][i] = 1;
	for(int t = 1; t <= len; t++){
		for(int i = 0; i < n; i++)
		for(int j = 0; j < n; j++) if(way[t - 1][i][j])
		for(int k = 0; k < n; k++)
			Add(way[t][i][k], mul(way[t - 1][i][j], a[j][k]));
		for(int i = 0; i < n; i++)
		for(int j = 0; j < n; j++) Add(ret[i & j][t], way[t][i][j]);
	}
	for(int t = 1; t <= len; t++){
		int vl = 0;
		for(int i = 0; i < n; i++) vl = add(mul(vl, 12345), ret[i][t]);
		f[t] = vl; BM::extend(t);
	}
	poly p = BM::cur;
	for(int i = 0; i < n; i++)
	for(int t = len + 1; t <= m; t++)
	for(int j = 1; j < p.size(); j++) Add(ret[i][t] , mul(ret[i][t - j], p[j]));
	for(int i = 0; i < n; i++){
		for(int j = i + 1; j < n; j++)
			for(int t = 1; t <= m; t++) if((i|j) == j) Add(ret[i][t], ret[j][t]);
		poly G; G.push_back(1);
		for(int j = 1; j <= m; j++) G.push_back(dec(0, ret[i][j]));
		ans[i] = NTT::inv(G);
	}
	int sum = 0;
	for(int t = 1; t <= m; t++){
		for(int i = n - 1; ~i; i--){
			for(int j = i + 1; j < n; j++)
				if((i | j) == j) Add(ans[i][t], Mod - ans[j][t]);
			sum ^= ans[i][t];
		} 
	} cout << sum; return 0;
}
评论
添加红包

请填写红包祝福语或标题

红包个数最小为10个

红包金额最低5元

当前余额3.43前往充值 >
需支付:10.00
成就一亿技术人!
领取后你会自动成为博主和红包主的粉丝 规则
hope_wisdom
发出的红包

打赏作者

FSYo

你的鼓励将是我创作的最大动力

¥1 ¥2 ¥4 ¥6 ¥10 ¥20
扫码支付:¥1
获取中
扫码支付

您的余额不足,请更换扫码支付或充值

打赏作者

实付
使用余额支付
点击重新获取
扫码支付
钱包余额 0

抵扣说明:

1.余额是钱包充值的虚拟货币,按照1:1的比例进行支付金额的抵扣。
2.余额无法直接购买下载,可以购买VIP、付费专栏及课程。

余额充值