bzoj4128 Matrix BSGS&hash

       顺便学了一下矩阵求逆（其实就是一个高斯消元），不过后来没有用。

       采用BSGS的思想，实际上令m=sqrt(p)，那么就是求最小的i*m-j使得：

      (A^m)^i=A^jB，其中A^j一共m个，可以预处理。

       因此关键是快速判断两个矩阵是否相等。可以将矩阵看成一个大的p进制数然后将其取模后得到的值设为其hash值，判断为O(1)；或者用另外一种方法，即将矩阵去乘一个(n*1)的矩阵，得到一个(n*1)的矩阵，然后就可以O(N)判断了。后面一种实现更简单，而且更快（用小号刷到rk2辣）。

       时间复杂度O(P^0.5N^2+N^3+PN)。

AC代码如下：

#include<iostream>
#include<cstdio>
#include<cstring>
#include<cstdlib>
#include<cmath>
#define N 75
using namespace std;

int n,mod,num[N];
struct matrix{
	int p[N][N];
	void init(){
		int i,j;
		for (i=1; i<=n; i++)
			for (j=1; j<=n; j++) scanf("%d",&p[i][j]);
	}
}a,b,c,hash[N<<1];
matrix tms(matrix x,matrix y){
	matrix z; int i,j,k;
	for (i=1; i<=n; i++)
		for (j=1; j<=n; j++){
			z.p[i][j]=0;
			for (k=1; k<=n; k++) z.p[i][j]=(z.p[i][j]+x.p[i][k]*y.p[k][j])%mod;
		}
	return z;
}
matrix mul(matrix x,matrix y){
	matrix z; int i,j;
	for (i=1; i<=n; i++)
		for (j=1; j<=n; j++){
			z.p[i][j]=0;
			z.p[i][1]=(z.p[i][1]+x.p[i][j]*y.p[j][1])%mod;
		}
	return z;
}
matrix ksm(matrix x,int y){
	matrix z=x;
	for (; y; y>>=1,x=tms(x,x)) if (y&1) z=tms(z,x);
	return z;
}
bool ok(matrix x,matrix y){
	int i; for (i=1; i<=n; i++) if (x.p[i][1]!=y.p[i][1]) return 0;
	return 1;
}
int main(){
	scanf("%d%d",&n,&mod); int m=(int)sqrt(mod)+1,i,j;
	srand(20160307); a.init(); b.init();
	for (i=1; i<=n; i++) c.p[i][1]=num[i]=rand()%mod;
	for (i=0; i<m; i++){
		hash[i]=mul(b,c); c=mul(a,c);
	}
	a=ksm(a,m-1);
	for (i=1; i<=n; i++) c.p[i][1]=num[i];
	for (i=1; i<=m; i++){
		c=mul(a,c);
		for (j=m-1; j>=0; j--)
			if (ok(c,hash[j])){ printf("%d\n",i*m-j); return 0; }
	}
	return 0;
}

by lych

2016.3.7