Matrix Power Series

链接:

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

---

题目：

Description

Given a n × n matrix A and a positive integer k, find the sum S=A1+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

---

题意：

  给你一个矩阵A，让你求S=∑ni=1Ai。

---

思路：

  这是矩阵乘法中关于等比矩阵的求法，不难建立转移矩阵：

[EOAA]

  其中E为单位矩阵，0为0矩阵，然后我们求这个转移矩阵的n次幂，所得到的矩阵的右上角就是我们所要求的S。

---

实现：

#include <iostream>#include <algorithm>#include <set>#include <string>#include <vector>#include <queue>#include <map>#include <stack>#include <iomanip>#include <functional>#include <cstdio>#include <cstring>#include <cmath>#include <cctype>#define il inline#define ll long long#define ull unsigned long longusing namespace std;const int maxn = 107, Matmatsize = 2;int Matrixsize = 2, mod = int(1e9)+7;struct Matrix {    int m[maxn][maxn];    Matrix(int i = 0) {        memset(m, 0, sizeof m);        if (i == 1) {            for (int I = 0; I < Matrixsize; I++) m[I][I] = 1;        }    }    Matrix operator * (const Matrix tmp) const {        Matrix ret;        long long x;        for(int i=0 ; i<Matrixsize ; i++) {            for(int j=0 ; j<Matrixsize ; j++) {                x=0;                for(int k=0 ; k<Matrixsize ; k++) {                    x+=((long long)m[i][k] * tmp.m[k][j]) % mod;                }                ret.m[i][j] = int(x % mod);            }        }        return ret;    }    Matrix operator+(const Matrix b) const {        Matrix a = *this;        for (int i = 0; i < Matrixsize; i++)            for (int j = 0; j < Matrixsize; j++)                a.m[i][j] = (a.m[i][j] + b.m[i][j]) % mod;        return a;    }    Matrix operator-(const Matrix b) const {        Matrix a = *this;        for (int i = 0; i < Matrixsize; i++)            for (int j = 0; j < Matrixsize; j++)                a.m[i][j] = (a.m[i][j] - b.m[i][j]) % mod;        return a;    }    Matrix operator-() const {        Matrix tmp = *this;        for (int i = 0; i < Matrixsize; i++)            for (int j = 0; j < Matrixsize; j++)                tmp.m[i][j] = -tmp.m[i][j];        return tmp;    }    Matrix qpow(long long n) {        Matrix ret = 1, tmp = *this;        while (n != 0) {            if (bool(n & 1)) ret = ret * tmp;            tmp = tmp * tmp;            n >>= 1;        }        return ret;    }};struct Matmat {    Matrix m[Matmatsize][Matmatsize];    Matmat (int i=0) {        for(int i=0 ; i<Matmatsize ; i++)            for(int j=0 ; j<Matmatsize ; j++)                m[i][j] = 0;        if(i == 1) for(int i=0 ; i<Matmatsize ; i++) m[i][i] = 1;    }    Matmat operator * (const Matmat tmp) const {        Matmat ret;        Matrix x;        for(int i=0 ; i<Matmatsize ; i++) {            for(int j=0 ; j<Matmatsize ; j++) {                x = 0;                for(int k=0 ; k<Matmatsize ; k++) {                    x = x + m[i][k] * tmp.m[k][j];                }                ret.m[i][j] = x;            }        }        return ret;    }    Matmat qpow(long long n) {        Matmat ret = 1, tmp = *this;        while (n != 0) {            if (bool(n & 1))                ret = ret * tmp;            tmp = tmp * tmp;            n >>= 1;        }        return ret;    }};int main() {#ifndef ONLINE_JUDGE    freopen("in.txt", "r", stdin);#endif    ios_base::sync_with_stdio(false);cin.tie(0);    int k;    while(cin >> Matrixsize >> k >> mod) {        Matrix E = 1, A = 0, ans;        for (int i = 0; i < Matrixsize; i++)            for (int j = 0; j < Matrixsize; j++) cin >> A.m[i][j];        Matmat tmp = 0;        tmp.m[0][0] = E, tmp.m[0][1] = tmp.m[1][1] = A;        ans = tmp.qpow(k).m[0][1];        for (int i = 0; i < Matrixsize; i++) {            for (int j = 0; j < Matrixsize - 1; j++) cout << ans.m[i][j] << ' ';            cout << ans.m[i][Matrixsize - 1] << '\n';        }    }    return 0;}