斐波那契数列

斐波那契数列（The Fibonacci sequence，OEIS A000045）的定义如下：

\[ F_0 = 0, F_1 = 1, F_n = F_{n-1} + F_{n-2} \]

该数列的前几项如下：

\[ 0, 1, 1, 2, 3, 5, 8, 13, 21, 34, 55, 89, \dots \]

卢卡斯数列

卢卡斯数列（The Lucas sequence，OEIS A000032）的定义如下：

\[ L_0 = 2, L_1 = 1, L_n = L_{n-1} + L_{n-2} \]

该数列的前几项如下：

\[ 2, 1, 3, 4, 7, 11, 18, 29, 47, 76, 123, 199, \dots \]

研究斐波那契数列，很多时候需要借助卢卡斯数列为工具．

斐波那契数列通项公式

第 \(n\) 个斐波那契数可以在 \(\Theta (n)\) 的时间内使用递推公式计算．但我们仍有更快速的方法计算．

解析解

解析解即公式解．我们有斐波那契数列的通项公式（Binet's Formula）：

\[ F_n = \frac{\left(\frac{1 + \sqrt{5}}{2}\right)^n - \left(\frac{1 - \sqrt{5}}{2}\right)^n}{\sqrt{5}} \]

这个公式可以很容易地用归纳法证明，当然也可以通过生成函数的概念推导，或者解一个方程得到．

当然你可能发现，这个公式分子的第二项总是小于 \(1\)，并且它以指数级的速度减小．因此我们可以把这个公式写成

\[ F_n = \left[\frac{\left(\frac{1 + \sqrt{5}}{2}\right)^n}{\sqrt{5}}\right] \]

这里的中括号表示取离它最近的整数．

这两个公式在计算的时候要求极高的精确度，因此在实践中很少用到．但是请不要忽视！结合模意义下二次剩余和逆元的概念，在 OI 中使用这个公式仍是有用的．

卢卡斯数列通项公式

我们有卢卡斯数列的通项公式：

\[ L_n = \left(\frac{1 + \sqrt{5}}{2}\right)^n + \left(\frac{1 - \sqrt{5}}{2}\right)^n \]

与斐波那契数列非常相似．事实上有：

\[ \frac{L_n + F_n\sqrt{5}}{2} = \left(\frac{1 + \sqrt{5}}{2}\right)^n \]

也就是说，\(L_n\) 和 \(F_n\) 恰好构成 \(\left(\frac{1 + \sqrt{5}}{2}\right)^n\) 二项式展开再合并同类项后的分子系数．也就是说，Pell 方程

\[ x^2-5y^2=-4 \]

的全体解，恰好是

\[ \frac{x_n + y_n\sqrt{5}}{2} = \frac{L_n + F_n\sqrt{5}}{2} \]

恰好是卢卡斯数列和斐波那契数列．因此有

\[ {L_n}^2-5{F_n}^2=-4 \]

矩阵形式

斐波那契数列的递推可以用矩阵乘法的形式表达：

\[ \begin{bmatrix}F_{n-1} & F_{n} \cr\end{bmatrix} = \begin{bmatrix}F_{n-2} & F_{n-1} \cr\end{bmatrix} \begin{bmatrix}0 & 1 \cr 1 & 1 \cr\end{bmatrix} \]

设 \(P = \begin{bmatrix}0 & 1 \cr 1 & 1 \cr\end{bmatrix}\)，我们得到

\[ \begin{bmatrix}F_n & F_{n+1} \cr\end{bmatrix} = \begin{bmatrix}F_0 & F_1 \cr\end{bmatrix} P^n \]

于是我们可以用矩阵乘法在 \(\Theta(\log n)\) 的时间内计算斐波那契数列．此外，前一节讲述的公式也可通过矩阵对角化的技巧来得到．

快速倍增法

使用上面的方法我们可以得到以下等式：

\[ \begin{aligned} F_{2k} &= F_k (2 F_{k+1} - F_{k}) \\ F_{2k+1} &= F_{k+1}^2 + F_{k}^2 \end{aligned} \]

于是可以通过这样的方法快速计算两个相邻的斐波那契数（常数比矩乘小）．代码如下，返回值是一个二元组 \((F_n,F_{n+1})\)．

pair<int, int> fib(int n) {
  if (n == 0) return {0, 1};
  auto p = fib(n >> 1);
  int c = p.first * (2 * p.second - p.first);
  int d = p.first * p.first + p.second * p.second;
  if (n & 1)
    return {d, c + d};
  else
    return {c, d};
}

性质

斐波那契数列拥有许多有趣的性质，这里列举出一部分简单的性质：

卡西尼性质（Cassini's identity）：\(F_{n-1} F_{n+1} - F_n^2 = (-1)^n\)．
附加性质：\(F_{n+k} = F_k F_{n+1} + F_{k-1} F_n\)．
取上一条性质中 \(k = n\)，我们得到 \(F_{2n} = F_n (F_{n+1} + F_{n-1})\)．
由上一条性质可以归纳证明，\(\forall k\in \mathbb{N},F_n|F_{nk}\)．
上述性质可逆，即 \(\forall F_a|F_b,a|b\)．
GCD 性质：\((F_m, F_n) = F_{(m, n)}\)．
以斐波那契数列相邻两项作为输入会使欧几里德算法达到最坏复杂度（具体参见维基 - 拉梅）．

斐波那契数列与卢卡斯数列的关系

不难发现，关于卢卡斯数列与斐波那契数列的等式，与三角函数公式具有很高的相似性．比如：

\[ \frac{L_n + F_n\sqrt{5}}{2} = \left(\frac{1 + \sqrt{5}}{2}\right)^n \]

与

\[ \cos nx + i\sin nx = \left(\cos x + i\sin x\right)^n \]

很像．以及

\[ {L_n}^2-5{F_n}^2=-4 \]

与

\[ \cos^2 x + \sin^2 x = 1 \]

很像．因此，卢卡斯数列与余弦函数很像，而斐波那契数列与正弦函数很像．比如，根据

\[ \left(\frac{1 + \sqrt{5}}{2}\right)^m\left(\frac{1 + \sqrt{5}}{2}\right)^n = \left(\frac{1 + \sqrt{5}}{2}\right)^{m+n} \]

可以得到两下标之和的等式：

\[ 2L_{m+n}=5F_mF_n+L_mL_n \]

\[ 2F_{m+n}=F_mL_n+L_mF_n \]

于是推论就有二倍下标的等式：

\[ L_{2n}={L_n}^2-2{\left(-1\right)}^n \]

\[ F_{2n}=F_nL_n \]

这也是一种快速倍增下标的办法．同样地，也可以仿照三角函数的公式，比如奇偶性、和差化积、积化和差、半角、万能代换等等，推理出更多有关卢卡斯数列与斐波那契数列的相应等式．

斐波那契编码

我们可以利用斐波那契数列为正整数编码．根据齐肯多夫定理，任何自然数 \(n\) 可以被唯一地表示成一些斐波那契数的和：

\[ N = F_{k_1} + F_{k_2} + \ldots + F_{k_r} \]

并且 \(k_1 \ge k_2 + 2,\ k_2 \ge k_3 + 2,\ \ldots,\ k_r \ge 2\)（即不能使用两个相邻的斐波那契数）

于是我们可以用 \(d_0 d_1 d_2 \dots d_s 1\) 的编码表示一个正整数，其中 \(d_i=1\) 则表示 \(F_{i+2}\) 被使用．编码末位我们强制给它加一个 1（这样会出现两个相邻的 1），表示这一串编码结束．举几个例子：

\[ \begin{aligned} 1 &=& 1 &=& F_2 &=& (11)_F \\ 2 &=& 2 &=& F_3 &=& (011)_F \\ 6 &=& 5 + 1 &=& F_5 + F_2 &=& (10011)_F \\ 8 &=& 8 &=& F_6 &=& (000011)_F \\ 9 &=& 8 + 1 &=& F_6 + F_2 &=& (100011)_F \\ 19 &=& 13 + 5 + 1 &=& F_7 + F_5 + F_2 &=& (1001011)_F \end{aligned} \]

给 \(n\) 编码的过程可以使用贪心算法解决：

从大到小枚举斐波那契数 \(F_i\)，直到 \(F_i\le n\)．
把 \(n\) 减掉 \(F_i\)，在编码的 \(i-2\) 的位置上放一个 1（编码从左到右以 0 为起点）．
如果 \(n\) 为正，回到步骤 1．
最后在编码末位添加一个 1，表示编码的结束位置．

解码过程同理，先删掉末位的 1，对于编码为 1 的位置 \(i\)（编码从左到右以 0 为起点），累加一个 \(F_{i+2}\) 到答案．最后的答案就是原数字．

模意义下周期性

对于模 \(m\) 意义下的斐波那契数列，可以容易地使用抽屉原理证明，该数列是有周期性的．由于斐波那契数每一项的计算都依赖于前两项的取值，所以需要用相邻斐波那契数组成的数对描述数列当且所处的状态．考虑模意义下前 \(m^2+1\) 个斐波那契数对：

\[ (F_0,\ F_1),\ (F_1,\ F_2),\ \ldots,\ (F_{m^2},\ F_{m^2 + 1}) \]

模 \(m\) 的剩余系大小为 \(m\)，这意味着至多只可能有 \(m^2\) 种互不相同的数对．因此，在前 \(m^2+1\) 个数对中必有两个相同的数对，于是从这两个数对可以往后生成相同的斐波那契数列．那么，斐波那契数列就是周期性的，且（最小正）周期不会超过 \(m^2\)．

Pisano 周期

模 \(m\) 意义下斐波那契数列的最小正周期被称为 Pisano 周期（Pisano period，皮萨诺周期，OEIS A001175）．本文中用 \(\pi(m)\) 表示模 \(m\) 的 Pisano 周期．

这一观察可以用于计算第 \(n\) 项斐波那契数模 \(m\) 的值．如果 \(n\) 非常大，就需要计算斐波那契数模 \(m\) 的周期．当然，只需要计算周期，不一定是最小正周期．

为此，有如下结论：

对于互素的模数 \(m_1,m_2\)，有 \(\pi(m_1m_2)=\operatorname{lcm}(\pi(m_1),\pi(m_2))\)．
对于素数 \(p\) 和正整数 \(e\)，有 \(\pi(p^{e})\mid p^{e-1}\pi(p)\)．
对于 \(m=2^e~(e\in\mathbf N_+)\)，有 \(\pi(m)=3\cdot 2^{e-1}\)．
对于 \(m=5^e~(e\in\mathbf N_+)\)，有 \(\pi(m)=4\cdot 5^e\)．
最后，对于素数 \(p\equiv\pm1\pmod{10}\)，有 \(\pi(p)\mid(p-1)\)；对于素数 \(p\equiv\pm3\pmod{10}\)，有 \(\pi(p)\mid 2(p+1)\)．

综合这些情形，可以说明：模 \(m\) 的 Pisano 周期不会超过 \(6m\)，等号当且仅当 \(m = 2\times 5^e~(e\in\mathbf N_+)\) 时取得．

利用上述结论，可以基于素因数分解算法，得到如下快速计算 Pisano 周期的方法：

参考代码

// Get a period of Fibonacci sequence mod m.
// Not necessarily be the exact Pisano period.
uint32_t calc_cycle_from_mod(uint32_t m) {
  uint32_t res = 1;
  for (auto pe : factorize(m)) {
    auto p = pe.first;
    auto e = pe.second;
    uint64_t cur = pow(p, e - 1);
    if (p == 2) {
      cur *= 3;
    } else if (p == 5) {
      cur *= 20;
    } else if (p % 5 == 1 || p % 5 == 4) {
      cur *= p - 1;
    } else {
      cur *= 2 * (p + 1);
    }
    res = lcm(res, cur);
  }
  return res;
}

这样得到的周期可能只是 Pisano 周期的一个倍数．要得到精确的 Pisano 周期，可以进一步考察该周期的因数；或者，可以直接通过 BSGS 算法以 \(O(\sqrt{m})\) 的时间复杂度计算．

证明

最后，本文简要证明上述关于 Pisano 周期的结论．值得说明的是，利用下文说明的方法，类似的结论可以推广到一般的二阶常系数线性齐次递推数列．尽管具体的常数有所差异，这些数列模 \(m\) 的 Pisano 周期都是 \(O(m)\) 的．

第一个观察是：利用中国剩余定理，可以将讨论限制在素数幂模的情形．设 \(m_1,m_2\) 是两个互素的模数．斐波那契数列在模 \(m_1\) 下的周期是 \(\pi(m_1)\) 及其倍数，在模 \(m_2\) 下的周期是 \(\pi(m_2)\) 及其倍数，所以它在模 \(m_1m_2\) 下的最小正周期则恰为 \(\pi(m_1)\) 和 \(\pi(m_2)\) 的最小公倍数．这就是前文的结论 1．

另一个观察是：模 \(m\) 下的 Pisano 周期，其实是最小的正整数 \(k\)，使得

\[ A^k = \begin{pmatrix} 1&1\\1&0 \end{pmatrix}^k \equiv I \pmod{m}. \]

也就是说，它其实是矩阵 \(A\) 在模 \(m\) 下¹的阶．

对于素数幂模 \(m=p^e\) 的情形，可以通过经典的升幂论证联系到相应的素数模的情形．设 \(k=\pi(p^e)\)，就存在二阶方阵 \(\Lambda\)，使得

\[ A^k = p^e\Lambda + I \]

成立．故而，由二项式定理可知

\[ A^{kp} = (p^e\Lambda + I)^p = I + \sum_{i=1}^p\binom{p}{i}(p^e\Lambda)^i \equiv I\pmod{p^{e+1}}. \]

因此，由阶的性质，有 \(\pi(p^{e+1})\mid kp = p\pi(p^e)\)．对 \(e\) 归纳可知，\(\pi(p^e)\mid p^{e-1}\pi(p)\) 总是成立．

对于素数模 \(p\) 的情形，本文讨论两种证明方式．

利用通项公式利用扩域

一种是利用斐波那契数列的通项公式：

\[ F_n = \dfrac{1}{\sqrt{5}}\left(\dfrac{1+\sqrt{5}}{2}\right)^n - \dfrac{1}{\sqrt{5}}\left(\dfrac{1-\sqrt{5}}{2}\right)^n. \]

将它用二项式定理展开，并消去根式项：

\[ F_n = \dfrac{1}{2^{n-1}}\sum_{i=0}^{\lfloor(n-1)/2\rfloor}\binom{n}{2i+1}5^i. \]

对于 \(p=2\)，这一表达式无法直接取模，但可以验证对应的 Pisano 周期为 \(\pi(2)=3\)．对于 \(p=5\)，有 \(F_n\equiv n\cdot 3^{n-1}\pmod{p}\)，可以直接验证对应的 Pisano 周期为 \(\pi(5)=20\)．对于剩余的奇素模数，可以分为两种情形：

如果 \(p\equiv 1,4\pmod{5}\)，就有

\[ \begin{aligned} F_{p} &\equiv \dfrac{1}{2^{p-1}}\binom{p}{p}5^{(p-1)/2} \equiv 1 \pmod{p},\\ F_{p+1} &\equiv \dfrac{1}{2^p}\left(\binom{p+1}{1} + \binom{p+1}{p}5^{(p-1)/2}\right) \equiv 1 \pmod{p}. \end{aligned} \]

化简过程中，利用了如下结论：由 Lucas 定理，对于 \(0 < k < p\) 都有 \(\dbinom{p}{k}\equiv 0\pmod{p}\)，而对于 \(1 < k < p\) 都有 \(\dbinom{p+1}{k}\equiv 0\pmod{p}\)；由 Fermat 小定理，有 \(2^{p-1}\equiv 5^{p-1}\equiv 1\pmod{p}\)；对于 \(p\equiv 1,4\pmod{5}\)，都有 \(p\) 是模 \(5\) 的二次剩余，利用二次互反律，也有 \(5\) 是模 \(p\) 的二次剩余，故而 \(5^{(p-1)/2} \equiv 1\pmod{p}\)．由此，有 \((F_p,F_{p+1}) \equiv (F_1,F_2) \pmod{p}\)，所以 \((p-1)\) 是模 \(p\) 的一个周期．所以，\(\pi(p)\mid(p-1)\)．
如果 \(p\equiv 2,3\pmod{5}\)，就有

\[ \begin{aligned} F_{2p} &\equiv \dfrac{1}{2^{2p-1}}\binom{2p}{p}5^{(p-1)/2} \equiv -1 \pmod{p},\\ F_{2p+1} &\equiv \dfrac{1}{2^{2p}}\left(\binom{2p+1}{1} + \binom{2p+1}{p}5^{(p-1)/2} + \binom{2p+1}{2p+1}5^p\right) \equiv -1\pmod{p}. \end{aligned} \]

化简过程中，利用了如下结论：由 Lucas 定理，对于 \(0 < k < p\) 和 \(p < k < 2p\) 都有 \(\dbinom{p}{k}\equiv 0\pmod{p}\)，以及 \(\dbinom{2p}{p}\equiv 2\pmod{p}\)，而对于 \(1 < k < p\) 和 \(p + 1 < k < 2p\) 都有 \(\dbinom{p}{k}\equiv 0\pmod{p}\)，以及 \(\dbinom{2p+1}{p}\equiv 2\pmod{p}\)；由 Fermat 小定理，有 \(2^{p-1}\equiv 5^{p-1}\equiv 1\pmod{p}\)；对于 \(p\equiv 2,3\pmod{5}\)，都有 \(p\) 是模 \(5\) 的二次非剩余，利用二次互反律，也有 \(5\) 是模 \(p\) 的二次非剩余，故而 \(5^{(p-1)/2} \equiv -1\pmod{p}\)．由此，有 \((F_{2p},F_{2p+1}) \equiv (F_{-2},F_{-1}) \pmod{p}\)，所以 \(2(p+1)\) 是模 \(p\) 的一个周期．所以，\(\pi(p)\mid 2(p+1)\)．

这就完成了证明．这一方法的局限性在于它高度依赖于斐波那契数列的通项公式，所以较难直接推广到一般的情形．

另一种证明方式则是试图直接计算矩阵 \(A=\begin{pmatrix}1&1\\1&0\end{pmatrix}\) 的阶．它的特征多项式是 \(f(x) = x^2-x-1\)，对应的判别式为 \(\Delta = 5\)．对于模 \(p=5\)，有 \(\Delta\equiv 0\pmod{5}\)，矩阵 \(A\) 有两个相同特征值 \(\lambda=3\)，且不能对角化，需要单独计算．对于模 \(p\equiv 1,4\pmod{5}\)，由二次互反律可知，判别式 \(\Delta=5\) 是模 \(p\) 的二次剩余，矩阵 \(A\) 在域 \(\mathbf F_p\) 内有两个相异特征值 \(\lambda_1\neq\lambda_2\)，矩阵 \(A\) 的阶就是 \(\operatorname{lcm}(\operatorname{ord}(\lambda_1),\operatorname{ord}(\lambda_2))\)，必然整除 \(|\mathbf F_p^\times|=p-1\)．对于模 \(p\equiv 2,3\pmod{5}\)，由二次互反律可知，判别式 \(\Delta=5\) 是模 \(p\) 的二次非剩余，矩阵 \(A\) 在域 \(\mathbf F_p\) 内没有特征值，而只有在扩域 \(\mathbf F_p[\sqrt{5}]\) 内才有两个相异特征值 \(\lambda_1\neq\lambda_2\)，由于 Frobenius 自同态 \(x\mapsto x^p\) 将两根交换，有 \(\lambda_2=\lambda_1^p\)，故而 \(\lambda_1^{p+1}=\lambda_2^{p+1}=\lambda_1\lambda_2=-1\)，亦即 \(\lambda_1^{2(p+1)}=\lambda_2^{2(p+1)}=1\)，由此，矩阵 \(A\) 的阶就是 \(\operatorname{lcm}(\operatorname{ord}(\lambda_1),\operatorname{ord}(\lambda_2))\)，必然整除 \(2(p+1)\)．这就得到了与前种方法一致的结论．

综上，对于不同的情形，相应地有：

\(\pi(2^e)=\dfrac{3}{2}\cdot 2^e,~\dfrac{1}{4}\pi(5^e)=5^e\)．
当 \(p\equiv\pm1\pmod{10}\) 时，\(\pi(p^e) \mid (p-1)p^{e-1}\)，所以 \(\pi(p^e)\le p^e\)．
当 \(p\equiv\pm3\pmod{10}\) 时，\(\dfrac{1}{4}\pi(p^e) \mid \dfrac{p+1}{2}p^{e-1}\)，所以 \(\dfrac{1}{4}\pi(p^e)\le p^e\)．

所以，利用结论 1，对于一般的模数 \(m=\prod_i p_i^{e_i}\)，有

\[ \begin{aligned} \pi(m)&=\operatorname{lcm}\{\pi(p_i^{e_i}):p_i\in\mathbf P\} \\ &\le \operatorname{lcm}\{\pi(p_i^{e_i}):p_i=2\text{ or }p_i\equiv\pm1~(\operatorname{mod}{10})\}\\ &\quad \cdot 4\cdot\operatorname{lcm}\{\pi(p_i^{e_i})/4:p_i=5\text{ or }p_i\equiv\pm3~(\operatorname{mod}{10})\}\\ &\le \prod\{\pi(p_i^{e_i}):p_i=2\text{ or }p_i\equiv\pm1~(\operatorname{mod}{10})\}\\ &\quad \cdot 4\cdot\prod\{\pi(p_i^{e_i})/4:p_i=5\text{ or }p_i\equiv\pm3~(\operatorname{mod}{10})\}\\ &\le \dfrac{3}{2}\cdot\prod\{p_i^{e_i}:p_i=2\text{ or }p_i\equiv\pm1~(\operatorname{mod}{10})\}\\ &\quad \cdot 4\cdot\prod\{p_i^{e_i}:p_i=5\text{ or }p_i\equiv\pm3~(\operatorname{mod}{10})\}\\ &= 6m. \end{aligned} \]

这就说明了斐波那契数列模 \(m\) 的 Pisano 周期总是不超过 \(6m\)，而且等号当且仅当在 \(m=2\cdot 5^e\) 处取得．

习题

参考文献与注释

本页面主要译自博文 Числа Фибоначчи 与其英文翻译版 Fibonacci Numbers．其中俄文版版权协议为 Public Domain + Leave a Link；英文版版权协议为 CC-BY-SA 4.0．内容有改动．

严格来说，它是矩阵 \(A\) 在一般线性群 \(GL_2(\mathbf Z_m)\) 中的阶． ↩

本页面最近更新：，更新历史
发现错误？想一起完善？在 GitHub 上编辑此页！
本页面贡献者：OI-wiki
本页面的全部内容在 CC BY-SA 4.0 和 SATA 协议之条款下提供，附加条款亦可能应用