【高科技】基于值域预处理的快速 GCD

· · 题解

功能简述

可在 \mathcal{O}(v)-\mathcal{O}(1) 的时间内多次询问 \gcd(n,m)(1\le n,m\le v)。

算法过程

证明： 采用数学归纳法，易得 $n=1$ 时成立。 对于正整数 $n$，设 $p$ 为其最小质因数，$a',b',c'$ 为 $\dfrac{n}{p}$ 的一合法分解（$a'\le b'\le c'$），则： - 若 $p\le\sqrt[4]{n}$，则因为 $a'\le\sqrt[3]{\dfrac{n}{p}}$，所以 $pa'\le p\sqrt[3]{\dfrac{n}{p}}=\sqrt[3]{np^2}\le\sqrt{n}$，所以 $pa',b',c'$ 为一符合条件的解； - 若 $p>\sqrt[4]{n}$，则若 $a'>p$，则 $n=pa'b'c'>p^4>n$，矛盾！故 $a'<p$，而又有 $p$ 为 $n$ 的最小质因数，所以 $a'=1$，故 $p,b,c$ 为一符合条件的解。 综上所述，原命题成立。 这也揭示了将 $n$ 分解为满足上述条件的 $a,b,c$ 的方法：线性筛后将 $n$ 最小素因子 $p$ 乘给 $\dfrac{n}{p}$ 分解成的三个数中最小的数即可。 所以当询问 $\gcd(n,m)$ 时，可以将 $n$ 分解为这样的 $abc$，每次求 $t=\gcd(a/b/c,m)$ 后，将最终结果乘上 $t$ 并将 $m$ 除以 $t$，即可得解。 具体地，若此时考虑到的数 $x$ 为素数，则求 $\gcd(x,m)$ 是简单的；否则必有 $x\le\sqrt{n}$，而 $\gcd(x,m)=\gcd(x,m\bmod x)$，所以只需预处理出 $\sqrt{v}$ 内任意两个数的 $\gcd$ 后 $\mathcal{O}(1)$ 查询即可。 时间复杂度：$\mathcal{O}(\sqrt{v}^2)=\mathcal{O}(v)$ 预处理，$\mathcal{O}(1)$ 查询。 ### 例题 $1$： - 给定数列 $a,b$，对每一个 $i$，求 $\sum\limits_{j=1}^{n}i^j\gcd(a_i,b_j)$。 - $n\le5000$，$a_i,b_i\le v=10^6$。 乍看一个推式子题，用如上高科技秒切。 时间复杂度：$\mathcal{O}(n^2+v) 
#include <bits/stdc++.h>

using namespace std;

const int N = 5005, M = 1e6 + 5, S = 1000, P = 998244353;

int _a[M], _b[M], _c[M];
int v[M], p[M], r;
int f[S + 1][S + 1];
int a[N], b[N];

void Init()
{
    _a[1] = _b[1] = _c[1] = 1;
    for (int i = 2; i <= M - 5; ++i) {
        if (!v[i]) {
            p[++r] = i;
            _a[i] = _b[i] = 1, _c[i] = i; 
        }
        int tp;
        for (int j = 1; j <= r && (tp = i * p[j]) <= M - 5; ++j) {
            v[tp] = 1;
            _a[tp] = _a[i] * p[j];
            _b[tp] = _b[i];
            _c[tp] = _c[i];
            if (_a[tp] > _b[tp]) {
                swap(_a[tp], _b[tp]);
                if (_b[tp] > _c[tp]) {
                    swap(_b[tp], _c[tp]);
                }
            }
            if (!(i % p[j])) {
                break;
            }
        }
    }

    for (int i = 1; i <= S; ++i) {
        f[0][i] = f[i][0] = i;
        for (int j = 1; j <= S; ++j) {
            f[i][j] = f[j % i][i];
        }
    }
    return;
}

int gcd(int x, int y)
{
    int A = 1, tp = f[_a[x]][y % _a[x]];
    A *= tp;
    y /= tp;
    tp = f[_b[x]][y % _b[x]];
    A *= tp;
    y /= tp;
    tp = (v[_c[x]] ? f[_c[x]][y % _c[x]] : (y % _c[x] ? 1 : _c[x]));
    A *= tp;
    y /= tp;
    return A;
}

signed main()
{
    Init();
    int n;
    scanf("%d", &n);
    for (int i = 1; i <= n; ++i) {
        scanf("%d", &a[i]);
    }
    for (int i = 1; i <= n; ++i) {
        scanf("%d", &b[i]);
    }
    for (int i = 1; i <= n; ++i) {
        int tp = 1, A = 0;
        for (int j = 1; j <= n; ++j) {
            tp = 1ll * tp * i % P;
            A = (A + 1ll * tp * gcd(a[i], b[j])) % P;
        }
        printf("%d\n", A);
    }
    return 0;
}