矩阵快速幂

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
namespace Matrix {
// Matrix mat(row, vec(col));
typedef vector<ll> vec;
typedef vector<vec> mat;
mat mul(mat& A, mat& B) {
mat C(A.size(), vec(B[0].size()));
for (int i = 0; i < A.size(); i++)
for (int k = 0; k < B.size(); k++)
if (A[i][k])
for (int j = 0; j < B[0].size(); j++)
C[i][j] = (C[i][j] + A[i][k] * B[k][j]) % mod;
return C;
}
mat Pow(mat A, ll n) {
mat B(A.size(), vec(A.size()));
for (int i = 0; i < A.size(); i++) B[i][i] = 1;
while (n) {
if (n & 1) B = mul(B, A);
A = mul(A, A);
n >>= 1;
}
return B;
}
} // namespace Matrix
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
struct Mat{
ll m[maxn][maxn];
};
Mat operator * (Mat a, Mat b){
Mat c;
memset(c.m,0,sizeof(c.m));
for(int i=0;i<n;i++){
for(int k=0;k<n;k++){
if(a.m[k][i] == 0) continue;
for(int p=0;p<n;p++){
if(b.m[i][p] == 0) continue;
c.m[k][p] = (c.m[k][p] + a.m[k][i] * b.m[i][p] % mod ) % mod;
}
}
}
return c;
}

Mat operator ^ (Mat a,int x){
Mat c;
for(int i=0;i<n;i++)
for(int k=0;k<n;k++)
c.m[i][k] = (i == k);
for(; x ; x >>= 1){
if(x & 1) c = c * a;
a = a * a;
}
return c;
}