P3704 [SDOI2017]数字表格
首先根据题意写出答案的表达式
[largeprod_{i=1}^nprod_{j=1}^mf_{gcd(i,j)}
]
按常规套路改为枚举 (d=gcd(i,j))
（不妨设 (nle m) ）
[largeprod_{d=1}^n{f_d}^{sum_{i=1}^nsum_{j=1}^m~[(i,j)=d]}
]
指数上的式子很熟悉了，单独拿出来推一下
[egin{aligned}
sum_{i=1}^nsum_{j=1}^m[(i,j)=d]
&=sum_{i=1}^{n/d}sum_{j=1}^{m/d}[(i,j)=1]\
&=sum_{i=1}^{n/d}sum_{j=1}^{m/d}sum_{tmid (i,j)}mu(t)\
&=sum_{t=1}^{n/d}mu(t)sum_{i=1}^{n/d}sum_{j=1}^{n/d}[t|i][t|j]\
&=sum_{t=1}^{n/d}mu(t)(n/dt)(m/dt)
end{aligned}]
代回原式
[largeprod_{d=1}^n{f_d}^{sum_{t=1}^{n/d}mu(t)(n/dt)(m/dt)}
]
设 (k=dt) ，改为枚举 (k,d)
[largeprod_{k=1}^nprod_{dmid k}{f_d}^{mu(k/d)(n/k)(m/k)}
]
[largeprod_{k=1}^nleft(prod_{dmid k}{f_d}^{mu(k/d)}
ight)^{(n/k)(m/k)}
]
[largeprod_{k=1}^nleft(prod_{dmid k}{f_{k/d}}^{mu(d)}
ight)^{(n/k)(m/k)}
]
括号里的东西可以 (O(nlog n)) 预处理，然后整除分块计算就可以了
时间复杂度 (O(nlog n+Tsqrt nlog n))
带 (log) 是因为用到了快速幂
到这里已经可以 AC 了
但还是介绍一下 (O(nloglog n+Tsqrt nlog n)) 的做法
可以理解为 dp ，设：
[large s_{i,n}=prod_{dmid n,d 只含前 i 种质因子}{f_{n/d}}^{mu(d)}
]
[large inv_{i,n}={s_{i,n}}^{-1}=prod_{dmid n,d 只含前 i 种质因子}{f_{n/d}}^{-mu(d)}
]
枚举质数 (p_i)
当 (p_i mid n) 时，显然 (s_{i,n}=s_{i-1,n},inv_{i,n}=inv_{i-1,n})
当 (p_imid n) 时，需要考虑 (d) 中含因子 (p_i) 的个数
(p_i mid d) 的情况对 (s_{i,n}) 的贡献即为 (s_{i-1,n})
({p_i}^2mid d) 时 (mu(d)=0) ，因此这时不产生贡献
若 (p_imid d) 且 ({p_i}^2 mid d) ，设 (d=p_id') ，则贡献为
[prod_{d'mid (n/p_i)}{f_{n/p_id'}}^{mu(p_id')}
]
(p_i) 与 (d') 互质，因此 (mu(p_id')=mu(p_i)mu(d')=-mu(d'))
代回之后发现贡献就是 (inv_{i-1,n/p_i})
故 (s_{i,n}=s_{i-1,n} imes inv_{i-1,n/p_i})
显然 (inv_{i,n}=inv_{i-1,n} imes s_{i-1,n/p_i})
综上，转移方程为
[s_{i,n}=egin{cases}s_{i-1,n}&(p_i
mid n)\ s_{i-1,n} imes invs_{i-1,n/p_i} & (p_imid n)end{cases}
]
[inv_{i,n}=egin{cases}inv_{i-1,n}&(p_i
mid n)\ inv_{i-1,n} imes s_{i-1,n/p_i} &(p_imid n)end{cases}
]
初始状态为 (s_{0,n}=f_n,inv_{0,n}={f_n}^{-1})
但是直接求逆元会导致复杂度升到 (O(nlog P)) ，就前功尽弃了
所以还要用一下线性求逆元的方法
#include<stdio.h>
const int N=1000010,P=1e9+7; int T,n,m,ans,t1,t2;
int t,cnt,suf,I,p[100000],s[N],inv[N],v[N],pre[N];
inline int min(int x,int y) { return x<y?x:y; }
inline int power(int x,int y) {
int s=1;
while (y) (y&1)&&(s=1ll*s*x%P),x=1ll*x*x%P,y>>=1;
return s;
}
int main() {
s[1]=pre[0]=pre[1]=inv[0]=suf=1;
for (int i=2; i<N; ++i) { //线性筛、初始化 s
v[i]||(p[++cnt]=i);
for (int j=1; j<=cnt&&(t=p[j]*i)<N; ++j) {
v[t]=1;
if (i%p[j]==0) break;
}
(s[i]=s[i-2]+s[i-1])>=P&&(s[i]-=P);
pre[i]=1ll*pre[i-1]*s[i]%P;
}
I=power(pre[N-1],P-2);
for (int i=N-1; i; --i) //初始化 inv
inv[i]=1ll*pre[i-1]*suf%P*I%P,
suf=1ll*suf*s[i]%P;
for (int i=1,j; i<=cnt; ++i) // dp
for (j=(N-1)/p[i],t=j*p[i]; j; --j,t-=p[i])
s[t]=1ll*s[t]*inv[j]%P,
inv[t]=1ll*inv[t]*s[j]%P;
for (int i=2; i<N; ++i) s[i]=1ll*s[i]*s[i-1]%P;
for (int i=2; i<N; ++i) inv[i]=1ll*inv[i]*inv[i-1]%P;
//存储 s 和 inv 的前缀积，便于整除分块计算
for (scanf("%d",&T); T; --T) {
scanf("%d%d",&n,&m),ans=1;
for (int l=1,r,mn=min(n,m); l<=mn; l=r+1) //整除分块
r=min(n/(t1=n/l),m/(t2=m/l)),
t=1ll*s[r]*inv[l-1]%P,
ans=1ll*ans*power(t,1ll*t1*t2%(P-1))%P;
printf("%d
",ans);
}
return 0;
}
P5518 [MtOI2019]幽灵乐团 / 莫比乌斯反演基础练习题
这个题要麻烦得多
全程大力推式子 虽然非常不基础 但确实是针对性很强的练习题
[largeprod_{i=1}^Aprod_{j=1}^Bprod_{k=1}^Cleft(dfrac{operatorname{lcm}(i,j)}{gcd(i,k)}
ight)^{f(type)}
]
等于
[largeprod_{i=1}^Aprod_{j=1}^Bprod_{k=1}^Cleft(dfrac{ij}{gcd(i,j)gcd(i,k)}
ight)^{f(type)}
]
可以转化为以下两个子问题：
[large M(A,B,C)=prod_{i=1}^Aprod_{j=1}^Bprod_{k=1}^Ci^{f(type)}
]
[large N(A,B,C)=prod_{i=1}^Aprod_{j=1}^Bprod_{k=1}^Cgcd(i,j)^{f(type)}
]
答案为
[large dfrac{M(A,B,C) imes M(B,A,C)}{N(A,B,C) imes N(A,C,B)}
]
下面按 (type) 的取值分为三种情况讨论
(Large extbf{1. type=0,f(type)=1})
[large M(A,B,C)=prod_{i=1}^Aprod_{j=1}^Bprod_{k=1}^Ci=prod_{i=1}^Ai^{BC}=(A!)^{BC}
]
预处理阶乘，快速幂回答询问，复杂度为 (O(n+Tlog n)) （ (n) 与 (A,B,C) 同级， (T) 为询问组数，下同）
[large N(A,B,C)=prod_{i=1}^Aprod_{j=1}^Bprod_{k=1}^Cgcd(i,j)=left(prod_{i=1}^Aprod_{j=1}^Bgcd(i,j)
ight)^C
]
括号内的式子和上题基本一致 过程不详细写了
设 (D=min(A,B)) ，则
[large N(A,B,C)=prod_{t=1}^Dleft(prod_{d mid t}d^{mu(t/d)}
ight)^{(A/t)(B/t)C}
]
依然是预处理括号内的式子 然后整除分块
时间复杂度 (O(nlog n+Tsqrt nlog n))
预处理也可以做到 (O(nlog log n)) ，但对本题而言基本没啥优化效果
(Large extbf{2. type=1,f(type)=i×j×k})
[large M(A,B,C)=prod_{i=1}^Aprod_{j=1}^Bprod_{k=1}^Ci^{i imes j imes k}=prod_{i=1}^Ai^{i(sum_{j=1}^Bj)(sum_{k=1}^Ck)}
]
令 (S(n)=sumlimits_{i=1}^ni=dfrac{n(n+1)}2)
[large M(A,B,C)=left(prod_{i=1}^Ai^i
ight)^{S(B)S(C)}
]
括号内的式子可以 (O(nlog n)) 预处理，然后快速幂回答询问
[large N(A,B,C)=prod_{i=1}^Aprod_{j=1}^Bprod_{k=1}^Cgcd(i,j)^{i imes j imes k}=left(prod_{i=1}^Aprod_{j=1}^Bgcd(i,j)^{i imes j}
ight)^{S(C)}
]
括号里的式子拿出来推一下
[large egin{aligned}prod_{i=1}^Aprod_{j=1}^Bgcd(i,j)^{i imes j}
&=prod_{d=1}^Dprod_{i=1}^Aprod_{j=1}^Bd^{i imes j imes [(i,j)=d]}
=prod_{d=1}^Dprod_{i=1}^{A/d}prod_{j=1}^{B/d}d^{id imes jd imes [(i,j)=1]}
=prod_{d=1}^Dd^{d^2sum_{i=1}^{A/d}sum_{j=1}^{B/d}ij[(i,j)=1]}\
&=prod_{d=1}^Dd^{d^2sum_{i=1}^{A/d}sum_{j=1}^{B/d}ijsum_{tmid (i,j)}~mu(t)}
=prod_{d=1}^Dd^{d^2sum_{t=1}^{D/d}mu(t)sum_{i=1}^{A/dt}itsum_{j=1}^{B/dt}jt}
=prod_{d=1}^Dd^{d^2sum_{t=1}^{D/d}mu(t)t^2S(A/dt)S(B/dt)}\
&=prod_{d=1}^Dprod_{t=1}^{D/d}d^{mu(t)(dt)^2S(A/dt)S(B/dt)}
=prod_{t'=1}^Dprod_{dmid t'}d^{mu(t'/d)t'^2S(A/t')S(B/t')}
=prod_{t'=1}^Dleft(prod_{dmid t'}d^{mu(t'/d)}
ight)^{t'^2S(A/t')S(B/t')}
end{aligned}]
推到这里就可以了，代回 (N(A,B,C)) 的表达式
[large N(A,B,C)=prod_{t=1}^Dleft(prod_{dmid t}d^{mu(t/d)}
ight)^{t^2S(A/t)S(B/t)S(C)}
]
括号里的式子预处理过了 然后还是整除分块
时间复杂度 (O(nlog n+Tsqrt nlog n))
(Large extbf{3. type=2,f(type)=gcd(i,j,k)})
[large M(A,B,C)=prod_{i=1}^Aprod_{j=1}^Bprod_{k=1}^C i^{gcd(i,j,k)}=prod_{i=1}^Ai^{sum_{j=1}^Bsum_{k=1}^Cgcd(i,j,k)}
]
指数可以考虑欧拉反演
[large sum_{j=1}^Bsum_{k=1}^Cgcd(i,j,k)
=sum_{j=1}^Bsum_{k=1}^Csum_{dmid gcd(i,j,k)}varphi(d)
=sum_{dmid i}varphi(d)sum_{j=1}^Bsum_{k=1}^C[d|j][d|k]
=sum_{dmid i}varphi(d)(B/d)(C/d)]
设 (E=min(A,B,C))
[large M(A,B,C)=prod_{i=1}^Ai^{sum_{dmid i}varphi(d)(B/d)(C/d)}
=prod_{d=1}^Eleft(prod_{i'=1}^{A/d}i'd
ight)^{varphi(d)(B/d)(C/d)}
=prod_{d=1}^Eleft((A/d)! imes d^{A/d}
ight)^{varphi(d)(B/d)(C/d)}]
[large M(A,B,C)=prod_{d=1}^E(A/d)!^{varphi(d)(B/d)(C/d)} imes prod_{d=1}^Ed^{varphi(d)(A/d)(B/d)(C/d)}
]
两部分都可以整除分块做
接下来是本题最难处理的式子了
[large N(A,B,C)=prod_{i=1}^Aprod_{j=1}^Bprod_{k=1}^Cgcd(i,j)^{gcd(i,j,k)}
]
改为枚举 (gcd(i,j))
[large N(A,B,C)=prod_{d=1}^Dd^{sum_{i=1}^Asum_{j=1}^B[(i,j)=d]sum_{k=1}^Cgcd(d,k)}
]
指数还是常规反演
[large sum_{i=1}^Asum_{j=1}^B[(i,j)=d]=sum_{i=1}^{A/d}sum_{j=1}^{B/d}sum_{tmid (i,j)}mu(t)=sum_{t=1}^{D/d}mu(t)(A/dt)(B/dt)
]
[large sum_{k=1}^Cgcd(d,k)=sum_{k=1}^Csum_{t'mid (d,k)}varphi(t')=sum_{t'mid d}varphi(t')(C/t')
]
[large N(A,B,C)=prod_{d=1}^Dd^{sum_{t=1}^{D/d}mu(t)(A/dt)(B/dt)sum_{t'mid d}varphi(t')(C/t')}
]
枚举因数看起来不对劲，换一下枚举方式
[large N(A,B,C)=prod_{t'=1}^Eprod_{d'=1}^{D/t'}(d't')^{varphi(t')(C/t')sum_{t=1}^{D/d't'}mu(t)(A/d't't)(B/d't't)}
]
式子有点丑 所以替换一次字母（
[large N(A,B,C)=prod_{i=1}^Eprod_{j=1}^{D/i}(ij)^{varphi(i)(C/i)sum_{k=1}^{D/ij}mu(k)(A/ijk)(B/ijk)}
]
似乎很难继续推了
但通过看官方题解我们了解到一个技巧，把 (i,j) 拆开分别算贡献
[large N(A,B,C)=prod_{i=1}^Eprod_{j=1}^{D/i}i^{varphi(i)(C/i)sum_{k=1}^{D/ij}mu(k)(A/ijk)(B/ijk)} imes prod_{i=1}^Eprod_{j=1}^{D/i}j^{varphi(i)(C/i)sum_{k=1}^{D/ij}mu(k)(A/ijk)(B/ijk)}
]
对于第一部分，枚举 (t=jk)
[largeprod_{i=1}^Ei^{varphi(i)(C/i)sum_{t=1}^{D/i}~(A/it)(B/it)sum_{kmid t}mu(k)}
]
[largeprod_{i=1}^Ei^{varphi(i)(C/i)sum_{t=1}^{D/i}~(A/it)(B/it)[t=1]}
]
发现只有 (t=1) 的情况能对答案产生贡献
[largeprod_{i=1}^Ei^{varphi(i)(A/i)(B/i)(C/i)}
]
而我们在 (M(A,B,C)) 中也算出过相同的式子：
[large M(A,B,C)=prod_{d=1}^E(A/d)!^{varphi(d)(B/d)(C/d)} imes prod_{d=1}^Ed^{varphi(d)(A/d)(B/d)(C/d)}
]
所以直接约掉就好了
最后来处理第二部分
[large prod_{i=1}^Eprod_{j=1}^{D/i}j^{varphi(i)(C/i)sum_{k=1}^{D/ij}mu(k)(A/ijk)(B/ijk)}
]
[large prod_{i=1}^Eleft(prod_{j=1}^{D/i}j^{sum_{k=1}^{D/ij}mu(k)(A/ijk)(B/ijk)}
ight)^{varphi(i)(C/i)}
]
同样枚举 (t=jk)
[large prod_{i=1}^Eleft(prod_{t=1}^{D/i}prod_{jmid t}j^{mu(t/j)(A/it)(B/it)}
ight)^{varphi(i)(C/i)}
]
简单整理一下
[large prod_{i=1}^Eleft(prod_{t=1}^{D/i}left(prod_{dmid t}d^{mu(t/d)}
ight)^{(A/it)(B/it)}
ight)^{varphi(i)(C/i)}
]
最内层括号中的式子预处理过了
可以两层整除分块做，时间复杂度 (O(n^{3/4}log n))
还有个 (O(n^{2/3}log n)) 的做法，虽然我没有写，而且也证不来复杂度，但还是提一下
之前我们算过
[largeprod_{i=1}^Aprod_{j=1}^Bgcd(i,j)=prod_{t=1}^Dleft(prod_{dmid T}d^{mu(t/d)}
ight)^{(A/t)(B/t)}
]
带入到现在的式子里就是
[largeprod_{i=1}^Eleft(prod_{j=1}^{A/i}prod_{k=1}^{B/i}gcd(j,k)
ight)^{varphi(i)(C/i)}
]
预处理 (A,Ble n^{2/3}) 的所有 (prodlimits_{i=1}^Aprodlimits_{j=1}^Bgcd(i,j))
(i) 需要一层整除分块
对于括号内的部分，当 (ige n^{1/3}) 时直接用预处理的值即可，当 (i<n^{1/3}) 时需要再套一层整除分块
时间复杂度证明可类比杜教筛 我都不会证
至此这道题的推式子部分终于结束了
写 (O(nlog n+Tn^{3/4}log n)) 或者 (O(n^{4/3}+Tn^{2/3}log n)) 的做法都能通过
细节很多 需要耐心调试
但很小的错误都能导致答案产生巨大误差 所以能过样例就基本就没问题了 调试难度不大
象征性地放一下代码（虽然丑到只有编译器看得懂
#include<stdio.h>
#define M0(A,B) power(fac[A],1ll*B*C%P2)
#define M1(A,B) power(pow[A],1ll*S1(B)*S1(C)%P2)
const int N=100001; int A,B,C,D,E,T,P,P2,t,cnt,ans,p[10000];
int s[N],inv[N],S[N],INV[N],fac[N],pow[N],phi[N],v[N];
inline int min(int x,int y) { return x<y?x:y; }
inline int S1(int x) { return (1ll*x*(x+1)>>1)%P2; }
inline int power(int x,int y) {
int ans=1;
while (y) (y&1)&&(ans=1ll*ans*x%P),x=1ll*x*x%P,y>>=1;
return ans;
}
void prework() {
phi[1]=fac[0]=pow[0]=S[0]=INV[0]=1;
s[0]=s[1]=inv[0]=inv[1]=1;
for (int i=2; i<N; ++i) {
v[i]||(p[++cnt]=i,phi[i]=i-1);
for (int j=1; j<=cnt&&(t=p[j]*i)<N; ++j) {
v[t]=1,phi[t]=phi[i]*(p[j]-1);
if (i%p[j]==0) { phi[t]+=phi[i]; break; }
}
(phi[i]+=phi[i-1])>=P2&&(phi[i]-=P2);
}
for (int i=2; i<N; ++i) s[i]=i,inv[i]=1ll*(P-P/i)*inv[P%i]%P;
for (int i=1; i<=cnt; ++i)
for (int j=(N-1)/p[i],t=j*p[i]; j; --j,t-=p[i])
s[t]=1ll*s[t]*inv[j]%P,inv[t]=1ll*inv[t]*s[j]%P;
for (int i=1; i<N; ++i)
S[i]=1ll*S[i-1]*power(s[i],1ll*i*i%P2)%P,
INV[i]=1ll*INV[i-1]*power(inv[i],1ll*i*i%P2)%P,
s[i]=1ll*s[i-1]*s[i]%P,inv[i]=1ll*inv[i-1]*inv[i]%P,
fac[i]=1ll*fac[i-1]*i%P,pow[i]=1ll*pow[i-1]*power(i,i)%P;
}
int N0(int A,int B,int C) {
int D=min(A,B),ans=1;
for (int l=1,r,t1,t2; l<=D; l=r+1)
r=min(A/(t1=A/l),B/(t2=B/l)),
ans=1ll*ans*power(1ll*s[r]*inv[l-1]%P,1ll*t1*t2%P2)%P;
return power(ans,P2-C);
}
int N1(int B,int C) {
D=min(A,B),ans=1;
for (int l=1,r,t,t1,t2; l<=D; l=r+1)
r=min(A/(t1=A/l),B/(t2=B/l)),
t=1ll*S1(t1)*S1(t2)%P2,
ans=1ll*ans*power(1ll*S[r]*INV[l-1]%P,t)%P;
return power(ans,P2-S1(C));
}
int M2(int A,int B) {
E=min(min(A,B),C),ans=1;
for (int l=1,r,t,t1,t2,t3; l<=E; l=r+1)
r=min(min(A/(t1=A/l),B/(t2=B/l)),C/(t3=C/l)),
t=1ll*t2*t3%P2*(phi[r]-phi[l-1]+P2)%P2,
ans=1ll*ans*power(fac[t1],t)%P;
return ans;
}
int N2(int B,int C) {
E=min(D=min(A,B),C),ans=1;
for (int l=1,r,t1,t2,t3; l<=E; l=r+1)
r=min(min(A/(t1=A/l),B/(t2=B/l)),C/(t3=C/l)),
ans=1ll*ans*N0(t1,t2,1ll*t3*(phi[r]-phi[l-1]+P2)%P2)%P;
return ans;
}
int main() {
scanf("%d%d",&T,&P),P2=P-1,prework();
while (T--)
scanf("%d%d%d",&A,&B,&C),t=1ll*M0(A,B)*M0(B,A)%P,
printf("%d ",1ll*t*N0(A,B,C)%P*N0(A,C,B)%P),
printf("%d ",1ll*M1(A,B)*M1(B,A)%P*N1(B,C)%P*N1(C,B)%P),
printf("%d
",1ll*M2(A,B)*M2(B,A)%P*N2(B,C)%P*N2(C,B)%P);
return 0;
}