矩阵幂和(SOJ 2919)
题意:计算$(\sum_{i=1}^{k}\mathbf{A}^{i})\%m$.
分析:$f_{k+1}=\mathbf{A}+\mathbf{A}^{2}+\cdots+\mathbf{A}^{k+1}=\mathbf{A}+\mathbf{A}(\mathbf{A}+\cdots+\mathbf{A}^{k})=\mathbf{A}+\mathbf{A}\cdot f_{k}$. 据此,我们可得递推式:
$\left(\begin{array}{c}f_{k+1}\\ \mathbf{A}\end{array}\right)=\left(\begin{array}{cc} \mathbf{A} &\mathbf{I}\\ \mathbf{0} &\mathbf{I} \end{array}\right)\left(\begin{array}{c} f_{k}\\ \mathbf{A}\end{array}\right)$
代码:(注意,这里我们计算矩阵的幂跟之前的方法有所不同,从二进制末位开始计算)
1 #include<iostream> 2 #include<string.h> 3 using namespace std; 4 int A[31][31]; 5 int m; 6 struct matrix 7 { 8 int P[62][62]; 9 }; 10 matrix matMul(matrix P1, matrix P2,int a,int b,int c) 11 { 12 matrix P3; 13 memset(P3.P, 0, sizeof(P3.P)); 14 int i, j, k; 15 for (i = 1; i <= a; i++) 16 for (j = 1; j <= c; j++) 17 for (k = 1; k <= b; k++) 18 P3.P[i][j] = (P3.P[i][j]+P1.P[i][k]*P2.P[k][j])%m; 19 return P3; 20 } 21 matrix matPow(matrix M, int n, int k) 22 { 23 matrix Res; 24 int i; 25 memset(Res.P, 0, sizeof(Res.P)); 26 for (i = 1; i <= 2 * n; i++) 27 Res.P[i][i] = 1; 28 while (k) 29 { 30 if (k & 1) 31 Res = matMul(Res,M,2*n,2*n,2*n); 32 k >>= 1; 33 M=matMul(M,M,2*n,2*n,2*n); 34 } 35 return Res; 36 } 37 int main() 38 { 39 int n, k; 40 int i, j; 41 matrix B; 42 while (scanf("%d%d%d", &n, &k, &m) == 3) 43 { 44 memset(B.P, 0, sizeof(B.P)); 45 for (i = 1; i <= n; i++) 46 for (j = 1; j <= n; j++) 47 { 48 scanf("%d", &A[i][j]); 49 B.P[i][j] = A[i][j]; 50 } 51 for (i = 1; i <= n; i++) 52 { 53 B.P[i][n + i] = 1; 54 B.P[n + i][n + i] = 1; 55 } 56 matrix Fin = matPow(B, n, k-1); 57 matrix Fin1, Fin2; 58 for (i = 1; i <= n; i++) 59 for (j = 1; j <= 2 * n; j++) 60 Fin1.P[i][j] = Fin.P[i][j]; 61 for (i = 1; i <= n; i++) 62 for (j = 1; j <= n;j++) 63 { 64 Fin2.P[i][j] = A[i][j]; 65 Fin2.P[n + i][j] = A[i][j]; 66 } 67 Fin = matMul(Fin1, Fin2, n, 2 * n, n); 68 for (i = 1; i <= n; i++) 69 { 70 for (j = 1; j < n; j++) 71 printf("%d ", Fin.P[i][j]); 72 printf("%d\n", Fin.P[i][n]); 73 } 74 } 75 return 0; 76 }