PKU 3233 經典矩陣乘法

S = A + A2 +
A3 + … + Ak. 的一個很好的求法是

構造這樣一個矩陣

A A

0 1

然後這個矩陣自乘K次即可，也就是矩陣套矩陣

 

#include <cstdio>#include <cstring>const int MAX = 65;int n, k, m, tn, mod;struct Mat{int mat[MAX][MAX];Mat() {memset(mat, 0, sizeof(mat));}void init() {for(int i = 0; i < n; i++) {for(int j = 0; j < n; j++)mat[i][j] = i == j;}}void print() {printf("**************\n");for(int i = 0; i <n; i++) {for(int j = 0; j < n; j++) {printf("%d ", mat[i][j]);}printf("\n");}printf("****************\n");}friend Mat operator *(Mat a, Mat b);friend Mat operator +(Mat a, Mat b);friend Mat operator ^(Mat a, int k);}E, A;int a[MAX][MAX];Mat operator +(Mat a, Mat b) {Mat c;for(int i = 0; i < n; i++)for(int j = 0; j < n; j++) {c.mat[i][j] = a.mat[i][j] + b.mat[i][j];if(c.mat[i][j] >= mod) c.mat[i][j] -= mod;}return c;}Mat operator *(Mat a, Mat b) {Mat ans;for(int i = 0; i < n; i++) {for(int j = 0; j < n; j++) {for(int k = 0; k <n; k++) {long long tmp = (long long)a.mat[i][k] * b.mat[k][j];if(tmp > mod) tmp %= mod;ans.mat[i][j] = ans.mat[i][j] + tmp;if(ans.mat[i][j] >= mod) ans.mat[i][j] -= mod;}}}return ans;}Mat operator ^(Mat a, int k) {Mat ans =  E;while(k) {if(k & 1) ans = ans * a;a = a * a, k >>= 1;}return ans;}void init() {              E.init();for(int i = 0; i < tn; i++) {    //A A for(int j = 0; j <tn; j++) {A.mat[i][j] = a[i][j];A.mat[i][j+tn] = a[i][j];}}for(int i = tn; i < n; i++) {   //0 1for(int j = tn; j < n; j++) {if(i == j) A.mat[i][j] = 1;}}}int main() {while(scanf("%d%d%d", &tn, &m, &mod) != EOF) {for(int i = 0; i < tn; i++) {for(int j = 0; j < tn; j++) {scanf("%d", &a[i][j]);}}n = 2 * tn;init();Mat ans = A ^(m);long long sum = 0;for(int i = 0; i < tn; i++)for(int j = tn; j < n; j++)printf(j == n-1?"%d\n":"%d ", ans.mat[i][j]);}return 0;}