题解 LOJ2075 「JSOI2016」位运算（状压DP，矩阵）

题目大意

题目链接

给定两个整数(n), (k)和一个01串(S)。我们设(R)是一个二进制数，它的二进制表示，就是(S)重复(k)次。请你选出(n)个不同的、小于(R)的非负整数（也就是值在([0,R-1])之间），使得它们的异或和为(0)。

数据范围：(3leq nleq 7,1leq kleq 10^5,1leq |S|leq 50)。

本题题解

假设我们已经知道了(R)（也就是把(S)重复(k)次大力展开），该怎么做？设选出的这(n)个数为(x_1,x_2,dots,x_n)。为了保证这(n)个数互不相同且都小于(R)，不妨设(R>x_1>x_2>dots >x_n)。不妨设(x_0=R)。

因为(n)很小，考虑状压DP。设(dp[i][ ext{mask}])表示考虑了(R)的前(i)位（从高到低）；( ext{mask})是一个长度为(n)的二进制数，对于第(j)个位置 ((1leq jleq n))，如果当前（只考虑数的前(i)位）(x_j=x_{j-1})，则( ext{mask})第(j)位为(1)，否则( ext{mask})第(j)位为(0)。转移时，枚举(x_1dots x_n)的第(i)位分别填什么，这样可以计算出新的( ext{mask}')。令(dp[i][ ext{mask}'] exttt{+=}dp[i-1][ ext{mask}])即可。

这个DP的时间复杂度是(O(|R|cdot (2^n)^2cdot n)=O(kcdot |S|cdot 2^{2n}cdot n))的，无法通过本题。

发现上面的这个DP，没有用到“(R)是由一个非常短的串(S)，重复(k)次得来的”，这一特殊条件。我们发现，在(S)的(k)次重复中，每一次的转移其实都是一样的。那么就容易想到做矩阵快速幂。

那么关键就是要求出转移矩阵：( ext{trans}[ ext{mask}_1][ ext{mask}_2])表示从(dp[xcdot|S|][ ext{mask}_1])转移到(dp[(x+1)cdot|S|][ ext{mask}_2])的系数 ((0leq x<k))。可以枚举( ext{mask}_1)，令(dp[0][ ext{mask}_1]=1)，然后对(S)做一遍上面的那个DP，就可以对所有( ext{mask}_2)求出( ext{trans}[ ext{mask}_1][ ext{mask}_2])了。这部分时间复杂度(O(|S|cdot 2^{3n}cdot n))。可以承受。

然后就对这个大小为(2^n imes 2^n)的矩阵做矩阵快速幂即可。这部分时间复杂度(O(2^{3n}log k))。

总时间复杂度(O(|S|cdot 2^{3n}cdot n+2^{3n}log k))。

参考代码：

//problem:LOJ2075
#include <bits/stdc++.h>
using namespace std;

#define pb push_back
#define mk make_pair
#define lob lower_bound
#define upb upper_bound
#define fi first
#define se second
#define SZ(x) ((int)(x).size())

typedef unsigned int uint;
typedef long long ll;
typedef unsigned long long ull;
typedef pair<int,int> pii;

const int MOD=1e9+7;
inline int mod1(int x){return x<MOD?x:x-MOD;}
inline int mod2(int x){return x<0?x+MOD:x;}
inline void add(int& x,int y){x=mod1(x+y);}
inline void sub(int& x,int y){x=mod2(x-y);}
inline int pow_mod(int x,int i){int y=1;while(i){if(i&1)y=(ll)y*x%MOD;x=(ll)x*x%MOD;i>>=1;}return y;}

const int MAXN=7,MAXK=1e5,MAXS=50;
const int SIZE2=1<<MAXN;
int n,K,len,bitcnt[SIZE2],dp[MAXS+5][SIZE2];
char s[MAXS+5];
struct Matrix{
	int a[SIZE2][SIZE2],sz;
	Matrix(){
		memset(a,0,sizeof(a));
		sz=0;
	}
};
Matrix operator*(const Matrix& X,const Matrix& Y){
	Matrix Z;
	assert(X.sz==Y.sz);
	Z.sz=X.sz;
	for(int i=0;i<=X.sz;++i){
		for(int j=0;j<=X.sz;++j){
			for(int k=0;k<=X.sz;++k){
				add(Z.a[i][j],(ll)X.a[i][k]*Y.a[k][j]%MOD);
			}
		}
	}
	return Z;
}
Matrix mat_pow(Matrix X,int i){
	Matrix Y;
	Y.sz=X.sz;
	for(int j=0;j<=X.sz;++j)Y.a[j][j]=1;
	while(i){
		if(i&1)Y=Y*X;
		X=X*X;
		i>>=1;
	}
	return Y;
}

int main() {
	cin>>n>>K;
	cin>>(s+1);len=strlen(s+1);
	int sz=(1<<n)-1;
	for(int i=1;i<=sz;++i)bitcnt[i]=bitcnt[i>>1]+(i&1);
	Matrix trans;
	trans.sz=sz;
	for(int st=0;st<=sz;++st){
		memset(dp,0,sizeof(dp));
		dp[0][st]=1;
		for(int i=1;i<=len;++i){
			for(int j=0;j<=sz;++j)if(dp[i-1][j]){
				for(int k=0;k<=sz;++k)if(bitcnt[k]%2==0){
					static int curs[MAXN+5];
					curs[0]=s[i]-'0';
					for(int l=1;l<=n;++l)curs[l]=((k>>(l-1))&1);
					bool fail=false;
					int newj=0;
					for(int l=1;l<=n;++l){
						if((j>>(l-1))&1){
							//之前是等于的
							if(curs[l]>curs[l-1]){fail=1;break;}
							if(curs[l]==curs[l-1])newj|=(1<<(l-1));
						}
					}
					if(fail)continue;
					add(dp[i][newj],dp[i-1][j]);
				}
			}
		}
		for(int ed=0;ed<=sz;++ed){
			trans.a[st][ed]=dp[len][ed];
		}
	}
	trans=mat_pow(trans,K);
	Matrix res;
	res.a[0][sz]=1;
	res.sz=sz;
	res=res*trans;
	cout<<res.a[0][0]<<endl;
	return 0;
}