题解 | #斐波那契公约数#

题目链接

PEEK151 斐波那契公约数

题目描述

给定经典的斐波那契数列 $f_i$ ，其定义如下： $f_i = \begin{cases}1,& 1 \le i \le 2; \\ f_{i-1}+f_{i-2},& i > 2. \end{cases}$

现给定两个正整数 $n$ 与 $m$ （ $1 \le n, m \le 10^9$ ），请你计算 $\gcd(f_n, f_m)$ 并对 $10^9 + 7$ 取模后输出结果。

解题思路

本题要求计算两个斐波那契数 $f_n$ 和 $f_m$ 的最大公约数（GCD）。由于 $n$ 和 $m$ 可能非常大，直接计算 $f_n$ 和 $f_m$ 是不可行的。

1. 关键性质

解决本题的核心在于利用斐波那契数列的一个重要数论性质：

\gcd(f_n, f_m) = f_{\gcd(n, m)}

这个定理指出，两个斐波那契数的最大公约数，等于它们下标的最大公约数所对应的那个斐波那契数。

这个优美的性质将问题极大地简化了。我们不再需要处理两个庞大的斐波那契数，而是可以先计算出下标的 $\gcd$ ，然后再求对应的斐波那契数。

2. 问题转化

根据上述性质，原问题 $\gcd(f_n, f_m) \pmod{10^9+7}$ 就等价于计算 $f_{\gcd(n, m)} \pmod{10^9+7}$ 。设 $g = \gcd(n, m)$ ，我们的新目标是计算 $f_g \pmod{10^9+7}$ 。

首先，我们可以用欧几里得算法（辗转相除法）快速地计算出 $g = \gcd(n, m)$ ，其时间复杂度为 $\mathcal{O}(\log(\min(n, m)))$ 。

3. 求解 $f_g$

计算出 $g$ 之后，问题就变成了“求解斐波那契数列第 $g$ 项模 $10^9+7$ 的值”。由于 $g$ 的值仍然可能非常大，因此需要使用矩阵快速幂算法。

这与 PEEK149 斐波那契数列 的解法完全一致：

构造转移矩阵: 斐波那契数列的递推关系可以表示为矩阵形式： $\begin{pmatrix} f_i \\ f_{i-1} \end{pmatrix} = \begin{pmatrix} 1 & 1 \\ 1 & 0 \end{pmatrix} \cdot \begin{pmatrix} f_{i-1} \\ f_{i-2} \end{pmatrix}$
矩阵快速幂: 我们可以通过计算转移矩阵的 $g-2$ 次幂，来从初始状态 $(f_2, f_1)$ 递推到 $(f_g, f_{g-1})$ ： $\begin{pmatrix} f_g \\ f_{g-1} \end{pmatrix} = \begin{pmatrix} 1 & 1 \\ 1 & 0 \end{pmatrix}^{g-2} \cdot \begin{pmatrix} f_2 \\ f_1 \end{pmatrix}$ 其中 $f_1=1, f_2=1$ 。这个幂次计算可以在 $\mathcal{O}(\log g)$ 时间内完成。

4. 算法流程总结

读入 $n, m$ 。
使用欧几里得算法计算 $g = \gcd(n, m)$ 。
如果 $g \le 2$ ，结果为 $1$ 。
否则，使用矩阵快速幂计算 $\begin{pmatrix} 1 & 1 \\ 1 & 0 \end{pmatrix}^{g-2} \pmod{10^9+7}$ ，得到结果矩阵 $B$ 。
最终答案为 $(B_{00} \cdot f_2 + B_{01} \cdot f_1) \pmod{10^9+7} = (B_{00} + B_{01}) \pmod{10^9+7}$ 。

代码

cpp
java
python

#include <iostream>
#include <vector>
#include <numeric>

using namespace std;

const long long MOD = 1000000007;

struct Matrix {
    long long mat[2][2];
    Matrix() {
        mat[0][0] = mat[0][1] = mat[1][0] = mat[1][1] = 0;
    }
};

Matrix multiply(Matrix A, Matrix B) {
    Matrix C;
    for (int i = 0; i < 2; i++) {
        for (int j = 0; j < 2; j++) {
            for (int k = 0; k < 2; k++) {
                C.mat[i][j] = (C.mat[i][j] + A.mat[i][k] * B.mat[k][j]) % MOD;
            }
        }
    }
    return C;
}

Matrix power(Matrix A, long long exp) {
    Matrix result;
    result.mat[0][0] = result.mat[1][1] = 1; // 单位矩阵
    while (exp > 0) {
        if (exp % 2 == 1) result = multiply(result, A);
        A = multiply(A, A);
        exp /= 2;
    }
    return result;
}

long long gcd(long long a, long long b) {
    return b == 0 ? a : gcd(b, a % b);
}

int main() {
    long long n, m;
    while (cin >> n >> m) {
        long long g = gcd(n, m);

        if (g <= 2) {
            cout << 1 << endl;
            continue;
        }

        Matrix T;
        T.mat[0][0] = 1; T.mat[0][1] = 1;
        T.mat[1][0] = 1; T.mat[1][1] = 0;

        T = power(T, g - 2);

        long long ans = (T.mat[0][0] * 1 + T.mat[0][1] * 1) % MOD;
        cout << ans << endl;
    }
    return 0;
}

import java.util.Scanner;

public class Main {
    static final long MOD = 1000000007;

    static long[][] multiply(long[][] a, long[][] b) {
        long[][] c = new long[2][2];
        for (int i = 0; i < 2; i++) {
            for (int j = 0; j < 2; j++) {
                for (int k = 0; k < 2; k++) {
                    c[i][j] = (c[i][j] + a[i][k] * b[k][j]) % MOD;
                }
            }
        }
        return c;
    }

    static long[][] power(long[][] a, long exp) {
        long[][] res = {{1, 0}, {0, 1}}; // 单位矩阵
        while (exp > 0) {
            if ((exp & 1) == 1) res = multiply(res, a);
            a = multiply(a, a);
            exp >>= 1;
        }
        return res;
    }

    static long gcd(long a, long b) {
        return b == 0 ? a : gcd(b, a % b);
    }

    public static void main(String[] args) {
        Scanner sc = new Scanner(System.in);
        while (sc.hasNextLong()) {
            long n = sc.nextLong();
            long m = sc.nextLong();

            long g = gcd(n, m);

            if (g <= 2) {
                System.out.println(1);
                continue;
            }

            long[][] T = {{1, 1}, {1, 0}};
            T = power(T, g - 2);

            long ans = (T[0][0] * 1 + T[0][1] * 1) % MOD;
            System.out.println(ans);
        }
    }
}

import sys

MOD = 1000000007

def gcd(a, b):
    return a if b == 0 else gcd(b, a % b)

def multiply(a, b):
    c = [[0, 0], [0, 0]]
    for i in range(2):
        for j in range(2):
            for k in range(2):
                c[i][j] = (c[i][j] + a[i][k] * b[k][j]) % MOD
    return c

def power(a, exp):
    res = [[1, 0], [0, 1]]  # 单位矩阵
    while exp > 0:
        if exp % 2 == 1:
            res = multiply(res, a)
        a = multiply(a, a)
        exp //= 2
    return res

for line in sys.stdin:
    if not line.strip():
        continue
    
    n, m = map(int, line.split())
    
    g = gcd(n, m)
    
    if g <= 2:
        print(1)
        continue
        
    T = [[1, 1], [1, 0]]
    T = power(T, g - 2)
    
    ans = (T[0][0] * 1 + T[0][1] * 1) % MOD
    print(ans)

算法及复杂度

算法: 欧几里得算法 + 矩阵快速幂
时间复杂度: 计算 $\gcd(n, m)$ 的时间复杂度为 $\mathcal{O}(\log(\min(n, m)))$ 。计算 $f_g$ 的矩阵快速幂时间复杂度为 $\mathcal{O}(\log g)$ 。因此，总时间复杂度为 $\mathcal{O}(\log n + \log m)$ 。
空间复杂度: 算法只需要常数个变量和矩阵来存储中间结果，因此空间复杂度为 $\mathcal{O}(1)$ 。