POJ3233 power of matrix(矩阵乘法+二进制优化)

传送门

题意:给定一个(n*n)的矩阵A,求(sum_{i=1}^{k}A^i),运算值对m取模.((n<=30,k<=10^9))

分析:如果题目是求(A^k),就是矩阵快速幂的模板题了.但现在我们要求的是(sum_{i=1}^{k}A^i),如果直接求时间复杂度是(O(kn^3)),T飞.因为本题最基础的肯定是矩阵乘法,n又这么小,所以时间复杂度(n^3)是不可避免的,所以我们要从k下手,因为(sqrt{k})都会超时,所以肯定是要想方设法优化到(logk).

对于(sum_{i=1}^{k}A^i)这种式子,它一定是可以由前面的(k比较小的)推出后面的(k比较大的),例如,假设我们已知原式(A^1+A^2+A^3),要求(A^4+A^5+A^6),则可以让原式乘上(A^3),不难发现如果再加上一个原式,就得到(A^1+A^2+A^3+A^4+A^5+A^6),这不就是我们想要的吗?

为了方便,我们设(B_x=A^1+A^2+...+A^x),即(B_3=A^1+A^2+A^3),则(B_6=B_3*A^3+B_3),我们要求的答案就是(B_k).

而且我们现在要把k优化到(logk),就要想到可以二进制优化.我们预处理求出(A^{2^0},A^{2^1},A^{2^2},A^{2^3},...,A^{2^{31}}),这同时也可以得到(B_{2^0},B_{2^1},B_{2^2},B_{2^3},...,B_{2^{31}}).当然了,数组不能开(2^{31})这么大,所以我们直接用A[i]表示(A^{2^i}),B[i]表示(B^{2^i}).

我们预处理完A,B数组后,考虑我们的答案(B_k),首先k一定能分解成几个2的整数次幂相加的形式,例如(k=7=2^0+2^1+2^2),这时我们想,要是(B_7=B_{1+2+4}=B_1+B_2+B_4)那该多好啊!

虽然这是不合法的,但它的确为我们提供了一种思路.结合上面那个(B_6=B_3*A^3+B_3)的例子,你能否推出(B_{x+y}=B_x+A^x*B_y)这个式子呢?(其中x和y都是2的整数次幂)

int n,m,K;
struct mat{
    int jz[35][35];
    void in(){
	for(int i=1;i<=n;i++)
	    for(int j=1;j<=n;j++)
			jz[i][j]=read()%m;
    }//读入
    void out(){
		for(int i=1;i<=n;i++){
	    	for(int j=1;j<=n;j++)
				printf("%d ",jz[i][j]%m);
	    	printf("
");
		}
    }//输出
//直接重载运算符'*'为矩阵乘法
    mat operator *(mat x){
		mat y;
		for(int i=1;i<=n;i++)
	    	for(int j=1;j<=n;j++)
				y.jz[i][j]=0;//记得初始化为0
		for(int i=1;i<=n;i++)
	    	for(int j=1;j<=n;j++)
				for(int k=1;k<=n;k++)
		    		y.jz[i][j]=y.jz[i][j]%m+(jz[i][k]*x.jz[k][j])%m;
		return y;
    }
//直接重载运算符'*'为矩阵加法
    mat operator +(mat x){
		for(int i=1;i<=n;i++)
	    	for(int j=1;j<=n;j++)
				x.jz[i][j]=(jz[i][j]+x.jz[i][j])%m;
		return x;
    }
}ans,a[32],b[32];
int main(){
    n=read();K=read();m=read();
    a[1].in();
//上面的分析中,说a[i]表示a[1<<i],但为了能够从
//下标1开始操作,代码中a[i]表示a[1<<(i-1)],b同理
//所以a[1]就是A^1
    b[1]=a[1];
    for(int i=1;i<=30;i++){
		a[i+1]=a[i]*a[i];//我们重载了乘号和加号
		b[i+1]=b[i]*a[i]+b[i];
    }//预处理a,b数组
    for(int i=1;i<=31;i++){//二进制分解K
		if(K&1)ans=ans*a[i]+b[i];
//K&1表示K分解成二进制数后的第i位是1,
//这里的ans可以理解为b[ans]
		K>>=1;
    }
    ans.out();
    return 0;
}