P5518-[MtOI2019]幽灵乐团【莫比乌斯反演,欧拉反演】
正題
題目鏈接:https://www.luogu.com.cn/problem/P5518
題目大意
TTT次給出A,B,CA,B,CA,B,C求以下三個式子
∏i=1A∏j=1B∏k=1Clcm(i,j)gcd(i,k)\prod_{i=1}^A\prod_{j=1}^B\prod_{k=1}^{C}\frac{lcm(i,j)}{gcd(i,k)}i=1∏A?j=1∏B?k=1∏C?gcd(i,k)lcm(i,j)?
∏i=1A∏j=1B∏k=1C(lcm(i,j)gcd(i,k))i×j×k\prod_{i=1}^A\prod_{j=1}^B\prod_{k=1}^{C}\left(\frac{lcm(i,j)}{gcd(i,k)}\right)^{i\times j\times k}i=1∏A?j=1∏B?k=1∏C?(gcd(i,k)lcm(i,j)?)i×j×k
∏i=1A∏j=1B∏k=1C(lcm(i,j)gcd(i,k))gcd(i,j,k)\prod_{i=1}^A\prod_{j=1}^B\prod_{k=1}^{C}\left(\frac{lcm(i,j)}{gcd(i,k)}\right)^{gcd(i,j,k)}i=1∏A?j=1∏B?k=1∏C?(gcd(i,k)lcm(i,j)?)gcd(i,j,k)
1≤T≤70,1≤A,B,C≤1051\leq T\leq 70,1\leq A,B,C\leq 10^51≤T≤70,1≤A,B,C≤105
解題思路
開始寫了個O(Tnlog?n)O(Tn\log n)O(Tnlogn)結果發現不能過,然后就多浪費了三個多小時
只需要用到兩個反演的式子
F(i)=∑i∣df(d)?f(i)=∑i∣dμ(di)F(d)F(i)=\sum_{i|d}f(d)\Rightarrow f(i)=\sum_{i|d}\mu(\fracze8trgl8bvbq{i})F(d)F(i)=i∣d∑?f(d)?f(i)=i∣d∑?μ(id?)F(d)
gcd(S)=∑?x∈S,d∣xφ(d)gcd(S)=\sum_{\forall x\in S,d|x}\varphi(d)gcd(S)=?x∈S,d∣x∑?φ(d)
然后因為推導過程出來冗長以外沒有太多難的部分,所以推薦自己手推到不會的再翻題解。
然后就開始吧,因為lcm(i,j)=ijgcd(i,j)lcm(i,j)=\frac{ij}{gcd(i,j)}lcm(i,j)=gcd(i,j)ij?所以問題可以化為兩個部分,求
(∏i=1A∏j=1B∏k=1Cij)f(type),(∏i=1A∏j=1B∏k=1C1gcd(i,j))f(type)\left(\prod_{i=1}^A\prod_{j=1}^B\prod_{k=1}^Cij\right)^{f(type)},\left(\prod_{i=1}^A\prod_{j=1}^B\prod_{k=1}^C\frac{1}{gcd(i,j)}\right)^{f(type)}(i=1∏A?j=1∏B?k=1∏C?ij)f(type),(i=1∏A?j=1∏B?k=1∏C?gcd(i,j)1?)f(type)
首先是第一個式子f(type)=1f(type)=1f(type)=1。
第一部分就是
∏i=1AiB×C×∏j=1BjA×C\prod_{i=1}^Ai^{B\times C}\times \prod_{j=1}^Bj^{A\times C}i=1∏A?iB×C×j=1∏B?jA×C
這個十分簡單,我們預處理階乘就可以做到O(log?P)O(\log P)O(logP)
然后第二部分考慮枚舉約數
∏d=11dC×∑i=1A∑j=1B[gcd(i,j)=d]\prod_{d=1}\frac{1}ze8trgl8bvbq^{C\times \sum_{i=1}^A\sum_{j=1}^B[gcd(i,j)=d]}d=1∏?d1?C×∑i=1A?∑j=1B?[gcd(i,j)=d]
然后莫反
∏d=11dC∑k=1μ(k)?Akd??Bkd?\prod_{d=1}\frac{1}ze8trgl8bvbq^{C\sum_{k=1}\mu(k)\lfloor\frac{A}{kd}\rfloor\lfloor\frac{B}{kd}\rfloor}d=1∏?d1?C∑k=1?μ(k)?kdA???kdB??
顯然的我們可以d,kd,kd,k都可以整除分塊,預處理一下逆元的前綴乘積就可以快速計算區間逆元乘積了。
第一個式子時間復雜度O(n34)O(n^\frac{3}{4})O(n43?)
然后第二個式子類似的,記S(n)=n×(n+1)2S(n)=\frac{n\times (n+1)}{2}S(n)=2n×(n+1)?那么有
∏i=1AiiS(B)S(C)×∏j=1BjjS(A)(C)\prod_{i=1}^Ai^{iS(B)S(C)}\times \prod_{j=1}^Bj^{jS(A)(C)}i=1∏A?iiS(B)S(C)×j=1∏B?jjS(A)(C)
維護一個iii^iii的前綴積即可。
第二部分
∏d=11dS(C)×∑i=1A∑j=1Bij[gcd(i,j)=d]?∏d=11dS(C)∑k=1μ(k)S(?Akd?)S(?Bkd?)\prod_{d=1}\frac{1}ze8trgl8bvbq^{S(C)\times \sum_{i=1}^A\sum_{j=1}^Bij[gcd(i,j)=d]}\Rightarrow \prod_{d=1}\frac{1}ze8trgl8bvbq^{S(C)\sum_{k=1}\mu(k)S(\lfloor\frac{A}{kd}\rfloor)S(\lfloor\frac{B}{kd}\rfloor)}d=1∏?d1?S(C)×∑i=1A?∑j=1B?ij[gcd(i,j)=d]?d=1∏?d1?S(C)∑k=1?μ(k)S(?kdA??)S(?kdB??)
需要提前處理μ(i)×i2\mu(i)\times i^2μ(i)×i2的前綴和,i2,1i2i^2,\frac{1}{i^2}i2,i21?的前綴積就好了
第二個式子時間復雜度O(n34)O(n^\frac{3}{4})O(n43?)
主要的難點在第三個式子
首先第一部分先只考慮iii
∏i=1A∏j=1B∏k=1Cigcd(i,j)\prod_{i=1}^A\prod_{j=1}^B\prod_{k=1}^Ci^{gcd(i,j)}i=1∏A?j=1∏B?k=1∏C?igcd(i,j)
考慮歐拉反演,枚舉約數
∏d=1(∏i?Ad?i×d)φ(d)?Bd??Cd?\prod_{d=1}\left(\prod_{i}^{\lfloor\frac{A}ze8trgl8bvbq\rfloor}i\times d\right)^{\varphi(d)\lfloor\frac{B}ze8trgl8bvbq\rfloor\lfloor\frac{C}ze8trgl8bvbq\rfloor}d=1∏????i∏?dA???i×d???φ(d)?dB???dC??
設fi=∏j=1iif_i=\prod_{j=1}^iifi?=∏j=1i?i那么有
∏d=1(f?Ad?d?Ad?)φ(d)?Bd??Cd?\prod_{d=1}\left(f_{\lfloor\frac{A}ze8trgl8bvbq\rfloor}d^{\lfloor\frac{A}ze8trgl8bvbq\rfloor}\right)^{\varphi(d)\lfloor\frac{B}ze8trgl8bvbq\rfloor\lfloor\frac{C}ze8trgl8bvbq\rfloor}d=1∏?(f?dA???d?dA??)φ(d)?dB???dC??
拆開來fff和ddd的部分
∏d=1(f?Ad?)φ(d)?Bd??Cd?×dφ(d)?Ad??Bd??Cd?\prod_{d=1}(f_{\lfloor\frac{A}ze8trgl8bvbq\rfloor})^{\varphi(d)\lfloor\frac{B}ze8trgl8bvbq\rfloor\lfloor\frac{C}ze8trgl8bvbq\rfloor}\times d^{\varphi(d)\lfloor\frac{A}ze8trgl8bvbq\rfloor\lfloor\frac{B}ze8trgl8bvbq\rfloor\lfloor\frac{C}ze8trgl8bvbq\rfloor}d=1∏?(f?dA???)φ(d)?dB???dC??×dφ(d)?dA???dB???dC??
預處理出fff數組,和ggg數組gi=∏j=1ijφ(j)g_i=\prod_{j=1}^i j^{\varphi(j)}gi?=∏j=1i?jφ(j),φ\varphiφ的前綴和就可以整除分塊搞了
然后是第二部分
∏i=1A∏j=1B∏k=1C1gcd(i,j)gcd(i,j,k)\prod_{i=1}^A\prod_{j=1}^B\prod_{k=1}^C\frac{1}{gcd(i,j)}^{gcd(i,j,k)}i=1∏A?j=1∏B?k=1∏C?gcd(i,j)1?gcd(i,j,k)
你可以試一下枚舉gcd(i,j)gcd(i,j)gcd(i,j),反正我推了好久都沒有對出來/kk
所以考慮枚舉gcd(i,j,k)gcd(i,j,k)gcd(i,j,k)的約數然后歐拉反演里面再莫反
∏d=1(∏k=11dk∑z=1μ(z)?Adkz??Bdkz?)φ(d)?Cd?\prod_{d=1}\left(\prod_{k=1}\frac{1}{dk}^{\sum_{z=1}\mu(z)\lfloor\frac{A}{dkz}\rfloor\lfloor\frac{B}{dkz}\rfloor}\right)^{\varphi(d)\lfloor\frac{C}ze8trgl8bvbq\rfloor}d=1∏?(k=1∏?dk1?∑z=1?μ(z)?dkzA???dkzB??)φ(d)?dC??
然后把1d\frac{1}ze8trgl8bvbqd1?和1k\frac{1}{k}k1?分開處理
∏d=1(∏k=11k∑z=1μ(z)?Adkz??Bdkz?)φ(d)?Cd?\prod_{d=1}\left(\prod_{k=1}\frac{1}{k}^{\sum_{z=1}\mu(z)\lfloor\frac{A}{dkz}\rfloor\lfloor\frac{B}{dkz}\rfloor}\right)^{\varphi(d)\lfloor\frac{C}ze8trgl8bvbq\rfloor}d=1∏?(k=1∏?k1?∑z=1?μ(z)?dkzA???dkzB??)φ(d)?dC??
×∏d=1(1d)φ(d)?Cd?∑k=1∑z=1μ(z)?Adkz??Bdkz?\times \prod_{d=1}\left(\frac{1}ze8trgl8bvbq\right)^{\varphi(d)\lfloor\frac{C}ze8trgl8bvbq\rfloor\sum_{k=1}\sum_{z=1}\mu(z)\lfloor\frac{A}{dkz}\rfloor\lfloor\frac{B}{dkz}\rfloor}×d=1∏?(d1?)φ(d)?dC??∑k=1?∑z=1?μ(z)?dkzA???dkzB??
(至于為什么要這么分開,會注意到第一個式子如果我們進行兩次整除分塊,那么對于每個1dk\frac{1}{dk}dk1?就會有兩個影響它的指數φ(d)\varphi(d)φ(d)和∑z=1μ(z)?Adkz??Bdkz?\sum_{z=1}\mu(z)\lfloor\frac{A}{dkz}\rfloor\lfloor\frac{B}{dkz}\rfloor∑z=1?μ(z)?dkzA???dkzB??所以很難處理,此時我們分開來就可以直接用區間的值來做了)
注意到這樣要做三次整除分塊,十分地慢,但是考慮后面那個式子∑z=1μ(z)?Adkz??Bdkz?\sum_{z=1}\mu(z)\lfloor\frac{A}{dkz}\rfloor\lfloor\frac{B}{dkz}\rfloor∑z=1?μ(z)?dkzA???dkzB??對于一個固定的dkdkdk是一個確定的值的,并且因為是整除分塊,所以不同的值不多我們可以考慮預處理這個東西,整除分塊一個dkdkdk然后里面再整除分塊計算就好了。
還有拆開后處理1d\frac{1}ze8trgl8bvbqd1?的式子需要預處理1dφ(d)\frac{1}ze8trgl8bvbq^{\varphi(d)}d1?φ(d)的前綴積。
然后就做完了,第三部分時間復雜度O(n+n34log?n)O(n+n^{\frac{3}{4}}\log n)O(n+n43?logn)
實際上其實最后的式子兩個部分可以約掉一些東西省一些常數,但是這樣做也能過但是得開int。
或者你可以用別的方法卡卡常反正我開int了。
code
因為中途long long改int所以代碼巨丑
#include<cstdio> #include<cstring> #include<algorithm> //#define int long long using namespace std; const int N=1e5+10; int T,P,Phi,A,B,C,ans,cnt;bool v[N]; int pri[N/10],mu[N],su[N],phi[N],f[N],g[N],h[N]; int inv[N],inc[N],inw[N],fac[N],gac[N],hac[N]; int power(int x,int b){int ans=1;b=(b%Phi+Phi)%Phi;while(b){if(b&1)ans=ans*1ll*x%P;x=x*1ll*x%P;b>>=1;}return ans; } int Sum(int n) {return n*1ll*(n+1)/2%Phi;} int Tum(int n) {return n*1ll*(n+1)*1ll*(2ll*1ll*n+1)/6ll%Phi;} void Prime(){mu[1]=phi[1]=1;for(int i=2;i<N;i++){if(!v[i])pri[++cnt]=i,mu[i]=-1,phi[i]=i-1;for(int j=1;j<=cnt&&i*1ll*pri[j]<N;j++){v[i*1ll*pri[j]]=1;if(i%pri[j]==0){phi[i*1ll*pri[j]]=phi[i]*1ll*pri[j];break;}mu[i*1ll*pri[j]]=-mu[i];phi[i*1ll*pri[j]]=phi[i]*1ll*(pri[j]-1);}}inc[1]=fac[0]=gac[0]=inw[0]=inv[0]=hac[0]=1;for(int i=2;i<N;i++)inc[i]=P-inc[P%i]*1ll*(P/i)%P;for(int i=1;i<N;i++)fac[i]=fac[i-1]*1ll*i%P;for(int i=1;i<N;i++)gac[i]=gac[i-1]*1ll*power(i,i)%P;for(int i=1;i<N;i++)hac[i]=hac[i-1]*1ll*power(i,i*1ll*i%Phi)%P;for(int i=1;i<N;i++)inw[i]=inw[i-1]*1ll*power(inc[i],i*1ll*i%Phi)%P;for(int i=1;i<N;i++)inv[i]=inv[i-1]*1ll*inc[i]%P;for(int i=1;i<N;i++)su[i]=su[i-1]+mu[i];g[0]=h[0]=1;for(int i=1;i<N;i++)g[i]=g[i-1]*1ll*power(inc[i],phi[i])%P;for(int i=1;i<N;i++)h[i]=h[i-1]*1ll*power(i,phi[i])%P;for(int i=1;i<N;i++)(phi[i]+=phi[i-1])%=Phi;return; } void qart1(int A,int B,int C){int n=min(A,B);for(int L=1,R;L<=n;L=R+1){R=min(A/(A/L),B/(B/L));int a=A/L,b=B/L,m=min(a,b),w=0;for(int l=1,r;l<=m;l=r+1){r=min(a/(a/l),b/(b/l));(w+=(su[r]-su[l-1])*(a/l)*1ll*(b/l)%Phi)%=Phi;}ans=ans*1ll*power(inv[R]*1ll*fac[L-1]%P,w*1ll*C%Phi)%P;}return; } void part1(){ans=power(fac[A],B*1ll*C%Phi)*1ll*power(fac[B],A*1ll*C%Phi)%P;qart1(A,B,C);qart1(A,C,B);printf("%d ",ans);return; } void qart2(int A,int B,int C){int n=min(A,B);for(int i=1;i<=n;i++)f[i]=(f[i-1]+mu[i]*1ll*i*1ll*i%Phi)%Phi;for(int L=1,R;L<=n;L=R+1){R=min(A/(A/L),B/(B/L));int a=A/L,b=B/L,m=min(a,b),w=0;for(int l=1,r;l<=m;l=r+1){r=min(a/(a/l),b/(b/l));(w+=Sum(a/l)*1ll*Sum(b/l)%Phi*1ll*(f[r]-f[l-1])%Phi)%=Phi;}w=w*1ll*Sum(C)%Phi;ans=ans*1ll*power(inw[R]*1ll*hac[L-1]%P,w)%P;}return; } void part2(){ans=power(gac[A],Sum(B)*1ll*Sum(C)%Phi)*1ll*power(gac[B],Sum(A)*1ll*Sum(C)%Phi)%P;qart2(A,B,C);qart2(A,C,B);printf("%d ",ans);return; } void qart3(int A,int B,int C){int n=min(A,B);for(int L=1,R;L<=n;L=R+1){R=min(A/(A/L),B/(B/L));int a=A/L,b=B/L,m=min(a,b),w=0;for(int l=1,r;l<=m;l=r+1){r=min(a/(a/l),b/(b/l));(w+=1ll*(a/l)*(b/l)*(su[r]-su[l-1])%Phi)%=Phi;}for(int i=L;i<=R;i++)f[i]=w;}n=min(n,C);for(int L=1,R;L<=n;L=R+1){R=min(A/(A/L),min(B/(B/L),C/(C/L)));int a=A/L,b=B/L,c=C/L,m=min(a,b),w=0,s=1;for(int l=1,r;l<=m;l=r+1){r=min(a/(a/l),b/(b/l));(w+=1ll*f[l*1ll*L]*(r-l+1)%Phi)%=Phi; s=1ll*s*power(1ll*inv[r]*fac[l-1]%P,f[l*1ll*L])%P;}ans=ans*1ll*power(1ll*g[R]*power(g[L-1],P-2)%P,1ll*w*(C/L)%Phi)%P;ans=ans*1ll*power(s,1ll*(phi[R]-phi[L-1])*(C/L)%Phi)%P;} // for(int i=1;i<=A;i++) // for(int j=1;j<=B;j++) // for(int k=1;k<=C;k++) // ans=ans*1ll*power(inc[__gcd(i,j)],__gcd(__gcd(i,j),k))%P;return; } void part3(){ans=1;int n=min(A,min(B,C));for(int L=1,R;L<=n;L=R+1){R=min(A/(A/L),min(B/(B/L),C/(C/L)));int a=A/L,b=B/L,c=C/L; // b=fac[b]*1ll*power(Sum(R)-Sum(L-1),b)%P;ans=1ll*ans*power(fac[a],1ll*(phi[R]-phi[L-1])*b*c%Phi)%P;ans=1ll*ans*power(fac[b],1ll*(phi[R]-phi[L-1])*a*c%Phi)%P;ans=1ll*ans*power(h[R]*1ll*power(h[L-1],P-2)%P,2ll*a*b*c%Phi)%P;} // int k=ans;qart3(A,B,C);qart3(A,C,B); // else ans=1ll*ans*ans%P*power(k,P-2)%P;printf("%d ",ans);return; } signed main() {scanf("%d%d",&T,&P);Phi=P-1;Prime();while(T--){scanf("%d%d%d",&A,&B,&C); // A=B=C=1e5;part1();part2();part3();putchar('\n');}return 0; }總結
以上是生活随笔為你收集整理的P5518-[MtOI2019]幽灵乐团【莫比乌斯反演,欧拉反演】的全部內容,希望文章能夠幫你解決所遇到的問題。
- 上一篇: 彩排的意思 何谓彩排
- 下一篇: 邢菲个人资料 邢菲个人资料简介