Matrix Power Series

Description

Time Limit: 3000MS Memory Limit: 131072K
Given a $n*n$ matrix $A$ and a positive integer $k$, find the sum $S = A + A^2 + A^3 + … + A^k$.

Input

The input contains exactly one test case. The first line of input contains three positive integers $n (n ≤ 30)$, $k (k ≤ 10^9)$ and $m (m < 10^4)$. 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

Solution & Code

法一：$k$的范围是$10^9$，考虑分治。
$S_i=A+A^2+...+A^i$
$S_{2i}=A^{i+1}+A^{i+2}+...+A^{2i}$
$S_{2i}=S_i*A^i$

#include <cstdio>
#include <cstring>
#include <iostream>
#include <algorithm>
using namespace std;

typedef long long ll;

const int maxr = 35;

int n, k;
ll m;

struct matrix{ll val[maxr][maxr];};
matrix U, V, Z, I;

matrix multiply(matrix A, matrix B){

    matrix C = Z;
    for(int i = 1; i <= n; ++i)
        for(int j = 1; j <= n; ++j)
            for(int k = 1; k <= n; ++k){
                C.val[i][j] += A.val[i][k] * B.val[k][j];
                C.val[i][j] %= m;
            }
    return C;
}

matrix add(matrix A, matrix B){

    matrix C = Z;
    for(int i = 1; i <= n; ++i)
        for(int j = 1; j <= n; ++j){
            C.val[i][j] = A.val[i][j] + B.val[i][j];
            C.val[i][j] %= m;
        }
    return C;
}

matrix powmod(matrix A, int x){

    matrix B = I;
    while(x){
        if(x & 1) B = multiply(B, A);
        x = x >> 1;
        A = multiply(A, A);
    }
    return B;
}

matrix calc(int x){

    if(x == 1) return U;
    int y = x >> 1;
    matrix L = calc(y);
    matrix R = multiply(L, powmod(U, y));
    if(y * 2 != x) return add(add(L, R), powmod(U, x));
        else return add(L, R);
}

int main(){

    cin >> n >> k >> m;
    for(int i = 1; i <= n; ++i){
        for(int j = 1; j <= n; ++j) cin >> U.val[i][j];
        I.val[i][i] = 1;
    }
    V = calc(k);
    for(int i = 1; i <= n; ++i){
        for(int j = 1; j <= n; ++j) printf("%lld ", V.val[i][j]);
        puts("");
    }

    return 0;
}

法二：如果求$s=a+a^2+...a^k$我们可以考虑用矩阵快速幂对其进行加速，求$S = A + A^2 + A^3 + … + A^k$时则可以构造一个元素为矩阵的矩阵，然后用矩阵快速幂加速。