洛谷 P5518 [MtOI2019]幽靈樂團 題解

Description

令：

\[f(\textit{type})=\begin{cases}1&\textit{type}=0\\ijk&\textit{type}=1\\\gcd(i,j,k)&\textit{type}=2\end{cases} \]

給定 \(A,B,C\)，求：

\[\textit{ans}_{\textit{type}}=\prod_{i=1}^A\prod_{j=1}^B\prod_{k=1}^C\left(\frac{\operatorname{lcm}(i,j)}{\gcd(i,k)}\right)^{f(\textit{type})} \]

多測，對給定質數 \(p\)

取模。

\(1\le A,B,C\le 10^5\)，\(10^7\le p\le 1.05\times10^9,p\in\mathbb{P}\)，\(T=70\)。

Solution

首先對原式化簡一下：

\[\begin{aligned} \textit{ans}_{\textit{type}}=&\prod_{i=1}^A\prod_{j=1}^B\prod_{k=1}^C\left(\frac{\operatorname{lcm}(i,j)}{\gcd(i,k)}\right)^{f(\textit{type})}\\ =&\prod_{i=1}^A\prod_{j=1}^B\prod_{k=1}^C\left(\frac{ij}{\gcd(i,j)\cdot\gcd(i,k)}\right)^{f(\textit{type})} \end{aligned} \]

於是原式可以分為兩個部分：

\[\begin{aligned} & \mathcal F(A,B,C) = \prod_{i=1}^A\prod_{j=1}^B\prod_{k=1}^Ci^{f(\textit{type})}\\ & \mathcal G(A,B,C)=\prod_{i=1}^A\prod_{j=1}^B\prod_{k=1}^C\gcd(i,j)^{f(\textit{type})} \end{aligned} \]

顯然有：

\[\textit{ans}_{\textit{type}}=\frac{\mathcal F(A,B,C)\cdot\mathcal F(B,A,C)}{\mathcal G(A,B,C)\cdot\mathcal G(A,C,B)} \]

接下來就開始推柿子吧！

---

\(\textit{type}=0\)

\[\begin{aligned} \mathcal F(A,B,C)=&\prod_{i=1}^A\prod_{j=1}^B\prod_{k=1}^Ci\\ =&\prod_{i=1}^A\prod_{j=1}^Bi^C\\ =&\prod_{i=1}^Ai^{B\cdot C}\\ =&(A!)^{B\cdot C} \end{aligned} \]

預處理階乘，每次計算快速冪即可。

\[\begin{aligned} \mathcal G(A,B,C)=&\prod_{i=1}^A\prod_{j=1}^B\prod_{k=1}^C\gcd(i,j)\\ =&\left(\prod_{i=1}^A\prod_{j=1}^B\gcd(i,j)\right)^C\\ \end{aligned} \]

然後開始莫反：

\[\begin{aligned} \prod_{i=1}^A\prod_{j=1}^B\gcd(i,j)=&\prod_{d=1}\prod_{i=1}^A\prod_{j=1}^Bd^{[\gcd(i,j)=d]}\\ =&\prod_{d=1}\prod_{i=1}^{\lfloor \frac Ad\rfloor}\prod_{j=1}^{\lfloor\frac Bd\rfloor}d^{[gcd(i,j)=1]}\\ =&\prod_{d=1}d^{\sum_{i=1}^{A/d}\sum_{j=1}^{B/d}[gcd(i,j)=1]} \end{aligned} \]

把指數提出來推一下：

\[\begin{aligned} \sum_{i=1}^{\lfloor\frac Ad\rfloor}\sum_{j=1}^{\lfloor\frac Bd\rfloor}[gcd(i,j)=1]=&\sum_{e=1}\sum_{i=1}^{\lfloor\frac Ad\rfloor}\sum_{j=1}^{\lfloor\frac Bd\rfloor}[e\mid i][e \mid j]\mu(e)\\ =&\sum_{e=1}\sum_{i=1}^{\lfloor \frac A{de}\rfloor}\sum_{j=1}^{\lfloor\frac B{de}\rfloor}\mu(e)\\ =&\sum_{e=1}\lfloor \frac A{de}\rfloor\lfloor\frac B{de}\rfloor\mu(e) \end{aligned} \]

帶回去就有：

\[\begin{aligned} \prod_{i=1}^A\prod_{j=1}^B\gcd(i,j)=&\prod_{d=1}d^{\sum_{e=1}\lfloor \frac A{de}\rfloor\lfloor\frac B{de}\rfloor\mu(e)}\\ =&\prod_{t=1}\prod_{d\mid t}d^{\lfloor\frac At\rfloor\lfloor\frac Bt\rfloor\mu(\frac td)}\\ =&\prod_{t=1}\left(\prod_{d\mid t}d^{\mu(\frac td)}\right)^{\lfloor\frac At\rfloor\lfloor\frac Bt\rfloor}\\ \mathcal G(A,B,C)=&\prod_{i=1}^A\prod_{j=1}^B\prod_{k=1}^C\gcd(i,j)\\ =&\left(\prod_{t=1}\left(\prod_{d\mid t}d^{\mu(\frac td)}\right)^{\lfloor\frac At\rfloor\lfloor\frac Bt\rfloor}\right)^C\\ \end{aligned} \]

最裡面括號內的東西可以用類似埃氏篩的方法預處理一下，即對於每個 \(d\)，將它的貢獻累乘到它的每個倍數上。預處理完成後就可以快樂的整除分塊了。

---

\(\textit{type}=1\)

定義：\(S(n)=\sum_{i=1}^ni=\frac{n(n+1)}2\)。

然後開始推柿子：

\[\begin{aligned} \mathcal F(A,B,C)=&\prod_{i=1}^A\prod_{j=1}^B\prod_{k=1}^Ci^{ijk}\\ =&\left(\prod_{i=1}^A\prod_{j=1}^Bi^{ij}\right)^{S(C)}\\ =&\left(\prod_{i=1}^Ai^i\right)^{S(B)\cdot S(C)} \end{aligned} \]

預處理 \(\prod_{i=1}^Ai^i\) 即可。

\[\begin{aligned} \mathcal G(A,B,C)=&\prod_{i=1}^A\prod_{j=1}^B\prod_{k=1}^C\gcd(i,j)^{ijk}\\ =&\left(\prod_{i=1}^A\prod_{j=1}^B\gcd(i,j)^{ij}\right)^{S(C)}\\ \prod_{i=1}^A\prod_{j=1}^B\gcd(i,j)^{ij}=&\prod_{d=1}\prod_{i=1}^A\prod_{j=1}^Bd^{[\gcd(i,j)=d]ij}\\ =&\prod_{d=1}\prod_{i=1}^{\lfloor\frac Ad\rfloor}\prod_{j=1}^{\lfloor\frac Bd\rfloor}d^{[\gcd(i,j)=1]ijd^2}\\ =&\prod_{d=1}d^{d^2\sum_{i=1}^{A/d}\sum_{j=1}^{B/d}[\gcd(i,j)=1]ij} \end{aligned} \]

把指數拎出來：

\[\begin{aligned} \sum_{i=1}^{\lfloor\frac Ad\rfloor}\sum_{j=1}^{\lfloor\frac Bd\rfloor}[\gcd(i,j)=1]ij=&\sum_{e=1}\sum_{i=1}^{\lfloor\frac Ad\rfloor}\sum_{j=1}^{\lfloor\frac Bd\rfloor}[e\mid i][e\mid j]\mu(e)ij\\ =&\sum_{e=1}\sum_{i=1}^{\lfloor\frac A{de}\rfloor}\sum_{j=1}^{\lfloor\frac B{de}\rfloor}e^2\mu(e)ij\\ =&\sum_{e=1}e^2\mu(e)\left(\sum_{i=1}^{\lfloor\frac A{de}\rfloor}\sum_{j=1}^{\lfloor\frac B{de}\rfloor}ij\right)\\ =&\sum_{e=1}e^2\mu(e)S(\lfloor\frac A{de}\rfloor)S(\lfloor\frac B{de}\rfloor) \end{aligned} \]

帶回去：

\[\begin{aligned} \prod_{i=1}^A\prod_{j=1}^B\gcd(i,j)^{ij}=&\prod_{d=1}d^{d^2\sum_{e=1}e^2\mu(e)S(\lfloor\frac A{de}\rfloor)S(\lfloor\frac B{de}\rfloor)}\\ =&\prod_{e=1}\prod_{d=1}d^{(de)^2\mu(e)S(\lfloor\frac A{de}\rfloor)S(\lfloor\frac B{de}\rfloor)}\\ =&\prod_{t=1}\prod_{d\mid t}d^{t^2\mu(\frac td)S(\lfloor\frac A{t}\rfloor)S(\lfloor\frac B{t}\rfloor)}\\ =&\prod_{t=1}\left(\prod_{d\mid t}d^{\mu(\frac td)}\right)^{t^2S(\lfloor\frac At\rfloor)S(\lfloor\frac Bt\rfloor)}\\ \mathcal G(A,B,C)=&\left(\prod_{i=1}^A\prod_{j=1}^B\gcd(i,j)^{ij}\right)^{S(C)}\\ =&\left(\left(\prod_{t=1}\left(\prod_{d\mid t}d^{\mu(\frac td)}\right)^{t^2}\right)^{S(\lfloor\frac At\rfloor)S(\lfloor\frac Bt\rfloor)}\right)^{S(C)} \end{aligned} \]

括號裡的東西其實就是第一問預處理的那個東西的 \(t^2\) 次方，預處理一下即可整除分塊。

---

\(\textit{type}=2\)

毒瘤。

\[\begin{aligned} \mathcal F(A,B,C)=&\prod_{i=1}^A\prod_{j=1}^B\prod_{k=1}^Ci^{\gcd(i,j,k)}\\ =&\prod_{i=1}^Ai^{\sum_{j=1}^B\sum_{k=1}^C\gcd(i,j,k)}\\ =&\prod_{d=1}\prod_{i=1}^Ai^{\sum_{j=1}^B\sum_{k=1}^C[\gcd(i,j,k)=d]d}\\ =&\prod_{d=1}\prod_{i=1}^{\lfloor\frac Ad\rfloor}(id)^{\sum_{j=1}^{B/d}\sum_{k=1}^{C/d}[\gcd(i,j,k)=1]d}\\ =&\prod_{d=1}\prod_{e=1}\prod_{i=1}^{\lfloor\frac A{de}\rfloor}(ide)^{\lfloor{\frac B{de}}\rfloor\lfloor{\frac C{de}}\rfloor\mu(e)d}\\ =&\prod_{t=1}\prod_{d\mid t}\prod_{i=1}^{\lfloor\frac At\rfloor}(it)^{\lfloor\frac Bt\rfloor\lfloor\frac Ct\rfloor\mu(\frac td)d}\\ =&\prod_{t=1}\left(\prod_{d\mid t}\left(\lfloor\frac At\rfloor!\cdot t^{\lfloor\frac At\rfloor}\right)^{\mu(\frac td)d}\right)^{\lfloor\frac Bt\rfloor\lfloor\frac Ct\rfloor}\\ =&\prod_{t=1}\left(\left(\lfloor\frac At\rfloor!\cdot t^{\lfloor\frac At\rfloor}\right)^{\sum_{d\mid t}\mu(\frac td)d}\right)^{\lfloor\frac Bt\rfloor\lfloor\frac Ct\rfloor}\\ =&\prod_{t=1}\left(\left(\lfloor\frac At\rfloor!\cdot t^{\lfloor\frac At\rfloor}\right)^{\varphi(t)}\right)^{\lfloor\frac Bt\rfloor\lfloor\frac Ct\rfloor}\\ \end{aligned} \]

這玩意有點複雜，把它拆成兩個部分：

\[\begin{aligned} & f_0(A,B,C)=\prod_{t=1}\left(\left(\lfloor\frac At\rfloor!\right)^{\lfloor\frac Bt\rfloor\lfloor\frac Ct\rfloor}\right)^{\varphi(t)}\\ & f_1(A,B,C)=\prod_{t=1}\left(t^{\varphi(t)}\right)^{\lfloor\frac At\rfloor\lfloor\frac Bt\rfloor\lfloor\frac Ct\rfloor} \end{aligned} \]

\(f_1\) 可以整除分塊了，把 \(f_2\) 先放著，看 \(\mathcal G\)：

\[\begin{aligned} \mathcal G(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}\prod_{i=1}^A\prod_{j=1}^B\prod_{k=1}^Cd^{[\gcd(i,j)=d]\gcd(d,k)}\\ =&\prod_{d=1}d^{\sum_{i=1}^{A/d}\sum_{j=1}^{B/d}[\gcd(i,j)=1]\sum_{k=1}^C\gcd(d,k)}\\ =&\prod_{d=1}d^{\sum_{e=1}\lfloor{\frac A{de}\rfloor}\lfloor\frac B{de}\rfloor\mu(e)\cdot\sum_{k=1}^C\gcd(d,k)}\\ =&\prod_{t=1}\left(\prod_{d\mid t}d^{\mu(\frac td)\cdot\sum_{k=1}^C\gcd(d,k)}\right)^{\lfloor\frac At\rfloor\lfloor\frac Bt\rfloor} \end{aligned} \]

指數上的東西：

\[\begin{aligned} \sum_{k=1}^C\gcd(d,k)=&\sum_{f\mid d}\sum_{k=1}^{\lfloor\frac Cf\rfloor}[\gcd(\frac df,k)=1]f\\ =&\sum_{f\mid d}\sum_{g\mid\frac df}\lfloor\frac C{fg}\rfloor\mu(g)f\\ =&\sum_{u\mid d}\varphi(u)\lfloor\frac Cu\rfloor \end{aligned} \]

然後帶回去：

\[\begin{aligned} \mathcal G(A,B,C)=&\prod_{t=1}\left(\prod_{d\mid t}d^{\mu(\frac td)\sum_{u\mid d}\varphi(u)\lfloor\frac Cu\rfloor}\right)^{\lfloor\frac At\rfloor\lfloor\frac Bt\rfloor}\\ \end{aligned} \]

這玩意不怎麼好預處理。考慮把底數的 \(d\) 拆成 \(u\) 和 \(\frac du\)：

\[\begin{aligned} & g_0(A,B,C)=\prod_{t=1}\left(\prod_{d\mid t}\prod_{u\mid d}\left(\frac du\right)^{\mu(\frac td)\varphi(u)\lfloor\frac Cu\rfloor}\right)^{\lfloor\frac At\rfloor\lfloor\frac Bt\rfloor}\\ & g_1(A,B,C)=\prod_{t=1}\left(\prod_{d\mid t}\prod_{u\mid d}u^{\mu(\frac td)\varphi(u)\lfloor\frac Cu\rfloor}\right)^{\lfloor\frac At\rfloor\lfloor\frac Bt\rfloor} \end{aligned} \]

分別推一下：

\[\begin{aligned} g_0(A,B,C)=&\prod_{t=1}\left(\prod_{d\mid t}\prod_{u\mid d}\left(\frac du\right)^{\mu(\frac td)\varphi(u)\lfloor\frac Cu\rfloor}\right)^{\lfloor\frac At\rfloor\lfloor\frac Bt\rfloor}\\ =&\prod_{t=1}\left(\prod_{u\mid t}\prod_{u\mid d\mid t}\left(\frac du\right)^{\mu(\frac td)\varphi(u)\lfloor\frac Cu\rfloor}\right)^{\lfloor\frac At\rfloor\lfloor\frac Bt\rfloor}\\ =&\prod_{t=1}\left(\prod_{u\mid t}\left(\prod_{d\mid \frac tu}d^{\mu(\frac td)}\right)^{\varphi(u)\lfloor\frac Cu\rfloor}\right)^{\lfloor\frac At\rfloor\lfloor\frac Bt\rfloor}\\ =&\prod_{u=1}\left(\prod_{t=1}\left(\prod_{d\mid t}d^{\mu(\frac td)}\right)^{\lfloor\frac A{tu}\rfloor\lfloor\frac B{tu}\rfloor}\right)^{\varphi(u)\lfloor\frac Cu\rfloor} \end{aligned} \]

括號裡的東西預處理過了，兩次整除分塊可以搞定。

\[\begin{aligned} g_1(A,B,C)=&\prod_{t=1}\left(\prod_{d\mid t}\prod_{u\mid d}u^{\mu(\frac td)\varphi(u)\lfloor\frac Cu\rfloor}\right)^{\lfloor\frac At\rfloor\lfloor\frac Bt\rfloor}\\ =&\prod_{u=1}\prod_{t=1}\left(\prod_{d\mid t}u^{\mu(\frac td)\varphi(u)\lfloor\frac Cu\rfloor}\right)^{\lfloor\frac A{tu}\rfloor\lfloor\frac B{tu}\rfloor}\\ =&\prod_{u=1}u^{\sum_{t=1}\sum_{d\mid t}\mu(\frac td)\varphi(u)\lfloor\frac Cu\rfloor\lfloor\frac A{tu}\rfloor\lfloor\frac B{tu}\rfloor}\\ =&\prod_{u=1}u^{\varphi(u)\lfloor\frac Cu\rfloor\sum_{t=1}\lfloor\frac A{tu}\rfloor\lfloor\frac B{tu}\rfloor\sum_{d\mid t}\mu(\frac td)}\\ =&\prod_{u=1}u^{\varphi(u)\lfloor\frac Cu\rfloor\sum_{t=1}\lfloor\frac A{tu}\rfloor\lfloor\frac B{tu}\rfloor[t=1]}\\ =&\prod_{u=1}u^{\varphi(u)\lfloor\frac Au\rfloor\lfloor\frac Bu\rfloor\lfloor\frac Cu\rfloor}\\ \end{aligned} \]

注意到，\(g_1\) 和之前的 \(f_1\) 其實是一樣的，而且它們的值不受 \(A,B,C\) 順序的影響，因此可以把它們約掉。

---

卡常小技巧（不過應該眾所周知）：把預處理的東西的逆元同時預處理出來可以大幅減小常數

code：

#include <cstdio>

const int N = 1e5,maxn = N + 5;

int T,p,phi,A,B,C,ta,tb,tc;
int pri[maxn],mu[maxn],Inv[maxn],f[maxn],Phi[maxn];
int fac[maxn],fd[maxn],facInv[maxn],fdInv[maxn];
int Fac[maxn],Fd[maxn],FdInv[maxn],ps[maxn];
bool isp[maxn];

inline int read() {
#define gc c = getchar()
    int d = 0,f = 0,gc;
    for(;c < 48 || c > 57;gc) f |= (c == '-');
    for(;c > 47 && c < 58;gc) d = (d << 1) + (d << 3) + (c ^ 48);
#undef gc
    return f ? -d : d;
}

inline int Mul(int a,int b,int mod) { return 1LL * a * b % mod; }
inline int Add(int a,int b,int mod) { a += b; return a > mod ? a - mod : a; }
inline int max(int a,int b) { return a > b ? a : b; }
inline int min(int a,int b) { return a < b ? a : b; }
inline int min(int a,int b,int c) { return min(min(a,b),c); }

inline int fpow(int a,int b,int mod) {
    int res = 1;
    for(;b;a = Mul(a,a,mod),b >>= 1) if(b & 1) res = Mul(res,a,mod);
    return res;
}

inline int GetInv(int a) { return a <= N ? Inv[a] : fpow(a,p - 2,p); }

namespace Type0 {
    inline int F(int A,int B,int C) { return fpow(fac[A],Mul(B,C,phi),p); }

    inline int G(int A,int B,int C) {
        int res,ans = 1;
        for(int r,l = 1;l <= A && l <= B;l = r + 1) {
            ta = A / l,tb = B / l; r = min(A / ta,B / tb);
            res = Mul(fd[r],fdInv[l - 1],p);
            ans = Mul(ans,fpow(res,Mul(ta,tb,phi),p),p);
        }
        return fpow(ans,C,p);
    }

    inline int solve() {
        int res = Mul(F(A,B,C),F(B,A,C),p);
        int inv = Mul(G(A,B,C),G(A,C,B),p);
        return Mul(res,GetInv(inv),p);
    }
}

namespace Type1 {
    inline int S(int n) { return 1LL * n * (n + 1) / 2 % phi; }

    inline int F(int A,int B,int C) { return fpow(Fac[A],Mul(S(B),S(C),phi),p); }

    inline int G(int A,int B,int C) {
        int res,ans = 1;
        for(int r,l = 1;l <= A && l <= B;l = r + 1) {
            ta = A / l,tb = B / l;
            r = min(A / ta,B / tb);
            res = Mul(Fd[r],FdInv[l - 1],p);
            ans = Mul(ans,fpow(res,Mul(S(ta),S(tb),phi),p),p);
        }
        return fpow(ans,S(C),p);
    }

    inline int solve() {
        int res = Mul(F(A,B,C),F(B,A,C),p);
        int inv = Mul(G(A,B,C),G(A,C,B),p);
        return Mul(res,GetInv(inv),p);
    }
}

namespace Type2 {
    inline int F(int A,int B,int C) { 
        int exp,res,ans = 1;
        for(int r,l = 1;l <= min(A,B,C);l = r + 1) {
            ta = A / l,tb = B / l,tc = C / l;
            r = min(A / ta,B / tb,C / tc);
            exp = Add(ps[r],phi - ps[l - 1],phi);
            res = fpow(fac[ta],Mul(tb,tc,phi),p);
            ans = Mul(ans,fpow(res,exp,p),p);
        }
        return ans;
    }

    inline int g0(int A,int B) {
        int res,ans = 1;
        for(int r,l = 1;l <= A && l <= B;l = r + 1) {
            ta = A / l,tb = B / l;
            r = min(A / ta,B / tb);
            res = Mul(fd[r],fdInv[l - 1],p);
            ans = Mul(ans,fpow(res,Mul(ta,tb,phi),p),p);
        }
        return ans;
    }

    inline int G(int A,int B,int C) {
        int res,ans = 1;
        for(int r,l = 1;l <= min(A,B,C);l = r + 1) {
            ta = A / l,tb = B / l,tc = C / l;
            r = min(A / ta,B / tb,C / tc);
            res = fpow(g0(ta,tb),tc,p);
            ans = Mul(ans,fpow(res,Add(ps[r],phi - ps[l - 1],phi),p),p);
        }
        return ans;
    }

    inline int solve() {
        int res = Mul(F(A,B,C),F(B,A,C),p);
        int inv = Mul(G(A,B,C),G(A,C,B),p);
        return Mul(res,GetInv(inv),p);
    }
}

inline void Init() {
    mu[1] = Phi[1] = f[1] = Inv[0] = Inv[1] = 1;
    fac[0] = fac[1] = Fac[0] = facInv[0] = 1;
    fd[0] = Fd[0] = fdInv[0] = FdInv[0] = 1;
    for(int i = 2;i <= N;i ++) {
        f[i] = 1,fac[i] = Mul(i,fac[i - 1],p);
        if(!isp[i]) pri[++ pri[0]] = i,mu[i] = -1,Phi[i] = i - 1;
        for(int j = 1;j <= pri[0] && i * pri[j] <= N;j ++) {
            isp[i * pri[j]] = true;
            if(!(i % pri[j])) {
                mu[i * pri[j]] = 0;
                Phi[i * pri[j]] = Phi[i] * pri[j];
                break;
            }
            mu[i * pri[j]] = -mu[i];
            Phi[i * pri[j]] = Phi[i] * Phi[pri[j]];
        }
    }
    facInv[N] = fpow(fac[N],p - 2,p);
    for(int i = N - 1;i;i --) {
        facInv[i] = Mul(facInv[i + 1],i + 1,p);
        Inv[i + 1] = Mul(facInv[i + 1],fac[i],p);
    }
    for(int i = 1;i <= N;i ++) {
        ps[i] = Add(Phi[i],ps[i - 1],phi);
        Fac[i] = Mul(Fac[i - 1],fpow(i,i,p),p);
        fd[i] = Mul(f[i],fd[i - 1],p);
        Fd[i] = Mul(fpow(f[i],Mul(i,i,phi),p),Fd[i - 1],p);
        fdInv[i] = GetInv(fd[i]),FdInv[i] = GetInv(Fd[i]);
        if(!mu[i]) continue;
        for(int j = 1;i * j <= N;j ++)
           f[i * j] = Mul(f[i * j],mu[i] == 1 ? j : Inv[j],p);
    }
}

int main() {
    T = read(),p = read();
    phi = p - 1; Init();
    for(;T;T --) {
        A = read(),B = read(),C = read();
        printf("%d %d %d\n",Type0::solve(),Type1::solve(),Type2::solve());
    }
    return 0;
}