技術標籤：迴圈斐波那契數列

考慮一個數列，它的第一項和第二項都是1，從第三項開始，每一項都是前面兩項的和。

1，1，2，3，5，8，13，21....

這個數列，就是數學中著名的斐波那契數列了。它被義大利數學家列昂納多·斐波那契(Leonardoda Fibonacci)在1202年提出，也因此得名。今天的問題，就是求這個數列的第N項是多少。

(解決這個問題，本文會提供三種方法。後面兩種需要用到快速冪，可從本公眾號演算法趣談—快速冪中獲取詳解)

我們把這個數列的定義數學化一下，口語化的定義是：“它的第一項和第二項都是1，從第三項開始，每一項都是前面兩項的和。”。而將其數學化，我們令該數列的第i項為f[i]，則把口語化表達翻譯過來就是：f[1]=1,f[2]=1,f[i]=f[i-1]+f[i-2](i>2)

像上面這樣的公式，我們稱其為遞推式。也就是通過前面能推出來後面的式子。凡是這種遞推式，都可以用迴圈來求出數列的第n項，程式碼應該比較好理解。

#include #include using namespace std;int n,f[100010];int main(){   cin>>n;  f[1]=1;  f[2]=1;  for(int i=3;i<=n;i++)    f[i]=f[i-1]+f[i-2];//遞推式  cout<  return 0;}

這種演算法，程式需要進行n次迴圈才可以算出來。他模擬的過程就相當於你一步一步推，f[2],f[3],f[4].....看起來十分費勁。那麼，有沒有更快的演算法？數學好的同學可能已經想出來了，我們是否可以把f[i]=f[i-1]+f[i-2]的通項公式給求出來，這樣就可以用一次計算求出來答案了。

首先，我們用待定係數法，設(其中α和β是兩個引數)：

把αan-1挪過去，就可以得到：

根據原本的遞推式，我們可以得知：

仔細觀察這個式子，如果α和β都是一個一元二次方程的根，那麼這個方程組豈不是像韋達定理？那我們構造一個一元二次方程，令α、β為這個方程的兩個根，且它們的關係符合上面的方程組。這樣的方程比較多，但無論如何它們的解都不會變。我們就構造一個簡單點的：

這個方程的解是：

也就是說：(α β可顛倒順序)

我們再把α和β帶回原來的遞推式(注意因為α和β可以交換順序，所以有兩個遞推式):

這樣，我們就把原本的遞推式轉化為了兩個等比數列，這樣，就可以先分別求出這兩個等比數列的通項公式。

最後，讓：

再化簡，就可以求出通項式了！：

好不容易求出來通項式，把它應用於計算機吧！

#include #include #include using namespace std;double n; double ksm(double a,int b){//快速冪  double ans=1;  while(b){    if(b%2==1)  ans*=a;    a*=a;    b=b/2;  }  return ans;}int main(){   cin>>n;  printf("%lf",(ksm((1+sqrt(5))/2,n)-ksm((1-sqrt(5))/2,n))/sqrt(5));  return 0;}

這裡面快速冪的時間複雜度是O(logn)，所以整個程式的時間複雜度也是O(logn).

我們來試驗一下程式碼：

十分正常，再換個大一點的數：

？？？為什麼好好的斐波那契數列，蹦出來一個小數？這是什麼情況？

回到程式仔細觀察，我們發現根號5在計算機中儲存的形式是小數，也就是計算機存的根號5根本不是根號5，只是近似於它的一個小數罷了。在n較小，比如等於10時，因為精度還夠，所以不會影響答案。但是！一旦n變大，計算機儲存的小數的精度就不夠了。因此，答案就會有所影響了！

那麼，現在怎麼辦，用這個方法，雖然優化了時間，卻消耗了精度啊。

所以，只能考慮別的做法了。觀察一個矩陣M，它長這個樣子：

如果把M平方，就會變成這個樣子：

而M的立方，是這個樣子：

仔細觀察，發現這個矩陣M的n次方後的矩陣中，M(1,1) 和M(2,1)分別對應著斐波那契數列的第n+1項和第n項。也就是，只用計算這個矩陣的n-1次方，它的第一行第一例的數就是斐波那契數列的第n項！

為什麼呢？因為在矩陣運算時：

當前矩陣的左上角記錄的是前一個矩陣左上和右上之和，當前矩陣的右上角記錄的是上一個矩陣的左上角。這樣就相當於完成了斐波那契的遞推。

這樣，我們把問題從求斐波那契數列的第n項變成了求這個矩陣的n-1次方。我們知道矩陣滿足乘法交換律，因此，我們可以用矩陣快速冪來解決此問題。

首先，我們要先寫出來矩陣的乘法。也就是目標矩陣的(i,j)為第一個矩陣的第i行和第二個矩陣的第j列每兩個數乘積之和。

inlinenodematri_mult(nodeA,nodeB){//其中node是一個含矩陣的結構體    node C;    for(int i=1;i<=2;i++)        for(int j=1;j<=2;j++)            C.a[i][j]=0;for(inti=1;i<=2;i++)        for(int j=1;j<=2;j++)for(intk=1;k<=2;k++)C.a[i][j]=C.a[i][j]+A.a[i][k]*B.a[k][j];    return C;}

然後，寫出快速冪。在這個快速冪中，不同的是把原先的乘法變成了矩陣乘法函式。而且，原本記錄冪的sum初始值是1，現在是一個對角線為1的單位矩陣，因為任何矩陣乘以這個單位矩陣都不會變，就像乘以1一樣。

inline ksm(node A,int b){    node sum;    for(int i=1;i<=2;i++)//構造單位矩陣        for(int j=1;j<2;j++)            if(i==j)  sum.a[i][j]=1;            else sum.a[i][j]=0;    while(b){//矩陣快速冪if(b&1)sum=matri_mult(sum,A);//把乘號替換成了矩陣乘法函式        A=matri_mult(A,A);        b/=2;    }    return sum;}

最後輸出矩陣n-2次方第一列兩個數的和就ok了!

完整程式碼如下：

#include #include #include #include #include #define MAXN 1010using namespace std;int n,k;struct node{    int a[110][110];};inline node matri_mult(node A,node B){//矩陣乘法    node C;    for(int i=1;i<=2;i++)        for(int j=1;j<=2;j++)            C.a[i][j]=0;    for(int i=1;i<=2;i++)        for(int j=1;j<=2;j++)            for(int k=1;k<=2;k++)        C.a[i][j]+=A.a[i][k]*B.a[k][j];    return C;}inline node ksm(node A,int b){//矩陣快速冪    node sum;    for(int i=1;i<=2;i++)        for(int j=1;j<=2;j++)            if(i==j)  sum.a[i][j]=1;            else sum.a[i][j]=0;    while(b){        if(b&1)  sum=matri_mult(sum,A);        A=matri_mult(A,A);        b/=2;    }    return sum;}node in;int main(){    cin>>n;in.a[1][1]=1,in.a[1][2]=1,in.a[2][1]=1,in.a[2][2]=0;//初始化矩陣    in=ksm(in,n-2);if(n<=2)cout<<1;//注意n<=2時矩陣快速冪不適用，需要特判。    else  cout<1][    return 0;}

總結一下，一共有三種求斐波那契數列第n項的方式。第一種是用遞推式，時間複雜度O(n)。第二種是用通項式，時間複雜度O(logn)，但是精度會不夠。第三種時矩陣快速冪，時間複雜度O(logn)，精度足夠。