浅谈矩阵快速幂,希望能够帮助到你!
现在给出一个实际的问题,要求2的十次方,很简单,直接for循环就行了。
可是现在如果说让求2的10000000000次方,又该怎么办呢?
这个时候我们可以考虑使用快速幂。
#include<algorithm>
#include<iostream>
#include<string.h>
#include<cstdio>
#define P pair<long long,long long>
#define ll long long
using namespace std;
ll x,y;
ll pow(ll x,ll y)
{
ll ans=1;
while(y)
{
if(y%2==1)
{
ans*=x;
}
x*=x;
y>>=1;
}
return ans;
}
int main()
{
cin>>x>>y;
printf("%lld\n",pow(x,y));
return 0;
}
矩阵快速幂分为两个部分:
1:奇数
2:偶数
从100直接到50,就相当于二分,所以能够大大降低时间复杂度。
现在把问题从快速幂转移到矩阵快速幂。
第一个问题是求斐波那契数:
在Fibonacci整数序列中,对于n≥2,F0 = 0,F1 = 1,Fn = Fn - 1 + Fn - 2.例如,Fibonacci序列的前十项是: 0,1,1,2,3,5,8,13,21,34 ...... 给定整数n,您的目标是计算Fn的最后4位数。 输入 输入测试文件将包含多个测试用例。每个测试用例由一行包含n(其中0≤n≤1,000,000,000)组成。文件结尾由包含数字-1的单行表示。 输出 对于每个测试用例,打印Fn的最后四位数字。如果Fnare的最后四位数全为零,则打印'0';否则,省略任何前导零(即,打印Fn mod 10000)。 样本输入 0 9 999999999 -1 样本输出 0 34 626
要求的是斐波那契数,我们知道对于斐波那契数的第n位(n>=3)满足f(n)=f(n-1)+f(n-2)。
可以推出A是:
也就是说,现在我们要求F(N),只需要对A进行N-1次方就行了。
具体的实现如下:
#include<algorithm>
#include<iostream>
#include<string.h>
#include<cstdio>
#include<cstring>
#define P pair<long long,long long>
#define ms(a,b) memset(a,b,sizeof(a))
#define ll long long
using namespace std;
const int mod=1e4;
//ll x,y;
//ll pow(ll x,ll y)
//{
// ll ans=1;
// while(y)
// {
// if(y%2==1)
// {
// ans*=x;
// }
// x*=x;
// y>>=1;
// }
// return ans;
//}
struct xiao
{
ll f[4][4];
} ans,mapp;
void init(xiao &ans)
{
for(int i=1; i<=2; ++i)
{
for(int j=1; j<=2; ++j)
{
if(i==j)
ans.f[i][j]=1;
else
ans.f[i][j]=0;
}
}
}
void power(xiao &x,xiao &y,xiao &z)
{
//z初始化,这里是存储相乘后矩阵信息的矩阵
ms(z.f,0);
for(int i=1; i<=2; ++i)
{
for(int j=1; j<=2; ++j)
{
if(x.f[i][j])
{
for(int k=1; k<=2; ++k)
{
z.f[i][k]+=x.f[i][j]*y.f[j][k];
z.f[i][k]%=mod;
}
}
}
}
/*
注意,这里把所有的信息都存储到z里面了,现在的情况只是说我有100毫升水在C杯子里面,
我往A中倒了50ml,现在还有50ml在c杯子里面,接下来就是把剩下的50ml水倒到b杯子里面。
实现这一步就是简单的赋值。具体实现在下面.
*/
}
ll judex(ll x)
{
init(ans);
/*要求第n位,对n-1位进行快速幂运算。
类比前面提到的快速幂的知识,我们首先定义一个单位矩阵,相当于前面的那个ans变量。
*/
xiao temp=mapp,t;
x--;
while(x)
{
if(x%2==1)
{
//类比前面的 ans*=2;
power(ans,temp,t);
ans=t;
}
//这里类比把2的100次方变成4的50次方。
power(temp,temp,t);
temp=t;
x/=2;
}
return ans.f[1][2]+ans.f[1][1];
}
int main()
{
mapp.f[1][1]=0;
mapp.f[1][2]=1;
mapp.f[2][1]=1;
mapp.f[2][2]=1;
ll n;//要求的是哪一位的
while(cin>>n)
{
if(n==-1)
break;
if(n==1)
{
printf("1\n");
continue;
}
else if(n==2 )
{
printf("1\n");
continue;
}
else if(n==0)
{
printf("0\n");
continue;
};
printf("%lld\n",judex(n)%mod);
}
return 0;
}
如果这道题解决了的话,接下来再看一道稍微复杂一点的题:
Recursive sequence
农夫约翰喜欢和他的N奶牛一起玩数学游戏。最近,它们被递归序列所吸引。在每个回合中,奶牛会排成一排,而约翰在黑板上写下两个正数a和b。然后,奶牛会逐一说出他们的身份证号码。第一头牛说第一个数字a,第二头说第二个数字b。在那之后,第i个牛说出第(i-2)个数,第(i-1)个数和i4的两倍之和。现在,您需要编写一个程序来计算第N头奶牛的数量,以检查John的奶牛是否可以做到正确。
输入
第一行输入包含一个整数t,即测试用例的数量。 t测试案例如下。 每个情况仅包含一行具有三个数字N,a和b,其中N,a,b <231,如上所述。
输出
对于每个测试用例,输出第N个牛的数量。这个数字可能非常大,所以你需要以模数2147493647输出它。
样本输入
2
3 1 2
4 1 10
样例输出
85
369
接下来,我们一起来一步步分析一下这道题:
首先,我们知道的是F(n)=2*f(n-2)+f(n-1)+n^4.
那么F(n+1)=2*f(n-1)+f(n)+(n+1)^4.
现在的思想还是类比的思想,我们在前面的斐波那契数的计算的时候知道,使用矩阵实现求得f(n),那么我们也可以使用矩阵实现从n^4到(n+1)^4的实现!
接下来请看下面的图:
也就是说,既然可以拆分成n^4,n^3,n^2,n^1,1这些项,那这些不就可以作为矩阵的内容吗?
对,那就可以实现了!接下来就是求A了。
因为矩阵里面有7个元素,所以我们构造的矩阵是7*7的矩阵。
举个例子,对于(n+1)^4来说,拆开就是
n^4+4*n^3+6*n^2+4*n+1.
按照这样的方法可以求出A
A=
构造算是矩阵快速幂的难点了,构造出来以后就好做了,接下来的一个难点是输出。
我们知道,输入的三个数分别为k,a,b,现在k如果是1的话,输出的是a,如果k是2的话,输出的是b,如果k是三的话,输出的就是2*a+b+1,这样就开始使用矩阵乘了。
return ((ans.f[1][0])*a+(ans.f[1][1])*b+(ans.f[1][2])*16+(ans.f[1][3])*8+(ans.f[1][4])*4+(ans.f[1][5])*2+(ans.f[1][6]));
附上完整的代码:
#include<algorithm>
#include<iostream>
#include<string.h>
#include<cstdio>
#include<cstring>
#define P pair<long long,long long>
#define ms(a,b) memset(a,b,sizeof(a))
#define ll long long
using namespace std;
const ll mod=2147493647;
ll N,a,b,k;
struct xiao
{
ll f[7][7];
} ans,mapp;
void init(xiao &ans)
{
for(int i=0; i<7; ++i)
{
for(int j=0; j<7; ++j)
{
if(i==j)
ans.f[i][j]=1;
else
ans.f[i][j]=0;
}
}
}
void power(xiao &x,xiao &y,xiao &z)
{
ms(z.f,0);
for(int i=0; i<7; ++i)
{
for(int j=0; j<7; ++j)
{
if(x.f[i][j])
{
for(int k=0; k<7; ++k)
{
z.f[i][k]+=(x.f[i][j]*y.f[j][k]);
z.f[i][k]%=mod;
}
}
}
}
}
ll solve(ll a,ll b,ll k)
{
init(ans);
ll number=k-2;
xiao temp=mapp,t;
while(number)
{
if(number%2==1)
{
power(ans,temp,t);
ans=t;
}
power(temp,temp,t);
temp=t;
number/=2;
}
return ((ans.f[1][0])*a+(ans.f[1][1])*b+(ans.f[1][2])*16+(ans.f[1][3])*8+(ans.f[1][4])*4+(ans.f[1][5])*2+(ans.f[1][6]));
}
int main()
{
mapp=
{
0,1,0,0,0,0,0,
2,1,1,4,6,4,1,
0,0,1,4,6,4,1,
0,0,0,1,3,3,1,
0,0,0,0,1,2,1,
0,0,0,0,0,1,1,
0,0,0,0,0,0,1
};
for(cin>>N; N--;)
{
cin>>k>>a>>b;
if(k==1)
{
printf("%lld\n",a);
}
else if(k==2)
{
printf("%lld\n",b);
}
else
printf("%lld\n",solve(a,b,k)%mod);
}
return 0;
}
最后针对矩阵的乘法进行补充:
最后得到的矩阵的形状是:
那么我们按照遍历一遍的方法去做的话,我们不能说每一个坐标点我们只计算一次,但我们可以肯定的是每一行我们只走一次,那么我们在走过这一行的时候就要把所有的数值都给计算出来。
对于上面的思想,我们用到的方法是矩阵乘法里面的把每一行都给提取出来进行计算(前面矩阵中的行)。