GRIS 的文章，没鸽。真的。星月妹妹什么时候骗过你了？而且我还准备开新坑讲一下 BEATLESS。
“回字有四样写法，你知道么？”
而且似乎不少人写过这玩意儿了。所以您也可以将本文作为茶余饭后的消遣，乐乐就好。
一、前言
今天在面试的时候被问到了这样一个问题：求斐波那契数列的第 $n$ 项，你能想到多少种算法？
自然，猴子都知道简单递归。再进一步，黑猩猩都能看出来递归可以用一个额外的表来速查（空间换时间）。有递归就自然要想能不能迭代。稍稍想一下的话可以知道是能迭代的。
然后又问，时间和空间复杂度是多少？
迭代简单，$T(n) = O(n)$，$S(n) = O(1)$。递归我就炸了（所以说我太水了），时间分析就满脸问号（用表的话……优化是优化了，但是和迭代也没差多少了，却还多了调用开销；简单递归分析卡壳了），没答上来，面试官也就没问空间复杂度。所以我当时只想到了两种；后来因为提到迭代的空间复杂度为 $O(1)$，面试官可能没太听清，说道：“（时间复杂度）$O(1)$ 的解法也是有的，不过你用了迭代只能是 $O(n)$ 了。”虽然和我的回答意思相差了一些，但是我还是惊诧了一下（没出声）：“这个 $O(1)$ 是怎么来的？”
所以下午我就搜索了一下。
……连个本科生（复杂度计算本科生都会）都比不过，我自裁罢。
二、准备
斐波那契数列 $f_n (n \in \mathbb{N})$ 的定义如下：
之后的代码以此为准。
为了简单起见，在之后的代码中使用 int，而不是更大的 long long，懒得用无符号数，也没有防溢出机制。
三、直接的代码解法
这部分的方法就是只用到代码、基础算法和优化，几乎不需要数学。
3.1 递归法
递归就是最简单直接的方法了，代码形式和数学形式接近，理解很好理解，大部分工作交给计算机完成。
3.1.1 朴素递归
朴素递归应该学过代码都会，代码就不需要多讲了。
int fib(int n) {
if (n < 1) {
return 0;
} else if (n <= 2) {
return 1;
} else {
return fib(n - 1) + fib(n - 2);
}
}
对应的调用如下：
int result = fib(15);
这个算法的空间复杂度比较好分析（虽然有点绕），是 $O(n)$。令初始调用为第0层，调用一次加深一层。可以知道，在某一层的所有调用结束之前，不会产生新的调用（即，栈深不会增加）。那么栈深最深的是什么时候呢？画一棵调用树就可以看出，是从 fib(n) 开始到 fib(1)（实际上可以算到 fib(2) 为止）的这条“单刀直入”的路径。或者可以这么想，每次 fib(n) 执行实际计算（求和）的时候，fib(n - 1) 已经计算完成了。所以空间复杂度是 $O(n)$（更准确地说是 $\Theta(n)$）。
那么，时间复杂度是多少呢？许多文章直接就说结论是 $O(2^n)$，但是也没给出推导过程，或许也是因为太简单了。观察调用树可以知道，一个非叶子节点最多有两个子节点，而树的深度是 $O(n)$；在每个节点的操作都是简单加法，时间为 $\Theta(1)$；所以总计就是 $O(2^n)$。
但是大多数人没想过，这是一个好的估算吗？我们知道，大 O 函数只是给出上界；给出比一个已知上界更高的上界也不算错。也就是说，我可以胡扯说复杂度是 $O(2^{2^n})$，你也不能说我错了，只能说这个上界“质量不好”。
再次观察调用树，可以很明显地发现不管树高度如何增长，总是有不可忽视的一部分节点“缺失”了。因此，很明显还能找出比 $O(2^n)$ 更紧的上界。
为了计算这个上界，我们不妨来观察函数调用的过程。设 fib(n) 的调用时间为 $t{n}$ ，我们可以知道：$t{n} = \Theta(1) + t{n-1} + t{n-2}$。那么，这不是和原来的斐波那契数列形式完全一样吗？回到原点了。
但是等一下！斐波那契数列是可以求出通项的（见 4.3 节）。提取“公因式” $\Theta(1)$（不是公因式，不过在复杂度分析中可以视作等量级；另外这个 $\Theta(1)$ 虽然上面是出现在递归式里，但是在通项中可以化为一个系数——这么讲不是很准确，不过意思就是这样）之后可以知道：
$\frac{1+\sqrt{5}}{2} \approx 1.618$，虽说看起来好像比2差不了多少，不过对于 $n$ 很大的时候差距就比较明显了。这个上界明显比 $O(2^n)$ 要好，也是能获得的最紧的上界。
3.1.2 尾递归
既然朴素递归开销这么大，有没有可能用无敌的尾递归想想办法呢？毕竟一些现代编译器是支持尾递归优化为循环（迭代）的。（注：其实尾递归并不是无敌的，有着严格的应用条件；这里只是玩梗。）
我们知道，尾递归的一个重要条件是对调用自身后返回的结果不可以进行处理（即，尾调用——对自身的调用必须是返回前的最后一条指令）。那么我们将对之前的两个状态的操作（相加）转移到参数中，不就成了尾递归了吗。代码如下：
int fib(int x1, int x2, int n) {
if (n < 1) {
return 0;
} else if (n <= 2) {
return 1;
} else if (n == 3) {
return x1 + x2;
} else {
return fib(x2, x1 + x2, n - 1);
}
}
对应的调用如下：
int result = fib(1, 1, 15);
这个的复杂度分析要考虑一下编译器。对于支持尾递归优化的编译器，这实际上跟下面的迭代是一样的，所以时间复杂度同为 $O(n)$，空间复杂度同为 $O(1)$。对于不支持的编译器，时间和空间的复杂度都很好分析，看参数 n 的变化就行——单次调用开销为 $\Theta(1)$，所以总复杂度同为 $O(n)$。
3.2 迭代法
不管是否想到了尾递归的方法，迭代都应该比较容易想到。如果想到了尾递归的话，很自然地可以进一步到迭代。
想法直截了当：既然每次只需要前两项的结果（状态），而 $f{n-3}$ 并不会被 $f_n$ 及之后的计算直接用到，那么为什么要保存其他状态呢？$f{n-1}$ 和 $f_{n-2}$ 已经很好地充当了“缓存”的功能。代码如下：
int fib(int n) {
if (n < 1) {
return 0;
} else if (n <= 2) {
return 1;
} else {
int n1 = 1, n2 = 1;
int tmp;
for (int i = 3; i < n; ++i) {
tmp = n1 + n2;
n1 = n2;
n2 = tmp;
}
return n2;
}
}
调用如下：
int result = fib(15);
显而易见，时间复杂度为 $O(n)$，空间复杂度为 $O(1)$（更准确地说是 $\Theta(1)$）。我觉得分析都不用写了。
四、含数学分析的解法
简单暴力的代码解法毕竟是太拿衣服。要往深了优化（或者给出无法继续优化的证据）还是得要数学。
4.1 数列分析
为了理解接下来的分析，首先要知道特征方程。数列的特征方程我记得是高中讲过的，可是后来一直都没怎么用过（高数里也是）；加上本来我数学就弱，所以基本就是忘光了。对不起了，宋老师。
斐波那契数列是一个二阶的递推数列，一般形式为 $x{n} = c{1}x{n-1} + c{2}x{n-2}$。假设能找到两个数 $p$ 和 $q$，使得 $x{n} - px{n-1} = q(x{n-1} - px{n-2})$，移项得 $x_n = (p+q)x{n-1} - pqx{n-2}$。所以有 $p + q = c_1$，$pq = -c_2$。不难看出 $p$ 和 $q$ 的地位是相等的，因此消去任意一个（以 $q$ 为例）可得 $p^2 - c{1}p - c_2 = 0$。这就是二阶递推数列的特征方程。求解该二次方程，两个实根（如果有）就是 $p$ 和 $q$。
代回原来的数列形式，若 $p \ne q$，则通项为 $x_n = \frac{x_2-qx_1}{p(p-q)}p^n + \frac{px_1-x_2}{q(p-q)}q^n$；若 $p=q$，则通项为 $x_n = (\frac{px_1-x_2}{p^2} + \frac{x_2-px_1}{p^2}n)p^n$。
将斐波那契数列代入（即 $c_1 = 1, c_2 = 1$）可以解得 $p = \frac{1+\sqrt{5}}{2}$，$q = \frac{1-\sqrt{5}}{2}$（$p$ 和 $q$ 可互换）。
4.2 改进的迭代：矩阵快速幂
（是搜索之后才知道有可以应用“矩阵快速幂”这么一类问题的。我没参加过 ACM，所以渣了。）
从 4.1 节可知，存在符合条件的 $p$ 和 $q$。于是，就可以将每次求值化为矩阵乘法。
对于 $k$ 阶线性递推数列 $xn = c{1}x{n-1} + c{2}x{n-2} + \cdots + c{k}x_{n-k}$，每次递推可以表示为：
其中大小为 $k \times k$ 的矩阵 $\mathbf{M}$ 就是……呃，好像没有标准称呼？就叫递推矩阵好了。由此可见，由于 $\mathbf{M}$ 在递推过程中不变，计算 $x_n$ 就被转化为了计算矩阵的幂的问题。这个矩阵受初始值选取和递推公式两方面的影响。
对于斐波那契数列，给定初始值 $[1, 1]$，可以得到
虽说是矩阵快速幂，但是对幂函数的计算优化其实和数的是一样的，只不过乘法的定义略有不同。
要计算一个数 $a$ 的 $n$ 次幂，最直接的想法是连乘，时间复杂度为 $O(n)$。但是只要见过质数筛这类的问题，或者是插值函数优化，就马上会意识到，其中许多次的乘法是多余的，结果是可以缓存的（甚至可以不用缓存）。那么，如何剔除这些多余的计算呢？
不是有对数函数嘛。大杀器。
考虑到实现的便利性，一般就取 2 为底了。将幂次 $n$ 分解为 2 的幂的和：$n = c{0}2^0 + c{1}2^1 + \cdots + c{m}2^{m}$，其中 $c_i \in {0, 1} (0 \le i \lt m), c_m = 1$。进一步展开到乘法，就可以看到，$a^n = a^{c{0}2^0} \times a^{c{1}2^1} \times \cdots \times a^{c{m}2^m}$。这样，从前往后“扫描”，由于 $a^{2^n} \times a^{2^n} = a^{2^{n+1}}$（废话），所以只需要保存最近的中间结果，就可以依次得到正确的乘积，并最终得到 $a^n$。代码如下：
// 简洁起见，仅处理 a > 0, n >= 0 的情况
int pown_fast(int a, int n) {
int result = 1;
int cur = a;
while (n > 0) {
if (n & 1) {
result = result * cur;
}
cur = cur * cur;
n = n >> 1;
}
return result;
}
可以见到，pown_fast() 的时间复杂度为 $O(\log n)$，空间复杂度为 $O(1)$（$\Theta(1)$）。对于矩阵的幂，由于斐波那契数列相关的递推矩阵大小是常数，所以时间和空间复杂度都和计算数的幂是一样的。再考虑初始值的输入等等，最终的时间复杂度为 $O(\log n) + \Theta(1) = O(\log n)$，空间复杂度为 $\Theta(1) + \Theta(1) = \Theta(1)$。
也有用类似思想，但不是快速幂的解法，而是采用新的递推公式，参见[1]。
4.3 通项公式
将 $p$ 和 $q$ 的值代入 4.1 节的通项公式可得： $f_n = \frac{1}{\sqrt{5}}((\frac{1+\sqrt{5}}{2})^n + (\frac{1-\sqrt{5}}{2})^n)$。
按照这个通项计算 $f_n$，理论上可以达到 $O(1)$ 的时间复杂度。
很美妙，不是吗？$O(1)$ 和 $O(1)$ 啊！但是简单地计算依然会出问题。这个公式中，多处涉及了（在最终结果前）无法约简的无理数。在现代的计算框架中，如果按照一般的思维，用浮点数来表示这些无理数，进行计算的话，一定有因浮点数的不精确性导致的误差。$n$ 越大，误差越大。且不说“四舍五入”这个行为在此正不正确，当误差超过 0.5 的时候，舍入都救不回来了。所以在实践中，这个方法是不可取的。
当然不是没有能保证精度的方法。采用符号代数就可以了。而且这个计算过程中的无理数只有 $\sqrt{5}$ 一个，而且 $(\sqrt{5})^2 = 5$ 还是个整数。不过计算过程就得根据二项式定理展开，然后累加并分类求和，所花的时间肯定会多于 $O(n)$，空间开销也不会小。
五、结论
以上我们见到了五种计算 $f_n$ 的方法：
- 朴素递归，$T(n) = O((\frac{1+\sqrt{5}}{2})^n)$，$S(n) = O(n)$；
- 尾递归，在有优化的情况下 $T(n) = O(n)$，$S(n) = O(1)$；
- 迭代，$T(n) = O(n)$，$S(n) = O(1)$；
- 快速幂，$T(n) = O(\log n)$，$S(n) = O(1)$；
- 通项，$T(n) = O(1)$，$S(n) = O(1)$。
在几种方法中，时间和空间复杂度都是从前往后越来越好。但是通项由于需要用浮点数表示无理数的值会导致误差所以基本上不能用。所以综合各种因素考虑，最优的算法是采用快速幂计算 $f_n$。
参考：