您的位置:首页 > 其它

矩阵快速幂 && 求解常系数齐次线性递推式

2019-07-18 10:34 141 查看
版权声明:本文为博主原创文章,遵循 CC 4.0 by 版权协议,转载请附上原文出处链接和本声明。 本文链接:https://blog.csdn.net/qq_43803508/article/details/95585003

先说矩阵A的n次幂(An)的快速幂的求法
因为矩阵乘满足结合律,所以矩阵的快速幂和前面普通快速幂原理都一样,主要在于矩阵乘法的运算

1、矩阵乘法

1.1 直接模拟矩阵乘(m行)

用结构体存矩阵(一个结构体可以直接赋值另一个结构体)

struct mat
{
LL a[m][m];
};
mat init=mat{0,0,0,0};
mat mul(mat a,mat b)
{
mat ans=init;
for(int i=0; i<m; i++)
for(int j=0; j<m; j++)
for(int k=0; k<m; k++)
ans.a[i][j]=(x.a[i][k]*y.a[k][j]%mod+ans.a[i][j])%mod;
return ans;
}

所以A,B两个矩阵相乘就可以调用函数mul===>mul(A,B)

1.2 重载*号,

C++的运算符重载

mat operator *(const mat &x,const mat &y)//重载运算符
{
mat ans=init;
for(int i=0; i<m; i++)
for(int j=0; j<m; j++)
for(int k=0; k<m; k++)
ans.a[i][j]=(x.a[i][k]*y.a[k][j]%mod+ans.a[i][j])%mod;
return ans;
}

现在A,B两个矩阵就可以直接相乘了===>A*B

2、求解常系数齐次线性递推式


以 Fibonacci 数列为例:
数列中 F0 = 0,F1 = 1, Fn = Fn−1 + Fn−2 ,它有两个单项式,我们就可以构造一个2行2列的矩阵,利用矩阵的乘法求Fn,我们先写出它的下一项
Fn = Fn−1 + Fn−2
Fn+1 = Fn + Fn−1

然后现在来填构造矩阵,因为Fn = 1* Fn−1 +1* Fn−2 ,所以第一列都是1

因为Fn-1 = 1* Fn−1 +0* Fn−2 ,所以第二列分别是1,0

然后我们用矩阵算F3
F3 = F2 + F1
F4 = F3 + F2

∣F 2 F 1 00∣\begin{vmatrix} F~2~ &amp; F~1~ \\ 0&amp;0 \\ \end{vmatrix}∣∣∣∣​F 2 0​F 1 0​∣∣∣∣​∣1110∣\begin{vmatrix} 1 &amp;1 \\ 1&amp;0 \\ \end{vmatrix}∣∣∣∣​11​10​∣∣∣∣​=∣F 3 F 2 00∣\begin{vmatrix} F~3~ &amp; F~2~ \\ 0&amp;0 \\ \end{vmatrix}∣∣∣∣​F 3 0​F 2 0​∣∣∣∣​
即:
∣1100∣\begin{vmatrix} 1 &amp;1 \\ 0&amp;0 \\ \end{vmatrix}∣∣∣∣< 8000 /span>​10​10​∣∣∣∣​∣1110∣\begin{vmatrix} 1 &amp;1 \\ 1&amp;0 \\ \end{vmatrix}∣∣∣∣​11​10​∣∣∣∣​=∣F 3 F 2 00∣\begin{vmatrix} F~3~ &amp; F~2~ \\ 0&amp;0 \\ \end{vmatrix}∣∣∣∣​F 3 0​F 2 0​∣∣∣∣​

所以Fn就等于矩阵∣F 2 F 1 00∣\begin{vmatrix} F~2~ &amp; F~1~ \\ 0&amp;0 \\ \end{vmatrix}∣∣∣∣​F 2 0​F 1 0​∣∣∣∣​乘以矩阵∣1110∣\begin{vmatrix} 1 &amp;1 \\ 1&amp;0 \\ \end{vmatrix}∣∣∣∣​11​10​∣∣∣∣​的n-2次幂

代码:poj 3070

#include<cstring>
#include<cstdio>
#include<iostream>
using namespace std;
typedef long long LL;
const int mod = 1e4;
struct mat
{
LL a[2][2];
};
mat init=mat{0,0,0,0};
mat init1=mat{1,0,0,1};//单位矩阵
mat operator *(const mat &x,const mat &y)//重载运算符
{
mat ans=init;
for(int i=0; i<2; i++)
for(int j=0; j<2; j++)
for(int k=0; k<2; k++)
ans.a[i][j]=(x.a[i][k]*y.a[k][j]%mod+ans.a[i][j])%mod;
return ans;
}
mat quick_pow(mat a,LL b)
{
mat ans=a;
while(b)
{
if(b&1)
ans=ans*a;
a=a*a;
b>>=1;
}
return ans;
}
int main()
{
int n;
mat a,b;
a=mat{1,1,0,0};
b=mat{1,1,1,0};
while(scanf("%d",&n),n!=-1)
{
if(n==0)
{
printf("0\n");
continue;
}
if(n==1||n==2)
{
printf("1\n");
continue;
}
mat ans=a*quick_pow(b,n-3);
printf("%lld\n",ans.a[0][0]);
}
}
内容来自用户分享和网络整理,不保证内容的准确性,如有侵权内容,可联系管理员处理 点击这里给我发消息
标签: