Pentium.Labs

System全家桶:https://zhuanlan.zhihu.com/c_1238468913098731520

导航

poj 3233 矩阵快速幂+YY

题意:给你矩阵A,求S=A+A^1+A^2+...+A^n

 

sol:直接把每一项解出来显然是不行的,也没必要。

我们可以YY一个矩阵:

其中1表示单位矩阵

 

然后容易得到:

可以看出这个分块矩阵的左下角那块就可以得到要求的解S

我们取这一块,再减去一个单位矩阵1即可。

 

 1 #include "iostream"
 2 #include "vector"
 3 #include "cstring"
 4 #include "cstdio"
 5 using namespace std;
 6 
 7 //typedef unsigned long int ULL;
 8 typedef vector<int> vec;
 9 typedef vector<vec> mat;
10 int n,m,k,P;
11 
12 mat mul(mat &A,mat &B)      //return A*B
13 {
14     mat C(A.size(),vec(B[0].size()));
15     for (int i=0;i<(int)A.size();i++)
16     {
17         for (int k=0;k<(int)B.size();k++)
18         {
19             for (int j=0;j<(int)B[0].size();j++)
20             {
21                 C[i][j]=(C[i][j]+A[i][k]*B[k][j])%P;
22             }
23         }
24     }
25     return C;
26 }
27 
28 mat m_pow(mat A,int m)      //return A^m
29 {
30     mat B(A.size(),vec(A.size()));
31     for (int i=0;i<(int)A.size();i++)
32         B[i][i]=1;
33     while (m>0)
34     {
35         if (m&1)    B=mul(B,A);
36         A=mul(A,A);
37         m>>=1;
38     }
39     return B;
40 }
41 
42 int main()
43 {
44     //freopen("in.txt","r",stdin);
45 
46     cin>>n>>k>>P;
47     mat A(2*n,vec(2*n));
48 
49     for (int i=0;i<n;i++)
50         for (int j=0;j<n;j++)
51             cin>>A[i][j];
52     for (int i=n;i<=2*n-1;i++)
53     {
54         A[i][i-n]=1;
55         A[i][i]=1;
56     }
57 
58     A=m_pow(A,k+1);
59 
60     for (int i=n;i<2*n;i++)
61     {
62         for (int j=0;j<n;j++)
63         {
64             if (i==j+n)
65                 A[i][j]--;
66             if (A[i][j]<0)  A[i][j]+=P;
67             cout<<A[i][j]<<" ";
68         }
69         cout<<endl;
70     }
71     return 0;
72 }

 

注意:挑战程序设计那本书P205上取了左上角那块,估计是写错了orz

 ----------------------------

本题还有一种方法,reference:  http://www.cnblogs.com/ACShiryu/archive/2011/08/09/poj3233.html

posted on 2014-10-30 21:02  Pentium.Labs  阅读(571)  评论(0编辑  收藏  举报



Pentium.Lab Since 1998