POJ 3040 浅谈矩阵快速幂优化LogN斐波拉契函数求解

世界真的很大

讲道理矩阵这个东西还完全不怎么会怎么办？

赶紧复习一波再找一道水题刷一刷

矩阵优化线性函数求解曾经也听说过，但今天不知怎么奇困无比然后公式推起来巨慢

但其实还是比较简单

看题先：

description：

求斐波拉契数列第n项

input：

The input test file will contain multiple test cases. Each test case consists of a single line containing n (where 0 ≤ n ≤ 1,000,000,000). The end-of-file is denoted by a single line containing the number −1.

output：

For each test case, print the last four digits of Fn. If the last four digits of Fn are all zeros, print ‘0’; otherwise, omit any leading zeros (i.e., print Fn mod 10000).

讲道理斐波拉契函数本来就可以O（n）求出来，看上去已经很优了

但是这道题的数据范围奇大无比，1e9

基于斐波拉契是一个线性函数，而且每一项的值只与前两项有关

而且是和的关系

只要碰上这样的函数，可以考虑由矩阵快速幂优化求解

由于答案与前两项有关，那么我们很自然想到搞一个只有两项的矩阵。一个1*2的矩阵a，a（1,1）是F（i），a（2,1）是F（i-1 )

这样就很轻松得到转移矩阵e

（ 1， 1 ）

（ 1， 0 ）

求F（n）就等于F（1）再乘以e^(n-1)，这一步用快速幂搞定

切记注意矩阵乘法的运算顺序，这是由于其没有交换律的缘故

完整代码：

#include<stdio.h>
#include<cstring>
using namespace std;
typedef long long dnt;

const int mod=10000;

struct Matrix
{
int x,y;
dnt a[5][5];
void reset()
{
memset(a,0,sizeof(a));
}
}e,F;

Matrix operator*(const Matrix &aa,const Matrix &bb)
{
Matrix rt;
rt.reset();
rt.x=aa.x,rt.y=bb.y;
for(int i=1;i<=rt.x;i++)
for(int j=1;j<=rt.y;j++)
for(int k=1;k<=aa.y;k++)
rt.a[i][j]=(rt.a[i][j]+aa.a[i][k]*bb.a[k][j])%mod;
return rt;
}

Matrix quickmub(Matrix a,int b)
{
Matrix rt;
rt.x=rt.y=a.x;
rt.reset();
for(int i=1;i<=rt.x;i++) rt.a[i][i]=1;
while(b)
{
if(b&1) rt=rt*a;
a=a*a;
b>>=1;
}
return rt;
}

int main()
{
F.x=2,F.y=1,F.reset();
F.a[1][1]=1,F.a[2][1]=0;
e.x=2,e.y=2,e.reset();
e.a[1][1]=1,e.a[1][2]=1,e.a[2][1]=1,e.a[2][2]=0;
while(1)
{
int n;
scanf("%d",&n);
if(n==-1) break ;
Matrix A;
A.reset();
if(n>0) A=quickmub(e,n-1)*F;
printf("%d\n",A.a[1][1]);
}
return 0;
}
/*
Whoso pulleth out this sword from this stone and anvil is duly born King of all England
*/

嗯，就是这样