标签：矩阵
|
Matrix Power Series
Description Given a n × n matrix A and a positive integer k, find the sum S = A + A2 + A3 + … + Ak. Input The input contains exactly one test case. The first line of input contains three positive integers n (n ≤ 30), k (k ≤ 109) and m (m < 104). Then follow n lines each containing n nonnegative integers below 32,768, giving A’s elements in row-major order. Output Output the elements of S modulo m in the same way as A is given. Sample Input 2 2 4
0 1
1 1
Sample Output 1 2
2 3
Source |
题意：给出一个n*n的矩阵A，求A+A^2+A^3+……+A^k mod m的结果是多少？
方法一：两次二分
Sk = A + A2 + A3 + … + Ak
=(1+A^(k/2))*(A + A^2 + A^3 + … + A^(k/2) )+{A^k}
=(1+A^(k/2))*(S(k/2) )+ {Ak} 当k为偶数时没有{Ak}
即
k%2==0: S[k]=F[k/2] (1+A[k/2]);
k%2==1: S[k]=F[k-1]+A[k];
#include <cstdio> #include <cstring> #include <algorithm> using namespace std; const int N = 32; struct Matrix { int mat[N][N]; Matrix() { memset(mat, 0, sizeof(mat)); for(int i = 0; i < N; i++) mat[i][i] = 1; } } E; int n, k, mod; Matrix Multi(Matrix a, Matrix b) { Matrix res; for(int i = 0; i < n; i++) { for(int j = 0; j < n; j++) { res.mat[i][j] = 0; for(int k = 0; k < n; k++) res.mat[i][j] = (res.mat[i][j] + a.mat[i][k] * b.mat[k][j]) % mod; } } return res; } Matrix Add(Matrix a, Matrix b) { Matrix res; for(int i = 0; i < n; i++) for(int j = 0; j < n; j++) res.mat[i][j] = (a.mat[i][j] + b.mat[i][j]) % mod; return res; } Matrix Pow(Matrix a, int n) { Matrix res; while(n) { if(n&1) res = Multi(res, a); a = Multi(a, a); n >>= 1; } return res; } Matrix Get_Ans(Matrix a, int k) { if(k == 1) return a; if(k&1) return Add(Pow(a, k), Get_Ans(a, k-1)); if(k % 2 == 0) { Matrix A = Get_Ans(a, k/2); Matrix B = Pow(a, k/2); Matrix C = Multi(A, B); return Add(C, A); } } int main() { Matrix A; while(~scanf("%d%d%d", &n, &k, &mod)) { for(int i = 0; i < n; i++) for(int j = 0; j < n; j++) { scanf("%d", &A.mat[i][j]); A.mat[i][j] %= mod; } Matrix ans = Get_Ans(A, k); for(int i = 0; i < n; i++) { for(int j = 0; j < n; j++) { if(j) printf(" "); printf("%d", ans.mat[i][j]); } printf("\n"); } } return 0; }方法二：
构造一个新矩阵B = | A E |
| 0 E |
则 B^(k+1) = | A^(k+1) A^k+A^(k-1)+……+A^2 + A + 1 |
| 0 E |
所以可以先直接用矩阵快速幂求出B^(k+1),然后用左上角的那一部分减去单位矩阵就是最后要求的矩阵。
#include <cstdio> #include <cstring> #include <algorithm> using namespace std; const int N = 61; int n, mod, k; struct Matrix { int mat[N][N]; Matrix() { memset(mat, 0, sizeof(mat)); for(int i = 0; i < N; i++) mat[i][i] = 1; } }; Matrix Multi(Matrix a, Matrix b) { Matrix res; memset(res.mat, 0, sizeof(res.mat)); for(int i = 0; i < n * 2; i++) { for(int k = 0; k < n * 2; k++) { if(a.mat[i][k]) { for(int j = 0; j < n * 2; j++) { res.mat[i][j] += a.mat[i][k] * b.mat[k][j]; res.mat[i][j] %= mod; } } } } return res; } Matrix Pow(Matrix x, int m) { Matrix res; while(m) { if(m&1) res = Multi(res, x); x = Multi(x, x); m >>= 1; } return res; } int main() { Matrix A; while(~scanf("%d%d%d",&n, &k, &mod)) { memset(A.mat, 0, sizeof(A.mat)); for(int i = 0; i < n; i++) { for(int j = 0; j < n; j++) { scanf("%d", &A.mat[i][j]); A.mat[i][j] %= mod; } A.mat[i][i+n] = 1; } for(int i = n; i < 2 * n; i++) A.mat[i][i] = 1; Matrix ans = Pow(A, k + 1); for(int i = 0; i < n; i++) { for(int j = n; j < n * 2; j++) { if(j > n) printf(" "); if(j - i == n) printf("%d", (ans.mat[i][j] - 1 + mod) % mod); else printf("%d", ans.mat[i][j]); } printf("\n"); } } return 0; }
POJ 3233 Matrix Power Series（矩阵快速幂）
标签：矩阵
原文地址：http://blog.csdn.net/lyhvoyage/article/details/39584771
