1. 程式人生 > >【POJ3233】Matrix Power Series

【POJ3233】Matrix Power Series

題面

給出矩陣A,求S = A + A2 + A3 + … + Ak.

分析

矩陣的乘方是可以通過快速冪很快的推出來,主要是相加的問題

但是別忘了,雖然沒有等比求和,但是矩陣是滿足結合律的

因此求出了A + A2 + A3 + … + Ak/2後只需要乘Ak/2,就可以得到A(k/2)+1 + A(k/2)+2 + A(k/2)+3 + … + Ak

類似地,一直二分下去。(上面忽略了k的奇偶,寫的時候要注意)

程式碼

#include<cstdio>
#include<cstring>
#include<algorithm>
#include<iostream>
using namespace std;
struct email
{
    int x[33][33];
}a,o;
int n,m,k;
email mul(email a,email b)
{
    email c;
    memset(c.x,0,sizeof(c.x));
    for(int i=1;i<=n;i++)
        
for(int j=1;j<=n;j++) for(int k=1;k<=n;k++) c.x[i][j]=(a.x[i][k]*b.x[k][j]+c.x[i][j])%m; return c; } email add(email a,email b) { int i,j; for(i=1;i<=n;i++) for(j=1;j<=n;j++) a.x[i][j]+=b.x[i][j],a.x[i][j]%=m; return
a; } email ksm(email a,int k) { if(k==1)return a; email b=ksm(a,k/2); b=mul(b,b); if(k&1)b=mul(b,a); return b; } email divide(email a,int k) { if(k==1)return a; if(k%2==0) { email b=divide(a,k/2); return mul(add(ksm(a,k/2),o),b); } else { email b=divide(a,(k)/2); return add(mul(add(ksm(a,k/2),o),b),ksm(a,k)); } } int main() { scanf("%d%d%d",&n,&k,&m); for(int i=1;i<=n;i++)o.x[i][i]=1; for(int i=1;i<=n;i++) for(int j=1;j<=n;j++) scanf("%d",&a.x[i][j]); email ans=divide(a,k); for(int i=1;i<=n;i++) { for(int j=1;j<=n;j++) printf("%d ",ans.x[i][j]); printf("\n"); } return 0; }