斐波那契数列

1. 定义与性质

我们先从最基础的定义开始。斐波那契数列 {Fn}\{F_n\} 通常定义如下：

F0=0F_0 = 0 F1=1F_1 = 1 $$F_n = F_{n-1} + F_{n-2} \quad (\text{对于 } n \ge 2) $$

数列的前几项是：0, 1, 1, 2, 3, 5, 8, 13, 21, 34, …\dots

这个定义简洁明了，描述了一个非常自然的增长过程。除了这个核心递推关系，斐波那契数列还有许多有趣的性质，比如：

  • 黄金分割比: 相邻两项的比值 FnFn−1\frac{F_n}{F_{n-1}} 随着 nn 的增大趋近于黄金分割比 ϕ=1+52≈1.618…\phi = \frac{1 + \sqrt{5}}{2} \approx 1.618\dots

  • 通项公式 (Binet 公式):

    Fn=ϕn−ψn5F_n = \frac{\phi^n - \psi^n}{\sqrt{5}}

    其中 ϕ=1+52\phi = \frac{1 + \sqrt{5}}{2} (黄金分割比) 和 ψ=1−52\psi = \frac{1 - \sqrt{5}}{2} 是特征方程 x2−x−1=0x^2 - x - 1 = 0 的两个根。这个公式在理论分析中有用，但在竞赛中直接计算涉及浮点数和精度问题，且效率不高，通常不用于求解具体数值，除非题目有特殊要求或需要推导性质。

  • 组合意义: FnF_n 可以表示用 1x1 和 1x2 的骨牌铺满 1x(n-1) 的地板的方案数，或者一个人爬 n-1 级楼梯，每次可以爬 1 级或 2 级的方案数（若定义 F1=1,F2=2F_1=1, F_2=2 则对应爬 nn 级楼梯）。这种组合意义是识别某些问题本质为斐波那契数列的关键。

  • 其他恒等式: 比如 Cassini 恒等式:

    Fn−1Fn+1−Fn2=(−1)nF_{n-1}F_{n+1} - F_n^2 = (-1)^n

    以及 和恒等式:

    Fm+n=Fm−1Fn+FmFn+1F_{m+n} = F_{m-1}F_n + F_m F_{n+1}

    这些恒等式有时能在特定问题中提供解题思路或简化计算。

我们今天的重点是计算 FnF_n，特别是当 nn 很大或者需要模运算时。

2. 朴素递归

最直观的实现方式就是直接根据定义写递归函数：

#include <bits/stdc++.h>

using namespace std;

// 朴素递归计算斐波那契数
long long fib_rec(int n) {
    // Base cases
    if (n == 0) return 0;
    if (n == 1) return 1;
    // Recursive step based on definition
    // F_n = F_{n-1} + F_{n-2}
    return fib_rec(n - 1) + fib_rec(n - 2);
}

int main() {
    ios::sync_with_stdio(0); cin.tie(0); cout.tie(0);

    int n;
    cin >> n;
    cout << fib_rec(n) << endl;

    return 0;
}
  • 代码解释: 这个 fib_rec 函数完全是定义 Fn=Fn−1+Fn−2F_n = F_{n-1} + F_{n-2} 的直接翻译，基础情况是 n=0n=0 和 n=1n=1。
  • 复杂度分析:
    • 时间复杂度: 呈指数级增长，大致为 O(ϕn)O(\phi^n) 或 O(2n)O(2^n)。为什么？计算 fib_rec(n) 需要计算 fib_rec(n-1) 和 fib_rec(n-2)，而这两者又会进一步递归。画出递归树你会发现，很多子问题（如 fib_rec(n-k)）被重复计算了无数次。对于 n=40n=40 左右，计算量就非常庞大了，在竞赛中几乎一定会超时（TLE）。
    • 空间复杂度: 递归深度为 O(n)O(n)，所以栈空间复杂度是 O(n)O(n)。

这个实现虽然简单，但效率极低，仅适用于 nn 非常小的情况，或者作为理解后续优化的起点。

动态规划优化：记忆化与迭代

递归方法的主要问题在于重复计算。解决重复计算的经典方法就是动态规划 (Dynamic Programming, DP)。对于斐波那契数列，DP 有两种主要形式：记忆化搜索（自顶向下）和迭代（自底向上）。

3. 记忆化搜索 (Memoization)

记忆化搜索本质上还是递归，但增加了一个“备忘录”（通常是数组或哈希表）来存储已经计算过的子问题的结果。每次计算 fib(k) 时，先检查备忘录里是否已有结果，如果有，直接返回；如果没有，才进行计算，并将结果存入备忘录后再返回。

#include <bits/stdc++.h>

using namespace std;

const int MAXN = 100005; // 根据题目 n 的范围调整
long long memo[MAXN];    // 备忘录数组，存储已计算的 F(k)
bool vis[MAXN];      // 标记数组，记录 memo[k] 是否已计算

// 记忆化搜索计算斐波那契数
long long fib_mem(int n) {
    // Base cases
    if (n == 0) return 0;
    if (n == 1) return 1;

    // 检查备忘录
    if (vis[n]) {
        return memo[n];
    }

    // 计算并存入备忘录
    // F_n = F_{n-1} + F_{n-2}
    memo[n] = fib_mem(n - 1) + fib_mem(n - 2);
    vis[n] = true; // 标记已计算

    // 注意：如果需要模运算，加法后要取模
    // memo[n] = (fib_mem(n - 1) + fib_mem(n - 2)) % MOD;

    return memo[n];
}

int main() {
    ios::sync_with_stdio(0); cin.tie(0); cout.tie(0);

    // 初始化备忘录 (可以只初始化 vis，或者用特殊值如 -1 标记 memo)
    // memset(vis, false, sizeof(vis));
    // 或者全局变量默认初始化为 0，可以利用 vis 区分

    int n;
    cin >> n;

    // 如果 n 很大，需要注意 long long 也会溢出，可能需要模运算
    // 如果需要模运算，需要在 fib_mem 函数内部进行
    cout << fib_mem(n) << endl;

    return 0;
}
  • 代码解释: 使用 memo 数组存储 FkF_k 的结果，vis 数组标记 memo[k] 是否有效。计算前先查 vis[n]，如果为 true 则直接返回 memo[n]。否则，递归计算 fib_mem(n-1) 和 fib_mem(n-2)，将结果相加（如果需要取模则在此处取模），存入 memo[n]，标记 vis[n] = true，然后返回结果。
  • 复杂度分析:
    • 时间复杂度: 每个状态 fib(k) (从 0 到 n) 只会被实际计算一次，因为计算后结果就被存储了。每次计算涉及常数次操作（加法、赋值、查表）。所以总时间复杂度是 O(n)O(n)。
    • 空间复杂度: 需要 O(n)O(n) 的额外空间存储备忘录 memo 和标记 vis，另外递归栈深度仍然是 O(n)O(n)。

记忆化搜索保留了递归的直观性，同时通过空间换时间，将复杂度降至线性，对于 nn 在 105∼10610^5 \sim 10^6 范围且不需要模运算（或模数不大，long long 能存下中间结果）的情况是可行的。

4. 迭代实现

迭代方法是自底向上计算斐波那契数。我们知道 F0F_0 和 F1F_1，然后可以计算 F2=F1+F0F_2 = F_1 + F_0，接着计算 F3=F2+F1F_3 = F_2 + F_1，依此类推，直到计算出 FnF_n。

#include <bits/stdc++.h>

using namespace std;

const int MAXN = 100005; // 同样根据 n 的范围调整
long long dp[MAXN];      // DP 数组，dp[i] 存储 F_i

// 迭代计算斐波那契数
long long fib_iter(int n) {
    if (n == 0) return 0;
    if (n == 1) return 1;

    dp[0] = 0;
    dp[1] = 1;

    for (int i = 2; i <= n; ++i) {
        // dp[i] = dp[i-1] + dp[i-2]
        dp[i] = dp[i - 1] + dp[i - 2];
        // 如果需要模运算
        // dp[i] = (dp[i - 1] + dp[i - 2]) % MOD;
    }
    return dp[n];
}

int main() {
    ios::sync_with_stdio(0); cin.tie(0); cout.tie(0);

    int n;
    cin >> n;

    // 同样，考虑 long long 溢出和模运算
    cout << fib_iter(n) << endl;

    return 0;
}
  • 代码解释: 初始化 dp[0] 和 dp[1]。然后通过一个循环，从 i=2i = 2 到 nn，利用 dp[i-1] 和 dp[i-2] 计算 dp[i]。最后返回 dp[n]。
  • 复杂度分析:
    • 时间复杂度: 一个简单的循环，从 2 到 n，每次循环内部是常数时间操作。总时间复杂度 O(n)O(n)。
    • 空间复杂度: 需要一个大小为 O(n)O(n) 的数组 dp 来存储 F0F_0 到 FnF_n。

迭代方法通常比记忆化搜索有更小的常数因子（没有函数调用开销），且空间复杂度分析更直接。

5. 空间优化：滚动数组

观察迭代方法，你会发现计算 dp[i]dp[i] 其实只需要 dp[i−1]dp[i-1] 和 dp[i−2]dp[i-2] 的值。我们并不需要存储从 dp[0]dp[0] 到 dp[i−3]dp[i-3] 的所有历史值。因此，可以用滚动数组的思想来优化空间。我们只需要维护最近的两个斐波那契数即可。

#include <bits/stdc++.h>

using namespace std;

// 迭代计算斐波那契数，空间优化 O(1)
long long fib_opt(int n) {
    if (n == 0) return 0;
    if (n == 1) return 1;

    long long a = 0; // 代表 F_{i-2}
    long long b = 1; // 代表 F_{i-1}
    long long cur;   // 代表 F_i

    for (int i = 2; i <= n; ++i) {
        // cur = F_{i-1} + F_{i-2}
        cur = a + b;
        // 如果需要模运算
        // cur = (a + b) % MOD;

        // 更新 a 和 b，为下一轮计算做准备
        a = b;      // a 更新为 F_{i-1}
        b = cur;    // b 更新为 F_i
    }
    // 循环结束后, b 存储的是 F_n
    return b;
}

int main() {
    ios::sync_with_stdio(0); cin.tie(0); cout.tie(0);

    int n;
    cin >> n;
    cout << fib_opt(n) << endl;

    return 0;
}
  • 代码解释: 使用三个变量 a, b, cur。a 存储 Fi−2F_{i-2}，b 存储 Fi−1F_{i-1}。在循环中计算 cur = a + b (即 FiF_i)，然后更新 a 为原来的 b (即 Fi−1F_{i-1})，b 更新为 cur (即 FiF_i)，这样就滚动到了下一个状态。
  • 复杂度分析:
    • 时间复杂度: 仍然是 O(n)O(n)，因为循环次数不变。
    • 空间复杂度: 只使用了常数个变量 (a, b, cur, i)，所以空间复杂度是 O(1)O(1)。

这种 O(n)O(n) 时间、O(1)O(1) 空间的迭代方法是计算单个斐波那契数（当 nn 不太大时）的标准且最高效的线性方法。对于 nn 在 106∼10710^6 \sim 10^7 级别，这个方法是首选。

矩阵快速幂：处理巨大 n 的利器

当 nn 非常大时，比如 nn 达到 10910^9 甚至 101810^{18}，O(n)O(n) 的时间复杂度就无法接受了。这时，我们需要一个更快的算法。斐波那契数列的递推关系是线性的，这提示我们可以使用矩阵快速幂 (Matrix Exponentiation) 将复杂度降到 O(log⁡n)O(\log n)。

6. 斐波那契数列与矩阵

我们的目标是找到一个 2x2 矩阵 MM，使得它能将状态向量 (FnFn−1)\begin{pmatrix} F_n \\ F_{n-1} \end{pmatrix} 转移到下一个状态向量 (Fn+1Fn)\begin{pmatrix} F_{n+1} \\ F_n \end{pmatrix}。即：

$$\begin{pmatrix} F_{n+1} \\ F_n \end{pmatrix} = \begin{pmatrix} ? & ? \\ ? & ? \end{pmatrix} \begin{pmatrix} F_n \\ F_{n-1} \end{pmatrix} $$

根据斐波那契定义：

  • Fn+1=1⋅Fn+1⋅Fn−1F_{n+1} = 1 \cdot F_n + 1 \cdot F_{n-1}
  • Fn=1⋅Fn+0⋅Fn−1F_n = 1 \cdot F_n + 0 \cdot F_{n-1}

将这个关系写成矩阵形式：

$$\begin{pmatrix} F_{n+1} \\ F_n \end{pmatrix} = \begin{pmatrix} 1 & 1 \\ 1 & 0 \end{pmatrix} \begin{pmatrix} F_n \\ F_{n-1} \end{pmatrix} $$

我们令这个转移矩阵

M=(1110)M = \begin{pmatrix} 1 & 1 \\ 1 & 0 \end{pmatrix}

那么，我们可以进行连续的转移：

$$\begin{pmatrix} F_n \\ F_{n-1} \end{pmatrix} = M \begin{pmatrix} F_{n-1} \\ F_{n-2} \end{pmatrix} = M^2 \begin{pmatrix} F_{n-2} \\ F_{n-3} \end{pmatrix} = \dots = M^{n-1} \begin{pmatrix} F_1 \\ F_0 \end{pmatrix} $$

注意这里的指数是 n−1n-1，因为我们是从 (F1F0)\begin{pmatrix} F_1 \\ F_0 \end{pmatrix} 开始转移 n−1n-1 次得到 (FnFn−1)\begin{pmatrix} F_n \\ F_{n-1} \end{pmatrix}。 由于 F1=1,F0=0F_1=1, F_0=0，我们有：

$$\begin{pmatrix} F_n \\ F_{n-1} \end{pmatrix} = M^{n-1} \begin{pmatrix} 1 \\ 0 \end{pmatrix} $$

如果我们计算出矩阵 $M^{n-1} = \begin{pmatrix} a & b \\ c & d \end{pmatrix}$，那么：

$$\begin{pmatrix} F_n \\ F_{n-1} \end{pmatrix} = \begin{pmatrix} a & b \\ c & d \end{pmatrix} \begin{pmatrix} 1 \\ 0 \end{pmatrix} = \begin{pmatrix} a \cdot 1 + b \cdot 0 \\ c \cdot 1 + d \cdot 0 \end{pmatrix} = \begin{pmatrix} a \\ c \end{pmatrix} $$

所以，FnF_n 就是 Mn−1M^{n-1} 这个矩阵左上角的元素 aa！（这里假设 n≥1n \ge 1。如果 n=0n=0，F0=0F_0=0 可以直接处理）。

  • 一个稍微不同的、更常用的形式: 使用状态向量 (Fn+1Fn)\begin{pmatrix} F_{n+1} \\ F_n \end{pmatrix} 和初始向量 (F1F0)\begin{pmatrix} F_1 \\ F_0 \end{pmatrix}。

    $$\begin{pmatrix} F_2 \\ F_1 \end{pmatrix} = M^1 \begin{pmatrix} F_1 \\ F_0 \end{pmatrix} = \begin{pmatrix} 1 & 1 \\ 1 & 0 \end{pmatrix} \begin{pmatrix} 1 \\ 0 \end{pmatrix} = \begin{pmatrix} 1 \\ 1 \end{pmatrix} $$

    推广可得：

    $$\begin{pmatrix} F_{n+1} \\ F_n \end{pmatrix} = M^n \begin{pmatrix} F_1 \\ F_0 \end{pmatrix} = M^n \begin{pmatrix} 1 \\ 0 \end{pmatrix} $$

    在这种情况下，计算 Mn=(abcd)M^n = \begin{pmatrix} a & b \\ c & d \end{pmatrix}，则：

    $$\begin{pmatrix} F_{n+1} \\ F_n \end{pmatrix} = \begin{pmatrix} a & b \\ c & d \end{pmatrix} \begin{pmatrix} 1 \\ 0 \end{pmatrix} = \begin{pmatrix} a \\ c \end{pmatrix} $$

    所以，FnF_n 是 MnM^n 矩阵的左下角元素 cc。 我们也可以用初始向量 $\begin{pmatrix} F_2 \\ F_1 \end{pmatrix} = \begin{pmatrix} 1 \\ 1 \end{pmatrix}$，然后计算 Mn−1M^{n-1}：

    $$\begin{pmatrix} F_{n+1} \\ F_n \end{pmatrix} = M^{n-1} \begin{pmatrix} F_2 \\ F_1 \end{pmatrix} = M^{n-1} \begin{pmatrix} 1 \\ 1 \end{pmatrix} $$

    如果 $M^{n-1} = \begin{pmatrix} a' & b' \\ c' & d' \end{pmatrix}$，则 Fn=c′⋅1+d′⋅1=c′+d′F_n = c' \cdot 1 + d' \cdot 1 = c' + d'。这个形式稍微复杂一点。

    最常用的结论：计算 $M^n = \begin{pmatrix} 1 & 1 \\ 1 & 0 \end{pmatrix}^n$，则 FnF_n 就是结果矩阵的 [1][0][1][0] 元素（0-based index，即左下角 cc）。 这个形式很方便，因为指数是 nn。

现在的问题转化为：如何快速计算矩阵 MM 的 nn 次幂？

7. 矩阵快速幂算法

计算一个数 xx 的 nn 次幂可以用快速幂 (Binary Exponentiation / Exponentiation by Squaring)，时间复杂度 O(log⁡n)O(\log n)。其原理是利用 nn 的二进制表示：

$$x^n = x^{\sum_{i=0}^k b_i 2^i} = \prod_{i=0}^k (x^{2^i})^{b_i} $$

其中 bi∈{0,1}b_i \in \{0, 1\} 是 nn 的二进制表示的第 ii 位。 我们可以通过反复平方 x→x2→x4→x8→…x \to x^2 \to x^4 \to x^8 \to \dots 来得到 x2ix^{2^i}，然后根据 nn 的二进制位 bib_i 决定是否将对应的 x2ix^{2^i} 乘入最终结果。

这个思想完全适用于矩阵。我们需要定义矩阵乘法，并注意矩阵乘法不满足交换律 A B ≠\ne B A。

#include <bits/stdc++.h>

using namespace std;

typedef long long ll; // 使用 long long 避免溢出

// 定义矩阵结构体
struct Mat {
    ll mat[2][2];
    // 构造函数
    Mat() { memset(mat, 0, sizeof(mat)); }
};

ll MOD = 1e9 + 7; // 假设模数为 10^9 + 7, 根据题目要求修改

// 矩阵乘法 (计算 C = A * B)
Mat mul(Mat a, Mat b) {
    Mat res; // 结果矩阵 C
    for (int i = 0; i < 2; ++i) {
        for (int j = 0; j < 2; ++j) {
            for (int k = 0; k < 2; ++k) {
                // C[i][j] = sum(A[i][k] * B[k][j])
                // 注意取模！
                res.mat[i][j] = (res.mat[i][j] + a.mat[i][k] * b.mat[k][j]) % MOD;
            }
        }
    }
    return res;
}

// 矩阵快速幂 (计算 base^p mod MOD)
Mat mat_pow(Mat base, ll p) {
    Mat res;
    // 初始化 res 为单位矩阵 I
    // I = [[1, 0], [0, 1]]
    res.mat[0][0] = res.mat[1][1] = 1;
    res.mat[0][1] = res.mat[1][0] = 0;

    while (p > 0) {
        // 如果 p 的当前最低位是 1 (p % 2 == 1)
        if (p & 1) {
            res = mul(res, base); // res = res * base^(2^k)
        }
        // base 自乘，准备下一位
        base = mul(base, base); // base = base^(2^k) -> base^(2^(k+1))
        // p 右移一位 (p = p / 2)
        p >>= 1;
    }
    return res;
}

// 使用矩阵快速幂计算 F(n) mod MOD
ll fib_mat(ll n) {
    if (n == 0) return 0;
    if (n == 1) return 1 % MOD; // 注意 n=1 时也要对 MOD 取模

    Mat base;
    // base = [[1, 1], [1, 0]]
    base.mat[0][0] = 1; base.mat[0][1] = 1;
    base.mat[1][0] = 1; base.mat[1][1] = 0;

    // 我们需要计算 base^n (对应 F_n = (M^n)[1][0] 的形式)
    Mat res_mat = mat_pow(base, n); // 计算 M^n mod MOD

    // F(n) 是结果矩阵的 [1][0] (左下角) 元素
    return res_mat.mat[1][0];
    // 如果用 F(n) = (M^(n-1))[0][0] 的形式，则调用 mat_pow(base, n-1) 并返回 res_mat.mat[0][0]
}

int main() {
    ios::sync_with_stdio(0); cin.tie(0); cout.tie(0);

    ll n;
    cin >> n;

    // 设置模数 (如果需要从输入读取)
    // cin >> MOD;

    cout << fib_mat(n) << endl;

    return 0;
}
  • 代码解释:
    • 定义 Mat 结构体存储 2x2 矩阵。
    • mul 函数实现两个 2x2 矩阵的乘法 C=A×BC = A \times B。注意每次加法和乘法后都要取模，即 (res.mat[i][j] + a.mat[i][k] * b.mat[k][j]) % MOD，防止中间结果溢出 long long 以及保证最终结果在模 MOD 下正确。
    • mat_pow 函数实现矩阵的快速幂 basep(modMOD)base^p \pmod{MOD}。res 初始化为单位矩阵（矩阵乘法的幺元）。循环体内，如果 p 的当前二进制位为 1，则将当前的 base (代表 M2kM^{2^k}) 乘入 res；然后 base 自乘变为 M2k+1M^{2^{k+1}}；p 右移一位。
    • fib_mat 函数处理 n=0,n=1n=0, n=1 的边界情况，然后构造基础转移矩阵 base，调用 mat_pow 计算 Mn(modMOD)M^n \pmod{MOD}，最后根据我们选择的递推形式返回结果矩阵的 [1][0] 元素作为 Fn(modMOD)F_n \pmod{MOD}。
  • 复杂度分析:
    • 时间复杂度: 矩阵乘法 mul 的复杂度是 O(k3)O(k^3)，其中 kk 是矩阵的维度。这里 k=2k=2，所以是 O(23)=O(8)O(2^3) = O(8)，即 O(1)O(1) 常数时间。矩阵快速幂 mat_pow 需要进行 O(log⁡p)O(\log p) 次矩阵乘法（pp 是指数，这里是 nn）。因此，总时间复杂度是 O(k3⋅log⁡n)=O(log⁡n)O(k^3 \cdot \log n) = \mathbf{O(\log n)}。
    • 空间复杂度: O(k2)=O(1)O(k^2) = O(1)，只使用了常数个矩阵变量和辅助变量。

矩阵快速幂是解决形如 Fn=a⋅Fn−1+b⋅Fn−2F_n = a \cdot F_{n-1} + b \cdot F_{n-2} 这类常系数线性齐次递推关系在 nn 很大时的标准武器。

模运算下的斐波那契数列

在算法竞赛中，斐波那契数列常常需要对一个模数 MM 取模，即计算 Fn(modM)F_n \pmod{M}。前面矩阵快速幂的代码已经考虑了模运算。但模运算下的斐波那契数列还有一些独特的性质。

8. 斐波那契数列模 M 的周期性 (Pisano Period)

将斐波那契数列的每一项都对 MM 取模，得到一个新的序列 {Fn(modM)}\{F_n \pmod{M}\}。这个序列是周期性的。这个周期被称为皮萨诺周期 (Pisano Period)，记作 π(M)\pi(M)。

例如，模 3： FnF_n: 0, 1, 1, 2, 3, 5, 8, 13, 21, 34, …\dots Fn(mod3)F_n \pmod 3: 0, 1, 1, 2, 0, 2, 2, 1, 0, 1, …\dots 序列是 0, 1, 1, 2, 0, 2, 2, 1, (0, 1, ...) 重复了。周期是从 (0,1)(0, 1) 这个状态对开始，到下一次出现 (0,1)(0, 1) 前结束。这里的周期长度是 8。即 π(3)=8\pi(3) = 8。

关键性质: 对于 n≥0n \ge 0，有

Fn(modM)=Fn(modπ(M))(modM)F_n \pmod{M} = F_{n \pmod{\pi(M)}} \pmod{M}

这意味着，如果我们需要计算 Fn(modM)F_n \pmod{M} 且 nn 非常巨大，我们可以先计算出皮萨诺周期 π(M)\pi(M)，然后计算 n′=n(modπ(M))n' = n \pmod{\pi(M)}，最后只需求 Fn′(modM)F_{n'} \pmod{M} 即可。由于 n′<π(M)n' < \pi(M)，计算 Fn′(modM)F_{n'} \pmod{M} 可以用 O(π(M))O(\pi(M)) 的迭代法或者 O(log⁡π(M))O(\log \pi(M)) 的矩阵快速幂。

如何求 π(M)\pi(M)? 没有已知的简单公式可以直接计算 π(M)\pi(M)。但我们可以通过模拟斐波那契数列模 MM 的生成过程来找到周期： 从 (F0(modM),F1(modM))=(0,1)(F_0 \pmod{M}, F_1 \pmod{M}) = (0, 1) 开始，不断计算下一项 (Fi+1(modM),Fi(modM))(F_{i+1} \pmod{M}, F_i \pmod{M})，直到再次遇到 (0,1)(0, 1) 这个状态对。由于模 MM 下的状态对 (a,b)(a, b) 只有 M×MM \times M 种可能（0≤a,b<M0 \le a, b < M），根据鸽巢原理，这个序列必然会出现重复的状态对，并且一旦 (0,1)(0, 1) 重复，整个序列就开始循环。周期长度 π(M)\pi(M) 满足 π(M)≤M2\pi(M) \le M^2。

#include <bits/stdc++.h>

using namespace std;

typedef long long ll;

// 计算斐波那契数列模 m 的皮萨诺周期 pi(m)
ll get_pisano_period(ll m) {
    if (m <= 1) return 1; // 特殊情况处理：pi(1)=1

    ll a = 0, b = 1;
    ll c; // F_{i+1} mod m
    // 周期从 (F_0 mod m, F_1 mod m) = (0, 1) 开始
    // 我们寻找下一个 (0, 1) 出现的位置
    for (ll i = 0; i < m * m; ++i) { // 周期长度理论上界是 m*m
        c = (a + b) % m;
        a = b;
        b = c;
        // 找到周期模式 (0, 1)
        // 注意 i=0 时 (a,b)=(1,1), i=1 时 (a,b)=(1,2) or (1,0) ...
        // 当 (a, b) 重新变为 (0, 1) 时，表示一个周期结束
        // 此时已经产生了 i+1 对 (从 (F_0,F_1) 到 (F_i, F_{i+1}))
        // 周期长度是 i+1
        if (a == 0 && b == 1) {
            return i + 1;
        }
    }
    // 理论上对于 m > 1，总能在 m*m 步内找到周期
    return m * m; // 实际上不会执行到这里，可以视为错误或未找到
}

// 迭代计算 F(n) mod m, 空间 O(1)
ll fib_mod(ll n, ll m) {
    if (m <= 1) return 0; // F(n) mod 1 = 0
    if (n <= 1) return n % m;

    ll a = 0, b = 1;
    ll cur;
    for (ll i = 2; i <= n; ++i) {
        cur = (a + b) % m;
        a = b;
        b = cur;
    }
    return b;
}

// (需要包含上面定义的 Mat 结构体, mul, mat_pow 函数)
// ll fib_mat(ll n, ll mod) { ... } // 使用矩阵快速幂计算 F(n) mod mod


int main() {
    ios::sync_with_stdio(0); cin.tie(0); cout.tie(0);

    ll n; // 可能非常大
    ll m; // 模数
    cin >> n >> m;

    if (m <= 1) {
        cout << 0 << endl;
        return 0;
    }

    ll period = get_pisano_period(m);
    ll n_reduced = n % period; // 使用周期性缩减 n

    // 计算 F(n_reduced) mod m
    // 方法1：使用 O(n_reduced) = O(period) 的迭代法
    cout << fib_mod(n_reduced, m) << endl;

    // 方法2：如果 period 仍然很大 (e.g., m 很大导致 m^2 很大),
    // 可以用 O(log n_reduced) = O(log period) 的矩阵快速幂
    // 注意：需要确保 MOD 在 fib_mat 中设置为当前的 m
    // MOD = m; // 临时设置全局 MOD 或将 m 作为参数传入 fib_mat
    // cout << fib_mat(n_reduced, m) << endl; // 假设 fib_mat 接受模数参数

    return 0;
}
  • 代码解释:
    • get_pisano_period(m) 函数通过迭代模拟 (Fi(modm),Fi+1(modm))(F_i \pmod m, F_{i+1} \pmod m) 的过程，寻找 (0,1)(0, 1) 的再次出现，返回周期长度。循环上界设为 m×mm \times m 保证能找到。
    • fib_mod(n, m) 是前面 O(1)O(1) 空间迭代计算 Fn(modm)F_n \pmod m 的函数。
    • main 函数先读入 nn 和 mm，处理 m≤1m \le 1 的边界。然后调用 get_pisano_period 得到周期 period。计算 n′=n(modperiod)n' = n \pmod{\text{period}}。最后调用 fib_mod (或者 fib_mat 如果 period\text{period} 仍然很大) 计算 Fn′(modm)F_{n'} \pmod m。
  • 复杂度分析:
    • get_pisano_period(m) 的时间复杂度是 O(π(M))O(\pi(M))，最坏情况 O(M2)O(M^2)。
    • 计算 Fn′(modm)F_{n'} \pmod m：
      • 如果用迭代法 fib_mod，时间复杂度 O(n′)=O(π(M))O(n') = O(\pi(M))，最坏 O(M2)O(M^2)。
      • 如果用矩阵快速幂 fib_mat，时间复杂度 O(log⁡n′)=O(log⁡π(M))O(\log n') = O(\log \pi(M))，最坏 O(log⁡(M2))=O(log⁡M)O(\log(M^2)) = O(\log M)。
    • 总复杂度主要取决于求周期和最后计算这两部分中较大的一个。通常是 O(π(M))O(\pi(M)) 或 O(M2)O(M^2)。
  • 适用场景: 当 nn 极其巨大 (例如 n≈10100n \approx 10^{100}，甚至更大，以字符串形式给出)，而模数 MM 相对较小 (例如 M≤106M \le 10^6) 时，皮萨诺周期是唯一可行的方法。你需要先对大数 nn 取模 π(M)\pi(M) (这需要大数模运算，或者如果 nn 是指数形式 aba^b，可以用费马小定理/欧拉定理结合指数循环节来求 n(modπ(M))n \pmod{\pi(M)})。

关于皮萨诺周期的一些补充:

  • π(M)\pi(M) 存在且 π(M)≤M2\pi(M) \le M^2。
  • 如果 M=pkM = p^k (pp 是素数)，有 π(pk)=pk−1π(p)\pi(p^k) = p^{k-1} \pi(p) 的关系（有一些例外，如 π(2)=3,π(2k)=3⋅2k−1\pi(2)=3, \pi(2^k)=3 \cdot 2^{k-1} for k≥3k \ge 3）。
  • 如果 M=p1a1p2a2…prarM = p_1^{a_1} p_2^{a_2} \dots p_r^{a_r} 是 MM 的标准素数分解，则 $\pi(M) = \text{lcm}(\pi(p_1^{a_1}), \pi(p_2^{a_2}), \dots, \pi(p_r^{a_r}))$。结合中国剩余定理，可以分解问题，但这通常比直接模拟找周期更复杂，竞赛中较少这样考。
  • π(10k)\pi(10^k) 的值有规律，可以预计算或查找 ($\pi(10)=60, \pi(100)=300, \pi(1000)=1500, \pi(10^k) = 15 \cdot 10^{k-1}$ for k≥3k \ge 3)。

斐波那契数列的推广与应用

斐波那契数列不仅仅是它本身，它的思想和解决方法（尤其是矩阵快速幂）可以推广到更广泛的线性递推问题。

9. 广义斐波那契数列与线性递推

形如

$$F_n = c_1 F_{n-1} + c_2 F_{n-2} + \dots + c_k F_{n-k} $$

的 kk 阶常系数线性齐次递推关系，都可以用矩阵快速幂在 O(k3⋅log⁡n)O(k^3 \cdot \log n) 的时间内求解 FnF_n（通常是模 MM 意义下）。

构造 k×kk \times k 转移矩阵: 状态向量通常选为 $\mathbf{v}_n = \begin{pmatrix} F_n \\ F_{n-1} \\ \vdots \\ F_{n-k+1} \end{pmatrix}$。 下一个状态是 $\mathbf{v}_{n+1} = \begin{pmatrix} F_{n+1} \\ F_n \\ \vdots \\ F_{n-k+2} \end{pmatrix}$。 我们需要找到矩阵 MM 使得 vn+1=M⋅vn\mathbf{v}_{n+1} = M \cdot \mathbf{v}_n。

根据递推关系 $F_{n+1} = c_1 F_n + c_2 F_{n-1} + \dots + c_k F_{n-k+1}$，以及 Fn=FnF_n = F_n, Fn−1=Fn−1F_{n-1} = F_{n-1}, …\dots, Fn−k+2=Fn−k+2F_{n-k+2} = F_{n-k+2}，可以构造出转移矩阵 MM：

$$M = \begin{pmatrix} c_1 & c_2 & \dots & c_{k-1} & c_k \\ 1 & 0 & \dots & 0 & 0 \\ 0 & 1 & \dots & 0 & 0 \\ \vdots & \vdots & \ddots & \vdots & \vdots \\ 0 & 0 & \dots & 1 & 0 \end{pmatrix} $$

这是一个 k×kk \times k 的矩阵。第一行是递推系数 c1c_1 到 ckc_k，下面是一个 (k−1)×(k−1)(k-1) \times (k-1) 的单位矩阵，嵌入在左上角，并向右移动一列，右侧和底行补 0。

然后，计算 MpM^p（pp 通常是 nn 或 n−k+1n-k+1 等，取决于初始状态向量的选择和 nn 的定义），再乘以初始状态向量 $\mathbf{v}_{\text{init}} = \begin{pmatrix} F_{k-1} \\ F_{k-2} \\ \vdots \\ F_0 \end{pmatrix}$（或其他合适的初始向量），就能得到目标状态向量 vn=Mp⋅vinit\mathbf{v}_n = M^p \cdot \mathbf{v}_{\text{init}}，从而求出 FnF_n。

带常数项的线性递推: Fn=c1Fn−1+⋯+ckFn−k+CF_n = c_1 F_{n-1} + \dots + c_k F_{n-k} + C 可以通过增广矩阵解决。构造 (k+1)×(k+1)(k+1) \times (k+1) 的矩阵。状态向量增加一个常数项 1： $\mathbf{v}'_n = \begin{pmatrix} F_n \\ F_{n-1} \\ \vdots \\ F_{n-k+1} \\ 1 \end{pmatrix}$。 转移矩阵 M′M' 类似：

$$M' = \begin{pmatrix} c_1 & c_2 & \dots & c_k & C \\ 1 & 0 & \dots & 0 & 0 \\ 0 & 1 & \dots & 0 & 0 \\ \vdots & \vdots & \ddots & \vdots & \vdots \\ 0 & 0 & \dots & 1 & 0 \\ 0 & 0 & \dots & 0 & 1 \end{pmatrix} $$

然后计算 (M′)p(M')^p 乘以初始状态向量 $\mathbf{v}'_{\text{init}} = \begin{pmatrix} F_{k-1} \\ \vdots \\ F_0 \\ 1 \end{pmatrix}$。

例题场景: 很多 DP 问题最后归结为一个线性递推关系，但 nn 的范围很大（如 101810^{18}），就需要用矩阵快速幂加速。例如，某些限制步数的路径计数、特定模式串的构造计数等。

10. 组合计数中的斐波那契

前面提到，FnF_n 可以解释为用 1x1 和 1x2 骨牌铺满 1x(n-1) 地板的方案数（假设 F1=1F_1=1 对应 1x0 地板， F2=1F_2=1 对应 1x1 地板）。

  • 铺 1x(n-1) 地板：
    • 最后一块是 1x1，则前面需要铺满 1x(n-2)，方案数 Fn−1F_{n-1} (对应 F(n−1)F(n-1) 项)
    • 最后一块是 1x2，则前面需要铺满 1x(n-3)，方案数 Fn−2F_{n-2} (对应 F(n−2)F(n-2) 项)
    • 总方案数 = Fn−1+Fn−2F_{n-1} + F_{n-2}。若定义 F1=1,F2=1F_1=1, F_2=1, 则 FnF_n 是铺 1x(n-1) 的方案数。
    • 若定义 f(1)=1f(1)=1 (铺 1x1), f(2)=2f(2)=2 (铺 1x2)，则 f(n)=f(n−1)+f(n−2)f(n) = f(n-1) + f(n-2) 是铺 1x n 地板的方案数。这 f(n)=Fn+1f(n) = F_{n+1} （如果 F0=0,F1=1F_0=0, F_1=1）。要看清题目定义和边界。

类似地，爬楼梯问题：每次爬 1 或 2 级，爬 nn 级楼梯的方案数。 令 f(n)f(n) 为爬 nn 级楼梯的方案数。

  • 最后一步是爬 1 级，之前在 n−1n-1 级，方案数 f(n−1)f(n-1)
  • 最后一步是爬 2 级，之前在 n−2n-2 级，方案数 f(n−2)f(n-2)
  • 总方案数 f(n)=f(n−1)+f(n−2)f(n) = f(n-1) + f(n-2)。
  • 初始条件：f(1)=1f(1)=1 (只能爬1级)，f(2)=2f(2)=2 (可以 1+1 或 2)。
  • 这个序列是 1, 2, 3, 5, 8, ... 即 f(n)=Fn+1f(n) = F_{n+1} (基于 F0=0,F1=1F_0=0, F_1=1)。

识别技巧: 当你发现一个问题的状态转移可以分解为来自前一个状态或前两个（或前 k 个）状态，且转移方式固定（系数是常数）时，就要警惕是不是斐波那契或类似的线性递推。

赛场实战策略与避坑指南

理解了各种方法后，如何在比赛中正确、快速地使用？

11. 如何选择合适的算法？

关键在于分析数据范围 nn 和模数 MM (如果有)。

  • nn 很小 (e.g., n≤90n \le 90 使得 FnF_n 不超过 long long):
    • O(n)O(n) 迭代 + O(1)O(1) 空间优化 (fib_opt) 是最佳选择。速度快，代码简单，不易出错。
    • 如果需要模运算，且 nn 在 106∼10710^6 \sim 10^7 内，这个方法也足够快。
  • nn 较大 (e.g., 105≤n≤10710^5 \le n \le 10^7)，通常需要模运算:
    • O(n)O(n) 迭代 + O(1)O(1) 空间 (fib_opt 带模运算) 仍然是首选。
    • O(n)O(n) DP 数组 (fib_iter 带模运算) 也可以，但空间占用稍大。
    • 记忆化搜索 (fib_mem 带模运算) 也行，但常数可能稍大，且有递归深度限制风险（虽然 C++ 默认栈深通常够用）。
  • nn 非常大 (e.g., n≥109n \ge 10^9, n≤1018n \le 10^{18}): 必须使用 O(log⁡n)O(\log n) 算法。
    • 矩阵快速幂 (fib_mat) 是标准解法。 需要熟练掌握矩阵乘法和快速幂模板，并正确处理模运算。复杂度 O(k3log⁡n)O(k^3 \log n)，对于斐波那契是 O(log⁡n)O(\log n)。
    • 如果题目只要求 Fn(modM)F_n \pmod{M} 且 MM 相对较小 (e.g., M≤106M \le 10^6)，可以考虑皮萨诺周期。先 O(M2)O(M^2) 或 O(π(M))O(\pi(M)) 求周期 π(M)\pi(M)，然后计算 Fn(modπ(M))(modM)F_{n \pmod{\pi(M)}} \pmod{M}。计算 Fn′F_{n'} 时，因为指数 n′=n(modπ(M))n' = n \pmod{\pi(M)} 变小了，可以用 O(π(M))O(\pi(M)) 迭代或 O(log⁡π(M))O(\log \pi(M)) 矩阵快速幂。当 nn 极其巨大 (远超 long long 范围，如以字符串给出) 时，这是唯一方法。
  • 问题本质是 kk 阶线性递推，nn 很大:
    • 构造 k×kk \times k 转移矩阵，使用 O(k3log⁡n)O(k^3 \log n) 的矩阵快速幂。关键是正确构造转移矩阵。

12. 常见陷阱与注意事项 (Pitfalls)

  1. 整数溢出 (Integer Overflow):

    • 斐波那契数列增长很快 (Fn≈ϕn5F_n \approx \frac{\phi^n}{\sqrt{5}})。即使 nn 只有 93 左右，FnF_n 也会超过 64 位整数 (unsigned long long 的最大值约为 1.8×10191.8 \times 10^{19})。
    • 解决方案:
      • 如果题目不需要模运算且 nn 较大，可能需要使用高精度计算 (模拟大数加法)。
      • 如果题目要求模 MM，必须在每次加法和乘法后都取模，尤其是在矩阵乘法中。 res.mat[i][j] = (res.mat[i][j] + a.mat[i][k] * b.mat[k][j]) % MOD; 这一步至关重要，不能只在最后结果取模，否则中间计算会溢出。
      • 确保使用 long long 存储中间结果和最终结果，即使模数 MM 只是 int 范围，因为两个 long long 相乘可能超过 long long，但在取模前如果结果保存在 long long 中是允许的（只要 M2M^2 不超 long long，或者乘法时小心处理）。更安全的做法是确保 a⋅b(modM)a \cdot b \pmod M 不会溢出，例如使用 (a % MOD * b % MOD) % MOD。

  2. 模运算错误 (Modulo Arithmetic Errors):

    • 负数取模: 在 C++ 中，(-a) % M 的结果可能是负数 (例如 -5 % 3 == -2)。如果需要保证结果在 [0,M−1][0, M-1] 范围内，标准做法是 (x % M + M) % M。虽然斐波那契数列非负，但在推广的线性递推或中间计算（如系数为负）中可能遇到负数，需要注意。对于斐波那契本身，只要 MM 是正数，(a % M + b % M) % M 不会产生负数。
    • 忘记取模: 在迭代、矩阵乘法等步骤中漏掉 % MOD。

  3. 边界条件错误 (Off-by-One Errors / Base Cases):

    • F0,F1F_0, F_1 的定义要明确（题目可能会有不同的起始定义 F0=a,F1=bF_0=a, F_1=b）。
    • 矩阵快速幂的指数是 nn 还是 n−1n-1？取决于你的状态转移矩阵和初始向量的选择。务必推导清楚或用小数据验证。例如，用 M=(1110)M = \begin{pmatrix} 1 & 1 \\ 1 & 0 \end{pmatrix}，若求 FnF_n，通常计算 MnM^n 作用于 (F1F0)\begin{pmatrix} F_1 \\ F_0 \end{pmatrix} 或 Mn−1M^{n-1} 作用于 (F2F1)\begin{pmatrix} F_2 \\ F_1 \end{pmatrix}。
    • 循环的起始和结束条件是否正确。比如迭代是从 i=2i=2 开始。

  4. 复杂度选择错误 (Complexity Misjudgment):

    • 看到 n=1018n=10^{18} 仍然使用 O(n)O(n) 算法导致 TLE。
    • 不必要地使用矩阵快速幂：如果 nn 只有 10510^5， O(n)O(n) 迭代比 O(log⁡n)O(\log n) 矩阵快速幂更快（因为矩阵乘法常数 k3=8k^3=8 较大）且更容易写对。

  5. 矩阵构造错误:

    • 对于广义线性递推，转移矩阵的系数 (cic_i)、1 的位置、0 的位置要仔细核对。特别是对于非齐次项 CC 的增广矩阵。

  6. 皮萨诺周期使用不当:

    • 求周期 get_pisano_period 算法写错，或者死循环 (忘记检查 a == 0 && b == 1)。
    • 计算 n′=n(modperiod)n' = n \pmod{\text{period}} 时，如果 nn 是超大数 (字符串)，需要用大数模运算。如果 n=abn = a^b 形式，可能需要指数循环节和欧拉定理/费马小定理。

总结

斐波那契数列，这个古老而优美的数学序列，在现代算法竞赛中依然扮演着重要角色。它不仅是一个常见的考点，更是学习和掌握以下核心算法思想的绝佳载体：

  • 递归与分治: 理解递归思想和子问题重复计算的弊端。
  • 动态规划: 学习记忆化搜索和迭代两种 DP 实现方式，以及状态压缩（滚动数组）优化空间。
  • 矩阵快速幂: 掌握处理线性递推在指数级 nn 下的高效计算方法 (O(k3log⁡n)O(k^3 \log n))，这是 ICPC 竞赛中的必备技能。
  • 模运算与数论: 理解模运算规则和周期性（皮萨诺周期 π(M)\pi(M)）在处理大数问题中的应用。
  • 组合意义: 培养从组合问题中识别斐波那契（或类似递推）模式的能力。

经典例题

P4994 终于结束的起点

题目描述

广为人知的斐波拉契数列 fib(n)\mathrm{fib}(n) 是这么计算的

$$\mathrm{fib}(n)=\begin{cases} 0,& n=0 \\ 1,& n=1 \\ \mathrm{fib}(n-1) + \mathrm{fib}(n-2),& n>1 \end{cases} $$

也就是 0,1,1,2,3,5,8,13⋯0, 1, 1, 2, 3, 5, 8, 13 \cdots，每一项都是前两项之和。

小 F 发现，如果把斐波拉契数列的每一项对任意大于 11 的正整数 MM 取模的时候，数列都会产生循环。

当然，小 F 很快就明白了，因为 (fib(n−1) mod M\mathrm{fib}(n - 1) \bmod M) 和 (fib(n−2) mod M)\mathrm{fib}(n - 2) \bmod M) 最多只有 M2M ^ 2 种取值，所以在 M2M ^ 2 次计算后一定出现过循环。

甚至更一般地，我们可以证明，无论取什么模数 MM，最终模 MM 下的斐波拉契数列都会是 0,1,⋯ ,0,1,⋯0, 1, \cdots, 0, 1, \cdots。

现在，给你一个模数 MM，请你求出最小的 n>0n > 0，使得 $\mathrm{fib}(n) \bmod M = 0, \mathrm{fib}(n + 1) \bmod M = 1$。

输入格式

输入一行一个正整数 MM。

输出格式

输出一行一个正整数 nn。

输入输出样例 #1

输入 #1

2

输出 #1

3

输入输出样例 #2

输入 #2

6

输出 #2

24

说明/提示

样例 1 解释

斐波拉契数列为 0,1,1,2,3,5,8,13,21,34,⋯0, 1, 1, 2, 3, 5, 8, 13, 21, 34, \cdots，在对 22 取模后结果为 0,1,1,0,1,1,0,1,1,0,⋯0, 1, 1, 0, 1, 1, 0, 1, 1, 0, \cdots。

我们可以发现，当 n=3n = 3 时，f(n) mod 2=0,f(n+1) mod 2=1f(n) \bmod 2= 0, f(n + 1) \bmod 2 = 1，也就是我们要求的 nn 的最小值。

数据范围

对于 30%30\% 的数据，M≤18M \leq 18；

对于 70%70\% 的数据，M≤2018M \leq 2018；

对于 100%100\% 的数据，2≤M≤706150=0xAC6662 \leq M \leq 706150=\verb!0xAC666!。

提示

如果你还不知道什么是取模 ( mod )(\bmod)，那我也很乐意告诉你，模运算是求整数除法得到的余数，也就是竖式除法最终「除不尽」的部分，也即

$$a \bmod M =k \iff a = bM + k\ (M > 0, 0 \leq k < M) $$

其中 a,b,ka, b, k 都是非负整数。

如果你使用 C / C++，你可以使用 % 来进行模运算。

如果你使用 Pascal，你可以使用 mod 来进行模运算。

参考代码

#include<bits/stdc++.h>
using namespace std;

int main() {
    ios::sync_with_stdio(0);
    cin.tie(0); cout.tie(0);

    int M;
    cin >> M;

    // 斐波那契数列初始化
    int a = 0, b = 1;
    int n = 0;

    // 从第2项开始计算斐波那契数列
    while (true) {
        n++;
        int c = (a + b) % M;  // 计算下一项
        a = b;  // 更新a
        b = c;  // 更新b
        
        // 检查是否满足条件：fib(n) % M == 0 和 fib(n+1) % M == 1
        if (a == 0 && b == 1) {
            cout << n << endl;
            break;
        }
    }

    return 0;
}

P4000 斐波那契数列

题目描述

大家都知道，斐波那契数列是满足如下性质的一个数列：

• f1=1f_1 = 1

• f2=1f_2 = 1

• fn=fn−1+fn−2f_n = f_{n-1} + f_{n-2} (n≥2n \geq 2 且 nn 为整数)

请你求出 fnmod  pf_n \mod p 的值。

输入格式

  • 第 1 行：一个整数 nn。

  • 第 2 行：一个整数 pp。

输出格式

  • 第 1 行：fnmod  pf_n \mod p 的值。

输入输出样例 #1

输入 #1

5
1000000007

输出 #1

5

输入输出样例 #2

输入 #2

10
1000000007

输出 #2

55

说明/提示

对于 100%100\% 的数据，n≤1030000000,p<231n \leq 10^{30000000}, p<2^{31}

参考代码

#include <bits/stdc++.h> 
using namespace std;

// 类型别名，方便书写
#define ll long long          // 使用 ll 作为 long long 整数的简写
#define ull unsigned long long // 使用 ull 作为 unsigned long long 整数的简写

// 全局变量和结构体
unordered_map < ull , ll > circ; // 哈希表，用于检测斐波那契数列模 MOD 的循环节（皮萨诺周期）
                                 // 使用类似 Pollard's Rho 的方法进行循环查找
                                 // 键 (Key): 状态表示，将 (F_k, F_{k+1}) 两个模 MOD 的值合并成一个 ull
                                 // 值 (Value): 首次遇到这个状态时的步数 'k'
ll len;                          // 存储检测到的循环节长度（皮萨诺周期）
int MOD;                         // 输入的模数 p
int MX = 1 << 18;                // 预计算块的大小 (2^18 = 262144)。用于类似 "Baby-step Giant-step" 的分块优化
mt19937_64 rnd(time(0));         // 64 位梅森旋转伪随机数生成器，使用当前时间作为种子

// 定义 2x2 矩阵结构体，用于斐波那契数列计算
struct matrix{
	ll arr[2][2]; // 存储矩阵元素的 2x2 数组

	// 默认构造函数：将矩阵初始化为零矩阵
	matrix(){memset(arr , 0 , sizeof(arr));}

	// 重载 [] 运算符，方便访问矩阵元素 (例如 matrix[row][col])
	ll* operator [](int x){return arr[x];}

	// 重载 * 运算符，实现矩阵乘法 (结果对 MOD 取模)
	friend matrix operator *(matrix p , matrix q){
		matrix x; // 结果矩阵，已由构造函数初始化为零
		// 标准矩阵乘法循环
		for(int i = 0 ; i < 2 ; ++i)
			for(int j = 0 ; j < 2 ; ++j)
				for(int k = 0 ; k < 2 ; ++k)
					// 累加乘积项
					// 注意：这里没有在循环内部取模，依赖于最后的结果取模
					// 这假设 MOD * MOD * 2 不会溢出 long long，对于 MOD < 2^31 是安全的
					x[i][k] += p[i][j] * q[j][k];

		// 将结果矩阵的每个元素对 MOD 取模
		for(int i = 0 ; i < 2 ; ++i)
			for(int j = 0 ; j < 2 ; ++j) x[i][j] %= MOD;

		return x; // 返回结果矩阵
	}
}G , T[2][1 << 18 | 1]; // G: 斐波那契转移矩阵 [[1, 1], [1, 0]] (虽然代码中没显式定义 G，但 T[0][1] 就是它)
                       // T: 预计算矩阵幂的数组
                       // T[0][i] = M^i mod MOD (M 是斐波那契转移矩阵)
                       // T[1][i] = (M^MX)^i mod MOD (用于大步计算)

signed main(){
	// 使用静态数组读取超长数字 n，避免动态分配的开销
	// 数组大小 300000003 足够容纳 10^30000000
	static char str[300000003];
	// 使用 scanf 读取，可能比 cin 更快
	scanf("%s %d" , str + 1 , &MOD); // 读取 n (存入 str+1 开始) 和模数 MOD

	// 初始化预计算矩阵数组 T
	// T[0][0] 和 T[1][0] 设置为单位矩阵 I = [[1, 0], [0, 1]]
	T[0][0][0][0] = T[0][0][1][1] = T[1][0][0][0] = T[1][0][1][1] = 1;

	// 初始化斐波那契转移矩阵 M = [[1, 1], [1, 0]]，并存储在 T[0][1] 中
	// T[0][1] 对应 M^1
	T[0][1][0][1] = T[0][1][1][0] = T[0][1][1][1] = 1;

	// 预计算 M^i mod MOD (小步)
	// T[0][i] = T[0][i-1] * T[0][1] = M^(i-1) * M = M^i
	for(int i = 2 ; i <= MX ; ++i) T[0][i] = T[0][i - 1] * T[0][1];

	// 预计算 (M^MX)^i mod MOD (大步)
	// T[1][1] 存储 M^MX
	T[1][1] = T[0][MX];
	// T[1][i] = T[1][i-1] * T[1][1] = (M^MX)^(i-1) * (M^MX) = (M^MX)^i
	for(int i = 2 ; i <= MX ; ++i) T[1][i] = T[1][i - 1] * T[1][1];

	// 使用类似 Pollard's Rho 的随机化方法查找循环节长度 (皮萨诺周期)
	while(1){
		// 生成一个随机步数 x (只取低 46 位，因为 18+18+10 < 64，但这里用了奇怪的位移 >> 28)
        // x 的范围理论上很大，但实际只需要在循环节内即可
        // (rnd() << 28 >> 28) 是一种有点奇怪的取低位随机数的方式，可能只是为了获得一个不太大的随机数
		ll x = (rnd() << 28 >> 28); // 取一个随机步数 k
		// 计算 M^x mod MOD
		// 通过预计算数组快速计算 M^x = M^(x % MX) * (M^MX)^(x / MX)
		// x & (MX - 1)  等价于 x % MX
		// x >> 18       等价于 x / MX (因为 MX = 2^18)
		matrix C = T[0][x & (MX - 1)] * T[1][x >> 18];

		// 将矩阵 C 的第一行 (F_{x+1}, F_x) 编码成一个 ull 作为状态标识符
		// val = (F_{x+1} << 32) | F_x
		ull val = ((1ull * C[0][0]) << 32) | C[0][1]; // C[0][0]=F_{x+1}, C[0][1]=F_x

		// 检查哈希表中是否已存在该状态
		if(circ.find(val) != circ.end()){
			// 如果找到，说明遇到了循环
			// 循环节长度 len = 当前步数 x - 第一次遇到该状态的步数 circ[val]
			len = abs(circ[val] - x); // 取绝对值，防止 x 恰好等于之前的 circ[val] (虽然概率极低)
			break; // 找到循环节，退出循环
		}
		// 如果没找到，将当前状态和步数存入哈希表
		circ[val] = x;
	} // 循环直到找到循环节

	// 计算 n mod len (n 是输入的大数，len 是循环节长度)
	ll sum = 0; // 用于存储 n mod len 的结果
	// 逐位读取大数 n 的字符串表示，计算模 len 的值
	for(int i = 1 ; str[i] ; ++i) sum = (sum * 10 + str[i] - '0') % len;

	// 计算最终结果 F_n mod MOD = F_{n mod len} mod MOD
	// 使用预计算的 T 数组快速计算 M^sum mod MOD
	// M^sum = M^(sum % MX) * (M^MX)^(sum / MX)
	// (T[0][sum & (MX - 1)] * T[1][sum >> 18]) 计算出 M^sum
	// 结果矩阵的第一行是 (F_{sum+1}, F_{sum})
	// 我们需要 F_sum，它在结果矩阵的 [0][1] 位置 (或 [1][0] 位置)
	cout << (T[0][sum & (MX - 1)] * T[1][sum >> 18])[0][1]; // 输出 F_sum mod MOD

	return 0; // 程序正常结束
}

P1962 斐波那契数列

题目描述

大家都知道，斐波那契数列是满足如下性质的一个数列：

$$F_n = \left\{\begin{aligned} 1 \space (n \le 2) \\ F_{n-1}+F_{n-2} \space (n\ge 3) \end{aligned}\right. $$

请你求出 Fn mod 109+7F_n \bmod 10^9 + 7 的值。

输入格式

一行一个正整数 nn

输出格式

输出一行一个整数表示答案。

输入输出样例 #1

输入 #1

5

输出 #1

5

输入输出样例 #2

输入 #2

10

输出 #2

55

说明/提示

【数据范围】
对于 60%60\% 的数据，1≤n≤921\le n \le 92；
对于 100%100\% 的数据，1≤n<2631\le n < 2^{63}。

参考代码

#include <bits/stdc++.h> 

using namespace std;
typedef long long ll;

// 定义模数常量，用于结果取模
const ll MOD = 1e9 + 7; // 即 1,000,000,007

// 定义 2x2 矩阵结构体
struct Matrix {
    ll m[2][2]; // 存储矩阵元素的二维数组
};

// 矩阵乘法函数，计算 a * b 并返回结果矩阵 res
// 参数 a 和 b 是常量引用，避免不必要的拷贝，提高效率
Matrix multiply(const Matrix& a, const Matrix& b) {
    Matrix res; // 用于存储结果的矩阵
    // 遍历结果矩阵的行
    for(int i = 0; i < 2; ++i) {
        // 遍历结果矩阵的列
        for(int j = 0; j < 2; ++j) {
            // 初始化结果矩阵当前元素为 0
            res.m[i][j] = 0;
            // 根据矩阵乘法定义计算 res.m[i][j] = sum(a.m[i][k] * b.m[k][j])
            for(int k = 0; k < 2; ++k) {
                // 累加乘积，注意每次乘法后取模，加法后也要取模，防止溢出
                res.m[i][j] = (res.m[i][j] + a.m[i][k] * b.m[k][j] % MOD) % MOD;
            }
        }
    }
    // 返回计算结果矩阵
    return res;
}

// 矩阵快速幂函数，计算 base 矩阵的 n 次幂 (base^n)
// 参数 base 是待求幂的矩阵，n 是指数
Matrix fast_pow(Matrix base, ll n) {
    Matrix res; // 用于存储最终结果的矩阵
    // 初始化 res 为单位矩阵 I = [[1, 0], [0, 1]]
    // 单位矩阵是矩阵乘法的幺元，相当于整数乘法中的 1
    res.m[0][0] = res.m[1][1] = 1;
    res.m[0][1] = res.m[1][0] = 0;

    // 使用二进制快速幂算法
    while(n > 0) { // 当指数 n 大于 0 时循环
        // 检查 n 的二进制表示的最低位是否为 1
        if(n & 1) {
            // 如果最低位为 1，说明当前幂次的 base 需要乘入结果 res
            res = multiply(res, base); // res = res * base^(2^k)
        }
        // base 自乘，为计算下一位做准备
        base = multiply(base, base); // base = base^(2^k) -> base^(2^(k+1))
        // 将 n 右移一位（相当于 n = n / 2），处理下一位二进制位
        n >>= 1;
    }
    // 返回最终结果矩阵 res = base^n
    return res;
}

// 计算斐波那契数列第 n 项的函数 (模 MOD)
// 注意：此实现基于 F_0=0, F_1=1, F_2=1, F_3=2...
ll Fibonacci(ll n) {
    // 定义斐波那契递推关系的转移矩阵
    // M = [[1, 1],
    //      [1, 0]]
    // 使得 [[F_{k+1}], [F_k]] = M * [[F_k], [F_{k-1}]]
    Matrix base;
    base.m[0][0] = 1; // M[0][0]
    base.m[0][1] = 1; // M[0][1]
    base.m[1][0] = 1; // M[1][0]
    base.m[1][1] = 0; // M[1][1]

    // 处理特殊/边界情况
    if (n == 0) return 0; // F_0 = 0
    if (n == 1) return 1; // F_1 = 1
    // n=2 的情况可以不特殊处理，由下面的公式计算

    // 根据矩阵乘法性质：
    // [[F_n], [F_{n-1}]] = M^(n-1) * [[F_1], [F_0]]
    // [[F_n], [F_{n-1}]] = M^(n-1) * [[1], [0]]
    // 设 M^(n-1) = [[a, b], [c, d]]
    // 则 F_n = a * 1 + b * 0 = a
    // 所以 F_n 是矩阵 M^(n-1) 的左上角元素 (m[0][0])
    // 注意：n >= 1
    return fast_pow(base, n - 1).m[0][0]; // 返回 M^(n-1) 的 [0][0] 元素
}

// 主函数
int main() {
    ll n; // 定义长整型变量 n 用于存储输入的斐波那契项数
    cin >> n; // 从标准输入读取 n
    // 调用 Fibonacci 函数计算第 n 项斐波那契数 (模 MOD) 并输出
    cout << Fibonacci(n) << endl;
    return 0; // 程序正常结束
}

P1939 矩阵加速（数列）

题目描述

已知一个数列 aa，它满足：

$$a_x= \begin{cases} 1 & x \in\{1,2,3\}\\ a_{x-1}+a_{x-3} & x \geq 4 \end{cases} $$

求 aa 数列的第 nn 项对 109+710^9+7 取余的值。

输入格式

第一行一个整数 TT，表示询问个数。

以下 TT 行，每行一个正整数 nn。

输出格式

每行输出一个非负整数表示答案。

输入输出样例 #1

输入 #1

3
6
8
10

输出 #1

4
9
19

说明/提示

  • 对于 30%30\% 的数据 n≤100n \leq 100；
  • 对于 60%60\% 的数据 n≤2×107n \leq2 \times 10^7；
  • 对于 100%100\% 的数据 1≤T≤1001 \leq T \leq 100，1≤n≤2×1091 \leq n \leq 2 \times 10^9。

参考代码

#include <bits/stdc++.h> 

using namespace std;

typedef long long ll;

// 定义模数常量
const ll MOD = 1e9 + 7;
// 定义矩阵维度
const int DIM = 3;

// 定义 DIM x DIM 矩阵结构体
struct Matrix {
    ll m[DIM][DIM];

    // 构造函数，默认初始化为零矩阵
    Matrix() {
        for (int i = 0; i < DIM; ++i) {
            for (int j = 0; j < DIM; ++j) {
                m[i][j] = 0;
            }
        }
    }
};

// 矩阵乘法函数 (计算 C = A * B), 结果模 MOD
Matrix multiply(const Matrix& a, const Matrix& b) {
    Matrix res; // 结果矩阵 C
    for (int i = 0; i < DIM; ++i) {
        for (int j = 0; j < DIM; ++j) {
            // 计算 C[i][j] = sum(A[i][k] * B[k][j]) for k from 0 to DIM-1
            for (int k = 0; k < DIM; ++k) {
                // 累加乘积，注意使用 long long 避免中间溢出，并在每次加法后取模
                res.m[i][j] = (res.m[i][j] + (a.m[i][k] * b.m[k][j]) % MOD) % MOD;
            }
        }
    }
    return res; // 返回结果矩阵
}

// 矩阵快速幂函数 (计算 base^p mod MOD)
Matrix fast_pow(Matrix base, ll p) {
    Matrix res;
    // 初始化 res 为单位矩阵
    for (int i = 0; i < DIM; ++i) {
        res.m[i][i] = 1;
    }

    // 二进制快速幂
    while (p > 0) {
        // 如果 p 的当前最低位是 1
        if (p & 1) {
            res = multiply(res, base); // res = res * base^(2^k)
        }
        // base 自乘，准备下一位
        base = multiply(base, base); // base = base^(2^k) -> base^(2^(k+1))
        // p 右移一位
        p >>= 1;
    }
    return res; // 返回 base^p
}

int main() {
    // 输入输出优化 (可选, 在数据量大或 T 大时有用)
    ios::sync_with_stdio(0);
    cin.tie(0);
    cout.tie(0);

    int T; // 测试用例数量
    cin >> T;
    while (T--) {
        ll n; // 数列项数
        cin >> n;

        // 处理 n <= 3 的基本情况
        if (n <= 3) {
            cout << 1 << "\n"; // a_1, a_2, a_3 都是 1
            continue; // 处理下一个测试用例
        }

        // 定义转移矩阵 M
        Matrix M;
        M.m[0][0] = 1; M.m[0][1] = 0; M.m[0][2] = 1;
        M.m[1][0] = 1; M.m[1][1] = 0; M.m[1][2] = 0;
        M.m[2][0] = 0; M.m[2][1] = 1; M.m[2][2] = 0;

        // 计算 M^(n-3)
        // 指数是 n-3，因为我们从 v_3 状态转移 n-3 次得到 v_n 状态
        Matrix Mn_3 = fast_pow(M, n - 3);

        // 计算最终结果 a_n
        // a_n = M^(n-3)[0][0]*a_3 + M^(n-3)[0][1]*a_2 + M^(n-3)[0][2]*a_1
        // 因为 a_1 = a_2 = a_3 = 1，所以 a_n = M^(n-3)[0][0] + M^(n-3)[0][1] + M^(n-3)[0][2]
        ll ans = (Mn_3.m[0][0] + Mn_3.m[0][1] + Mn_3.m[0][2]) % MOD;

        // 确保结果非负 (虽然在此特定问题中系数和初始值非负，不太可能出现负数，但好习惯)
        ans = (ans % MOD + MOD) % MOD;

        cout << ans << "\n"; // 输出结果
    }
    return 0; 
}
  1. 评测队列
  2. 服务状态
  3. 兼容模式
  4. 主题
  2. 陕ICP备2025061476号
  3. Worker 0, 252ms
  4. Powered Dong Dong