Description
给定 (n imes n) 矩阵 (mathbf{A}) 和正整数 (k),求和 (mathbf{T}=mathbf{A}+mathbf{A}^2+mathbf{A}^3+cdots+mathbf{A}^k)。矩阵元素对 (m) 取模。
(n leq 30,kleq 10^9,m<10^4)
Solution
像一些数列的前 (n) 项和很容易想到用矩阵快速幂实现。
而对于等比矩阵的前 (n) 项和我们同样设法用矩阵实现。
我们记矩阵 [mathbf S=egin{bmatrix}mathbf A &mathbf E\0&mathbf E\end{bmatrix}]
其中 (mathbf E) 为单位矩阵。
可知 ({mathbf S}^2=egin{bmatrix}mathbf A^2 &mathbf A+mathbf E\0&mathbf E\end{bmatrix},{mathbf S}^3=egin{bmatrix}mathbf A^3 &mathbf A^2+mathbf A+mathbf E\0&mathbf E\end{bmatrix},cdots,{mathbf S}^{k+1}=egin{bmatrix}mathbf A^{k+1} &mathbf A^k+mathbf A^{k-1}+cdots+mathbf A^2+mathbf A+mathbf E\0&mathbf E\end{bmatrix})
故答案为 ({mathbf S}^{k+1}) 的右上角矩阵减去一个单位矩阵。
Code
#include <bits/stdc++.h>
#define ll long long
using namespace std;
const int N = 65;
int n, k, m;
struct mat {
int a[N][N];
mat() {memset(a, 0, sizeof(a)); }
mat operator * (const mat &b) const {
mat ans;
for (int i = 1; i <= (n<<1); i++)
for (int j = 1; j <= (n<<1); j++)
for (int k = 1; k <= (n<<1); k++)
(ans.a[i][j] += (a[i][k]*b.a[k][j])) %= m;
return ans;
}
} S, E;
int main() {
scanf("%d%d%d", &n, &k, &m);
for (int i = 1; i <= n; i++) {
for (int j = 1; j <= n; j++)
scanf("%d", &S.a[i][j]), S.a[i][j] %= m;
S.a[i][i+n] = S.a[i+n][i+n] = 1;
}
E = S;
while (k) {
if (k&1) S = S*E;
k >>= 1, E = E*E;
}
for (int i = 1; i <= n; i++) S.a[i][i+n]--;
for (int i = 1; i <= n; i++)
for (int j = 1; j <= n; j++)
printf("%d%c", (S.a[i][j+n]+m)%m, "
"[j == n]);
return 0;
}