矩阵快速幂学习总结
前言
我们都知道,斐波那契数列的递推公式为 \(a_{i}\)=\(a_{i-1}+a_{i-2}\) , 我们要求其第n位值的时候,时间复杂度是O(n)的,但是当 n 达到\(10^{8}\)甚至以上级别的时候,就会TLE,那该怎么办呢?
矩阵快速幂
前置芝士 1.快速幂
洛谷 P1226 【模板】快速幂
所谓快速幂,其实就是乘法结合律和二进制的结合产物
我们都知道,一个数可以根据其二进制位分成: \(2^{p_0}\)+\(2^{p_1}\)+…+\(2^{p_k}\) , 如23=1+2+4+16
乘法结合律就是a * b * c = a * (b * c)
那么\(a^n\)=a * a * a * … * a = a * (a * a)* … *(a * a * … * a)
这样计算的速度就从O(n)降到了O($\log_{2}{n} $)
C++快速幂模板:
点击查看代码
int qpow(int a,int b,int mod)
{
	int res=1;
	while(b)
	{
		if(b&1)res=res*a%mod;
		a=a*a%mod;
		b>>=1;
	}
	return res;
}
前置芝士 2.矩阵乘法

例如:\(\begin{bmatrix} 1&2 \\ 3&4\end{bmatrix}\) $\times $ \(\begin{bmatrix} 1 \\ 1\end{bmatrix}\) = \(\begin{bmatrix} 3 \\ 7\end{bmatrix}\)
C++矩阵乘法模板:
点击查看代码
#define vvi vector<vector<int>>
vvi mul(vvi a,vvi b,int mod){
	int n=a.size(),m=a[0].size(),k=b[0].size();
	vvi c(n,vector<int>(k,0));
	for(int i=0;i<n;i++)	
		for(int j=0;j<k;j++)
			for(int u=0;u<m;u++)
				c[i][j]=(c[i][j]+a[i][u]*b[u][j])%mod;
	return c;
}
实现矩阵快速幂
将快速幂模板改一改就可以了,注意需要先构造单位矩阵:\(\begin{bmatrix}  1&0&0 \\  0&1&0 \\ 0&0&1 \end{bmatrix}\)
C++矩阵快速幂模板:
点击查看代码
#define vvi vector<vector<int>>
vvi mul(vvi a,vvi b,int mod){
	int n=a.size(),m=a[0].size(),k=b[0].size();
	vvi c(n,vector<int>(k,0));
	for(int i=0;i<n;i++)	
		for(int j=0;j<k;j++)
			for(int u=0;u<m;u++)
				c[i][j]=(c[i][j]+a[i][u]*b[u][j])%mod;
	return c;
}
vvi qpow(vvi a,int b,int mod){
	int n=(int)a.size();
	vector<vector<int>>c(n,vector<int>(n,0));
	for(int i=0;i<n;i++)c[i][i]=1;//构造单位矩阵
	while(b){
		if(b&1)c=mul(c,a,mod);
		a=mul(a,a,mod);
		b>>=1;
	}
	return c;
}
矩阵快速幂专项练习题单
P1962 斐波那契数列
我们构造这样的矩阵:\(\begin{bmatrix}  1&1 \\  1&0  \end{bmatrix}\) $\times $ \(\begin{bmatrix}  a_n \\  a_{n-1}  \end{bmatrix}\) = \(\begin{bmatrix}  a_{n+1} \\  a_{n}  \end{bmatrix}\) 。
怎么想到这个构造的呢?我觉得可以根据递推公式,观察第n项是由前面哪些项得来的,根据矩阵运算法则反推即可 。
利用矩阵快速幂的代码计算出其迭代n-1次后的新矩阵,这样第1行第1列的元素就是\(a_n\)了(第1行第2列的元素就是\(a_{n-1}\))
点击查看代码
#include <bits/stdc++.h>
#define vvi vector<vector<int>>
#define int long long
#define endl "\n"
using namespace std;
const int N=1e6+10,M=2010,mod=1e9+7,INF=0x3f3f3f3f;
vvi mul(vvi a,vvi b,int mod){
	int n=a.size(),m=a[0].size(),k=b[0].size();
	vvi c(n,vector<int>(k,0));
	for(int i=0;i<n;i++)	
		for(int j=0;j<k;j++)
			for(int u=0;u<m;u++)
				c[i][j]=(c[i][j]+a[i][u]*b[u][j])%mod;
	return c;
}
vvi qpow(vvi a,int b,int mod){
	int n=(int)a.size();
	vector<vector<int>>c(n,vector<int>(n,0));
	for(int i=0;i<n;i++)c[i][i]=1;//构造单位矩阵
	while(b){
		if(b&1)c=mul(c,a,mod);
		a=mul(a,a,mod);
		b>>=1;
	}
	return c;
}
void solve(){
	int n;
	cin>>n;
	vector<vector<int>>a={{1,1},{1,0}};
	vector<vector<int>>b=qpow(a,n-1,mod);
	cout<<b[0][0]<<endl;
}
signed main(){
    ios::sync_with_stdio(0);cin.tie(0),cout.tie(0);
    int _=1;
    //cin>>_;
    while(_--) solve();
    return 0;
}
P1349 广义斐波那契数列
考虑构造:\(\begin{bmatrix} p&q \\ 1&0 \end{bmatrix}\) $\times $ \(\begin{bmatrix} a_n \\ a_{n-1} \end{bmatrix}\) = \(\begin{bmatrix} a_{n+1} \\ a_{n} \end{bmatrix}\)
点击查看代码
#include <bits/stdc++.h>
#define int long long
#define pb push_back
#define vvi vector<vector<int>>
#define lowbit(x) (x&(-x))
#define fi first
#define se second
#define endl "\n"
using namespace std;
const int N=5e5+10,M=2010,mod=1e9+7,INF=0x3f3f3f3f;
const int inf=0x3f3f3f3f3f3f3f3f;
vvi mul(vvi a,vvi b,int mod){//a的列数应与b的行数相等
	//a为n行m列的矩阵,b为m行d列的矩阵,相乘后获得n行d列的矩阵
	int n=a.size(),m=a[0].size(),d=b[0].size();
	vvi c(n,vector<int>(d,0));
	for(int i=0;i<n;i++)
		for(int j=0;j<d;j++)
			for(int k=0;k<m;k++)
				c[i][j]=(c[i][j]+a[i][k]*b[k][j]%mod)%mod;
	return c;
}
vvi qpow(vvi a,int n,int m){
	vvi c(2,vector<int>(2,0));
	c[0][0]=c[1][1]=1;
	while(n){
		if(n&1)c=mul(c,a,m);
		a=mul(a,a,m);
		n>>=1;
	}
	return c;
}
void solve(){
	int p,q,a1,a2,n,m;
	cin>>p>>q>>a1>>a2>>n>>m;
	if(n==1)cout<<a1%m<<endl;
	else if(n==2)cout<<a2%m<<endl;
	else {
		vvi a={{p%m,q%m},{1,0}};
		vvi b={{a2%m},{a1%m}};
		a=qpow(a,n-2,m);
		a=mul(a,b,m);
		cout<<a[0][0]<<endl;
	}
	
}
signed main(){
    ios::sync_with_stdio(0);cin.tie(0),cout.tie(0);
    int _=1;
    //cin>>_;
    while(_--) solve();
    return 0;
}
P1939 矩阵加速(数列)
考虑构造\(\begin{bmatrix} 1&0&1 \\ 1&0&0 \\ 0&1&0 \end{bmatrix}\) $\times $ \(\begin{bmatrix} a_{n} \\ a_{n-1} \\ a_{n-2} \end{bmatrix}\) = \(\begin{bmatrix} a_{n+1} \\ a_{n} \\ a_{n-1} \end{bmatrix}\)
点击查看代码
#include <bits/stdc++.h>
#define int long long
#define pb push_back
#define fi first
#define se second
#define endl "\n"
using namespace std;
const int N=1e6+10,M=2010,mod=1e9+7,INF=0x3f3f3f3f;
int n,m,k;
vector<vector<int>>mul(vector<vector<int>>a,vector<vector<int>>b){
	vector<vector<int>>c(3,vector<int>(3,0));
	for(int i=0;i<3;i++)
		for(int j=0;j<3;j++)
			for(int k=0;k<3;k++)
				c[i][j]=(c[i][j]+a[i][k]*b[k][j]%mod)%mod;
	return c;
}
vector<vector<int>>qpow(vector<vector<int>>a,int b){
	vector<vector<int>>c(3,vector<int>(3,0));
	for(int i=0;i<3;i++)c[i][i]=1;
	while(b){
		if(b&1)c=mul(c,a);
		a=mul(a,a);
		b>>=1;
	}
	return c;
}
void solve(){
	cin>>n;
	vector<vector<int>>a={{1,0,1},{1,0,0},{0,1,0}};
	vector<vector<int>>b=qpow(a,n);
	cout<<b[1][0]<<endl;
}
signed main(){
    ios::sync_with_stdio(0);cin.tie(0),cout.tie(0);
    int _=1;
    cin>>_;
    while(_--) solve();
    return 0;
}
U583138 多重前缀和
考虑构造下三角矩阵\(\begin{bmatrix} 1&0&0&…&0 \\ 1&1&0&…&0 \\ … \\ 1&1&1&…&1 \end{bmatrix}\) $\times $ \(\begin{bmatrix} a_{1} \\ a_{2} \\ … \\ a_{n} \end{bmatrix}\) = \(\begin{bmatrix} a_{2} \\ a_{3} \\ … \\ a_{n+1} \end{bmatrix}\)
点击查看代码
#include <bits/stdc++.h>
#define int long long
#define pb push_back
#define vvi vector<vector<int>>
#define lowbit(x) (x&(-x))
#define fi first
#define se second
#define endl "\n"
using namespace std;
const int N=5e5+10,M=2010,mod=1e9+7,INF=0x3f3f3f3f;
const int inf=0x3f3f3f3f3f3f3f3f;
vvi mul(vvi a,vvi b){
	int n=a.size(),m=a[0].size(),k=b[0].size();
	vvi c(n,vector<int>(k,0));
	for(int i=0;i<n;i++)
		for(int j=0;j<k;j++)
			for(int u=0;u<m;u++)
				c[i][j]=(c[i][j]+a[i][u]*b[u][j]%mod)%mod;
	return c;
}
vvi qpow(vvi a,int m){
	int n=a.size();
	vvi c(n,vector<int>(n,0));
	for(int i=0;i<n;i++)
			c[i][i]=1;
			
	while(m){
		if(m&1)c=mul(c,a);
		a=mul(a,a);
		m>>=1;
	}
	return c;
}
void solve(){
	int n,m;
	cin>>n>>m;
	vvi a(n,vector<int>(1,0));
	for(int i=0;i<n;i++)
		cin>>a[i][0];
	
	vvi b(n,vector<int>(n,0));
	for(int i=0;i<n;i++)
		for(int j=0;j<=i;j++)
			b[i][j]=1;
	b=qpow(b,m);
	
	vvi c=mul(b,a);
	for(int i=0;i<n;i++)
		cout<<c[i][0]<<" ";
}
signed main(){
    ios::sync_with_stdio(0);cin.tie(0),cout.tie(0);
    int _=1;
    //cin>>_;
    while(_--) solve();
    return 0;
}
P2044 [NOI2012] 随机数生成器
考虑构造矩阵:\(\begin{bmatrix} a&c \\ 0&1 \end{bmatrix}\) $\times $ \(\begin{bmatrix} x_0 \\ 1 \end{bmatrix}\) = \(\begin{bmatrix} x_1 \\ 1 \end{bmatrix}\)
点击查看代码
#include <bits/stdc++.h>
#define int __int128
#define LL long long
#define pb push_back
#define fi first
#define se second
#define endl "\n"
using namespace std;
const int N=1e6+10,M=2010,mod=998244353,INF=0x3f3f3f3f;
#define vvi vector<vector<int>>
vvi mul(vvi a,vvi b,int mod){
	int n=a.size(),m=a[0].size(),k=b[0].size();
	vvi c(n,vector<int>(k,0));
	for(int i=0;i<n;i++)	
		for(int j=0;j<k;j++)
			for(int u=0;u<m;u++)
				c[i][j]=(c[i][j]+a[i][u]*b[u][j])%mod;
	return c;
}
vvi qpow(vvi a,int b,int mod){
	int n=(int)a.size();
	vector<vector<int>>c(n,vector<int>(n,0));
	for(int i=0;i<n;i++)c[i][i]=1;//构造单位矩阵
	while(b){
		if(b&1)c=mul(c,a,mod);
		a=mul(a,a,mod);
		b>>=1;
	}
	return c;
}
void solve(){
	LL m,a,c,x0,n,g;
	cin>>m>>a>>c>>x0>>n>>g;
	vvi p={{x0%m},{1}};
	vvi pp={{a,c%m},{0,1}};
	pp=qpow(pp,n,m);
	pp=mul(pp,p,m);
	cout<<(LL)pp[0][0]%g<<endl;
}
signed main(){
    ios::sync_with_stdio(0);cin.tie(0),cout.tie(0);
    int _=1;
    //cin>>_;
    while(_--) solve();
    return 0;
}
 
 
        
 
                
            
         浙公网安备 33010602011771号
浙公网安备 33010602011771号