poj 3233 Matrix Power Series(midhard)

http://poj.org/problem?id=3233

  矩阵快速幂的题,求ΣA^i(1≤i≤k)的矩阵并输出。将它二分,用快速幂的方法来O(nlogn)的复杂度完成计算。

代码如下:

View Code
  1 #include <cstdio>
  2 #include <cstring>
  3 #include <cstdlib>
  4 #include <algorithm>
  5 
  6 using namespace std;
  7 typedef int ll;
  8 const int matSize = 31;
  9 const int stdMod = 1000000007;
 10 int calSize = matSize;
 11 int mod = stdMod;
 12 
 13 struct Matrix {
 14     ll val[matSize][matSize];
 15 
 16     Matrix(bool Init = false) {
 17         for (int i = 0; i < calSize; i++) {
 18             for (int j = 0; j < calSize; j++) {
 19                 val[i][j] = 0;
 20             }
 21             if (Init) val[i][i] = 1;
 22         }
 23     }
 24 
 25     void print() {
 26         for (int i = 0; i < calSize; i++) {
 27             for (int j = 0; j < calSize; j++) {
 28                 if (j) putchar(' ');
 29                 printf("%d", val[i][j]);
 30             }
 31             puts("");
 32         }
 33     }
 34 } Base;
 35 
 36 Matrix operator * (Matrix &_a, Matrix &_b) {
 37     Matrix ret = Matrix();
 38 
 39     for (int i = 0; i < calSize; i++) {
 40         for (int k = 0; k < calSize; k++) {
 41             if (_a.val[i][k]) {
 42                 for (int j = 0; j < calSize; j++) {
 43                     ret.val[i][j] += _a.val[i][k] * _b.val[k][j];
 44                     ret.val[i][j] %= mod;
 45                 }
 46             }
 47         }
 48     }
 49 
 50     return ret;
 51 }
 52 
 53 Matrix operator ^ (Matrix &_a, ll _p) {
 54     Matrix ret = Matrix(true);
 55     Matrix __a = _a;
 56 
 57     while (_p) {
 58         if (_p & 1) {
 59             ret = ret * __a;
 60         }
 61         __a = __a * __a;
 62         _p >>= 1;
 63     }
 64 
 65     return ret;
 66 }
 67 
 68 Matrix operator + (Matrix &_a, Matrix &_b) {
 69     Matrix ret;
 70 
 71     for (int i = 0; i < calSize; i++) {
 72         for (int j = 0; j < calSize; j++) {
 73             ret.val[i][j] = (_a.val[i][j] + _b.val[i][j]) % mod;
 74         }
 75     }
 76     //puts("~~~~~~");
 77     //ret.print();
 78 
 79     return ret;
 80 }
 81 
 82 void deal(int k) {
 83     Matrix tmp = Base, ans = Matrix(), ep = Base, one = Matrix(true), as;
 84 
 85     while (k) {
 86 //        puts("tmp");
 87 //        tmp.print();
 88         if (k & 1) {
 89             ans = ans * ep;
 90             ans = ans + tmp;
 91         }
 92         as = one + ep;
 93         tmp = tmp * as;
 94         ep = ep * ep;
 95         k >>= 1;
 96 //        puts("!!!");
 97 //        ans.print();
 98 //        puts("~~~");
 99     }
100 
101     ans.print();
102 }
103 
104 
105 int main() {
106     int k;
107 
108 //    freopen("in", "r", stdin);
109     while (~scanf("%d%d%d", &calSize, &k, &mod)) {
110         for (int i = 0; i < calSize; i++) {
111             for (int j = 0; j < calSize; j++) {
112                 scanf("%d", &Base.val[i][j]);
113                 Base.val[i][j] %= mod;
114             }
115         }
116         deal(k);
117     }
118 
119     return 0;
120 }

——written by Lyon

原文地址:https://www.cnblogs.com/LyonLys/p/poj_3233_Lyon.html