题目大意
\(T\) 组数据，每组数据给定 \(A,B,C\)，求：
\[\prod_{i=1}^A\prod_{j=1}^B\prod_{k=1}^C\Big(\frac{\text{lcm}(i,j)}{\gcd(i,k)}\Big)^{f(type)}\bmod p
\]
其中，\(type \in \{0,1,2\}\)，\(f(0)=1,f(1)=i\times j\times k,f(2)=\gcd(i,j,k)\)。
思路分析
神经污染题，三合一大礼包。
首先简单化一下式子：
\[\prod_{i=1}^A\prod_{j=1}^B\prod_{k=1}^C\Big(\frac{\text{lcm}(i,j)}{\gcd(i,k)}\Big)^{f(type)}=\prod_{i=1}^A\prod_{j=1}^B\prod_{k=1}^C\Big(\frac{ij}{\gcd(i,j)\gcd(i,k)}\Big)^{f(type)}=\Bigg(\frac{\Big(\prod_{i=1}^A\prod_{j=1}^B\prod_{k=1}^Ci\Big)\Big(\prod_{i=1}^A\prod_{j=1}^B\prod_{k=1}^Cj\Big)}{\Big(\prod_{i=1}^A\prod_{j=1}^B\prod_{k=1}^C\gcd(i,j)\Big)\Big(\prod_{i=1}^A\prod_{j=1}^B\prod_{k=1}^C\gcd(i,k)\Big)}{\Bigg)}^{f(type)}
\]
设三元函数 \(f_1(a,b,c)=\prod_{i=1}^a\prod_{j=1}^b\prod_{k=1}^ci^{f(type)}\)，\(f_2(a,b,c)=\prod_{i=1}^a\prod_{j=1}^b\prod_{k=1}^c\gcd(i,j)^{f(type)}\)，则原式可化为
\[\frac{f_1(A,B,C)\times f_1(B,A,C)}{f_2(A,B,C)\times f_2(A,C,B)}
\]
现在的问题就变成了如何计算 \(f_1\) 和 \(f_2\)。
- \(type=0\)
此时 \(f_1(a,b,c)=\prod_{i=1}^a\prod_{j=1}^b\prod_{k=1}^ci=(a!)^{bc}\)，预处理阶乘后直接用快速幂计算即可。
考虑 \(f_2\)，有
\[\begin{aligned}
f_2(a,b,c)&=\prod_{i=1}^a\prod_{j=1}^b\prod_{k=1}^c\gcd(i,j)\\
&=\Big(\prod_{i=1}^a\prod_{j=1}^b\gcd(i,j)\Big)^c
\end{aligned}\]
括号里是一个经典式子，简单化一下：
\[\begin{aligned}
\prod_{i=1}^a\prod_{j=1}^b\gcd(i,j)&=\prod_{d=1}^{a}d^{\sum_{i=1}^{a/d}\sum_{j=1}^{b/d}[\gcd(i,j)=1]}\\
&=\prod_{d=1}^ad^{\sum_{x=1}^{a/d}\mu(x)\lfloor\frac{a}{xd}\rfloor\lfloor\frac{b}{xd}\rfloor}\\
&=\prod_{T=1}^a\Bigg(\prod_{d|T}d^{\mu(\frac{T}{d})}\Bigg)^{\lfloor\frac{a}{T}\rfloor\lfloor\frac{b}{T}\rfloor}
\end{aligned}\]
中间部分可以直接 \(O(n\ln n)\) 枚举倍数预处理出来，然后再整除分块就可以做到 \(O(\sqrt n\log n)\) 回答询问。
- \(type=1\)
先考虑 \(f_1\)，此时有：
\[f_1(a,b,c)=\prod_{i=1}^a\prod_{j=1}^b\prod_{k=1}^ci^{ijk}=\prod_{i=1}^ai^{i\sum_{j=1}^bj\sum_{k=1}^ck}
\]
设 \(S(x)=\frac{x(x+1)}{2}\)，则
\[f_1(a,b,c)=\prod_{i=1}^ai^{iS(b)S(c)}=\Bigg(\prod_{i=1}^ai^i\Bigg)^{S(b)S(c)}
\]
括号里面的部分可以暴力 \(O(n\log n)\) 预处理出来，\(O(\log n)\) 回答询问。
再看 \(f_2\)，此时有：
\[\begin{aligned}
f_2(a,b,c)&=\prod_{i=1}^a\prod_{j=1}^b\prod_{k=1}^c\gcd(i,j)^{ijk}\\
&=\Bigg(\prod_{i=1}^a\prod_{j=1}^b\gcd(i,j)^{ij}\Bigg)^{S(c)}\\
&=\Bigg(\prod_{d=1}^nd^{d^2\sum_{i=1}^{a/d}\sum_{j=1}^{b/d}ij[\gcd(i,j)=1]}\Bigg)^{S(c)}\\
&=\Bigg(\prod_{d=1}^nd^{d^2\sum_{x=1}^{a/d}\mu(x)x^2S(\lfloor\frac{a}{xd}\rfloor)S(\lfloor\frac{b}{xd}\rfloor)}\Bigg)^{S(c)}\\
&=\Bigg(\prod_{T=1}^a\Bigg(\prod_{d|T}d^{\mu(\frac{T}{d})T^2}\Bigg)^{S(\lfloor\frac{a}{T}\rfloor)S(\lfloor\frac{b}{T}\rfloor)}\Bigg)^{S(c)}
\end{aligned}\]
中间的部分也可以 \(O(n\ln n)\) 预处理出来，依旧是 \(O(\sqrt n\log n)\) 回答单次询问。
- \(type=2\)
先看 \(f_1\)：
\[\begin{aligned}
f_1(a,b,c)&=\prod_{i=1}^a\prod_{j=1}^b\prod_{k=1}^ci^{\gcd(i,j,k)}\\
&=\prod_{d=1}^a\prod_{i=1}^{a/d}(id)^{d\sum_{j=1}^{b/d}\sum_{k=1}^{c/d}[\gcd(i,j,k)=1]}\\
&=\prod_{d=1}^a\prod_{x=1}^a\prod_{i=1}^{a/xd}(ixd)^{\mu(x)\lfloor\frac{a}{xd}\rfloor\lfloor\frac{b}{xd}\rfloor}\\
&=\prod_{T=1}^a\Bigg(\Big\lfloor\frac{a}{T}\Big\rfloor !T^{\lfloor\frac{a}{T}\rfloor}\Bigg)^{\varphi(T)\lfloor\frac{b}{T}\rfloor\lfloor\frac{c}{T}\rfloor}\\
&=\Bigg(\prod_{T=1}^a(T^{\varphi(T)})^{\lfloor\frac{a}{T}\rfloor\lfloor\frac{b}{T}\rfloor\lfloor\frac{c}{T}\rfloor}\Bigg)\times \Bigg(\prod_{T=1}^a(\Big\lfloor\frac{a}{T}\Big\rfloor !)^{\varphi(T)\lfloor\frac{b}{T}\rfloor\lfloor\frac{c}{T}\rfloor}\Bigg)
\end{aligned}\]
先别急着看能不能做（做是肯定能做的，预处理 \(T^{\varphi(T)}\) 和 \(\varphi\) 的前缀和就可以了），推最后一个式子（也是最难推的式子）：
\[\begin{aligned}
f_2(a,b,c)&=\prod_{i=1}^a\prod_{j=1}^b\prod_{k=1}^c\gcd(i,j)^{\gcd(i,j,k)}\\
&=\prod_{d=1}^a\prod_{i=1}^{a/d}\prod_{j=1}^{b/d}(d\times \gcd(i,j))^{d\sum_{k=1}^{c/d}[\gcd(i,j,k)=1]}\\
&=\prod_{d=1}^a\prod_{t=1}^a\prod_{i=1}^{a/td}\prod_{j=1}^{b/td}(td\times \gcd(i,j))^{d\mu(t)\lfloor\frac{c}{td}\rfloor}\\
&=\Bigg(\prod_{d=1}^a\prod_{t=1}^a(td)^{\mu(t)d\lfloor\frac{a}{td}\rfloor\lfloor\frac{b}{td}\rfloor\lfloor\frac{c}{td}\rfloor}\Bigg)\times \Bigg(\prod_{d=1}^a\prod_{t=1}^a\prod_{i=1}^{a/td}\prod_{j=1}^{b/td}\gcd(i,j)^{d\mu(t)\lfloor\frac{c}{td}\rfloor}\Bigg)\\
&=\Bigg(\prod_{T=1}^a(T^{\varphi(T)})^{\lfloor\frac{a}{T}\rfloor\lfloor\frac{b}{T}\rfloor\lfloor\frac{c}{T}\rfloor}\Bigg)\times \Bigg(\prod_{d=1}^a\prod_{t=1}^a\prod_{i=1}^{a/td}\prod_{j=1}^{b/td}\gcd(i,j)^{d\mu(t)\lfloor\frac{c}{td}\rfloor}\Bigg)
\end{aligned}\]
将 \(f_1\) 和 \(f_2\) 同时消去左边部分，得
\[f_1'=\prod_{T=1}^a(\Big\lfloor\frac{a}{T}\Big\rfloor !)^{\varphi(T)\lfloor\frac{b}{T}\rfloor\lfloor\frac{c}{T}\rfloor}
\]
\[f_2'=\prod_{d=1}^a\prod_{t=1}^a\prod_{i=1}^{a/td}\prod_{j=1}^{b/td}\gcd(i,j)^{d\mu(t)\lfloor\frac{c}{td}\rfloor}
\]
此时 \(f_1'\) 容易计算，继续推 \(f_2\)：
\[\begin{aligned}
f_2'(a,b,c)&=\prod_{d=1}^a\prod_{t=1}^a\prod_{d'=1}^a\prod_{t'=1}^ad'^{\mu(t')\mu(t)d\lfloor\frac{c}{td}\rfloor\lfloor\frac{a}{dd'tt'}\rfloor\lfloor\frac{b}{dd'tt'}\rfloor}\\
&=\prod_{T=1}^a\Bigg(\prod_{T'=1}^a\Bigg(\prod_{d|T'}d^{\mu(\frac{T'}{d})}\Bigg)^{\lfloor\frac{a}{TT'}\rfloor\lfloor\frac{b}{TT'}\rfloor}\Bigg)^{\varphi(T)\lfloor\frac{c}{T}\rfloor}
\end{aligned}\]
对比中间部分和 \(type=0\) 时 \(f_2\) 的中间部分的形式，得
\[f_2'(a,b,c)=\prod_{T=1}^a\Bigg(\prod_{i=1}^{a/T}\prod_{j=1}^{b/T}\gcd(i,j)\Bigg)^{\varphi(T)\lfloor\frac{c}{T}\rfloor}
\]
中间部分按照类似计算 \(type=0\) 时 \(f_2\) 的计算方式计算，再套一个整除分块，可以做到 \(O(n^{\frac{3}{4}}\log n)\) 回答单次询问。
故总时间复杂度是 \(O(n\log n+Tn^{\frac{3}{4}}\log n)\)。
思路总结
- 预处理：
阶乘，欧拉函数，欧拉函数前缀和，莫比乌斯函数，\(g_1(n)=\prod_{d|n}d^{\mu(\frac{n}{d})}\)，\(g_2(n)=\prod_{i=1}^ni^i\)，\(g_3(n)=\prod_{d|n}d^{\mu(\frac{n}{d})n^2}\)，\(g_1\) 的前缀积和前缀积的逆元，\(g_3\) 的前缀积和前缀积的逆元。
- \(type=0\)
分子用阶乘和快速幂算，分母用整除分块和 \(g_1\) 的前缀积算。
- \(type=1\)
算分子用到 \(g_2\) 和快速幂，算分母的整除分块用到 \(g_3\) 的前缀积。
- \(type=2\)
用到阶乘和欧拉函数的前缀和，整除分块算分子。
分母用到整除分块套整除分块，其中外层用到欧拉函数的前缀和，内层与 \(type=0\) 时的整除分块相同。
代码
（这绝对是我到目前为止打过的最长的数论题的代码了）
#include <iostream>
#include <cstring>
#include <cmath>
#include <cstdio>
#include <algorithm>
using namespace std;
const int N=100100,L=100000;
typedef long long ll;
#define S2(n) ((1ll*n*(n+1)/2)%mod1)
int mod,mod1,T,A,B,C,cnt;
int prime[N],v[N],mu[N],phi[N];
ll fac[N],dmuTd[N],ii[N],Sphi[N];
ll PdmuTd[N],invPdmuTd[N];
ll PdmuTdTT[N],invPdmuTdTT[N];
ll q_pow(ll a,ll b){
ll res=1;
while(b){
if(b&1) res=res*a%mod;
a=a*a%mod;
b>>=1;
}
return res;
}
void init(){
phi[1]=mu[1]=fac[0]=fac[1]=1;
Sphi[1]=dmuTd[1]=ii[0]=1;
PdmuTd[0]=invPdmuTd[0]=1;
PdmuTdTT[0]=invPdmuTdTT[0]=1;
for(int i=2;i<=L;i++){
if(!v[i]){
prime[++cnt]=i;
mu[i]=-1;
phi[i]=i-1;
}
for(int j=1;j<=cnt&&i*prime[j]<=L;j++){
v[i*prime[j]]=1;
if(i%prime[j]==0){
phi[i*prime[j]]=phi[i]*prime[j];
break;
}
phi[i*prime[j]]=phi[i]*phi[prime[j]];
mu[i*prime[j]]=-mu[i];
}
dmuTd[i]=1;
fac[i]=fac[i-1]*i%mod;
Sphi[i]=Sphi[i-1]+phi[i];
}
for(int i=1;i<=L;i++){
ii[i]=ii[i-1]*q_pow(i,i)%mod;
for(int j=i;j<=L;j+=i){
if(mu[i]==1) dmuTd[j]=dmuTd[j]*(j/i)%mod;
if(mu[i]==-1) dmuTd[j]=dmuTd[j]*q_pow(j/i,mod-2)%mod;
}
}
for(int i=1;i<=L;i++){
PdmuTd[i]=PdmuTd[i-1]*dmuTd[i]%mod;
invPdmuTd[i]=q_pow(PdmuTd[i],mod-2);
PdmuTdTT[i]=PdmuTdTT[i-1]*q_pow(dmuTd[i],1ll*i*i%mod1)%mod;
invPdmuTdTT[i]=q_pow(PdmuTdTT[i],mod-2);
}
}
ll TYPE0_solve(int A,int B){
int M=min(A,B);
ll res=1;
for(int l=1,r;l<=M;l=r+1){
int Al=A/l,Bl=B/l;
r=min(A/Al,B/Bl);
res=res*q_pow(PdmuTd[r]*invPdmuTd[l-1]%mod,1ll*Al*Bl%mod1)%mod;
}
return res;
}
ll TYPE0_calc(){
ll res1=q_pow(fac[A],1ll*B*C%mod1);
ll res2=q_pow(fac[B],1ll*A*C%mod1);
ll res3=q_pow(TYPE0_solve(A,B),C);
ll res4=q_pow(TYPE0_solve(A,C),B);
ll res5=res1*res2%mod;
ll res6=res3*res4%mod;
return res5*q_pow(res6,mod-2)%mod;
}
ll TYPE1_solve(int A,int B){
int M=min(A,B);
ll res=1;
for(int l=1,r;l<=M;l=r+1){
int Al=A/l,Bl=B/l;
r=min(A/Al,B/Bl);
res=res*q_pow(PdmuTdTT[r]*invPdmuTdTT[l-1]%mod,S2(Al)*S2(Bl)%mod1)%mod;
}
return res;
}
ll TYPE1_calc(){
ll res1=q_pow(ii[A],S2(B)*S2(C)%mod1);
ll res2=q_pow(ii[B],S2(A)*S2(C)%mod1);
ll res3=q_pow(TYPE1_solve(A,B),S2(C));
ll res4=q_pow(TYPE1_solve(A,C),S2(B));
ll res5=res1*res2%mod;
ll res6=res3*res4%mod;
return res5*q_pow(res6,mod-2)%mod;
}
ll TYPE2_solve1(int A,int B,int C){
int M=min({A,B,C});
ll res=1;
for(int l=1,r;l<=M;l=r+1){
int Al=A/l,Bl=B/l,Cl=C/l;
r=min({A/Al,B/Bl,C/Cl});
res=res*q_pow(fac[Al],(1ll*Bl*Cl%(mod-1))*(Sphi[r]-Sphi[l-1])%mod1)%mod;
}
return res;
}
ll TYPE2_solve2(int A,int B,int C){
int M=min({A,B,C});
ll res=1;
for(int l=1,r;l<=M;l=r+1){
int Al=A/l,Bl=B/l,Cl=C/l;
r=min({A/Al,B/Bl,C/Cl});
res=res*q_pow(TYPE0_solve(Al,Bl),Cl*(Sphi[r]-Sphi[l-1])%mod1)%mod;
}
return res;
}
ll TYPE2_calc(){
ll res1=TYPE2_solve1(A,B,C);
ll res2=TYPE2_solve1(B,A,C);
ll res3=TYPE2_solve2(A,B,C);
ll res4=TYPE2_solve2(A,C,B);
ll res5=res1*res2%mod;
ll res6=res3*res4%mod;
return res5*q_pow(res6,mod-2)%mod;
}
signed main(){
scanf("%d%d",&T,&mod);mod1=mod-1;
init();
while(T--){
scanf("%d%d%d",&A,&B,&C);
ll ans0=TYPE0_calc();
ll ans1=TYPE1_calc();
ll ans2=TYPE2_calc();
printf("%lld %lld %lld\n",ans0,ans1,ans2);
}
return 0;
}