Matrix Power Series 题解

描述

Given a n × n matrix A and a positive integer k, find the sum S = A + A2 + A3 + … + Ak.

输入

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 the elements of S modulo m in the same way as A is given.

样例输入

2 2 4
0 1
1 1
 

样例输出

1 2
2 3
 

---

题意是对一个n*n的矩阵，求上述式子的结果。

一般来说用等比数列求和来计算，但是这里是矩阵没法求1-A和1-A^k，所以用分治的方式去求，其中需要求矩阵的乘方，用矩阵快速幂求解。

类似于一般的快速幂求法，在求矩阵快速幂里的ans数组要先赋值成单位矩阵，其他的和快速幂一样，因为不能直接返回数组，所以这个用一个结构体表示矩阵（一开始形参不小心写了二维数组就错了，救命）。

然后用分治的方式求等比数列的和。令等比数列从1到n的和为f[n]，分为两种情况：
        n为奇数：f[n] = f[n >> 1] + a ^ n;

        n为偶数：f[n] = f[n >> 1]

然后递归处理就可以得到答案了

#include<bits/stdc++.h>

using namespace std;
typedef long long LL;
typedef pair<int, int> PII;
const int N = 50;
int MOD;
LL n, m, k;
struct Martix {
	LL a[N][N];
	Martix(){memset(a, 0, sizeof a);}
};

Martix a, ans;
// int h[N], e[N], ne[N], idx;
// void add(int a, int b) {
// 	e[idx] = b, ne[idx] = h[a], h[a] = idx ++;
// }

Martix mul(Martix a, Martix b) {
    Martix c;
    for(int i = 0;i < n;i ++) {
        for(int j = 0;j < n;j ++) {
            for(int k = 0;k < n;k ++)
                c.a[i][j] = (c.a[i][j] + a.a[i][k] * b.a[k][j]) % MOD;
        }
    }

	return c;
}

Martix add(Martix a, Martix b) {
    Martix c;
    for(int i = 0;i < n;i ++) {
        for(int j = 0;j < n;j ++) {
    	   c.a[i][j] = (a.a[i][j] + b.a[i][j]) % MOD;
        }
    }
	return c;
}

Martix ksm(Martix a, LL b) {
	Martix res;
	for(int i = 0;i < n;i ++)
		res.a[i][i] = 1;
	while(b) {
		if(b & 1) res = mul(res, a);
		a = mul(a, a);
		b >>= 1;
	}
	
	return res;
}

Martix solve(int x) {
	if(x == 1) return a;
	Martix t = solve(x >> 1);
	if(x & 1) {
		return add(ksm(a, x), add(t, mul(t, ksm(a, x >> 1))));
	}
	return add(t, mul(t, ksm(a, x >> 1)));
}

int main() { 
	ios::sync_with_stdio(false);
	cin.tie(0), cout.tie(0);
	cin >> n >> k >> MOD;
	for(int i = 0;i < n;i ++)
		for(int j = 0;j < n;j ++)
			cin >> a.a[i][j];

	ans = solve(k);
	for(int i = 0;i < n;i ++){
		for(int j = 0;j < n;j ++){
			if(j) cout << ' ';
			cout << ans.a[i][j] % MOD;
		}
		cout << '\n';
	}		
	
	return 0;
}