The 2021 CCPC Guangzhou Onsite - K. Magus Night 题解
title: The 2021 CCPC Guangzhou Onsite - K. Magus Night 题解
date: 2021-11-18 18:19:00
tags: [inclusion-exclusion, number theory, math, team training]
mathjax: true
题意
Marisa Kirisame 有 \(n\) 个魔法珠，他可以进行一次充能，使得每个珠子 \(s_i\) 随机充能到 \([1,m]\) 的一个值，如果这些珠子的魔法值的 \(gcd\leq p\) 且 \(lcm\geq q\) 那么这算一次成功的充能，带来 \(\prod_{i=0}^n s_i\) 的贡献，否则贡献为 \(0\) ，记贡献期望值为 \(E\) ，要求输出 \(E \cdot m^n \ mod\ \ 998\ 244\ 353\) 的结果
思路
容斥套容斥套容斥，三层容斥，赛上没想到这么多层，也没想到gcd限制+lcm限制下是可以在2s内跑出来的，麻了
划分
先来个容斥：由于 \(lcm\geq p\) 的条件不好求贡献，考虑转换为 \(lcm
记（写法可能不标准）
对数组 \(S = \{s_1,s_2,s_3, \dots ,s_n\}\)
\[\begin{aligned} Z &= \{S\ |\ 1\leq s_i\leq m\} \\ A &= \{Z\ |\ gcd(S) \leq q\ \land lcm(S)\geq p\} \\ B &= \{Z\ |\ gcd(S) > q \} \\ C &= \{Z\ |\ gcd(S) \leq q\ \land lcm(S)Z
\[\begin{aligned} &(1+2+3+\dots +m)^n \\ =&\big(\frac{1}{2}\cdot m \cdot (m+1)\big)^n \end{aligned} \]复杂度 \(O(log\ n)\)
B
如果你不会求这个……那么你应该学习 luogu P4450 后再来看这篇题解
总体来说，就是枚举 \(s_i\) 为 \(1,2,3,\dots ,m\) 的倍数后乘上一个对应的莫比乌斯函数
贡献即为（ \(k\) 表示 \(gcd=k\) 时的贡献）
\[\begin{aligned} & \sum_{k=1}^m \sum _{i=1} ^{\lfloor \frac{m}{k} \rfloor} \mu(i) \cdot (i+2i+3i+ \dots + {\lfloor \frac{m}{i\cdot k} \rfloor}\cdot i)^n \\ =& \sum_{k=1}^m \sum _{i=1} ^{\lfloor \frac{m}{k} \rfloor} \mu(i) \cdot \big(\frac{1}{2}\cdot i \cdot {\lfloor \frac{m}{i\cdot k} \rfloor} ({\lfloor \frac{m}{i\cdot k} \rfloor}+1)\big) ^ n \end{aligned} \]复杂度\(O(m\cdot log\ m \cdot log\ n)\)
C
前方容斥套容斥警告
记 \(f(j)\) 表示 \(lcm(S)=j\) 时的贡献
\[\begin{aligned} f(j) &= \sum_{d|j}\mu(\frac{j}{d}) \cdot \Big(\sum d\Big)^n \\ \end{aligned} \]后面那个 \(\sum d\) 可以预处理一下
所以 \(C\) 的贡献是
\[\begin{aligned} & \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 \Big(\sum d\Big)^n \end{aligned} \]复杂度O(能过)
复杂度我也分析不清楚，但是不会超过\(O(m\cdot log\ m \cdot log^2\ p)\) ，预处理 \(\sum d\) 下
代码
#include
using namespace std;
using i64 = long long;
const int N = 2e5+6;
vector 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) 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 >> m >> p >> q;
for(int i=0;i