矩阵快速幂1

最新推荐文章于 2025-02-20 15:41:56 发布

原创最新推荐文章于 2025-02-20 15:41:56 发布 · 330 阅读

0 ·

CC 4.0 BY-SA版权

文章标签：

#线性代数 #算法 #矩阵

数学同时被 2 个专栏收录

11 篇文章

订阅专栏

矩阵快速幂

1 篇文章

订阅专栏

文章目录

关于快速幂

可参考以下博客：
快速幂算法

关于构造矩阵

可参考以下两篇博客：
矩阵构造方法
 根据递推公式构造系数矩阵用于快速幂

以斐波那契数列为例 $f[n]=f[n-1]+f[n-2]\quad (f[1]=f[2]=1)$
我们可得
$f[n]=1\times f[n-1]+1\times f[n-2]$
$f[n-1]=1\times f[n-1]+0\times f[n-2]$
进行矩阵构造
$\left[ \begin{matrix} f[n] \\ f[n-1] \\ \end{matrix} \right]=\left[ \begin{matrix} 1 &1 \\ 1 & 0 \\ \end{matrix} \right]\times \left[ \begin{matrix} f[n-1] \\ f[n-2] \\ \end{matrix} \right]$
$\times$ 左面矩阵为 $f [n], f [n - 1]$ 等式右侧系数，右侧为等号右侧的项，我们可使用以上式子进行递推，设 $\times$ 左边矩阵为 $A$ ，则最后可得 $f[n]=A^{n-2}\left[ \begin{matrix} f[2] \\ f[1] \\ \end{matrix} \right]$
设 $A^{n-2}$ 为 $B$ ，则 $f[n]=B_{1\,1}+B_{1\,2}$
矩阵快速幂中，最核心的是计算出基础矩阵A。要满足K阶常系数线性递推关系。
以下为两道例题：

P3390 【模板】矩阵快速幂

题目链接

P3390 【模板】矩阵快速幂

代码

#include <iostream>
#include <cstdio>
#define ll long long
using namespace std;

const int maxn=105,p=1e9+7 ;
int n;
ll k ;
struct ymatrix
{
	ll a[maxn][maxn] ;
}A ;

inline ll read()
{
	ll x=0,f=1;char c=getchar();
	while(c>'9'||c<'0') {if(c=='-') f=-1; c=getchar();}
	while(c<='9'&&c>='0') {x=x*10+c-'0'; c=getchar();}
	return f*x; 
}

ymatrix mul ( ymatrix a , ymatrix b )
{
	ymatrix ans ;
	for ( int i = 1 ; i <= n ; ++i )
		for ( int j = 1 ; j <= n ; ++j)
		{
			ans.a[i][j] = 0 ;
			for ( int k = 1 ; k <= n ; ++k )
			{
				ans.a[i][j] += a.a[i][k] * b.a[k][j] ;
				ans.a[i][j] %= p ;
			}
		}
	return ans ;
}

ymatrix quickpow( ymatrix a , ll b )
{
	ymatrix ans , base = a ;
	for ( int i = 1 ; i <= n ; ++i )
		for ( int j = 1 ; j <= n ; ++j )
			ans.a[i][j] = 0 ;
	for ( int i = 1 ; i <= n ; ++i ) ans.a[i][i] = 1 ;
	while ( b )
	{
		if ( b & 1 ) ans = mul ( ans , base ) ;
		base = mul ( base , base ) ;
		b >>= 1 ; 
	}
	return ans ;
}

int main()
{
	scanf("%d",&n) ; k = read() ;
	for ( int i = 1 ; i <= n ; ++i )
		for ( int j = 1 ; j <= n ; ++j )
			A.a[i][j] = read() ;
	ymatrix ans = quickpow(A,k) ;
	for ( int i = 1 ; i <= n ; ++i )
	{
		for ( int j = 1 ; j <= n ; ++j )
			printf("%lld ",ans.a[i][j]) ;
		printf("\n") ;
	}
	return 0 ;
}

P1962 斐波那契数列

题目链接

P1962 斐波那契数列

代码

#include <iostream>
#include <cstdio>
#define ll long long
using namespace std;

const int maxn=5,p=1e9+7 ;
ll n;
struct ymatrix
{
	ll a[maxn][maxn] ;
}A ;

inline ll read()
{
	ll x=0,f=1;char c=getchar();
	while(c>'9'||c<'0') {if(c=='-') f=-1; c=getchar();}
	while(c<='9'&&c>='0') {x=x*10+c-'0'; c=getchar();}
	return f*x; 
}

ymatrix mul ( ymatrix a , ymatrix b )
{
	ymatrix ans ;
	for ( int i = 1 ; i <= 2 ; ++i )
		for ( int j = 1 ; j <= 2 ; ++j)
		{
			ans.a[i][j] = 0 ;
			for ( int k = 1 ; k <= 2 ; ++k )
			{
				ans.a[i][j] += a.a[i][k] * b.a[k][j] ;
				ans.a[i][j] %= p ;
			}
		}
	return ans ;
}

ymatrix quickpow( ymatrix a , ll b )
{
	ymatrix ans , base = a ;
	for ( int i = 1 ; i <= 2 ; ++i )
		for ( int j = 1 ; j <= 2 ; ++j )
			ans.a[i][j] = 0 ;
	for ( int i = 1 ; i <= 2 ; ++i ) ans.a[i][i] = 1 ;
	while ( b )
	{
		if ( b & 1 ) ans = mul ( ans , base ) ;
		base = mul ( base , base ) ;
		b >>= 1 ; 
	}
	return ans ;
}

int main()
{
	n = read() ;
	if( n == 1 || n == 2 ) printf("1\n") ;
	else
	{
		A.a[1][1] = A.a[1][2] = A.a[2][1] = 1 ;
		A.a[2][2] = 0 ;
		ymatrix ans = quickpow( A , n-2 ) ;
		printf("%lld\n",(ans.a[1][1] + ans.a[1][2]) % p) ;
	}
	return 0 ;
}