程序算法藝術與實踐:遞歸策略之Fibonacci數列 · 程序算法藝術與實踐

**背景：** 假定你有一雄一雌一對剛出生的兔子，它們在長到一個月大小時開始交配，在第二月結束時，雌兔子產下另一對兔子，過了一個月后它們也開始繁殖，如此這般持續下去。每只雌兔在開始繁殖時每月都產下一對兔子，假定沒有兔子死亡，在一年后總共會有多少對兔子？在一月底，最初的一對兔子交配，但是還只有1對兔子；在二月底，雌兔產下一對兔子，共有2對兔子；在三月底，最老的雌兔產下第二對兔子，共有3對兔子；在四月底，最老的雌兔產下第三對兔子，兩個月前生的雌兔產下一對兔子，共有5對兔子；……如此這般計算下去，兔子對數分別是：1, 1, 2, 3, 5, 8, 13, 21, 34, 55,89, 144, ...看出規律了嗎？從第3個數目開始，每個數目都是前面兩個數目之和。這就是著名的斐波那契（Fibonacci）數列。 **數學表示：** Fibonacci數列的數學表達式就是： **F(n) = F(n-1) + F(n-2)** **F(1) = 1 ,F(2) = 1** **遞歸程序1：** **Fibonacci數列可以用很直觀的二叉遞歸程序來寫，用C++語言的描述如下：** ~~~ long fib1(int n){ if (n <= 2){ return 1;} else{ return fib1(n-1) + fib1(n-2);} } ~~~ 看上去程序的遞歸使用很恰當，可是在用VC2010的環境下測試n=37的時候用了大約3s，而n=45的時候基本下樓打完飯也看不到結果……顯然這種遞歸的效率太低了！！**遞歸效率分析：**例如，用下面一個測試函數： ~~~ long fib1(int n, int* arr){ arr[n]++; if (n <= 2){ return 1;} else{return fib1(n-1, arr) + fib1(n-2, arr);} } ~~~ 這時，可以得到每個fib(i)被計算的次數： fib(10) = 1???? fib(9) = 1????? fib(8) = 2????? fib(7) = 3 fib(6) = 5????? fib(5) = 8????? fib(4) = 13??? fib(3) = 21 fib(2) = 34?? fib(1) = 55??? fib(0) = 34 可見，計算次數呈反向的Fibonacci數列，這顯然造成了大量重復計算。我們令T(N)為函數fib(n)的運行時間，當N>=2的時候我們分析可知： ***T(N) = T(N-1) + T(N-2) + 2*** 而fib(n) = fib(n-1) + fib(n-2)，所以有T(N) >= fib(n)，歸納法證明可得： ***fib(N) < (5/3)^N*** ***當N>4時，fib（N）>= (3/2)^N*** ***標準寫法：![](https://box.kancloud.cn/612faf49da5e301191db5540aea7e863_187x67.jpg) *** 顯然這個***O（(3/2)^N）***是以指數增長的算法，基本上是最壞的情況。 **其實，這違反了遞歸的一個規則：合成效益法則。*合成效益法則（Compound interest rule）：在求解一個問題的同一實例的時候，切勿在不同的遞歸調用中做重復性的工作。***所以在上面的代碼中調用fib(N-1)的時候實際上同時計算了fib(N-2)。這種小的重復計算在遞歸過程中就會產生巨大的運行時間。 **遞歸程序2：** 用一叉遞歸程序就可以得到近似線性的效率，用C++語言的描述如下： ~~~ long fib(int n, long a, long b, int count){ if (count == n) return b; return fib(n, b, a+b, ++count); } long fib2(int n){ return fib(n, 0, 1, 1); } ~~~ **這種方法雖然是遞歸了，但是并不直觀，而且效率上相比下面的迭代循環并沒有優勢。** **迭代解法：** Fibonacci數列用迭代程序來寫也很容易，用C++語言的描述如下： ~~~ //也可以用數組將每次計算的f(n)存儲下來，用來下次計算用（空間換時間） long fib3 (int n){ long x = 0, y = 1; for (int j = 1; j < n; j++) { y = x + y; x = y - x; } return y; } ~~~ **這時程序的效率顯然為*O（N）*，N = 45的時候<1s就能得到結果。** 我們將數列寫成： Fibonacci[0] = 0，Fibonacci[1] = 1 Fibonacci[n] = Fibonacci[n-1] + Fibonacci[n-2] (n >= 2) 可以將它寫成矩陣乘法形式： ![](https://box.kancloud.cn/4fc5153d5dc4c40ddcd3edf8f918afdb_324x71.jpg) 將右邊連續的展開就得到： ![](https://box.kancloud.cn/73f038caa63d50f550de9118e4881801_328x58.jpg) 下面就是要用O(log(n))的算法計算：![](https://box.kancloud.cn/c586cea186feaa1ae902b410ff6bd650_67x51.jpg) 顯然用二分法來求，結合一些面向對象的概念，C++代碼如下： ~~~ class Matrix{ public: long matr[2][2]; Matrix(const Matrix&rhs); Matrix(long a, long b, long c, long d); Matrix& operator=(const Matrix&); friend Matrix operator*(const Matrix& lhs, const Matrix& rhs){ Matrix ret(0,0,0,0); ret.matr[0][0] = lhs.matr[0][0]*rhs.matr[0][0] + lhs.matr[0][1]*rhs.matr[1][0]; ret.matr[0][1] = lhs.matr[0][0]*rhs.matr[0][1] + lhs.matr[0][1]*rhs.matr[1][1]; ret.matr[1][0] = lhs.matr[1][0]*rhs.matr[0][0] + lhs.matr[1][1]*rhs.matr[1][0]; ret.matr[1][1] = lhs.matr[1][0]*rhs.matr[0][1] + lhs.matr[1][1]*rhs.matr[1][1]; return ret; } }; Matrix::Matrix(long a, long b, long c, long d){ this->matr[0][0] = a; this->matr[0][1] = b; this->matr[1][0] = c; this->matr[1][1] = d; } Matrix::Matrix(const Matrix &rhs){ this->matr[0][0] = rhs.matr[0][0]; this->matr[0][1] = rhs.matr[0][1]; this->matr[1][0] = rhs.matr[1][0]; this->matr[1][1] = rhs.matr[1][1]; } Matrix& Matrix::operator =(const Matrix &rhs){ this->matr[0][0] = rhs.matr[0][0]; this->matr[0][1] = rhs.matr[0][1]; this->matr[1][0] = rhs.matr[1][0]; this->matr[1][1] = rhs.matr[1][1]; return *this; } Matrix power(const Matrix& m, int n){ if (n == 1) return m; if (n%2 == 0) return power(m*m, n/2); else return power(m*m, n/2) * m; } long fib4 (int n){ Matrix matrix0(1, 1, 1, 0); matrix0 = power(matrix0, n-1); return matrix0.matr[0][0]; } ~~~ **這時程序的效率為*O（log(N)）*。公式解法：在*O（1）*的時間就能求得到F(n)了：** ***![](https://box.kancloud.cn/10ad360359a6abd122f88ca89eabbde6_260x64.jpg) *** **注意：其中[x]表示取距離x最近的整數。** **用C++寫的代碼如下：** ~~~ long fib5(int n){ double z = sqrt(5.0); double x = (1 + z)/2; double y = (1 - z)/2; return (pow(x, n) - pow(y, n))/z + 0.5; } ~~~ **這個與數學庫實現開方和乘方本身效率有關的，我想應該還是在O(log(n))的效率。** **總結：** 上面給出了5中求解斐波那契數列的方法，用測試程序主函數如下： ~~~ int main(){ cout << fib1(45) << endl; cout << fib2(45) << endl; cout << fib3(45) << endl; cout << fib4(45) << endl; cout << fib5(45) << endl; return 0; } ~~~ 函數fib1會等待好久，其它的都能很快得出結果，并且相同為：1134903170。而后面兩種只有在n = 1000000000的時候會顯示出優勢。由于我的程序都沒有涉及到高精度，所以要是求大數據的話，可以通過取模來獲得結果的后4位來測試效率與正確性。另外斐波那契數列在實際工作中應該用的很少，尤其是當數據n很大的時候（例如：1000000000），所以綜合考慮基本普通的非遞歸O(n)方法就很好了，沒有必要用矩陣乘法。程序全部源碼 ~~~ class Matrix{ public: long matr[2][2]; Matrix(const Matrix&rhs); Matrix(long a, long b, long c, long d); Matrix& operator=(const Matrix&); friend Matrix operator*(const Matrix& lhs, const Matrix& rhs) { Matrix ret(0,0,0,0); ret.matr[0][0] = lhs.matr[0][0]*rhs.matr[0][0] + lhs.matr[0][1]*rhs.matr[1][0]; ret.matr[0][1] = lhs.matr[0][0]*rhs.matr[0][1] + lhs.matr[0][1]*rhs.matr[1][1]; ret.matr[1][0] = lhs.matr[1][0]*rhs.matr[0][0] + lhs.matr[1][1]*rhs.matr[1][0]; ret.matr[1][1] = lhs.matr[1][0]*rhs.matr[0][1] + lhs.matr[1][1]*rhs.matr[1][1]; return ret; } }; Matrix::Matrix(long a, long b, long c, long d){ this->matr[0][0] = a; this->matr[0][1] = b; this->matr[1][0] = c; this->matr[1][1] = d; } Matrix::Matrix(const Matrix &rhs){ this->matr[0][0] = rhs.matr[0][0]; this->matr[0][1] = rhs.matr[0][1]; this->matr[1][0] = rhs.matr[1][0]; this->matr[1][1] = rhs.matr[1][1]; } Matrix& Matrix::operator =(const Matrix &rhs){ this->matr[0][0] = rhs.matr[0][0]; this->matr[0][1] = rhs.matr[0][1]; this->matr[1][0] = rhs.matr[1][0]; this->matr[1][1] = rhs.matr[1][1]; return *this; } Matrix power(const Matrix& m, int n){ if (n == 1) return m; if (n%2 == 0) return power(m*m, n/2); else return power(m*m, n/2) * m; } //普通遞歸 long fib1(int n){ if (n <= 2) { return 1; } else { return fib1(n-1) + fib1(n-2); } } long fib(int n, long a, long b, int count){ if (count == n) return b; return fib(n, b, a+b, ++count); } //一叉遞歸 long fib2(int n){ return fib(n, 0, 1, 1); } //非遞歸方法O(n) long fib3 (int n){ long x = 0, y = 1; for (int j = 1; j < n; j++) { y = x + y; x = y - x; } return y; } //矩陣乘法O(log(n)) long fib4 (int n){ Matrix matrix0(1, 1, 1, 0); matrix0 = power(matrix0, n-1); return matrix0.matr[0][0]; } //公式法O(1) long fib5(int n){ double z = sqrt(5.0); double x = (1 + z)/2; double y = (1 - z)/2; return (pow(x, n) - pow(y, n))/z + 0.5; } int main(){ //n = 45時候fib1()很慢 int n = 10; cout << fib1(n) << endl; cout << fib2(n) << endl; cout << fib3(n) << endl; cout << fib4(n) << endl; cout << fib5(n) << endl; return 0; } ~~~ 關于[程序算法藝術與實踐](http://blog.csdn.net/column/details/tac-programalgrithm.html)更多討論與交流，敬請關注本博客和新浪微博[songzi_tea](http://weibo.com/songzitea).