The 2021 CCPC Guangzhou Onsite - K. Magus Night 题解

题意

Marisa Kirisame 有 $n$ 个魔法珠，他可以进行一次充能，使得每个珠子 $s_i$ 随机充能到 $[1,m]$ 的一个值，如果这些珠子的魔法值的 $gcd\le p$ 且 $lcm\ge q$ 那么这算一次成功的充能，带来 $\prod _{i=0}^n{s}_i$ 的贡献，否则贡献为 $0$ ，记贡献期望值为 $E$ ，要求输出 $E\cdot m^{n}mod998244353$ 的结果

思路

容斥套容斥套容斥，三层容斥，赛上没想到这么多层，也没想到gcd限制+lcm限制下是可以在2s内跑出来的，麻了

---

划分

先来个容斥：由于 $lcm\ge p$ 的条件不好求贡献，考虑转换为 $lcm<p$ 后枚举 $lcm$（因为枚举 $lcm\ge p$ 的话，要注意到 $lcm$ 会超过 $m$ ）

记（写法可能不标准）

对数组 $S=\{s_1,s_2,s_3,\dots ,s_n\}$

$$\begin{array}{rl}Z& =\{S|1\le s_i\le m\}\\ A& =\{Z|gcd(S)\le q\wedge lcm(S)\ge p\}\\ B& =\{Z|gcd(S)>q\}\\ C& =\{Z|gcd(S)\le q\wedge lcm(S)<p\}\end{array}$$

故

$$A=Z-B-C$$

Z

$$\begin{array}{rl}& (1+2+3+\cdots +m{)}^n\\ =& (\frac{1}{2}\cdot m\cdot (m+1){)}^n\end{array}$$

复杂度 $O(logn)$

B

如果你不会求这个……那么你应该学习 luogu P4450 后再来看这篇题解

总体来说，就是枚举 $s_i$ 为 $1,2,3,\dots ,m$ 的倍数后乘上一个对应的莫比乌斯函数

贡献即为（ $k$ 表示 $gcd=k$ 时的贡献）

$$\begin{array}{rl}& \sum _{k=q+1}^m\sum _{i=1}^{\lfloor \frac{m}{k}\rfloor }\mu (i)\cdot (i+2i+3i+\cdots +\lfloor \frac{m}{i\cdot k}\rfloor \cdot i{)}^n\\ =& \sum _{k=q+1}^m\sum _{i=1}^{\lfloor \frac{m}{k}\rfloor }\mu (i)\cdot (\frac{1}{2}\cdot ik\cdot \lfloor \frac{m}{i\cdot k}\rfloor (\lfloor \frac{m}{i\cdot k}\rfloor +1){)}^n\end{array}$$

复杂度$O(m\cdot logm\cdot logn)$

C

前方容斥套容斥警告

记 $f(j)$ 表示 $lcm(S)=j$ 时的贡献

$$\begin{array}{rl}f(j)& =\sum _{d|j}\mu (\frac{j}{d})\cdot (\sum d{)}^n\end{array}$$

后面那个 $\sum d$ 可以预处理一下

所以 $C$ 的贡献是

$$\begin{array}{rl}& \sum _{k=1}^m\sum _{i=1}^{\lfloor \frac{m}{k}\rfloor }\mu (i)\cdot \sum _{j=1}^{\lfloor \frac{p-1}{i\cdot k}\rfloor }f(j)\cdot (i\cdot k{)}^n\\ =& \sum _{k=1}^m\sum _{i=1}^{\lfloor \frac{m}{k}\rfloor }\mu (i)\cdot (i\cdot k{)}^n\cdot \sum _{j=1}^{\lfloor \frac{p-1}{i\cdot k}\rfloor }\sum _{d|j}\mu (\frac{j}{d})\cdot (\sum d{)}^n\end{array}$$

复杂度O(能过)

复杂度我也分析不清楚，但是不会超过$O(m\cdot logm\cdot lo{g}^{2}p)$ ，预处理 $\sum d$ 下

代码

#include<bits/stdc++.h>
using namespace std;
using i64 = long long;

const int N = 2e5+6;
vector<int> divs[N];
i64 dsum[N] = {0};
int prms[N],cnt=0,mob[N];
bool notP[N] = {false};
const int MOD = 998244353;

i64 qpow(i64 b,int p) {
    i64 ret=1;
    for(;p;p>>=1,b=b*b%MOD) if(p&1) ret = ret*b%MOD;
    return ret;
}

void get_mob() {
    mob[1] = 1;
    for(int i=2;i<N;++i) {
        if(!notP[i]) {
            prms[cnt] = i;
            ++cnt;
            mob[i] = -1;
        }

        for(int j=0;j<cnt;++j) {
            if(i*prms[j]>=N) break;

            mob[i*prms[j]] = -mob[i];
            notP[i*prms[j]] = true;

            if(i%prms[j]==0) {
                mob[i*prms[j]] = 0;
                break;
            }
        }
    }
}

void norm(i64 &x) {
    while(x >= MOD) x -= MOD;
    while(x < 0) x += MOD;
}

void get_d() {
    for(int i=1;i<N;++i) {
        for(int j=i;j<N;j+=i) {
            divs[j].push_back(i);
            dsum[j] += i;
            norm(dsum[j]);
        }
    }
}

int main(int argc, char const *argv[])
{
    ios_base::sync_with_stdio(false);
    cin.tie(nullptr); cout.tie(nullptr);

    get_mob();
    get_d();

    i64 n;
    int m,p,q;
    cin >> n >> m >> p >> q;

    for(int i=0;i<N;++i) {
        dsum[i] = qpow(dsum[i],n);
    }

    i64 ans2=0,ans3=0,ans4=0;

    // 2
    for(int k=q+1;k<=m;++k) {
        int M = m/k;
        for(int i=1;i<=M;++i) {
            ans2 += mob[i] * qpow(1ll*(M/i)*(M/i+1)/2%MOD*i%MOD*k%MOD, n);
            norm(ans2);
        }
    }

    // cout << ans2 << endl;

    // 3
    for(int k=1;k<=q;++k) {
        int M = m/k;
        for(int i=1;i<=M;++i) {
            i64 tmp = 0;
            int P = (p-1)/i/k;
            for(int j=1;j<=P;++j) {
                for(int d : divs[j]) {
                    tmp += mob[j/d] * dsum[d] % MOD;
                    norm(tmp);
                }
            }
            ans3 += mob[i] * tmp * qpow(i*k%MOD,n) % MOD;
            norm(ans3);
        }
    }

    // cout << ans3 << endl;

    // 4
    ans4 = qpow(1ll*m*(m+1)/2%MOD,n);

    i64 ans = ans4 - ans2 - ans3;
    norm(ans);
    cout << ans;

    return 0;
}