文章目录
- 基础知识
- 1. 矩阵结构
- 2. 重载 * 运算符
- 3. 矩阵快速幂
- 例题1: 矩阵幂求和
- 例题2: 矩阵快速幂
- 例题3: 斐波那契数列
- 例题4: 矩阵加速
- 例题5: 广义斐波那契
- 例题6: 斐波那契公约数
- 例题7: 这个勇者明明超强却过分慎重
基础知识
1. 矩阵结构
struct Matrix {
int g[N][N];
// 矩阵初始化。type 为 true 则初始化为 E,type 为 false 则初始化为 O。
void init(bool type) {
for (int i = 1; i <= n; i++)
for (int j = 1; j <= n; j++)
if(i != j) g[i][j] = 0;
else g[i][j] = type == true ? 1 : 0;
}
};
2. 重载 * 运算符
Matrix operator*(const Matrix &o) const {
Matrix t;
t.init(false); // 0 矩阵
for (int i = 1; i <= n; i++) // 三重循环
for (int j = 1; j <= n; j++)
for (int k = 1; k <= n; k++)
t.g[i][j] = (t.g[i][j] + o.g[i][k] * g[k][j]) % mod;
return t;
}
3. 矩阵快速幂
// 计算a^k
Matrix operator^(Matrix a, int k) { // 重载矩阵快速幂
Matrix res;
res.init(true);
while (k) {
if (k & 1) res = res * a;
a = a * a;
k >>= 1;
}
return res;
}
例题1: 矩阵幂求和
矩阵幂求和
S = A + A 2 + A 3 + . . . + A k S=A+A^2+A^3+...+A^k \\\\ S=A+A2+A3+...+Ak
- 推导如下:
1. 当 k 为偶数:
S = A + . . . + A k / 2 + A k / 2 + 1 + . . . + A k = A + . . . + A k / 2 + A k / 2 ∗ ( A + . . . + A k / 2 ) = ( A k / 2 + E ) ∗ S k / 2 \begin{align} S&=A+...+A^{k/2}+A^{k/2+1}+...+A^k\\ &=A+...+A^{k/2}+A^{k/2}*(A+...+A^{k/2}) \\ &=(A^{k/2}+E)*S_{k/2} \end{align} S=A+...+Ak/2+Ak/2+1+...+Ak=A+...+Ak/2+Ak/2∗(A+...+Ak/2)=(Ak/2+E)∗Sk/2
2. 当 k 为奇数
令
k
=
2
∗
n
+
1
,
2
∗
n
=
k
−
1
令 \ k = 2*n+1, \ \ \ 2*n = k - 1 \\
令 k=2∗n+1, 2∗n=k−1
S
=
A
+
A
2
+
A
3
+
.
.
.
+
A
2
n
+
1
=
A
+
A
∗
(
A
+
.
.
.
+
A
2
n
)
=
A
+
A
∗
S
2
n
=
A
+
A
∗
S
k
−
1
\begin{align} S&=A+A^2+A^3+...+A^{2n+1} \\ &=A+A*(A+...+A^{2n}) \\ &=A + A*S_{2n}\\ &=A+A*S_{k-1} \end{align}
S=A+A2+A3+...+A2n+1=A+A∗(A+...+A2n)=A+A∗S2n=A+A∗Sk−1
代码如下:
#include <iostream>
using namespace std;
// https://www.acwing.com/solution/content/15850/
typedef long long ll;
const int N = 110;
int n, mod, k;
// const int mod = 1e9 + 7;
struct Matrix {
int g[N][N];
// 矩阵初始化。type 为 true 则初始化为 E,type 为 false 则初始化为 O。
void init(bool type) {
for (int i = 1; i <= n; i++)
for (int j = 1; j <= n; j++)
if(i != j) g[i][j] = 0; // 两种矩阵的共同点:不在 i = j 对角线上的数皆为 0
else g[i][j] = type == true ? 1 : 0; // 不同点:在 i = j 对角线上,E 为 1,O 为 0
}
Matrix operator*(const Matrix &o) const {
Matrix t;
t.init(false);
for (int i = 1; i <= n; i++)
for (int j = 1; j <= n; j++)
for (int k = 1; k <= n; k++)
t.g[i][j] = (t.g[i][j] + o.g[i][k] * g[k][j]) % mod;
return t;
}
Matrix operator+(const Matrix &o) const // 重载矩阵加法
{
Matrix t;
for (int i = 1; i <= n; i++)
for (int j = 1; j <= n; j++)
t.g[i][j] = (g[i][j] + o.g[i][j]) % mod;
return t;
}
void print() { // 将该矩阵输出
for (int i = 1; i <= n; i++) {
for (int j = 1; j <= n; j++)
printf("%d ", g[i][j]);
puts("");
}
}
void read() { // 读入该矩阵
for (int i = 1; i <= n; i++)
for (int j = 1; j <= n; j++) {
scanf("%d", &g[i][j]);
g[i][j] %= mod;
}
}
};
// 计算a^k
Matrix operator^(Matrix a, int k) { // 重载矩阵快速幂。由于要用到乘法,所以在结构体外重载。
Matrix res;
res.init(true);
while (k) {
if (k & 1) res = res * a;
a = a * a;
k >>= 1;
}
return res;
}
Matrix E, v; // E 即 E 矩阵,v 为读入矩阵
// 计算S=A+A**2+A**3+...+A**k
Matrix S(int k) {
if (k == 1) return v; // 如果 k 为 1,那么返回 v,终止递归。
if (k & 1) return v + v * S(k - 1); // 如果 k 是奇数,那么返回上述 A + A * S(k - 1)
return (E + (v ^ k >> 1)) * S(k >> 1); // 否则返回上述 (E + A ^ (k / 2)) * S(k / 2)
}
int main() {
cin >> n >> k >> mod; // scanf
E.init(true); // 初始化矩阵E
v.read(); // 读入矩阵v
S(k).print(); // 计算S(k)并输出
return 0;
}
例题2: 矩阵快速幂
矩阵快速幂
- 框架如上题,PS: long long
scanf("%lld %lld", &n, &k);
E.init(true);
v.read();
(v^k).print();
例题3: 斐波那契数列
f n = { 1 , ( n < = 2 ) f n − 1 + f n − 2 , ( n > = 3 ) f_n=\begin{cases} 1, & (n<=2) \\ f_{n-1} + f_{n-2}, & (n>=3)\\ \end{cases} fn={1,fn−1+fn−2,(n<=2)(n>=3)
- 分析如下:
通过递推式,推出
1
×
2
∗
2
×
2
=
=
1
×
2
的方阵(
2
×
2
)
f
n
=
f
n
−
1
∗
1
+
f
n
−
2
∗
1
f
n
−
1
=
f
n
−
1
∗
1
+
f
n
−
2
∗
0
推出(
2
×
2
)矩阵:
∣
1
1
1
0
∣
通过递推式, 推出 1 × 2 \ * 2 × 2 == 1 × 2 的方阵(2×2)\\ \begin{align} f_n &= f_{n-1} * 1 + f_{n-2}* 1 \\ f_{n-1}&=f_{n-1}*1+f_{n-2}* 0 \\ \end{align} \\ 推出(2×2)矩阵: \begin{vmatrix} \mathbf{1} & \mathbf{1} \\ \mathbf{1} & \mathbf{0} \\ \end{vmatrix} \\
通过递推式,推出1×2 ∗2×2==1×2的方阵(2×2)fnfn−1=fn−1∗1+fn−2∗1=fn−1∗1+fn−2∗0推出(2×2)矩阵:
1110
最终结果:
[
f
n
f
n
−
1
]
=
∣
1
1
1
0
∣
n
−
1
∗
[
f
1
f
0
]
最终结果: [f_n \ f_{n-1}] = \begin{vmatrix} \mathbf{1} & \mathbf{1} \\ \mathbf{1} & \mathbf{0} \\ \end{vmatrix} ^{n-1} * [f_1 \ f_0]
最终结果:[fn fn−1]=
1110
n−1∗[f1 f0]
#include <iostream>
#include <cstring>
using namespace std;
typedef long long ll;
const int mod = 1e9 + 7;
const int N = 3;
struct Matrix {
ll g[N][N];
};
Matrix E, v; // E是单位阵(也是最终答案矩阵(快速幂)), v是推导出的初始的2×2 矩阵
void init() {
memset(E.g, 0, sizeof(E));
for (int i = 1; i <= 2; i++) E.g[i][i] = 1;
memset(v.g, 0, sizeof(v.g));
v.g[1][1] = v.g[1][2] = v.g[2][1] = 1;
}
Matrix mul(Matrix a, Matrix b) {
Matrix t;
memset(t.g, 0, sizeof(t));
for (int i = 1; i <= 2; i++)
for (int j = 1; j <= 2; j++)
for (int k = 1; k <= 2; k++)
t.g[i][j] = (t.g[i][j] + a.g[i][k] * b.g[k][j]) % mod;
return t;
}
void qmi(ll k) {
while(k) {
if(k & 1) E = mul(E, v);
v = mul(v, v);
k >>= 1;
}
}
void print(Matrix a) {
for (int i = 1; i <= 2; i++) {
for (int j = 1; j <= 2; j++)
cout << a.g[i][j] << " ";
cout << endl;
}
}
int main() {
ll n; cin >> n;
init();
qmi(n - 1);
cout << E.g[1][1] << endl;
return 0;
}
例题4: 矩阵加速
矩阵加速
#include <iostream>
#include <cstring>
using namespace std;
typedef long long ll;
const int mod = 1e9 + 7;
const int N = 4;
int T, n;
struct Matrix {
ll g[N][N];
};
Matrix E, v;
void init() {
memset(E.g, 0, sizeof(E.g));
for (int i = 1; i <= 3; i++) E.g[i][i] = 1;
memset(v.g, 0, sizeof(v.g));
v.g[1][1] = v.g[1][3] = v.g[2][1] = v.g[3][2] = 1;
}
Matrix mul(Matrix a, Matrix b) {
Matrix t;
memset(t.g, 0, sizeof(t.g));
for (int i = 1; i <= 3; i++)
for (int j = 1; j <= 3; j++)
for (int k = 1; k <= 3; k++)
t.g[i][j] = (t.g[i][j] + a.g[i][k] * b.g[k][j]) % mod;
return t;
}
void qmi(int k) {
init();
while(k) {
if(k & 1) E = mul(E, v);
v = mul(v, v);
k >>= 1;
}
}
void print(Matrix a) {
for (int i = 1; i <= 3; i++) {
for (int j = 1; j <= 3; j++)
cout << a.g[i][j] << " ";
cout << endl;
}
}
int main() {
cin >> T;
while(T--) {
cin >> n;
if(n <= 3) {
cout << 1 << endl;
continue;
}
qmi(n);
cout << E.g[2][1] << endl;
}
return 0;
}
例题5: 广义斐波那契
广义斐波那契
广义斐波那契数列是指形如: a n = p × a N − 1 + q × a n − 2 的数列 广义斐波那契数列是指形如:a_n = p × a_{N-1} + q × a_{n-2} \ 的数列 广义斐波那契数列是指形如:an=p×aN−1+q×an−2 的数列
- 推导如下:
目标矩阵: [ F n F n − 1 ] ⬇ F n = p × F n − 1 + q × F n − 2 ⬇ [ F n − 1 F N − 2 ] × ∣ p 1 q 0 ∣ ⬇ [ ( p × F n − 1 + q × F n − 2 ) ( F n − 1 × 1 + F N − 2 × 0 ) ] ⬇ [ F n F n − 1 ] 目标矩阵: [F_n \ F_{n-1}] \\ ⬇ \\ F_n=p×F_{n-1}+q×F_{n-2} \\ ⬇ \\ [F_{n-1} \ F_{N-2}] × \begin{vmatrix} p & 1 \\ q & 0 \\ \end{vmatrix} \\ ⬇ \\ [(p×F_{n-1}+q×F_{n-2}) \ (F_{n-1} × 1 + F_{N-2} ×0)] \\ ⬇ \\ [F_{n} \ F_{n-1}] 目标矩阵:[Fn Fn−1]⬇Fn=p×Fn−1+q×Fn−2⬇[Fn−1 FN−2]× pq10 ⬇[(p×Fn−1+q×Fn−2) (Fn−1×1+FN−2×0)]⬇[Fn Fn−1]
void init() {
memset(E.g, 0, sizeof(E.g));
E.g[1][1] = a2, E.g[1][2] = a1;
memset(v.g, 0, sizeof(v.g));
v.g[1][1] = p;
v.g[1][2] = 1;
v.g[2][1] = q;
}
例题6: 斐波那契公约数
斐波那契公约数
- 结论: g c d ( F n , F m ) = F ( g c d ( n , m ) ) gcd(F_n, \ F_m) = F(gcd(n, \ m)) gcd(Fn, Fm)=F(gcd(n, m))
#include <iostream>
#include <cstring>
#include <algorithm>
using namespace std;
typedef long long ll;
const int mod = 1e8;
const int N = 3;
struct Matrix {
ll g[N][N];
};
Matrix E, v; // E是单位阵(也是最终答案矩阵(快速幂)), v是推导出的初始的2×2 矩阵
void init() {
memset(E.g, 0, sizeof(E));
for (int i = 1; i <= 2; i++) E.g[i][i] = 1;
memset(v.g, 0, sizeof(v.g));
v.g[1][1] = v.g[1][2] = v.g[2][1] = 1;
}
Matrix mul(Matrix a, Matrix b) {
Matrix t;
memset(t.g, 0, sizeof(t));
for (int i = 1; i <= 2; i++)
for (int j = 1; j <= 2; j++)
for (int k = 1; k <= 2; k++)
t.g[i][j] = (t.g[i][j] + a.g[i][k] * b.g[k][j]) % mod;
return t;
}
void qmi(ll k) {
init();
while(k) {
if(k & 1) E = mul(E, v);
v = mul(v, v);
k >>= 1;
}
}
void print(Matrix a) {
for (int i = 1; i <= 2; i++) {
for (int j = 1; j <= 2; j++)
cout << a.g[i][j] << " ";
cout << endl;
}
}
int main() {
int n, m; cin >> n >> m;
n = __gcd(n, m);
if(n <= 2) {
cout << 1 << endl;
return 0;
}
qmi(n - 1);
// print(E);
cout << E.g[1][1] << endl;
return 0;
}
例题7: 这个勇者明明超强却过分慎重
这个勇者明明超强却过分慎重
#include <bits/stdc++.h>
#pragma GCC optimize(2)
#pragma GCC optimize(3)
using namespace std;
#define js ios::sync_with_stdio(false);cin.tie(0); cout.tie(0)
typedef long long ll;
inline ll read() { ll s = 0, w = 1; char ch = getchar(); while (ch < 48 || ch > 57) { if (ch == '-') w = -1; ch = getchar(); } while (ch >= 48 && ch <= 57) s = (s << 1) + (s << 3) + (ch ^ 48), ch = getchar(); return s * w; }
const int mod = 666666;
const int N = 2; //矩阵规模
const int M = 2;
struct Node {
ll matrix[N][M];//结构体中的矩阵
Node() {//默认构造函数
memset(matrix, 0, sizeof(matrix));
}
void one() {//单位矩阵
for (int i = 0; i < N; ++i)
matrix[i][i] = 1;
}
Node operator*(Node other) {//矩阵运算重载*运算符
Node ans;//记录乘积
for (int i = 0; i < N; i++)
for (int j = 0; j < M; j++)
for (int k = 0; k < N; k++) {
ans.matrix[i][j] += matrix[i][k] * other.matrix[k][j];
ans.matrix[i][j] %= mod;
}
return ans;
}
};
Node power(Node a, ll b) {//快速幂求a的b次方
Node res;
res.one();//单位矩阵
while (b) {
if (b & 1)
res = res * a;
a = a * a;
b >>= 1;
}
return res;
}
int main() {
Node st, p;
// [f[1], f[2]] = [4, 233];
st.matrix[0][0] = 4;
st.matrix[0][1] = 233;
// f(n) = 4 * f(n - 1) + 3 * f(n - 2)
// 0 3
// 1 4
p.matrix[0][1] = 3;
p.matrix[1][0] = 1;
p.matrix[1][1] = 4;
int n, m;
while(~scanf("%d %d", &n, &m))
printf("%d\n", m - (st * power(p, n - 1)).matrix[0][0]);
return 0;
}
ps: 矩阵拓展