题解 P3704 【[SDOI2017]数字表格】

· · 题解

在博客园食用更佳：https://www.cnblogs.com/PinkRabbit/p/10011541.html。

题意简述：

求 \prod_{i=1}^{N}\prod_{j=1}^{M}F_{\gcd(i,j)}\bmod mod ，其中 F_{i} 是斐波那契数列的第 i 项， mod=10^9+7 。

共有 T 组数据。

题解：

喜闻乐见的推式子时间。

不失一般性，假设 N\le M 。

\begin{aligned}&=\prod_{i=1}^{N}\prod_{j=1}^{M}F_{\gcd(i,j)} \\ &=\prod_{k=1}^{N}{F_{k}}^{\left(\sum_{i=1}^{N}\;\sum_{j=1}^{M}\;\left[\gcd(i,j)=k\right]\right)}\end{aligned}

右上角的指数部分是老套路了。

\begin{aligned}&= \sum_{i=1}^{N}\sum_{j=1}^{M}\left[\gcd(i,j)=k\right]\\&= \sum_{i=1}^{\left\lfloor\frac{N}{k}\right\rfloor}\sum_{j=1}^{\left\lfloor\frac{M}{k}\right\rfloor}\left[\gcd(i,j)=1\right]\\&= \sum_{d=1}^{\left\lfloor\frac{N}{k}\right\rfloor}\mu(d)\left\lfloor\frac{N}{kd}\right\rfloor\left\lfloor\frac{M}{kd}\right\rfloor\end{aligned}

所以

\begin{aligned} &= \prod_{k=1}^{N}{F_{k}}^{\left(\sum_{d=1}^{\left\lfloor\frac{N}{k}\right\rfloor}\mu(d)\left\lfloor\frac{N}{kd}\right\rfloor\left\lfloor\frac{M}{kd}\right\rfloor\right)}\\ &= \prod_{T=1}^{N}\left(\prod_{k|T}{F_{k}}^{\mu(\frac{T}{k})}\right)^{\left\lfloor\frac{N}{T}\right\rfloor\left\lfloor\frac{M}{T}\right\rfloor} \end{aligned}

令 f(n)=\prod_{d|n}{F_{d}}^{\mu(\frac{n}{d})} 。

则

=\prod_{T=1}^{N}{f(T)}^{\left\lfloor\frac{N}{T}\right\rfloor\left\lfloor\frac{M}{T}\right\rfloor}

外层数论分块求出。内层的 f(T) 直接暴力预处理，每个数直接乘到它的倍数中，复杂度 \Theta(n\log n) 。

注意实现的时候的时间复杂度，我因为实现多了快速幂的一个 \log 被卡了。

正确的时间复杂度应该是 \Theta(N(\log N+\log mod)+T\sqrt{N}\log mod) 。

#include<cstdio>
#include<algorithm>
using namespace std;

#define mod 1000000007
#define LL long long

int Pow(int b, LL e) {
    if (e < 0) e += mod - 1;
    int a = 1;
    for (; e; b = (LL)b * b % mod, e >>= 1)
        if (e & 1) a = (LL)a * b % mod;
    return a;
}

bool ip[1000005];
int p[80005], pc;
int mu[1000005];
int f[1000005], fr[1000005];

void init() {

    ip[1] = 1;
    mu[1] = 1;

    for (int i = 2; i <= 1000000; ++i) {
        if (!ip[i]) {
            p[++pc] = i;
            mu[i] = -1;
        }
        for (int j = 1; j <= pc && (LL)p[j] * i <= 1000000; ++j) {
            register int k = p[j] * i;
            ip[k] = 1;
            if (i % p[j]) mu[k] = -mu[i];
            else break;
        }
    }

    for (int i = 1; i <= 1000000; ++i)
        f[i] = 1, fr[i] = 1;

    int A = 1, B = 0;
    for (int i = 1; i <= 1000000; ++i) {
        B = (A + B) % mod;
        A = (B - A + mod) % mod;
        int G[3] = {Pow(B, -1), 1, B};
        for (int j = i, k = 1; j <= 1000000; j += i, ++k) {
            f[j] = (LL)f[j] * G[mu[k] + 1] % mod,
            fr[j] = (LL)fr[j] * G[1 - mu[k]] % mod;
        }
    }

    f[0] = fr[0] = 1;
    for (int i = 1; i <= 1000000; ++i)
        f[i] = (LL)f[i - 1] * f[i] % mod,
        fr[i] = (LL)fr[i - 1] * fr[i] % mod;
}

int main() {
    init();
    int T;
    scanf("%d", &T);
    while (T--) {
        int N, M;
        scanf("%d%d", &N, &M);
        if (N > M) swap(N, M);
        int A = 1;
        for (int i = 1, j; i <= N; i = j + 1) {
            j = min(N / (N / i), M / (M / i));
            A = (LL)A * Pow((LL)f[j] * fr[i - 1] % mod, (LL)(N / i) * (M / i)) % mod;
        }
        printf("%d\n", A);
    }
    return 0;
}