Luogu P3704 [SDOI2017]数字表格
题意
\(T\) 组询问，每组询问给定 \(n,m\)，求 \(\prod\limits_{i=1}^{n}\prod\limits_{j=1}^{m}F_{\gcd(i,j)}\) 对 \(10^9+7\) 取模的值，\(F\) 是 Fibonacci 数列。
\(\texttt{Data Range:}1\leq T\leq 10^3,1\leq n,m\leq 10^6\)
题解
\(\texttt{Y}\)\(\texttt{LWang}\) 是神！
考虑套路地枚举一发 \(\gcd\)，有
\[\prod\limits_{i=1}^{n}\prod\limits_{j=1}^{m}F_{\gcd(i,j)}=\prod\limits_{d=1}^{n}\prod\limits_{i=1}^{n}\prod\limits_{j=1}^{m}F_{\gcd(i,j)}[\gcd(i,j)=d]
\]
再考虑套路地把 \(\gcd\) 提出来
\[\prod\limits_{i=1}^{n}\prod\limits_{j=1}^{m}F_{\gcd(i,j)}=\prod\limits_{d=1}^{n}\prod\limits_{i=1}^{\left\lfloor\frac{n}{d}\right\rfloor}\prod\limits_{j=1}^{\left\lfloor\frac{m}{d}\right\rfloor}F_{d}[\gcd(i,j)=1]
\]
然后我们可以将这个式子写成 \(F_d\) 的幂的形式
\[\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}^{\left\lfloor\frac{n}{d}\right\rfloor}\sum\limits_{j=1}^{\left\lfloor\frac{m}{d}\right\rfloor}[\gcd(i,j)=1]}
\]
注意到上面其实就是一个套路反演
\[\prod\limits_{i=1}^{n}\prod\limits_{j=1}^{m}F_{\gcd(i,j)}=\prod\limits_{d=1}^{n}F_d^{\sum\limits_{k=1}^{\left\lfloor\frac{n}{d}\right\rfloor}\mu(k)\left\lfloor\frac{n}{dk}\right\rfloor\left\lfloor\frac{m}{dk}\right\rfloor}
\]
然后给你写回来
\[\prod\limits_{i=1}^{n}\prod\limits_{j=1}^{m}F_{\gcd(i,j)}=\prod\limits_{d=1}^{n}\prod\limits_{k=1}^{\left\lfloor\frac{n}{d}\right\rfloor}F_d^{\mu(k)\left\lfloor\frac{n}{dk}\right\rfloor\left\lfloor\frac{m}{dk}\right\rfloor}
\]
现在做代换 \(T=dk\)
\[\prod\limits_{i=1}^{n}\prod\limits_{j=1}^{m}F_{\gcd(i,j)}=\prod\limits_{T=1}^{n}\prod\limits_{d\mid T}F_d^{\mu\left(\frac{T}{d}\right)\left\lfloor\frac{n}{T}\right\rfloor\left\lfloor\frac{m}{T}\right\rfloor}
\]
见证奇迹的时刻到了
\[\prod\limits_{i=1}^{n}\prod\limits_{j=1}^{m}F_{\gcd(i,j)}=\prod\limits_{T=1}^{n}\left(\prod\limits_{d\mid T}F_d^{\mu\left(\frac{T}{d}\right)}\right)^{\left\lfloor\frac{n}{T}\right\rfloor\left\lfloor\frac{m}{T}\right\rfloor}
\]
括号里面的东西可以预处理出来，括号外面整除分块即可，由于可能要求很多次逆元，所以建议写个线性求逆。
代码
#include<bits/stdc++.h>
using namespace std;
typedef int ll;
typedef long long int li;
const ll MAXN=1e6+51,MOD=1e9+7;
ll test,ptot,n,m,x;
ll np[MAXN],prime[MAXN],mu[MAXN],fib[MAXN],f[MAXN],prod[MAXN];
ll pr[MAXN],invf[MAXN],invp[MAXN];
inline ll read()
{
register ll num=0,neg=1;
register char ch=getchar();
while(!isdigit(ch)&&ch!='-')
{
ch=getchar();
}
if(ch=='-')
{
neg=-1;
ch=getchar();
}
while(isdigit(ch))
{
num=(num<<3)+(num<<1)+(ch-'0');
ch=getchar();
}
return num*neg;
}
inline ll qpow(ll base,ll exponent)
{
ll res=1;
while(exponent)
{
if(exponent&1)
{
res=(li)res*base%MOD;
}
base=(li)base*base%MOD,exponent>>=1;
}
return res;
}
inline void sieve(ll limit)
{
np[1]=mu[1]=f[1]=fib[1]=1;
for(register int i=2;i<=limit;i++)
{
fib[i]=(fib[i-1]+fib[i-2])%MOD,f[i]=1;
if(!np[i])
{
prime[++ptot]=i,mu[i]=-1;
}
for(register int j=1;i*prime[j]<=limit;j++)
{
np[i*prime[j]]=1;
if(i%prime[j]==0)
{
mu[i*prime[j]]=0;
break;
}
mu[i*prime[j]]=-mu[i];
}
}
}
inline ll calc(ll n,ll m)
{
ll res=1,cur;
for(register int l=1,r;l<=min(n,m);l=r+1)
{
r=min(n/(n/l),m/(m/l)),cur=1,cur=(li)prod[r]%MOD*invp[l-1]%MOD;
res=(li)res*qpow(cur,(li)(n/l)*(m/l)%(MOD-1))%MOD;
}
return res;
}
int main()
{
test=read(),sieve(1e6+10),prod[0]=pr[0]=1;
for(register int i=1;i<=1e6;i++)
{
pr[i]=(li)pr[i-1]*fib[i]%MOD;
}
x=qpow(pr[1000000],MOD-2),invf[1]=1;
for(register int i=1e6-1;i;i--)
{
invf[i+1]=(li)pr[i]*x%MOD,x=(li)x*fib[i+1]%MOD;
}
for(register int i=1;i<=1e6;i++)
{
for(register int j=1;i*j<=1e6;j++)
{
f[i*j]=(li)f[i*j]*(mu[j]==1?fib[i]:mu[j]==-1?invf[i]:1)%MOD;
}
}
for(register int i=1;i<=1e6;i++)
{
prod[i]=(li)prod[i-1]*f[i]%MOD,pr[i]=(li)pr[i-1]*prod[i]%MOD;
}
x=qpow(pr[1000000],MOD-2),invp[0]=invp[1]=1;
for(register int i=1e6-1;i;i--)
{
invp[i+1]=(li)pr[i]*x%MOD,x=(li)x*prod[i+1]%MOD;
}
for(register int i=0;i<test;i++)
{
n=read(),m=read(),printf("%d\n",calc(n,m));
}
}
