#莫比乌斯反演，乘法逆元，快速幂，整除分块#JZOJ 100006 洛谷 3704 bzoj 4816 数字表格

题目

求 ∏ i = 1 n ∏ j = 1 m F g c d ( i , j ) \prod_{i=1}^n\prod_{j=1}^mF_{gcd(i,j)} i=1∏n​j=1∏m​Fgcd(i,j)​,F是斐波那契数列

分析

使 n ≤ m n\leq m n≤m

= ∏ d = 1 n F d ∑ i = 1 n ∑ j = 1 m [ g c d ( i , j ) = = d ] \large=\prod_{d=1}^n{F_d}^{\sum_{i=1}^n\sum_{j=1}^m[gcd(i,j)==d]} =d=1∏n​Fd​∑i=1n​∑j=1m​[gcd(i,j)==d]

重点是指数

∑ i = 1 n ∑ j = 1 m [ g c d ( i , j ) = = d ] = ∑ i = 1 ⌊ n d ⌋ ∑ j = 1 ⌊ m d ⌋ [ g c d ( i , j ) = = 1 ] \large\sum_{i=1}^n\sum_{j=1}^m[gcd(i,j)==d]=\sum_{i=1}^{\lfloor\frac{n}{d}\rfloor}\sum_{j=1}^{\lfloor\frac{m}{d}\rfloor}[gcd(i,j)==1] i=1∑n​j=1∑m​[gcd(i,j)==d]=i=1∑⌊dn​⌋​j=1∑⌊dm​⌋​[gcd(i,j)==1]

= ∑ i = 1 ⌊ n d ⌋ ∑ j = 1 ⌊ m d ⌋ ∑ k ∣ i , k ∣ j μ ( k ) \large=\sum_{i=1}^{\lfloor\frac{n}{d}\rfloor}\sum_{j=1}^{\lfloor\frac{m}{d}\rfloor}\sum_{k|i,k|j}\mu(k) =i=1∑⌊dn​⌋​j=1∑⌊dm​⌋​k∣i,k∣j∑​μ(k)

= ∑ k = 1 ⌊ n d ⌋ μ ( k ) ⌊ n k d ⌋ ⌊ m k d ⌋ \large=\sum_{k=1}^{\lfloor\frac{n}{d}\rfloor}\mu(k)\lfloor\frac{n}{kd}\rfloor\lfloor\frac{m}{kd}\rfloor =k=1∑⌊dn​⌋​μ(k)⌊kdn​⌋⌊kdm​⌋

所以 = ∏ d = 1 n F d ∑ k = 1 ⌊ n d ⌋ μ ( k ) ⌊ n k d ⌋ ⌊ m k d ⌋ \large=\prod_{d=1}^n{F_d}^{\sum_{k=1}^{\lfloor\frac{n}{d}\rfloor}\mu(k)\lfloor\frac{n}{kd}\rfloor\lfloor\frac{m}{kd}\rfloor} =d=1∏n​Fd​∑k=1⌊dn​⌋​μ(k)⌊kdn​⌋⌊kdm​⌋

令 i = k d i=kd i=kd

∏ i = 1 n ( ∏ k ∣ i F k μ ( i k ) ) ⌊ n i ⌋ ⌊ m i ⌋ \large\prod_{i=1}^n(\prod_{k|i}F_{k}^{\mu(\frac{i}{k})})^{\lfloor\frac{n}{i}\rfloor\lfloor\frac{m}{i}\rfloor} i=1∏n​(k∣i∏​Fkμ(ki​)​)⌊in​⌋⌊im​⌋

时间复杂度 O ( n log ⁡ n ( log ⁡ n + log ⁡ m o d ) + q n log ⁡ m o d ) O(n\log n(\log n+\log mod)+q\sqrt n\log mod) O(nlogn(logn+logmod)+qn

​logmod)

代码

#include <cstdio>
#include <cctype>
#include <cstring>
#define rr register
#define min(a,b) ((a)<(b)?(a):(b))
using namespace std;
const int mod=1000000007,N=1000101;
int mu[N],dp[N],inv[N],prime[80001],cnt; bool v[N];
inline signed iut(){
    rr int ans=0; rr char c=getchar();
    while (!isdigit(c)) c=getchar();
    while (isdigit(c)) ans=(ans<<3)+(ans<<1)+(c^48),c=getchar();
    return ans;
}
inline signed ksm(int x,int y){
    rr int ans=1;
    for (;y;y>>=1,x=1ll*x*x%mod)
        if (y&1) ans=1ll*ans*x%mod;
    return ans;
}
inline void init(){
    dp[0]=dp[1]=inv[1]=inv[0]=mu[1]=1;
    for (rr int i=2;i<=N;++i){
        dp[i]=inv[i]=1;
        if (!v[i]) prime[++cnt]=i,mu[i]=-1;
        for (rr int j=1;prime[j]*i<=N;++j){
            v[i*prime[j]]=1;
            if (i%prime[j]) mu[i*prime[j]]=-mu[i];
                else break;
        }
    }
    for (rr int i=1,a=1,b=0;i<=N;++i){
        b=b+a-mod*(b+a>=mod),a=b-a+mod*(b<a);
        rr int t[3]={ksm(b,mod-2),1,b};
        for (rr int j=i,k=1;j<=N;j+=i,++k)
            dp[j]=1ll*dp[j]*t[mu[k]+1]%mod,
                inv[j]=1ll*inv[j]*t[1-mu[k]]%mod;
    }
    for (rr int i=1;i<=N;++i) dp[i]=1ll*dp[i-1]*dp[i]%mod,inv[i]=1ll*inv[i-1]*inv[i]%mod;
}
void print(int ans){
    if (ans>9) print(ans/10);
    putchar(ans%10+48);
}
signed main(){
    init();
    for (rr int t=iut();t;--t){
        rr int n=iut(),m=iut();
        if (n>m) n^=m,m^=n,n^=m;
        rr int ans=1;
        for (rr int l=1,r;l<=n;l=r+1){
            r=min(n/(n/l),m/(m/l));
            ans=1ll*ans*ksm(1ll*dp[r]*inv[l-1]%mod,1ll*(n/l)*(m/l)%(mod-1))%mod;
        }
        print(ans); putchar(10);
    }
    return 0;
}