#斐波那契数列#
题目描述
斐波那契数列是由如下递推式定义的数列
F
0
=
0
F_0=0
F0=0
F
1
=
1
F_1=1
F1=1
F
n
+
1
=
F
n
+
1
+
F
n
F_{n+1}=F_{n+1}+F_n
Fn+1=Fn+1+Fn
求这个数列第n项的值对
1
0
4
10^4
104取余后的结果。
限制条件
⋅
\huge ·
⋅
0
⩽
n
⩽
1
0
16
0\leqslant n\leqslant10^{16}
0⩽n⩽1016
我也是刚刚才搞懂了矩阵乘法(如果你不知道什么是矩阵乘法的话,右转百度百科),于是来应用一下新知识,如有表述不到位的地方请见谅。
下面进入正文
首先,我们先介绍一下对于斐波那契数列如何求解。把斐波那契数列的递推式表示成矩阵就得到下面的式子
(
F
n
+
2
F
n
+
1
)
=
{ \begin{pmatrix} F_{n+2}\\ F_{n+1} \end{pmatrix} }=
(Fn+2Fn+1)=
(
1
1
1
0
)
{ \begin{pmatrix} 1 & 1\\ 1 & 0 \end{pmatrix} }
(1110)
(
F
n
+
1
F
n
)
{ \begin{pmatrix} F_{n+1}\\ F_n \end{pmatrix} }
(Fn+1Fn)
我们发现式子里有个固定的矩阵 ( 1 1 1 0 ) { \begin{pmatrix} 1&1\\ 1&0 \end{pmatrix} } (1110)
记这个矩阵为A,则有
(
F
n
+
1
F
n
)
=
{ \begin{pmatrix} F_{n+1}\\ F_n \end{pmatrix} }=
(Fn+1Fn)=
A
n
A^n
An
(
F
1
F
0
)
=
{ \begin{pmatrix} F_1\\ F_0\\ \end{pmatrix} }=
(F1F0)=
A
n
A^n
An
(
1
0
)
{ \begin{pmatrix} 1\\ 0 \end{pmatrix} }
(10)
因此只要求出 A n A^n An就可以求出 F n F_n Fn了。关于 A n A^n An的计算我们可以采用类似快速幂的算法,在 O ( l o g n ) O(logn) O(logn)时间里求出第n项的值。
用结构体定义矩阵 ( a 0 a 1 a 2 a 3 ) { \begin{pmatrix} a0&a1\\ a2&a3 \end{pmatrix} } (a0a2a1a3),并重载矩阵乘法。
struct Matrix{
int a[4];
Matrix operator *(const Matrix &tmp)const{
Matrix res;
res.a[0]=(a[0]*tmp.a[0]+a[1]*tmp.a[2])%P;
res.a[1]=(a[0]*tmp.a[1]+a[1]*tmp.a[3])%P;
res.a[2]=(a[2]*tmp.a[0]+a[3]*tmp.a[2])%P;
res.a[3]=(a[2]*tmp.a[1]+a[3]*tmp.a[3])%P;
return res;
}
Matrix(){memset(a,0,sizeof(a));}
}A[Max];
计算前要先初始化
void Init(){
A[0].a[0]=1;
A[0].a[1]=1;
A[0].a[2]=1;
for(int i=1;i<Max;i++)
A[i]=A[i-1]*A[i-1];
}
求答案时要注意初始的res为 ( 1 0 0 1 ) { \begin{pmatrix} 1&0\\ 0&1 \end{pmatrix} } (1001)
int main(){
Init();
Rd(n);
Matrix res;
res.a[0]=1;
res.a[3]=1;
for(int i=Max-1;i>=0;i--)
if((n>>i)&1)res=res*A[i];
printf("%d\n",res.a[2]);
return 0;
}
完整的代码如下:
#include<stdio.h>
#include<string.h>
#include<algorithm>
#include<iostream>
#define P 10000
#define M 1000005
#define Max 25
using namespace std;
void Rd(int &res){
char c;res=0;int k=1;
while(c=getchar(),c<48&&c!='-');
if(c=='-'){k=-1;c='0';}
do{
res=(res<<3)+(res<<1)+(c^48);
}while(c=getchar(),c>=48);
res*=k;
}
struct Matrix{
int a[4];
Matrix operator *(const Matrix &tmp)const{
Matrix res;
res.a[0]=(a[0]*tmp.a[0]+a[1]*tmp.a[2])%P;
res.a[1]=(a[0]*tmp.a[1]+a[1]*tmp.a[3])%P;
res.a[2]=(a[2]*tmp.a[0]+a[3]*tmp.a[2])%P;
res.a[3]=(a[2]*tmp.a[1]+a[3]*tmp.a[3])%P;
return res;
}
Matrix(){memset(a,0,sizeof(a));}
}A[Max];
void Init(){
A[0].a[0]=1;
A[0].a[1]=1;
A[0].a[2]=1;
for(int i=1;i<Max;i++)
A[i]=A[i-1]*A[i-1];
}
int n;
int main(){
Init();
Rd(n);
Matrix res;
res.a[0]=1;
res.a[3]=1;
for(int i=Max-1;i>=0;i--)
if((n>>i)&1)res=res*A[i];
printf("%d\n",res.a[2]);
return 0;
}
第一次写博客,谢谢各位支持。O(∩_∩)O