一道究极恶心的毒瘤六合一题，式子推了我满满两面 A4 纸……
首先我们可以将式子拆成：
[ans=prodlimits_{i=1}^Aprodlimits_{j=1}^Bprodlimits_{k=1}^C(dfrac{ij}{gcd(i,j)gcd(i,j)})^{f(type)}
]
也就是说我们需要算出以下四项式子的值：
[prodlimits_{i=1}^Aprodlimits_{j=1}^Bprodlimits_{k=1}^Ci^{f(type)}
]
[prodlimits_{i=1}^Aprodlimits_{j=1}^Bprodlimits_{k=1}^Cj^{f(type)}
]
[prodlimits_{i=1}^Aprodlimits_{j=1}^Bprodlimits_{k=1}^C(dfrac{1}{gcd(i,j)})^{f(type)}
]
[prodlimits_{i=1}^Aprodlimits_{j=1}^Bprodlimits_{k=1}^C(dfrac{1}{gcd(i,k)})^{f(type)}
]
显然前两项与后两项是等价的，因此我们只需算出：
[f_1(type)=prodlimits_{i=1}^Aprodlimits_{j=1}^Bprodlimits_{k=1}^Ci^{f(type)}
]
与
[f_2(type)=prodlimits_{i=1}^Aprodlimits_{j=1}^Bprodlimits_{k=1}^C(gcd(i,j))^{f(type)}
]
即可求出答案。
考虑对 (type) 进行分类讨论，首先是 (type=0)，那么
[f_1(0)=prodlimits_{i=1}^Aprodlimits_{j=1}^Bprodlimits_{k=1}^Ci
]
考虑每个 (i) 的贡献，稍微想想即可得到：
[f_1(0)=(A!)^{BC}
]
然后是 (f_2(0))，套路地枚举 (gcd(i,j)) 有
[egin{aligned}
f_2(0)&=prodlimits_{i=1}^Aprodlimits_{j=1}^Bprodlimits_{k=1}^Cgcd(i,j)\
&=(prodlimits_{d=1}^{min(A,B)}d^{sumlimits_{i=1}^Asumlimits_{j=1}^B[gcd(i,j)=d]})^C\
&=(prodlimits_{d=1}^{min(A,B)}d^{sumlimits_{i=1}^{A/d}sumlimits_{j=1}^{B/d}[gcd(i,j)=1]})^C\
&=(prodlimits_{d=1}^{min(A,B)}d^{sumlimits_{p}mu(p)lfloorfrac{A}{dp}
floorlfloorfrac{B}{dp}
floor})^C\
end{aligned}
]
最右边那个 (k) 次方显然可以忽略掉不管它，最后求个快速幂即可。考虑对里面的 (dp) 进行二维的整除分块，那么答案的式子可以化为：
[egin{aligned}
f_2(0)=(prodlimits_{dp}g_1(dp)^{lfloorfrac{A}{dp}
floorlfloorfrac{B}{dp}
floor})^C
end{aligned}
]
其中
[g_1(x)=prodlimits_{d·p=x}d^{mu(p)}
]
考虑整除分块，对于一段区间 ([L,R])，满足 (forall xin[L,R]) 均有 (lfloordfrac{A}{x} floor=lfloordfrac{A}{L} floor,lfloordfrac{B}{x} floor=lfloordfrac{B}{L} floor)，我们这样计算它们的贡献：
[prodlimits_{i=L}^R(g_1(i)^{lfloorfrac{A}{i}
floorlfloorfrac{B}{i}
floor})^C=(prodlimits_{i=L}^Rg_1(i))^{lfloorfrac{A}{i}
floorlfloorfrac{B}{i}
floor·C}
]
预处理前缀积即可 (mathcal O(1)) 计算，时间复杂度 (mathcal O(log nsqrt{n}))。
接下来是 (type=1) 的情况，个人感觉与 (type=0) 的情况大差不差，毕竟指数上都只与 (i,j,k) 本身而不涉及到它们的 (gcd) 之类，只不过指数上枚举变量的次数稍微高了一点点，导致情况较于 (type=0) 略有一点点繁琐。
首先是 (f_1(1))
[f_1(1)=prodlimits_{i=1}^Aprodlimits_{j=1}^Bprodlimits_{k=1}^Ci^{ijk}
]
还是考虑每个 (i) 的贡献被累计了多少次：
[egin{aligned}
f_1(0)&=(A!)^{sumlimits_{j=1}^Bsumlimits_{k=1}^Cjk}\
&=(A!)^{frac{B(B+1)}{2}·frac{C(C+1)}{2}}
end{aligned}
]
一波快速幂带走。
其次是 (f_2(1))
[egin{aligned}
f_2(1)&=prodlimits_{i=1}^Aprodlimits_{j=1}^Bprodlimits_{k=1}^Cgcd(i,j)^{ijk}\
&=(prodlimits_{d=1}^{min(A,B)}d^{sumlimits_{i=1}^Asumlimits_{j=1}^Bij[gcd(i,j)=d]})^{sumlimits_{k=1}^Ck}\
&=(prodlimits_{d=1}^{min(A,B)}d^{sumlimits_{i=1}^{A/d}sumlimits_{j=1}^{B/d}ijd^2[gcd(i,j)=1]})^{frac{C(C+1)}{2}}\
&=(prodlimits_{d=1}^{min(A,B)}d^{sumlimits_{p}mu(p)d^2p^2·s(lfloorfrac{A}{dp}
floor)s(lfloorfrac{B}{dp}
floor)})^{frac{C(C+1)}{2}}\
end{aligned}
]
其中 (s(i)=dfrac{i(i+1)}{2})
那么我们还是枚举 (dp)，按照 (f_2(0)) 的套路设一个 (g_2(x))，定义如下：
[g_2(x)=prodlimits_{d·p=x}d^{d^2p^2mu(p)}
]
那么考虑对 (dp) 进行整除分块，那么
[f_2(1)=(prodlimits_{dp=1}^{min(A,B)}g_2(dp)^{s(lfloorfrac{A}{dp}
floor)s(lfloorfrac{B}{dp}
floor)})^{frac{C(C+1)}{2}}
]
预处理 (g_2(dp)) 的前缀积然后对 (dp) 整除分块即可。
最后是 (type=2) 的情况
先考虑 (f_1(2))
[f_1(2)=prodlimits_{i=1}^Aprodlimits_{j=1}^Bprodlimits_{k=1}^Ci^{gcd(i,j,k)}
]
上来先把 (gcd) 反演掉：
[egin{aligned}
f_1(2)&=prodlimits_{d=1}^{min(A,B,C)}prodlimits_{i=1}^{lfloorfrac{A}{d}
floor}(id)^{d·sumlimits_{j=1}^{lfloorfrac{B}{d}
floor}sumlimits_{k=1}^{lfloorfrac{C}{d}
floor}[gcd(i,j,k)=1]}\
&=prodlimits_{d=1}^{min(A,B,C)}prodlimits_{i=1}^{lfloorfrac{A}{d}
floor}(id)^{d·sumlimits_{pmid i}mu(p)lfloorfrac{B}{dp}
floorlfloorfrac{C}{dp}
floor}
end{aligned}
]
把 (p) 提到外面来
[egin{aligned}
f_1(2)&=prodlimits_{d=1}^{min(A,B,C)}prodlimits_{p}(prodlimits_{iin[1,lfloorfrac{A}{d}
floor]&pmid i}(id)^{d·mu(p)})^{lfloorfrac{B}{dp}
floorlfloorfrac{C}{dp}
floor}
end{aligned}
]
然后按照套路枚举 (dfrac{i}{p})
[egin{aligned}
f_1(2)&=prodlimits_{d=1}^{min(A,B,C)}prodlimits_{p}(prodlimits_{i=1}^{lfloorfrac{A}{dp}
floor}(idp)^{d·mu(p)})^{lfloorfrac{B}{dp}
floorlfloorfrac{C}{dp}
floor}
end{aligned}
]
然后枚举 (dp)，根据 (mu*i=varphi) 可知 (sumlimits_{d·p=x}d·mu(p)=varphi(x))，于是
[f_1(2)=prodlimits_{x=1}^{min(A,B,C)}(prodlimits_{i=1}^{lfloorfrac{A}{x}
floor}(ix)^{varphi(x)})^{lfloorfrac{B}{x}
floorlfloorfrac{C}{x}
floor}
]
外面的东西显然整除分块一下就好了，里面的东西
[prodlimits_{i=1}^{lfloorfrac{A}{x}
floor}(ix)^{varphi(x)}
]
显然等于
[egin{aligned}
&prodlimits_{i=1}^{lfloorfrac{A}{x}
floor}i^{varphi(x)}·x^{varphi(x)}\
=&((lfloordfrac{A}{x}
floor)!)^{varphi(x)}·(x^{varphi(x)})^{lfloorfrac{A}{x}
floor}
end{aligned}
]
然后套路地预处理 (g_3=x^{varphi(x)}) 的前缀积，以及 (varphi(x)) 的前缀和即可在整除分块的过程中 (mathcal O(1)) 求出式子的值，注意 (varphi(x)) 的前缀和应 (mod(MOD-1)) 而不是 (mod MOD)，因为 (varphi(x)) 的前缀和作用在指数上。
时间复杂度 (sqrt{n}log n)
然后是最精神污染的一个式子：
[f_2(2)=prodlimits_{i=1}^Aprodlimits_{j=1}^Bprodlimits_{k=1}^C(gcd(i,j))^{gcd(i,j,k)}
]
按照 P1587 的套路，碰到两个 (gcd) 咱们最好不要莽，要一个个反演，因此考虑先反演下面这个 (gcd)：
[egin{aligned}
f_2(2)&=prodlimits_{d=1}^{min(A,B)}d^{sumlimits_{i=1}^{lfloorfrac{A}{d}
floor}sumlimits_{j=1}^{lfloorfrac{B}{d}
floor}[gcd(i,j)=1]sumlimits_{k=1}^Cgcd(d,k)}\
&=prodlimits_{d=1}^{min(A,B)}d^{sumlimits_{p}^{lfloorfrac{min(A,B)}{d}
floor}mu(p)lfloorfrac{A}{dp}
floorlfloorfrac{B}{dp}
floorsumlimits_{k=1}^Cgcd(d,k)}
end{aligned}
]
发现后面有个 (sumlimits_{k=1}^Cgcd(d,k))，考虑对这个东西单独推个式子：
[egin{aligned}
g_4(d,C)&=sumlimits_{k=1}^Cgcd(d,k)\
&=sumlimits_{xmid d}x·sumlimits_{i=1}^{lfloorfrac{C}{x}
floor}[iperpdfrac{d}{x}]\
&=sumlimits_{xmid d}x·sumlimits_{ymidfrac{d}{x}}mu(y)lfloordfrac{C}{xy}
floor
end{aligned}
]
套路地枚举 (xy=T) 可得：
[g_4(d,C)=sumlimits_{Tmid d}lfloordfrac{C}{T}
floorsumlimits_{xmid T}xmu(dfrac{T}{x})
]
这已经是本题中第二次看到这个式子了：
[sumlimits_{xmid T}xmu(dfrac{T}{x})=varphi(T)
]
于是
[g_4(d,C)=sumlimits_{Tmid d}lfloordfrac{C}{T}
floorvarphi(T)
]
带回去
[f_2(2)=prodlimits_{d=1}^{min(A,B)}d^{sumlimits_{p}^{lfloorfrac{min(A,B)}{d}
floor}mu(p)lfloorfrac{A}{dp}
floorlfloorfrac{B}{dp}
floorsumlimits_{Tmid d}lfloorfrac{C}{T}
floorvarphi(T)}
]
按照这里总结出来的套路，看到先枚举 (i) 再枚举 (jmid i) 的求和/积式我们可以考虑交换求和/积的顺序，先枚举 (j) 再枚举 (i)，这样会出现下取整，就可以整除分块了。
因此考虑先枚举 (T) 再枚举 (d)，有：
[f_2(2)=prodlimits_{T=1}^{min(A,B)}prodlimits_{d=1}^{lfloorfrac{min(A,B)}{T}
floor}(dT)^{sumlimits_{p}mu(p)lfloorfrac{A}{dTp}
floorlfloorfrac{B}{dTp}
floorlfloorfrac{C}{T}
floorvarphi(T)}
]
然后考虑拆开来：
[f_2(2)=prodlimits_{T=1}^{min(A,B)}(T^{varphi(T)lfloorfrac{C}{T}
floor})^{sumlimits_{d}sumlimits_{p}mu(p)lfloorfrac{A}{dTp}
floorlfloorfrac{B}{dTp}
floor}·prodlimits_{d=1}^{lfloorfrac{min(A,B)}{T}
floor}(d^{sumlimits_{p}mu(p)lfloorfrac{A}{dTp}
floorlfloorfrac{B}{dTp}
floor})^{varphi(T)lfloorfrac{C}{T}
floor}
]
先考虑前面的式子：
[(T^{varphi(T)lfloorfrac{C}{T}
floor})^{sumlimits_{d}sumlimits_{p}mu(p)lfloorfrac{A}{dTp}
floorlfloorfrac{B}{dTp}
floor}
]
考虑枚举 (dp=x)，那么：
[egin{aligned}
ext{原式}&=(T^{varphi(T)lfloorfrac{C}{T}
floor})^{sumlimits_{x}sumlimits_{pmid x}mu(p)lfloorfrac{A}{Tx}
floorlfloorfrac{B}{Tx}
floor}\
&=(T^{varphi(T)lfloorfrac{C}{T}
floor})^{epsilon(x)lfloorfrac{A}{Tx}
floorlfloorfrac{B}{Tx}
floor}\
&=(T^{varphi(T)lfloorfrac{C}{T}
floor})^{lfloorfrac{A}{T}
floorlfloorfrac{B}{T}
floor}
end{aligned}
]
对 (T) 整除分块一下，预处理 (varphi(T)) 的前缀和即可 (mathcal O(1)) 算出。
然后是后面的式子（胜利就在眼前！）
[prodlimits_{T=1}^{min(A,B)}prodlimits_{d=1}^{lfloorfrac{min(A,B)}{T}
floor}(d^{sumlimits_{p}mu(p)lfloorfrac{A}{dTp}
floorlfloorfrac{B}{dTp}
floor})^{varphi(T)lfloorfrac{C}{T}
floor}
]
还是对 (T) 整除分块，然后枚举 (dp=x)，那么上面的式子可以写成：
[prodlimits_{T=1}^{min(A,B)}prodlimits_{x=1}^{lfloorfrac{min(A,B)}{T}
floor}((prodlimits_{dmid x}d^{mu(p)})^{lfloorfrac{A}{Tx}
floorlfloorfrac{B}{Tx}
floor})^{varphi(T)lfloorfrac{C}{T}
floor}
]
发现最里面的括号的东西就是在求 (f_2(0)) 时用到的 (g_1(x))，那么我们再套一个对 (x) 的整除分块即可。
根据整除分块里再套一个整除分块复杂度是 (sumlimits_{x,exists k,s.t.lfloorfrac{n}{k} floor=x}sqrt{x}=n^{0.75}) 可知这一部分复杂度为 (n^{0.75}log n)
于是这题就做完了，时间复杂度 (nlog n+Tn^{0.75}log n)。
关于此题的常数问题，由于取模运算较多，可以使用快速取模优化常数，具体可见 chenxia25 神仙的这篇博客。
const int MAXV=1e5;
ll mod;
int pr[MAXV/5+5],prcnt=0,vis[MAXV+5],mu[MAXV+5],smu[MAXV+5],phi[MAXV+5];
int fac[MAXV+5],prd[MAXV+5],prd_inv[MAXV+5],prd_sq[MAXV+5],prd_sq_inv[MAXV+5];
int pre_ii[MAXV+5],inv_pre_ii[MAXV+5],pre_mul[MAXV+5],inv_pre_mul[MAXV+5];
int inv[MAXV+5],prd_phi[MAXV+5],prd_phi_inv[MAXV+5];
int sphi[MAXV+5];
__int128_t _base1=1,_base2=1;
inline int mol1(__int128_t x){return x-mod*(_base1*x>>64);}
inline int mol2(__int128_t x){return x-(mod-1)*(_base2*x>>64);}
int qpow(int x,int e){
if(e<0) e+=mod-1;int ret=1;
for(;e;e>>=1,x=mol1(1ll*x*x)) if(e&1) ret=mol1(1ll*ret*x);
return ret;
}
int work(int x,int y){return (!y)?1:((~y)?x:inv[x]);}
void sieve(int n){
for(int i=(inv[0]=inv[1]=1)+1;i<=n;i++) inv[i]=mol1(1ll*inv[mod%i]*(mod-mod/i));
mu[1]=phi[1]=1;
for(int i=(fac[0]=1);i<=n;i++) fac[i]=mol1(1ll*fac[i-1]*i);
for(int i=2;i<=n;i++){
if(!vis[i]) mu[i]=-1,pr[++prcnt]=i,phi[i]=i-1;
for(int j=1;j<=prcnt&&pr[j]*i<=n;j++){
vis[pr[j]*i]=1;
if(i%pr[j]==0){phi[i*pr[j]]=phi[i]*pr[j];break;}
mu[i*pr[j]]=-mu[i];phi[i*pr[j]]=phi[i]*phi[pr[j]];
}
}
for(int i=1;i<=n;i++) smu[i]=smu[i-1]+mu[i];
for(int i=0;i<=n;i++) prd[i]=prd_sq[i]=1;
for(int i=1;i<=n;i++) for(int j=i;j<=n;j+=i) prd[j]=mol1(1ll*prd[j]*work(i,mu[j/i]));
for(int i=1;i<=n;i++) prd[i]=mol1(1ll*prd[i-1]*prd[i]);
for(int i=0;i<=n;i++) prd_inv[i]=qpow(prd[i],-1);
for(int i=1;i<=n;i++) for(int j=i;j<=n;j+=i)
prd_sq[j]=mol1(1ll*prd_sq[j]*
qpow(qpow(i,mol2(1ll*i*i)),mol2(1ll*mu[j/i]*(j/i)*(j/i))));
for(int i=1;i<=n;i++) prd_sq[i]=mol1(1ll*prd_sq[i-1]*prd_sq[i]);
for(int i=0;i<=n;i++) prd_sq_inv[i]=qpow(prd_sq[i],-1);
pre_ii[0]=1;for(int i=1;i<=n;i++) pre_ii[i]=mol1(1ll*pre_ii[i-1]*qpow(i,i));
for(int i=0;i<=n;i++) inv_pre_ii[i]=qpow(pre_ii[i],-1);
pre_mul[0]=1;for(int i=1;i<=n;i++) pre_mul[i]=mol1(1ll*pre_mul[i-1]*work(i,mu[i]));
for(int i=0;i<=n;i++) inv_pre_mul[i]=qpow(pre_mul[i],-1);
prd_phi[0]=1;for(int i=1;i<=n;i++) prd_phi[i]=mol1(1ll*prd_phi[i-1]*qpow(i,phi[i]));
for(int i=0;i<=n;i++) prd_phi_inv[i]=qpow(prd_phi[i],-1);
for(int i=1;i<=n;i++) sphi[i]=mol2(sphi[i-1]+phi[i]);
}
int calc1(int x,int y){
int res=1;
for(int l=1,r;l<=min(x,y);l=r+1){
r=min(x/(x/l),y/(y/l));
res=1ll*res*qpow(1ll*prd[r]*prd_inv[l-1]%mod,1ll*(x/l)*(y/l)%(mod-1))%mod;
}
return res;
}
int solve1(int a,int b,int c){
int res=1ll*qpow(calc1(a,b),-c)*qpow(calc1(a,c),-b)%mod;
res=1ll*res*qpow(fac[a],1ll*b*c%(mod-1))%mod;
res=1ll*res*qpow(fac[b],1ll*a*c%(mod-1))%mod;
return res;
}
ll get(int x){return mol2(1ll*x*(x+1)/2);}
int calc2(int x,int y){
int res=1;
for(int l=1,r;l<=min(x,y);l=r+1){
r=min(x/(x/l),y/(y/l));
res=1ll*res*qpow(1ll*prd_sq[r]*prd_sq_inv[l-1]%mod,
1ll*get(x/l)*get(y/l)%(mod-1))%mod;
}
return res;
}
int solve2(int a,int b,int c){
int res=1;
res=1ll*res*qpow(pre_ii[a],(1ll*b*(b+1)>>1)%(mod-1))%mod;
res=1ll*res*qpow(pre_ii[b],(1ll*a*(a+1)>>1)%(mod-1))%mod;
res=1ll*res*qpow(calc2(a,b),mod-2)%mod;
res=qpow(res,(1ll*c*(c+1)>>1)%(mod-1));
res=1ll*res*qpow(calc2(a,c),-(1ll*b*(b+1)>>1)%(mod-1))%mod;
return res;
}
int calc3(int x,int y,int z){
int res=1;
for(int l=1,r;l<=x;l=r+1){
r=1e9;
if(x/l) chkmin(r,x/(x/l));
if(y/l) chkmin(r,y/(y/l));
if(z/l) chkmin(r,z/(z/l));
int mul=qpow(1ll*prd_phi[r]*prd_phi_inv[l-1]%mod,x/l);
mul=1ll*mul*qpow(fac[x/l],sphi[r]-sphi[l-1])%mod;
res=1ll*res*qpow(mul,1ll*(y/l)*(z/l)%(mod-1))%mod;
} return res;
}
int calc4(int x,int y,int z){
int res=1;
for(int l=1,r;l<=min(x,y);l=r+1){
r=1e9;
if(x/l) chkmin(r,x/(x/l));
if(y/l) chkmin(r,y/(y/l));
if(z/l) chkmin(r,z/(z/l));
res=1ll*res*qpow(1ll*prd_phi[r]*prd_phi_inv[l-1]%mod,1ll*(x/l)*(y/l)*(z/l)%(mod-1))%mod;
int X=x/l,Y=y/l,Z=z/l,sm=1ll*Z*(sphi[r]-sphi[l-1])%(mod-1);
for(int L=1,R;L<=min(X,Y);L=R+1){
R=1e9;
if(X/L) chkmin(R,X/(X/L));
if(Y/L) chkmin(R,Y/(Y/L));
res=1ll*res*qpow(1ll*prd[R]*prd_inv[L-1]%mod,1ll*(X/L)*(Y/L)*sm%(mod-1))%mod;
}
}
return res;
}
int solve3(int x,int y,int z){
return 1ll*calc3(x,y,z)*calc3(y,x,z)%mod*
qpow(calc4(x,y,z),-1)%mod*qpow(calc4(x,z,y),-1)%mod;
}
int main(){
int qu;scanf("%d%lld",&qu,&mod);
_base1=(_base1<<64)/mod;_base2=(_base2<<64)/(mod-1);
sieve(MAXV);
while(qu--){
int a,b,c;scanf("%d%d%d",&a,&b,&c);
printf("%d %d %d
",solve1(a,b,c),solve2(a,b,c),solve3(a,b,c));
}
return 0;
}