luoguP3704 [SDOI2017]数字表格

题意

默认\(n\leqslant m\)
所求即为：\(\prod\limits_{i=1}^n\prod\limits_{j=1}^mf[\gcd(i,j)]\)
枚举\(\gcd(i,j)\)变为：
\(\prod\limits_{k=1}^{n}f(k)^{\sum\limits_{i=1}^n\sum\limits_{j=1}^m[\gcd(i,j)=k]}\)
上面那个是莫比乌斯反演套路形式：
\(\sum\limits_{i=1}^n\sum\limits_{j=1}^m[\gcd(i,j)=k]\)
\(\sum\limits_{i=1}^{\frac{n}{k}}\sum\limits_{j=1}^{\frac{m}{k}}[\gcd(i,j)=1]\)
\(\sum\limits_{i=1}^{\frac{n}{k}}\sum\limits_{j=1}^{\frac{m}{k}}\sum\limits_{x|\gcd(i,j)}\mu(x)\)
\(\sum\limits_{x=1}^{\frac{n}{k}}\mu(x)\sum\limits_{i=1}^{\frac{n}{k}}\sum\limits_{j=1}^{\frac{m}{k}}[x|\gcd(i,j)]\)
\(\sum\limits_{x=1}^{\frac{n}{k}}\mu(x)\sum\limits_{i=1}^{\frac{n}{k*x}}\sum\limits_{j=1}^{\frac{m}{k*x}}1\)
\(\sum\limits_{x=1}^{\frac{n}{k}}\mu(x)\frac{n}{k*x}\frac{m}{k*x}\)
代回原式：
\(\prod\limits_{k=1}^{n}f(k)^{\sum\limits_{x=1}^{\frac{n}{k}}\mu(x)\frac{n}{k*x}\frac{m}{k*x}}\)
设\(T=k*x\)，转而枚举\(T\)：
\(\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}}\)
显然指数部分可以除法分块，考虑如何求\(\prod\limits_{d|T}f(d)^{\mu(\frac{T}{d})}\)：
设\(g(T)=\prod\limits_{d|T}f(d)^{\mu(\frac{T}{d})}\)
在算到\(f(d)\)时乘到\(g(T)\)即可。
答案即为：
\(\prod\limits_{T=1}^{n}g(T)^{\frac{n}{T}\frac{m}{T}}\)

code:

#include<bits/stdc++.h>
using namespace std;
typedef long long ll;
const int maxn=1e6+10;
const ll mod=1e9+7;
int T,n,m;
ll ans=1;
ll f[maxn],g[maxn],mu[maxn],invf[maxn];
bool vis[maxn];
vector<int>prime;
inline ll power(ll x,ll k,ll mod)
{
    ll res=1;
    while(k)
    {
        if(k&1)res=res*x%mod;
        x=x*x%mod;k>>=1;
    }
    return res;
}
inline void pre_work(int n)
{
    g[0]=g[1]=1;
    vis[1]=1;mu[1]=1;
    for(int i=2;i<=n;i++)
    {
        g[i]=1;
        if(!vis[i])prime.push_back(i),mu[i]=-1;
        for(unsigned int j=0;j<prime.size()&&i*prime[j]<=n;j++)
        {
            vis[i*prime[j]]=1;
            if(i%prime[j]==0)break;
            mu[i*prime[j]]=-mu[i];  
        }   
    }
    f[0]=0,f[1]=1;invf[1]=1;
    for(int i=2;i<=n;i++)f[i]=(f[i-1]+f[i-2])%mod,invf[i]=power(f[i],mod-2,mod);
    for(int i=1;i<=n;i++)
    {
        if(!mu[i])continue;
        for(int j=i;j<=n;j+=i)
            g[j]=g[j]*(mu[i]==1?f[j/i]:invf[j/i])%mod;
    }
    for(int i=2;i<=n;i++)g[i]=g[i]*g[i-1]%mod;
}
int main()
{
    pre_work(1000000);
    scanf("%d",&T);
    while(T--)
    {
        ans=1;
        scanf("%d%d",&n,&m);
        if(n>m)swap(n,m);
        for(int l=1,r;l<=n;l=r+1)
        {
            r=min(n/(n/l),m/(m/l));
            ans=ans*power(g[r]*power(g[l-1],mod-2,mod)%mod,1ll*(n/l)*(m/l)%(mod-1),mod)%mod;
        }
        printf("%lld\n",(ans+mod)%mod);
    }
    return 0;
}

原文地址：https://www.cnblogs.com/nofind/p/11947576.html