傳送門

做莫比烏斯反演題顯著提高了我的\(\LaTeX\)水平

推式子(默認\(N \leq M\)，分數下取整，會省略大部分過程)

\(\begin{align*} \prod\limits_{i=1}^N \prod\limits_{j=1}^M f[gcd(i,j)] & = \prod\limits_{d=1}^N f[d]^{\sum\limits_{i=1}^\frac{N}{d} \sum\limits_{j=1}^\frac{M}{d}[gcd(i,j)==1]} \\ & = \prod\limits_{d = 1}^N f[d]^{\sum\limits_{p=1}^\frac{N}{d} \mu(p) \frac{N}{dp} \frac{M}{dp}} \end{align*}\)

推到這裏可以\(O(TN)\)地通過兩個數論分塊得出答案，可以獲得70pts

當然這還不夠，按照老套路枚舉\(dp\)繼續推式子

\(\begin{align*} \prod\limits_{d = 1}^N f[d]^{\sum\limits_{p=1}^\frac{N}{d} \mu(p) \frac{N}{dp} \frac{M}{dp}} & = \prod\limits_{T=1} ^ N \prod\limits_{d | T} f[d] ^ {\mu (\frac{T}{d}) \frac{N}{T} \frac{M}{T}} \\ & = \prod\limits_{T=1} ^ N \prod\limits_{d | T} (f[d] ^ {\mu (\frac{T}{d})} ) ^ { \frac{N}{T} \frac{M}{T}} \end{align*}\)

\(\frac{N}{T} \frac{M}{T}\)可以數論分塊，所以要求\(f[d]^{\mu(\frac{T}{d})}\)的前綴積

當你以為要線性篩的時候……\(N \leq 10^6\)直接枚舉倍數乘進去就行了……

總復雜度\(O(NlogN + T\sqrt{N})\)

#include<bits/stdc++.h>
//This code is written by Itst
using namespace std;

inline int read(){
    int a = 0;
    char c = getchar();
    bool f = 0;
    while(!isdigit(c)){
        if(c == ‘-‘)
            f = 1;
        c = getchar();
    }
    while(isdigit(c)){
        a = (a << 3) + (a << 1) + (c ^ ‘0‘);
        c = getchar();
    }
    return f ? -a : a;
}

const int MOD = 1e9 + 7 , MAXN = 1e6 + 7;
bool nprime[MAXN];
int fib[MAXN] , inv[MAXN] , mu[MAXN] , times[MAXN] , prime[MAXN];
int N , M , cnt;

inline int poww(long long a , int b){
    int times = 1;
    while(b){
        if(b & 1)
            times = times * a % MOD;
        a = a * a % MOD;
        b >>= 1;
    }
    return times;
}

void init(){
    fib[1] = inv[1] = times[1] = times[0] = mu[1] = 1;
    for(int i = 2 ; i <= N ; ++i){
        times[i] = 1;
        fib[i] = (fib[i - 1] + fib[i - 2]) % MOD;
        inv[i] = poww(fib[i] , MOD - 2);
    }
    for(int i = 2 ; i <= N ; ++i){
        if(!nprime[i]){
            prime[++cnt] = i;
            mu[i] = -1;
        }
        for(int j = 1 ; j <= cnt && i * prime[j] <= N ; ++j){
            nprime[i * prime[j]] = 1;
            if(i % prime[j] == 0)
                break;
            mu[i * prime[j]] = -1 * mu[i];
        }
    }
    for(int i = 1 ; i <= N ; ++i)
        for(int j = 1 ; j * i <= N ; ++j)
            if(mu[j])
                times[j * i] = 1ll * times[j * i] * (mu[j] > 0 ? fib[i] : inv[i]) % MOD;
    for(int i = 2 ; i <= N ; ++i)
        times[i] = 1ll * times[i - 1] * times[i] % MOD;
    for(int i = 0 ; i <= N ; ++i)
        inv[i] = poww(times[i] , MOD - 2);
}

int main(){
    N = 1e6;
    init();
    for(int T = read() ; T ; --T){
        N = read();
        M = read();
        if(N > M)
            swap(N , M);
        int ans = 1;
        for(int i = 1 , pi ; i <= N ; i = pi + 1){
            pi = min(N / (N / i) , M / (M / i));
            ans = 1ll * ans * poww(1ll * times[pi] * inv[i - 1] % MOD , 1ll * (N / i) * (M / i) % (MOD - 1)) % MOD;
        }
        printf("%d\n" , ans);
    }
    return 0;
}
 

BZOJ4816 SDOI2017 數字表格 莫比烏斯反演