POJ-3233+3070（矩阵加、乘、乘方）（）

题意：已知一个n*n的矩阵A，和一个正整数k，求S = A + A² + A³ + … + A^k。

思路：矩阵快速幂。首先我们知道 A^x 可以用矩阵快速幂求出来（具体可见poj 3070）。其次可以对k进行二分，每次将规模减半，分k为奇偶两种情况，如当k = 6和k = 7时有：

      k = 6 有： S(6) = (1 + A^3) * (A + A^2 + A^3) = (1 + A^3) * S(3)。

      k = 7 有： S(7) = A + (A + A^4) * (A + A^2 + A^3) = A + (A + A^4) * S(3)。

ps：对矩阵定义成结构体Matrix，求S时用递归，程序会比较直观，好写一点。当然定义成数组，然后再进行一些预处理，效率会更高些。

3070是一个缩水版，不说它了。

#define MAX 35int MOD = 10000;int M, N;struct matrix{int v[MAX][MAX];} t;matrix matrix_add(matrix a, matrix b)//a+b{matrix c;int i, j;for (i = 0; i < M; ++i) {for (j = 0; j < N; ++j) {c.v[i][j] = a.v[i][j] + b.v[i][j];c.v[i][j] %= MOD;}}return c;}matrix matrix_mul(matrix a, matrix b)//a*b{int i, j, k;matrix c;for (i = 0; i < M; ++i) {for (j = 0; j < N; ++j) {c.v[i][j] = 0;for (k = 0; k < N; ++k) {c.v[i][j] += a.v[i][k] * b.v[k][j];c.v[i][j] %= MOD;}}}return c;}matrix matrix_pow(matrix m, int k)//m^k{int i, j;if (k == 0) {memset(m.v, 0, sizeof(m.v));for (i = 0; i < M; ++i)m.v[i][i] = 1;return m;}if (k == 1) return m;matrix ans = matrix_pow(m, k / 2);if (k % 2 == 0)return matrix_mul(ans, ans);elsereturn matrix_mul(matrix_mul(ans, ans), m);}matrix while_matrix_pow(matrix m, int k)//m^k非递归版本{matrix tmp = m;matrix result = matrix_pow(m, 0);//得到单位矩阵while (k) {if ((k & 1) == 1)result = matrix_mul(result, tmp);tmp = matrix_mul(tmp, tmp);k >>= 1;}return result;}matrix matrix_cal(matrix a, int k)//求S(k) = A + A2 + A3 + … + Ak{if (k == 1) return a;matrix x = while_matrix_pow(a, (k + 1) / 2);matrix y = matrix_cal(a, k / 2);if (k % 2 == 0) {return matrix_mul(matrix_add(matrix_pow(a, 0), x), y);} else {return matrix_add(a, matrix_mul(matrix_add(a, x), y));}}int main(){int a, b, c;while (scanf("%d%d%d", &a, &b, &c) == 3) {M = a;N = a;MOD = c;matrix t;int i, j;for (i = 0; i < M; ++i) {for (j = 0; j < N; ++j) {scanf("%d", &t.v[i][j]);}}matrix ans = matrix_cal(t, b);for (i = 0; i < M; ++i) {for (j = 0; j < N - 1; ++j) {printf("%d ", ans.v[i][j]);}printf("%d\n", ans.v[i][j]);}}    return 0;}