http://poj.org/problem?id=3233
1 #include <iostream> 2 3 using namespace std; 4 5 const int N = 50; 6 int n, k, mod; 7 8 // 定义结构体 9 struct mat { 10 int m[N][N]; 11 mat () { 12 memset(m, 0, sizeof m); 13 for (int i = 0; i < N; i ++) { 14 m[i][i] = 1; 15 } 16 } 17 }; 18 19 // 矩阵乘法 20 mat multi(mat a, mat b) { 21 mat res; 22 for (int i = 0; i < n; i ++) { 23 for (int j = 0; j < n; j ++ ) { 24 res.m[i][j] = 0; 25 for (int k = 0; k < n; k ++) { 26 res.m[i][j] = (res.m[i][j] + a.m[i][k] * b.m[k][j]) % mod; 27 } 28 } 29 } 30 return res; 31 } 32 33 // 矩阵相加 34 mat add(mat a, mat b) { 35 mat res; 36 for (int i = 0; i < n; i ++) { 37 for (int j = 0; j < n; j++) { 38 res.m[i][j] = (a.m[i][j] + b.m[i][j]) % mod; 39 } 40 } 41 return res; 42 } 43 44 // 矩阵快速幂 45 mat fastpow(mat a, int k) { 46 mat res; 47 while (k) { 48 if (k & 1) res = multi(res, a); 49 a = multi(a, a); 50 k >>= 1; 51 } 52 return res; 53 } 54 55 // 矩阵求和 56 // sum(k) = (E+A^{k/2}) * sum(k/2); // 偶数 57 // sum(k) = (E+A^{k-1}/2) * sum((k - 1)/2) + A^k; // 奇数 58 mat sum(mat a, int k) { 59 mat tmp; 60 if (k == 1) return a; 61 if (k & 1) { 62 mat t = multi(sum(a, (k-1)/2),add(tmp, fastpow(a, (k-1) / 2))); 63 return add(t, fastpow(a, k)); 64 } else { 65 return multi(sum(a, k/2),add(tmp, fastpow(a, k / 2))); 66 } 67 } 68 69 int main() { 70 cin >> n >> k >> mod; 71 72 // 读矩阵 73 mat a; 74 for (int i = 0; i < n; i++) { 75 for (int j = 0; j < n; j ++) { 76 cin >> a.m[i][j]; 77 } 78 } 79 80 mat res = sum(a, k); 81 82 // 输出矩阵 83 for (int i = 0; i < n; i++) { 84 for (int j = 0; j < n; j ++) { 85 cout << res.m[i][j] << " "; 86 } 87 cout << endl; 88 } 89 90 91 return 0; 92 }