【题意】：给出一个n * n矩阵A和一个正整数k,求
S =
A +
A2 +
A3 + … +
Ak
【题解】：以下是matrix67的题解：
这道题两次二分，相当经典。首先我们知道，A^i可以二分求出。然后我们需要对整个题目的数据规模k进行二分。比如，当k=6时，有：
A + A^2 + A^3 + A^4 + A^5 + A^6 =(A + A^2 + A^3) + A^3*(A + A^2 + A^3)
应用这个式子后，规模k减小了一半。我们二分求出A^3后再递归地计算A + A^2 + A^3，即可得到原问题的答案。
【代码】：
1 #include "iostream"
2 #include "cstring"
3 #include "cstdio"
4 using namespace std;
5 #define maxn 35
6 int n, k, m;
7 struct Mat {
8 int val[maxn][maxn];
9 void unit() {//单位矩阵
10 for(int i = 0; i < maxn; i++) val[i][i] = 1;
11 }
12 void zero() {//零矩阵
13 memset(val, 0, sizeof(val));
14 }
15 }x;
16
17 Mat operator *(const Mat &a, const Mat &b) {//矩阵乘法
18 Mat tmp;
19 tmp.zero();
20 for(int k = 1; k <= n; k++) {
21 for(int i = 1; i <= n; i++)
22 if(a.val[i][k])
23 for(int j = 1; j <= n; j++) {
24 tmp.val[i][j] += a.val[i][k] * b.val[k][j];
25 if(tmp.val[i][j] >= m) tmp.val[i][j] %= m;
26 }
27 }
28 return tmp;
29 }
30
31 Mat operator ^(Mat x, int n) {//矩阵快速幂
32 Mat tmp;
33 tmp.zero();
34 tmp.unit();
35 while(n) {
36 if(n & 1) tmp = tmp * x;
37 x = x * x;
38 n >>= 1;
39 }
40 return tmp;
41 }
42
43 Mat operator +(const Mat &a, const Mat &b) {//矩阵加法
44 Mat tmp;
45 for(int i = 1; i <= n; i++)
46 for(int j = 1; j <= n; j++)
47 tmp.val[i][j] = (a.val[i][j] + b.val[i][j]) % m;
48 return tmp;
49 }
50
51 Mat sum(int k) {
52 if(k == 1) return x;
53 else {
54 Mat tmp = sum(k >> 1);
55 if(k & 1) {
56 Mat tmp2 = x ^ ((k >> 1) + 1);
57 return tmp + tmp2 + tmp * tmp2;
58 } else {
59 Mat tmp2 = x ^ (k >> 1);
60 return tmp + tmp * tmp2;
61 }
62 }
63 }
64
65 int main() {
66 while(~scanf("%d%d%d", &n, &k, &m)) {
67 for(int i = 1; i <= n; i++) {
68 for(int j = 1; j <= n; j++) {
69 scanf("%d", &x.val[i][j]);
70 x.val[i][j] %= m;
71 }
72 }
73 Mat ans = sum(k);
74 for(int i = 1; i <= n; i++) {
75 for(int j = 1; j < n; j++)
76 printf("%d ", ans.val[i][j]);
77 printf("%d\n", ans.val[i][n]);
78 }
79 break;
80 }
81 return 0;
82 }