1. 程式人生 > >斐波拉契

斐波拉契

names while 圖片 16px signed AR define event 有一個

在OI簡單數論中斐波拉契是常常出現的東西 是什麽
斐波那契數列,又稱黃金分割數列、因數學家列昂納多·斐波那契以兔子繁殖為例子而引入,故又稱為“兔子數列”,指的是這樣一個數列:1、1、2、3、5、8、13、21、34、……在數學上,斐波納契數列以如下被以遞歸的方法定義:F(0)=1,F(1)=1, F(n)=F(n-1)+F(n-2)(n>=2,n∈N*)在現代物理、準晶體結構、化學等領域,斐波納契數列都有直接的應用. (百度百科)
模板
1,遞歸
技術分享圖片
1 inline int Fib(int n)
2 { 3 if (n == 1 || n == 2) 4 { 5 return 1; 6 } 7 return Fib(n - 1) + Fib(n - 2); 8 }
View Code


復雜度分析
由於對於每一個1都是最後一層遞歸返回上來的故會遞歸F(n)次 , 由於斐波拉契數列是隨著指數上升的 故復雜度約為O(2^n) 2,循環
技術分享圖片
1 inline int Fib(int n)
2 {
3     F[1] = F[2] = 1;
4     for(int i = 3 ; i <= n ; ++ i)
5 F[i] = F[i - 1] + F[i - 2]; 6 return F[n]; 7 }
View Code

復雜度分析 :O(n) 3,矩陣乘法優化
講數列放在矩陣乘法中會發現 $$\left[\begin{matrix}F[n - 1]\\F[n - 2]\end{matrix}\right] * \left[\begin{matrix}1 & 1\\1 & 0\end{matrix}\right] = \left[\begin{matrix}F[n]\\F[n - 1]\end{matrix}\right]$$ 即如果求F[n]
$$\left[\begin{matrix}F[n]\\F[n - 1]\end{matrix}\right] = \left[\begin{matrix}F[n - 1]\\F[n - 2]\end{matrix}\right] * \left[\begin{matrix}1 & 1\\1 & 0\end{matrix}\right] = \left[\begin{matrix}F[n - 2]\\F[n - 3]\end{matrix}\right] * \left[\begin{matrix}1 & 1\\1 & 0\end{matrix}\right] ^ 2 = .... = \left[\begin{matrix}F[1] \\ F[2]\end{matrix}\right] * \left[\begin{matrix}1 & 1\\1 & 0\end{matrix}\right]^{n-1}$$ 將$$\left[\begin{matrix}F[1] \\ F[2]\end{matrix}\right]$$
拓展到2 * 2即
$$\left[\begin{matrix}1 & 0 \\ 0 & 1\end{matrix}\right]$$
技術分享圖片
 1 #include <iostream>
 2 #include <cstdio>
 3 #include <cstring>
 4 #define inc(i) (++ i)
 5 #define int long long
 6 using namespace std;
 7 const int Mod = 1000000000 + 7;
 8 struct JX
 9 {
10     int a[3][3];
11     friend JX operator * (JX a , JX b)
12     {
13         JX c;
14         memset(c.a , 0 , sizeof(c.a));
15         for(int i = 1 ; i <= 2 ; inc(i))
16             for(int j = 1 ; j <= 2 ; inc(j))
17                 for(int k = 1 ; k <= 2 ; inc(k))
18                     c.a[i][j] = (c.a[i][j] + a.a[i][k] * b.a[k][j]) % Mod;
19         return c;
20     }
21 }Ans , A;
22 int n;
23 inline void KSM()
24 {
25     int p = n - 2;
26     while(p)
27     {
28         if(p & 1) Ans = Ans * A;
29         p >>= 1 , A = A * A;
30     }
31 }
32 signed main()
33 {
34     scanf("%ld" , &n);
35     if(n <= 2)  putchar(49);
36     else
37     {
38         Ans.a[1][1] = A.a[1][1] = 1;//直接跳到n = 2時的矩陣
39         Ans.a[1][2] = A.a[1][2] = 1;
40         Ans.a[2][1] = A.a[2][1] = 1;
41         KSM();
42         printf("%lld" , Ans.a[1][1]);
43     }
44     return 0;
45 }
View Code

4 , 算不上什麽優化的優化 2(斐波拉契數列膜p的循環節)
轉自[Fib數模n的循環節](https://blog.csdn.net/acdreamers/article/details/10983813)(代碼原創) 我們知道Fibonacci數列,現在我們來求一個Fib數模n的循環節的長度. 對於一個正整數n,我們求Fib數模n的循環節的長度的方法如下: (1)把n素因子分解,即
(2)分別計算Fib數模每個的循環節長度,假設長度分別是
(3)那麽Fib模n的循環節長度 從上面三個步驟看來,貌似最困難的是第二步,那麽我們如何求Fib模的循環節長度呢?

這裏有一個優美的定理:
>Fib數模的最小循環節長度等於,其中表示Fib數模素數的最小循環節長度。可以看出我們現在最重要的就是求


對於求我們利用如下定理:

如果5是模的二次剩余,那麽循環節的的長度是的因子,否則,循環節的長度是的因子。 順便說一句,對於小於等於5的素數,我們直接特殊判斷,loop(2)=3,loop(3)=8,loop(5)=20。 那麽我們可以先求出所有的因子,然後用矩陣快速冪來一個一個判斷,這樣時間復雜度不會很大。
模板代碼:(巨醜 見諒)
技術分享圖片
1#include<iostream>2#include<cstdio>3#include<cctype>4#include<cstring>5#include<cmath>6#include<algorithm>7#defineintlonglong8#defineinc(i)(++i)9#definedec(i)(--i)10usingnamespacestd;11constintN=100000+7;12intn,Prime[N+7],tot,Pri[N],Cnt[N],cnt,TOT;13intFactor[N],P;14boolNo_Prime[N+7];15structJX16{17inta[4][4];18friendJXoperator*(JXa,JXb)19{20JXc;21memset(c.a,0,sizeof(c.a));22for(inti=1;i<=2;inc(i))23for(intj=1;j<=2;inc(j))24for(intk=1;k<=2;inc(k))25c.a[i][j]=(c.a[i][j]+a.a[i][k]*b.a[k][j]%P)%P;26returnc;27}28}A,X;29inlineJXKSM_2(JXa,intq,intP)30{//矩陣優化斐波拉契數列31JXAns=X,x=a;32while(q)33{34if(q&1)Ans=Ans*x;35q>>=1,x=x*x;36}37returnAns;38}39inlinevoidGet_Prime()40{//得出質數便於之後分解質因數41No_Prime[1]=1;42for(inti=2;i<=N;inc(i))43{44if(!No_Prime[i])Prime[inc(tot)]=i;45for(intj=1;j<=tot&&i*Prime[j]<=N;inc(j))46{47No_Prime[i*Prime[j]]=1;48if(i%Prime[j]==0)break;49}50}51}52inlinevoidResolve(intn,intPri[],intCnt[])53{//分解質因數,pri存質因子,cnt存次數54cnt=1;55intt=sqrt(n),TOT;56for(inti=1;Prime[i]<=t;inc(i))57if(n%Prime[i]==0)58{59TOT=0;60Pri[cnt]=Prime[i];61while(n%Prime[i]==0)62{63n/=Prime[i];64inc(TOT);65}66Cnt[cnt]=TOT;67inc(cnt);68}69if(n^1)70{71Pri[cnt]=n,Cnt[cnt]=1;72inc(cnt);73}74dec(cnt);75}7677inlineintKSM_1(intx,intq,intP)78{//普通的快速冪79intAns=1;80x%=P;81while(q)82{83if(q&1)Ans=Ans*x%P;84q>>=1,x=x*x%P;85}86returnAns;87}88inlinevoidWork(intn)89{//求n的因數90TOT=0;91intt=sqrt(n);92for(inti=1;i<=t;inc(i))93{94if(n%i==0)95{96if(i*i==n)Factor[inc(TOT)]=i;97else98{99Factor[inc(TOT)]=i;100Factor[inc(TOT)]=n/i;101}102}103}104}105inlineintGcd(inta,intb)106{107returna%b==0?b:Gcd(b,a%b);108}109inlineintFind_Loop(intn)110{111Resolve(n,Pri,Cnt);112intAns=1,loop;113for(inti=1;i<=cnt;inc(i))114{115loop=1;116P=Pri[i];117switch(P)118{119case2:loop=3;break;120case3:loop=8;break;121case5:loop=20;break;122default:123{124if(KSM_1(5,(P-1)>>1,P)==1)125Work(P-1);126else127Work(2*(P+1));128sort(Factor,Factor+TOT+1);129for(intj=1;j<=TOT;inc(j))130{131JXB=KSM_2(A,Factor[j]-1,P);132intx=(B.a[1][1]+B.a[1][2])%P;133inty=(B.a[2][1]+B.a[2][2])%P;134if(x==1&&y==0)135{136loop=Factor[j];137break;138}139}140}break;141}142loop*=KSM_1(P,Cnt[i]-1,9223372036854775807);143if(loop<Ans)swap(loop,Ans);144Ans=Ans/Gcd(loop,Ans)*loop;145}146returnAns;147}148inlineintMod(char*a,intb)//高精度a除以低精度b149{150intd=0,l=strlen(a);151for(inti=0;i<l;inc(i))d=(d*10+(a[i]-0))%b;//求出余數152returnd;153}154charpp[30000007];155signedmain()156{157A.a[1][1]=A.a[1][2]=A.a[2][1]=1;158X.a[1][1]=X.a[2][2]=1;159Get_Prime();160scanf("%s%lld",pp,&n);161intloop=Find_Loop(n);162intAns=Mod(pp,loop);163if(!Ans)Ans=loop;164if(Ans<=0)putchar(48);165else166if(Ans<=2)printf("%lld",1ll%n);167else168{169A.a[1][1]=A.a[1][2]=A.a[2][1]=1;170X.a[1][1]=X.a[2][2]=1;171P=n;172X=KSM_2(A,Ans-1,n);173printf("%lld",X.a[1][1]);174}175return0;176}
Ugly Code

斐波拉契