算法 · 深入理解 Fibonacci 数列计算及黄金分割

心宿二_Sco
2,339 阅读10分钟

在自然界中，有一串数字我们时常能看到，它的前几项序列是：

0, 1, 1, 2, 3, 5, 8, 13, 21, 34, 55, 89, 144, 233...

如果用数学语言描述，它的递归形式是这样的：

{F0=0F1=1Fn=Fn−1+Fn−2(n≥2)\begin{cases} F_0 = 0 \\ F_1 = 1 \\ F_n = F_{n-1}+F_{n-2} &(n\geq2) \end{cases}

也就是说，除了 第0项 和 第1项，之后的数字是由这个数紧邻的前两个数字，相加得到的。

如果用几何图形表示，你可能见过这张图：

image.png

这就是斐波那契数，它形成的数列就被称为斐波那契数列。

名称源自1170年出生的意大利比萨数学家列奥纳多·斐波那契：Leonardo Fibonacci。

（歪个楼：列奥纳多·迪·塞尔·皮耶罗·达·芬奇，出生于1452年）

求解和利用斐波那契数及数列，是现代计算机科学与传统古典艺术、现代新媒体艺术设计的热门交叉点，包括但不限于图形学、性能优化、先锋艺术、影视创作等等分支领域。

其中最为著名的，就是黄金螺旋分割线。所以此篇追根溯源，深入到算法和数学的最底层面，理解斐波那契数列的计算方式。

① Fibonacci_1(n)：定义式二分递归计算

根据上述定义式，fib() 函数的伪代码描述为（默认 n 为非负数）：

fib(n)={n(ifn≤1)fib(n−1)+fib(n−2)(ifn≥2)fib(n)= \begin{cases} n & (if &n\leq1)\\ fib(n-1)+fib(n-2) & (if &n\geq2) \end{cases}
// C++
int fib1(int n){
    return (2 > n) ? n : fib1(n - 1) + fib1(n - 2);
}

显然，这种算法浅显易懂，直接根据定义式编码即可。

但弊端也很显而易见，二分递归的形式，其计算效率非常低下。

假设，计算 ***fib(n)***所需的时间为 T(n)。那么根据此算法思路，计算 ***fib(n-1)***则相应的需要 ***T(n-1)***的时间，直到最后花费一个单位的时间，将其累加。

故，***T(n)***的递推式为：

T(n)={1(ifn≤1)T(n−1)+T(n−2)+1(ifn≥2)T(n)= \begin{cases} 1 & (if &n\leq1)\\ T(n-1)+T(n-2) + 1 & (if &n\geq2) \end{cases}

仔细观察，就会发现 ***T(n)***与 ***fib(n)***长得非常像。

那么，令 S(n) = [T(n)+1]/2，

S(n)=[T(n)+1]/2={(1+1)/2[T(n−1)+T(n−2)+1+1]/2={(1+1)/2{[2∗S(n−1)−1]+[2∗S(n−2)−1]+1+1]}/2={(1+1)/2[2∗S(n−1)+2∗S(n−2)]/2={1(ifn≤1)S(n−1)+S(n−2)(ifn≥2)\begin{aligned} S(n)&= [T(n)+1]/2 \\ &= \begin{cases} (1+1)/2 \\ [T(n-1)+T(n-2)+1+1]/2 \\ \end{cases} \\ & = \begin{cases} (1+1)/2 \\ \{[2*S(n-1)-1] + [2*S(n-2)-1] +1 +1 ]\}/2 \end{cases} \\ & = \begin{cases} (1+1)/2 \\ [2*S(n-1) + 2*S(n-2)]/2 \\ \end{cases} \\ &= \begin{cases} 1 & (if &n\leq1)\\ S(n-1)+S(n-2) & (if &n\geq2) \end{cases} \\ \end{aligned}

由此可以发现，除了起始项目外，***S(n)***与 fib(n) 在形式上一致。

而且：

{S(0)=(T(0)+1)/2=1=fib(1)S(1)=(T(1)+1)/2=1=fib(2)\begin{cases} S(0) = (T(0)+1)/2 = 1 = fib(1) \\ S(1) = (T(1)+1)/2 = 1 = fib(2) \\ \end{cases}

也就是说，S(n)=fib(n+1)

进而：

S(n)=fib(n+1)=(ϕn+1−ϕ^n+1)/5, ϕ=(1+5)/2, ϕ^=(1−5)/2T(n)=2∗S(n)−1=2∗fib(n+1)−1=O(ϕn+1)≈O(2n)\begin{aligned} & S(n)=fib(n+1)=(\phi^{n+1}-\hat{\phi}^{n+1})/\sqrt{5},\ \phi = (1+\sqrt{5})/2,\ \hat{\phi}=(1-\sqrt{5})/2 \\ & T(n) = 2*S(n)-1 = 2*fib(n+1)-1=O(\phi^{n+1})\approx O(2^n) \\ \end{aligned}

也就是说，二分递归的算法时间复杂度达到了恐怖的指数量级，而这个结果是远远不能被接受的（关于 φ 值是如何计算得出的稍后详解）。

如果说指数级复杂度，是给这个算法做了定性分析，那么我们定量举个例子。

最普通的 PC 运算性能大致为 1GHz，意味着每秒可以做 10^9 次浮点预算，也就是十亿次。

当我们计算第 100 项斐波那契数时，大约需要 2^100 次运算，那么大概耗时为 2^100/10^9 ，这是多少呢？

粗略估算，2^10 也就是 1024，约等于 10^3，也就是 1000。

那么：

2100=(210)10≈(103)102100/109≈(103)10/109=1021\begin{aligned} & 2^{100} = (2^{10})^{10} \approx (10^3)^{10} \\ & 2^{100}/10^9 \approx (10^3)^{10}/10^9 = 10^{21} \\ \end{aligned}

问题来了，10^21 秒是多少：

{1 day=24hr∗60min∗60sec=86400sec≈25∗4000=105sec100 year≈100∗365 day≈3∗104day=3∗109sec\begin{cases} 1 \ day = 24hr * 60min*60sec = 86400 sec \approx 25*4000 = 10^5 sec \\ 100 \ year \approx 100*365\ day \approx 3 * 10^4 day = 3* 10^9 sec \end{cases}

我们可以粗略计算，按照上式，三个世纪大概为 10^10 秒。

则，10^21 秒相当于 10^21/10^10 = 10^11 个“三个世纪”。

也就是三千亿个世纪

也就是三十万亿年

现有可观测宇宙的寿命，学界普遍认为在140亿年。可以毫不夸张的说，如果按照定义式的方式交给计算机做运算，大概我们现存的宇宙都重新诞生了上千次了。

而这漫长的一切，只不过是为了求解 第100项 斐波那契数 而已。

所以，这种二分递归，层层往复计算的算法，显然是不能被接受的。

② Fibonacci_2(n)：线性递归计算

从二分递归的算法中，不难看出问题所在。

每一次计算新的值时，需要重新回到它的前一项，然后层层往回走，直到最开始的起始项。然后再重新开始从起始项一直加到当前项。

在这个过程中，有大量相同的项目被无数次重复计算，这也是为什么时间复杂度如此之高，计算资源如此被大量浪费的原因。

顺着这种思路，如果说我可以让程序拥有一种记忆力，将已经计算过的项目保存下来，在需要的时候直接提取，而不是返到起点再重头计算一遍。这样是不是就可以节省大量的重复劳动了？

这种开辟存储空间作为计算辅助的做法，就是非常典型的“空间换时间”的思维。

在这个算法中，单独开辟一个长度为 n 的数组空间，程序运算到第几项，就将这一项的结果存储在空间中。

从而当需要当先项数据时，只需直接调取即可。计算第 n 项时，程序只需要从起点跑一趟，就可以累算得到当前项目值。这种方式也更贴近我们人类计算这个数列的方法。

// C++
int fib2(int n, int& prev){
    if (0 == n){
        prev = 1;
        return 0;
    }
    else {
        int prevPrev;
        prev = fib2(n - 1, prevPrev);
        return prevPrev + prev;
    }
}

这种算法的时间复杂度为 O(n)，同时，由于辅助的存储空间与项目所在数列的位置（距离/长度/坐标...）成正比，所以，空间复杂度同样也为 O(n)。

如果这里我们要用微机计算 第100项 斐波那契数，理论上只需要 100/10^9 秒的时间，远远小于一秒钟。这种算法上的优化效果，是远高于第一种二分递归算法的。

③ Fibonacci_3(n)：迭代计算

然而，仔细观察，线性递归算法还有它的问题所在。

对于计算机来说，“空间”是一种限量的资源，数字小一点还好，如果当数据量变大，空间消耗会成为不得不面对的现实。

而在上述开辟的数组空间里，其实真正用到的，仅仅为最后两项数字。也就是斐波那契数的当前结果和前第一项（前第二项，可以用当前结果减去第一项计算得到）。

对于利用过的小项目结果，再往后其实程序就不会再访问它，那么它存在的实际意义就消失了。

如果程序可以自动销毁之前的结果，让整个辅助空间保持在一定的大小范围内，同计算项的增长而动态跟随，始终保存最新的两个项目，避免浪费。那么空间消耗，将不再与计算规模挂钩。换言之，空间资源的开支，成为常数级，也就是说，空间复杂度将来到极致的 O(1)。

顺着这种思路继续优化下去。

我们其实不需要另外开一个程序去实时紧盯不用的计算结果并删除它，换一个角度想，我们只要在固定的空间内反复写入新结果，覆盖旧数据，就可以实现同样的效果。

事实上，我们只需要两个变量即可。

一个存储当前项，一个存储前一项。

// C++
int fib3(int n){
    int f = 1, g = 0;
    while (0 < n--) { g += f; f = g - f; }
    return g;
}

程序中，g 为最终结果，f 为前一项。

没有高昂的数组空间开支，也没有复杂的资源管理。

只需要两个变量，仿佛二龙戏珠盘旋向上一样，从 0 一直爬升计算到当前项即可。

算法的时间复杂度同样依然为 O(n)，而空间复杂度，已经被优化到了极致的 O(1)，也有称之为【就地算法】。

可以看到，程序非常的简洁，效率也足够的优秀。短短几个字符，凝结了人类顶级的智慧，不得不感叹算法的魅力。

④ Fibonacci_4(n)：通项公式求解

但是！

还没有结束！

人类追求极致智慧的脚步是不会停下来的。

空间复杂度已经被优化到了极致的 O(1)，那么评价算法优劣的指标剩下了另一项，也就是时间复杂度。

到目前为止，任意项斐波那契数求解计算的时间复杂度还依然是 O(n)。

接下来要做的，就是将它进一步优化，直到完美的 O(1)。

你可能会困惑，依托前 n-1 和 n-2 项的斐波那契数，无论如何都无法跳脱从起点到第 n 项的这个遍历过程。

是的，所以如果想实现最为极致的优化。

那么，欢迎来到数学的世界。  

注意到了吗，我们最开始，对斐波那契数的定义，是通过它的自然形式，自然而然推出来的、形象的递归公式。

生动的描述方式固然便于理解，但是难免有些原始。我们接下来要做的，是通过递归描述，推导出计算任意项斐波那契数的通项公式。   也就是说，只有一个变量，和一堆系数及数学方法，不经过逐个推导，直接用一个方程算出结果。

为了便于理解，我们用高中生都会的初等数学知识推导：

目前已知：

a1=1a2=1an=an−1+an−2\begin{aligned} & a_1 = 1 \\ & a_2 = 1 \\ & a_n = a_{n-1} + a_{n-2} \\ \end{aligned}

设：

an+α∗an−1=β∗(an−1+α∗an−2)a_n + \alpha*a_{n-1} = \beta*(a_{n-1}+\alpha*a_{n-2})

即：

an=(β−α)∗an−1+α∗β∗an−2a_n = (\beta-\alpha)*a_{n-1}+\alpha*\beta*a_{n-2}

根据：

{an=an−1+an−2an=(β−α)∗an−1+α∗β∗an−2\begin{cases} a_n = a_{n-1} + a_{n-2} \\ a_n = (\beta-\alpha)*a_{n-1}+\alpha*\beta*a_{n-2} \\ \end{cases}

可以求得：

{β−α=1α∗β=1=>{β=1+αα∗(1+α)=1=>α2+α−1=0\begin{cases} \beta - \alpha = 1 \\ \alpha * \beta = 1 \\ \end{cases} => \begin{cases} \beta =1 + \alpha \\ \alpha*(1+\alpha)=1 \\ \end{cases} => \alpha^2+\alpha-1=0
α1,2=−b±b2−4ac2a=−1±1+42=−1±52\begin{aligned} \alpha_{1,2} &= \frac{-b\pm\sqrt{b^2-4ac}} {2a} \\ &= \frac{-1\pm\sqrt{1+4}}{2} \\ & = \frac{-1\pm\sqrt{5}}{2} \end{aligned}

这里可以取正数，那么带回公式，可求得：

{α=5−12β=5+12\begin{cases} \alpha=\frac{\sqrt{5}-1}{2} \\ \beta=\frac{\sqrt{5}+1}{2} \\ \end{cases}

因为等比数列

Tn=an+α∗an−1Tn−1=an−1+α∗an−2Tn=β∗Tn−1\begin{aligned} & T_n = a_{n}+\alpha*a_{n-1} \\ & T_{n-1} = a_{n-1}+\alpha*a_{n-2} \\ & T_{n}=\beta*T_{n-1} \end{aligned}

所以：

an+1+α∗an=(a2+α∗a1)∗βn−1=(1+α)∗βn−1=βn\begin{aligned} a_{n+1}+\alpha*a_{n}&=(a_2+\alpha*a_1)*\beta^{n-1} \\ &= (1+\alpha)*\beta^{n-1} \\ &= \beta^n \end{aligned}

等式两边同时除去 β 的 n+1 次，可得：

an+1+α∗anβn+1=βnβn+1an+1βn+1+αβ∗anβn=1β\begin{aligned} & \frac{a_{n+1}+\alpha*a_{n}}{\beta^{n+1}}=\frac{\beta^{n}}{\beta^{n+1}} \\ & \frac{a_{n+1}}{\beta^{n+1}}+\frac{\alpha}{\beta}*\frac{a_n}{\beta^n}=\frac{1}{\beta} \\ \end{aligned}

令：

bn=anβnb_n=\frac{a_n}{\beta^n}

则代入上式：

bn+1+αβbn=1βb_{n+1}+\frac{\alpha}{\beta}b_n=\frac{1}{\beta}

设：

bn+1+λ=−αβ(bn+λ)b_{n+1}+\lambda = - \frac{\alpha}{\beta}(b_n+\lambda)

则：

bn+1+αβbn=−αβλ−λb_{n+1}+\frac{\alpha}{\beta}b_n=-\frac{\alpha}{\beta}\lambda-\lambda

即：

1β=−αβλ−λ\frac{1}{\beta} =-\frac{\alpha}{\beta}\lambda-\lambda

得：

λ=−1α+β\lambda = -\frac{1}{\alpha+\beta}

所以：

{bn+λ=(−αβ)n−1(b1+λ)b1=α1β=1βbn=anβnα=5−12β=5+12\begin{cases} b_n+\lambda=(-\frac{\alpha}{\beta})^{n-1}(b_1+\lambda) \\ b_1 = \frac{\alpha_1}{\beta} = \frac{1}{\beta} \\ b_n = \frac{a_n}{\beta^n} \\ \alpha = \frac{\sqrt{5}-1}{2} \\ \beta = \frac{\sqrt{5}+1}{2} \\ \end{cases}

代入可得：

anβn+λ=(−αβ)n−1∗(1β+λ)an−15∗βn=−1∗(−α)n−1∗β∗(1β−15)\begin{aligned} & \frac{a_n}{\beta^n}+\lambda = (-\frac{\alpha}{\beta})^{n-1}*(\frac{1}{\beta}+\lambda) \\ & a_n-\frac{1}{\sqrt{5}}*\beta^{n}=-1*(-\alpha)^{n-1}*\beta*(\frac{1}{\beta}-\frac{1}{\sqrt{5}}) \\ \end{aligned}

进一步求解：

an=−1∗(−α)n−1∗(1−β5)+15∗βn=−1∗(−(5−1)2)n−1∗(12)∗(2∗55−5+15)+15∗(5+12)n=15∗[(1+52)n−(1−52)n]\begin{aligned} a_n &= -1*(-\alpha)^{n-1}*(1-\frac{\beta}{\sqrt{5}})+\frac{1}{\sqrt{5}}*\beta^n \\ & =-1*(\frac{-(\sqrt{5}-1)}{2})^{n-1}*(\frac{1}{2})*(\frac{2*\sqrt{5}}{\sqrt{5}}-\frac{\sqrt{5}+1}{\sqrt{5}})+\frac{1}{\sqrt{5}}*(\frac{\sqrt{5}+1}{2})^n \\ &= \frac{1}{\sqrt{5}}*[(\frac{1+\sqrt{5}}{2})^n-(\frac{1-\sqrt{5}}{2})^n] \end{aligned}

最终可得，斐波那契数列的通项公式为：

fib(n)=(ϕn−ϕ^n)/5, ϕ=(1+5)/2, ϕ^=(1−5)/2fib(n)=(\phi^{n}-\hat{\phi}^{n})/\sqrt{5},\ \phi = (1+\sqrt{5})/2,\ \hat{\phi}=(1-\sqrt{5})/2

这也就是 解法① 中提到的时间复杂度为什么是指数级。

当然，求解通项公式的方法还有很多，这里就不赘述了，详情可以自己查询。

有了通项公式，程序就十分容易编写：

// C++
int fib4(int n)
{
    float var1 = (1 + sqrt(5)) / 2;
    float var2 = (1 - sqrt(5)) / 2;
    float var3 = sqrt(5) / 5;
    return var3 * (pow(var1, n) - pow(var2, n));
}

据查 C++ pow 函数可达 O(1)，那么很自然的，时间复杂度为 O(1)，空间复杂度也为 O(1)。

即便使用快速幂运算，时间复杂度依然可以降至 O(logn)

⑤ 黄金分割

我们经常听说的黄金分割，是一种比例。

定义为：

a+ba=ab=ϕ (a>b>0)\frac{a+b}{a}=\frac{a}{b}=\phi \ (a>b>0)

这个比值取小数点后三位为 1.618。

image.png

是不是感觉这个形式很面熟？

德国天文学家开普勒发现，斐波那契数列的前后两项之比:

11,12,23,35,58,813,1321,2134,...\frac{1}{1}, \frac{1}{2}, \frac{2}{3}, \frac{3}{5}, \frac{5}{8}, \frac{8}{13}, \frac{13}{21}, \frac{21}{34}, ...

在 n 趋于无穷大的时候，结果就是黄金分割比例了。

lim⁡n→∞fibn+1fibn=12(1+5)=ϕ≈1.618...\lim_{n \rightarrow \infty}{\frac{fib_{n+1}}{fib_{n}}} = \frac{1}{2}(1+\sqrt5)=\phi \approx 1.618...

image.png

黄金分割的应用非常广泛，非常著名的一个示意图就是下面蒙娜丽莎的这张：

image.png

斐波那契数列最早被人们观察到是在自然界中，花瓣的生长分布、向日葵花盘中葵花籽的分布、鹦鹉螺螺纹的变化，甚至银河系旋臂的走势，其实它就在我们身边。

作为自然界神奇的表现，斐波那契数列抽象描述了这一普遍而又充满神秘魅力的现象。

之所以人们普遍认为黄金分割是美的，多半是源自它抽象出的自然界的真实。

当然远远不止这些，

比如 UI 设计中界面的划分占比，页面排版的美观，logo 设计中的美感，建筑设计里的秩序。

甚至在计算机科学领域，还有专门的【斐波那契查找】这种独特的查找算法。

这里不再多述应用场景了。  

此篇为带大家深入理解斐波那契数列背后的数理知识和计算方式，当我们从最底层了解它的原理，表面的众多纷繁，也自然有了必然的联系。也为我们的进一步创新，提供更为扎实的理论和思路。

【此文原创，欢迎转发，禁止搬运】

avatar
文章
阅读
粉丝