poj 3233 Matrix Power Series

分类：分治 难度：1.5

 

题意：给定矩阵a，求s=a^1+…+a^k，结果mod m。矩阵乘法，两次分治，计算a^k和a^1+…+a^k。即a^k=(a*a)^(k/2), a^1+...+a^k=a^1+...a^(k/2) + a^(k/2) * [ a^1+...a^(k/2) ]。

#include<cstdio>
#include<cstring>

int n,k,m;
const int N = 32;

typedef struct
{
int d

;
}ma;
ma e;

ma madd(ma x, ma y)
{
int i,j;
ma z;
for(i=0; i<n; i++)
{
for(j=0; j<n; j++)
{
z.d[i][j] = x.d[i][j] + y.d[i][j];
z.d[i][j] %= m;
}
}
return z;
}

ma mmult(ma x, ma y)
{
int i,j,l;
int ans

;
ma z;
for(i=0; i<n; i++)
{
for(j=0; j<n; j++)
{
z.d[i][j] = 0;
for(l=0; l<n; l++)
{
z.d[i][j] += x.d[i][l]*y.d[l][j];
z.d[i][j] %= m;
}
}
}
return z;
}

ma mpow(ma a,int k)
{
if(k==0)
return e;
if(k==1)
return a;
if(k%2)
return mmult(mpow(a,k-1),a);
return mpow(mmult(a,a),k/2);
}

void print(ma x)
{
for(int i=0; i<n; i++)
{
for(int j=0; j<n; j++)
printf("%d ",x.d[i][j]);
printf("\n");
}
}

ma doit(ma a, int k)
{
if(k<=1)
return a;
if(k%2)
return madd(doit(a,k-1),mpow(a,k));
ma temp = doit(a,k/2);
return madd(temp,mmult(temp,mpow(a,k/2)));
}

int main()
{
int i,j;
ma a;
ma s;

e.d[0][0] = e.d[1][1] = 1;
e.d[0][1] = e.d[1][0] = 0;
scanf("%d%d%d",&n,&k,&m);

for(i=0; i<n; i++)
{
for(j=0; j<n; j++)
{
scanf("%d",&a.d[i][j]);
s.d[i][j] = 0;
}
}

s = doit(a,k);
print(s);
getchar();getchar();
}