小D的Lemon
题意
已知
\[g(x)=\left\{\begin{matrix}
1&,x=1\\
\sum_{i=1}^qk_i&,otherwise
\end{matrix}\right.
\]
1&,x=1\\
\sum_{i=1}^qk_i&,otherwise
\end{matrix}\right.
\]
求
\[\prod_{i=1}^n\prod_{j=1}^mg(\gcd(i,j))
\]
\]
说明
多组数据，\(T\le 1000,n,m\le 250000\)
式子太久没推都推不好了...
\[\begin{aligned}
\prod_{i=1}^n\prod_{j=1}^m g(\gcd(i,j))=&\prod_{d=1}^ng(d)^{\sum\limits_{i=1}^n\sum\limits_{j=1}^m[\gcd(i,j)=d]}\\
=&\prod_{d=1}^{\min(n,m)}g(d)^{\sum\limits_{k=1}^{\min(\lfloor\frac{n}{d}\rfloor,\rfloor\frac{m}{d}\rfloor)}\mu(k)\lfloor\frac{n}{dk}\rfloor\lfloor\frac{m}{dk}\rfloor}\\
=&\prod_{T=1}^{\min(n,m)}(\prod_{d|T}g(d)^{\mu(\frac{T}{d})})^{\lfloor\frac{n}{T}\rfloor\lfloor\frac{m}{T}\rfloor}
\end{aligned}
\]
\prod_{i=1}^n\prod_{j=1}^m g(\gcd(i,j))=&\prod_{d=1}^ng(d)^{\sum\limits_{i=1}^n\sum\limits_{j=1}^m[\gcd(i,j)=d]}\\
=&\prod_{d=1}^{\min(n,m)}g(d)^{\sum\limits_{k=1}^{\min(\lfloor\frac{n}{d}\rfloor,\rfloor\frac{m}{d}\rfloor)}\mu(k)\lfloor\frac{n}{dk}\rfloor\lfloor\frac{m}{dk}\rfloor}\\
=&\prod_{T=1}^{\min(n,m)}(\prod_{d|T}g(d)^{\mu(\frac{T}{d})})^{\lfloor\frac{n}{T}\rfloor\lfloor\frac{m}{T}\rfloor}
\end{aligned}
\]
把括号里面的预处理出来就可以了
复杂度\(O(n\sqrt n+T\log n\sqrt n)\)或者\(O(n\ln n+T\log n\sqrt n)\)
Code:
#include <cstdio>
const int N=250000;
const int mod=1e9+7;
int min(int x,int y){return x<y?x:y;}
inline int add(int x,int y){return x+y>=mod?x+y-mod:x+y;}
#define mul(x,y) (1ll*(x)*(y)%mod)
int qp(int d,int k){int f=1;while(k){if(k&1)f=mul(f,d);d=mul(d,d),k>>=1;}return f;}
int g[N+10],mu[N+10],ispri[N+10],pri[N+10],yuu[N+10],yuuinv[N+10],cnt,inv[21];
int cal(int x,int y)
{
if(y==1) return x;
if(y==0) return 1;
return inv[x];
}
void init()
{
g[1]=mu[1]=1;
for(int i=2;i<=N;i++)
{
if(!ispri[i])
{
pri[++cnt]=i;
g[i]=1;
mu[i]=-1;
}
for(int j=1;j<=cnt&&pri[j]*i<=N;j++)
{
int x=pri[j]*i;
ispri[x]=1;
g[x]=g[i]+1;
if(i%pri[j]) mu[x]=-mu[i];
else break;
}
}
for(int i=1;i<=20;i++) inv[i]=qp(i,mod-2);
yuu[0]=yuuinv[0]=1;
for(int j,i=1;i<=N;i++)
{
int bee=1;
for(j=1;j*j<i;j++)
{
if(i%j) continue;
bee=mul(bee,cal(g[j],mu[i/j]));
bee=mul(bee,cal(g[i/j],mu[j]));
}
if(j*j==i) bee=mul(bee,cal(g[j],mu[i/j]));
yuu[i]=mul(yuu[i-1],bee);
}
for(int i=1;i<=N;i++) yuuinv[i]=qp(yuu[i],mod-2);
}
int main()
{
init();
int T,n,m;
scanf("%d",&T);
while(T--)
{
scanf("%d%d",&n,&m);
int ans=1;
for(int l=1,r;l<=min(n,m);l=r+1)
{
r=min(n/(n/l),m/(m/l));
ans=mul(ans,qp(mul(yuu[r],yuuinv[l-1]),1ll*(n/l)*(m/l)%(mod-1)));
}
printf("%d\n",ans);
}
return 0;
}
2019.2.16