神秘题

· · 题解

神力 / Drunkards

给定一个长度为 n 的序列 a，其中 a_i \in \{-1, 1\}。初始位置为 0。

在 i \in [1, n] 的每个时刻，原计划执行的移动操作是使得当前位置增加 a_i。但是，有 p 的概率在第 i 个时刻不执行移动操作，即位置不变；有 1 - p 的概率执行移动操作。各个时刻的操作是否执行相互独立。

对于所有可能的终止位置 i \in [-n, n]，求出最终经过过该位置的概率，结果对 10^9+7 取模。

尝试 dp？设 f(i,j) 表示时刻 i 到达 j 的概率，你发现求出这个东西后似乎也没有什么用（一个位置可能会被重复经过），如果尝试容斥，是可行的，但时间复杂度 O(n^3)\sim O(n^5) 爆炸。

既然容斥不行，考虑钦定计入答案规则来去重，假如我们只让第一次到达 j 计入答案，那么设 g(i,j) 表示时刻 i 第一次到达 j 的概率，转移形如：

g(i,j)=f(i,j)-\sum_{k<i} f(k,j) P(k+1,i,0)

其中 P(x,y,z) 表示时刻 x 在任意位置 a，在时刻 y 到达 a+z 的概率。那么这个可以 O(n^3) dp 得到，时间复杂度爆炸。

尝试换一个计入答案规则，让最后一次到达 j 计入贡献，那么 j 的答案就是：

\sum_{i} f(i,j) P(i + 1)

其中 P(i) 表示从时刻 i+1 在任意位置 a，且 i+1 选择了移动操作，到时刻 n，(i+1,n] 时刻区间内没有到达 a 的概率。

这个怎么算呢？如果 a_{i+1}=1，那么假如我们称在时刻 k 停下为 a_k=0 的话，串 [i+1,n] 的前缀和都应该 >0（a_{i+1}=-1 类似，故不赘述）。

一个直观的想法是设 h(i,j) 表示以 i 开头的后缀，最小前缀和为 j 的概率，那么实际上已经可以转移了：

\begin{aligned} h(i,\min(0, j))&\gets p h(i+1, j)\\ h(i, a_i+\min(0, j))&\gets (1-p) h(i+1, j) \end{aligned}

最后可以直接反推出答案：

P(i)=(1-p)\sum_{k+1>0} h(i+1,k)

特别地，P(n+1)=1。时间复杂度 O(n^2)，可能小小卡常。

namespace solution {
    const int N = 5005, mod = 1e9 + 7;
    using mint = GF<int, mod>;
    int n, p_, a[N];
    template<class T, int N> struct neg_arr {
        T a[N << 1];
        T& operator[](int i){ return a[i + N]; }
    };
    neg_arr<mint, N> f[N], h_min[N], h_max[N];
    mint p, q, P[N];
    void solve(){
        cin >> n >> p_, p = mint(p_) / 100, q = mint(1) - p;
        for(int i = 1; i <= n; i++) cin >> a[i];
        f[0][0] = 1;
        for(int i = 0; i < n; i++){
            for(int j = -n; j <= n; j++){
                if(!f[i][j].v) continue;
                f[i + 1][j + a[i + 1]] += q * f[i][j], f[i + 1][j] += p * f[i][j];
            }
        }
        h_min[n + 1][0] = h_max[n + 1][0] = 1;
        for(int i = n; i; i--){
            for(int j = -n; j <= n; j++){
                if(h_min[i + 1][j].v){
                    h_min[i][std::min(0, j)] += p * h_min[i + 1][j];
                    h_min[i][a[i] + std::min(0, j)] += q * h_min[i + 1][j];
                }
                if(h_max[i + 1][j].v){
                    h_max[i][std::max(0, j)] += p * h_max[i + 1][j];
                    h_max[i][a[i] + std::max(0, j)] += q * h_max[i + 1][j];
                }
            }
        }
        P[n + 1] = 1;
        for(int i = 1; i <= n; i++){
            if(a[i] == 1) for(int k = 0; k <= n; k++) P[i] += h_min[i + 1][k];
            else for(int k = -n; k <= 0; k++) P[i] += h_max[i + 1][k];
            P[i] *= q;
        }
        for(int i = -n; i <= n; i++){
            mint ans = 0;
            for(int j = 0; j <= n; j++) ans += f[j][i] * P[j + 1];
            std::cout << ans << ' ';
        }
        std::cout << '\n';
    }
}

开关灯 / Flipping Game

有 n 盏灯，初始时灯 i 是否亮着记为 a_i，目标状态灯 i 是否亮着记作 b_i。

选择其中的三盏灯构成无序不重三元组 (x,y,z)，有 N=\binom{n}{3} 种选法，你需要从这 N 个三元组种不重复地选出若干个，然后对于每一个选出的三元组 (x,y,z)，按下灯 x,y,z 的开关。若可以将初始状态转换为目标状态，则称这个三元组集合是好的，你需要计算大小为 m 的好的三元组集合的数量。答案对 10007 取模。

先不考虑不重复选三元组的限制，一个最朴素的做法，设 f(i,S) 表示现在进行了 i 次操作，集合 S 内的灯被按了奇数次开关的方案数，状态数直接到了 O(m2^n) 不可承受。不过你发现大小相同的两个集合 S,T，f(i,S)=f(i,T)，因为它们完全对称。这样就可以将状态削减到 f(i,j)，其中 j=|S|。

转移的话可以枚举这一步的三元组怎么选：

时间复杂度 O(nm)，但是直接求出这个 f 显然没有用处，因为选出的三元组可能重复，尝试整体容斥又非常困难，尝试单步容斥，即每次转移的时候都尝试减去不合法的情况。

对于 f(i+1,j)，如果有不合法的，只能是 i+1 和 [1,i] 中的某一次重复了，假设是 x。等价于 i-1 次就翻转完了 j 盏灯，接下来做了两次相同的操作，方案数就是 f(i-1,j)i(\binom{n}{3}-(i-1))，还挺直观的。最后答案就是：

\frac{f(m,\text{difference})}{m! \binom{n}{\text{difference}}}

意义就是除掉顺序 m!，以及钦定选择翻转的就是哪些不同的位置 \binom{n}{\text{difference}}（因为你计算了所有同构的贡献，但是你只需要计算等价类中的一个元素）。于是总时间复杂度也可以保持 O(nm)。

namespace solution {
    const int N = 1005, mod = 10007;
    using mint = GF<int, mod>;
    int n, m;
    mint fact[N], ifact[N], f[N][N];
    std::string from, to;
    mint binom(int n, int m){ return (n < m || m < 0) ? 0 : fact[n] * ifact[m] * ifact[n - m]; }
    bool solve(){
        cin >> n >> m >> from >> to; int _ = std::max(n, m);
        if(!_) return false;
        fact[0] = 1; for(int i = 1; i <= _; i++) fact[i] = fact[i - 1] * i;
        ifact[_] = ~fact[_]; for(int i = _; i; i--) ifact[i - 1] = ifact[i] * i;
        f[0][0] = 1;
        for(int i = 0; i < m; i++){
            for(int j = 0; j <= n; j++) f[i + 1][j] = 0;
            for(int j = 0; j <= n; j++){
                if(j >= 3) f[i + 1][j - 3] += f[i][j] * binom(j, 3);
                if(j >= 1) f[i + 1][j - 1] += f[i][j] * binom(j, 2) * (n - j);
                f[i + 1][j + 1] += f[i][j] * j * binom(n - j, 2);
                f[i + 1][j + 3] += f[i][j] * binom(n - j, 3);
            }
            if(!i) continue;
            for(int j = 0; j <= n; j++) f[i + 1][j] -= f[i - 1][j] * i * (binom(n, 3) - (i - 1));
        }
        int cnt = 0; for(int i = 0; i < n; i++) cnt += from[i] != to[i];
        std::cout << f[m][cnt] * (~binom(n, cnt)) * ifact[m] << '\n';
        return true;
    }
}

排列问题 / Tiket

给定三个长度为 n 的排列 a,b,c，计算：

\sum_{1\leq i,j\leq n}[a_i<a_j][b_i<b_j][c_i<c_j]

看起来是一个经典的三维偏序问题，直接用经典的树套树或者 cdq 分治做，时间复杂度 O(n\log^2 n) 难以接受。

我们记 S(a,b,\cdots,\overline{A},\overline{B},\cdots) 表示在 a,b,\cdots 这些排列上满足 \square_i<\square_j，在 A,B 则不满足的二元组 (i,j) 集合。那么我们要求的就是 |S(a,b,c)|。

我们虽然暂时不会 O(n\log^2 n) 求这个东西，但是通过排序 + 树状数组，很容易求出 |S(a,b)|,|S(a,c)|,|S(b,c)|（这一步我相信你肯定会）。

考虑 |S(a,b)|+|S(a,c)|+|S(b,c)| 这个 naive 的想法，它算的是什么？对于 S(a,b,c) 中的元素，它算了三次，对于 S(a,b,\overline{c}) 之类的元素，它算了一次，其余的没有计算，那么我们最好是减去 \Delta=|S(a,b,c)|+|S(a,b,\overline{c})|+|S(a,c,\overline{b})|,|S(b,c,\overline{a})|，最后除以 2 就是答案。

惊人的事实： \Delta=\frac{1}{2}n(n-1)，为什么？发现若 (i,j)\in S(a,b,c)，那么 (j,i)\in S(\overline{a},\overline{b},\overline{c})，所以 |S(a,b,c)|=|S(\overline{a},\overline{b},\overline{c})|，同理 |S(a,b,\overline{c})|=|S(c,\overline{a},\overline{b})|。而这些 S 的情形正好覆盖了所有 i\neq j 的 (i,j) 对，总方案数为 n(n-1)，因此 \Delta 就是总方案数除以 2。证毕。

时间复杂度 O(n\log n) 常数较大但足以通过。

namespace solution {
    using i64 = int64_t;
    const int N = 2e6 + 5;
    int n, a[N], b[N], c[N], tmp[N];
    i64 seedA, seedB, seedC;
    i64 rand(i64& seed){ return seed = ((seed * 19260817) ^ 233333) & ((1 << 24) - 1); }
    struct BIT {
        int t[N];
        void clear(){ std::fill(t + 1, t + n + 1, 0); }
        void update(int p, int v){ for(; p <= n; p += p & -p) t[p] += v; }
        int query(int p){ int res = 0; for(; p; p -= p & -p) res += t[p]; return res; }
    } bit;
    i64 paritial2d(int* a, int* b){
        for(int i = 1; i <= n; i++) tmp[i] = i;
        std::sort(tmp + 1, tmp + n + 1, [=](int x, int y){ return a[x] < a[y]; });
        bit.clear(); i64 res = 0;
        for(int i = 1; i <= n; i++) res += bit.query(b[tmp[i]]), bit.update(b[tmp[i]], 1);
        return res;
    }
    void solve(){
        cin >> n >> seedA >> seedB >> seedC;
        for(int i = 1; i <= n; i++) a[i] = b[i] = c[i] = i;
        for(int i = 1; i <= n; i++) std::swap(a[i], a[rand(seedA) % i + 1]);
        for(int i = 1; i <= n; i++) std::swap(b[i], b[rand(seedB) % i + 1]);
        for(int i = 1; i <= n; i++) std::swap(c[i], c[rand(seedC) % i + 1]);
        i64 ans = paritial2d(a, b) + paritial2d(b, c) + paritial2d(a, c);
        ans -= 1ll * n * (n - 1) / 2;
        std::cout << (ans / 2) << '\n';
    }
}

数数 / Fancy Arrays

给定 m 和 q，q 次询问有多少个长度为 n 的数列 A 满足所有元素都是 m 的因子且任意相邻两数不互质，答案对 998244353 取模。

一个朴素的 dp 是记 f(i,j) 表示考虑前缀 [1,i]，且 A_i=j 的方案数，那么状态数达到了 O(nd(m))，总转移时间复杂度 O(nd^2(m)) 难以接受。

由于发现每一次转移都形如 f(i+1,j)\gets f(i,k)，其中 j\mid m, k\mid m,\gcd(j,k)\neq 1。那么这个转移可以轻易地使用矩阵乘法写出来，然后就可以用矩阵快速幂优化到 O(d^3(m)\log n) 同样难以接受。因为 d(m)\sim \exp\left((\ln 2+o(1)) \frac{\ln m}{\ln \ln m}\right) 很大（S. Wigert, 1907 / S. Ramanujan, 1915。在 m\leq 10^{18} 的规模可以达到 103680）。尝试削减状态数。

首先可以观察到我们不需要记录具体的因数，事实上只需要记录 A_i 有哪些质因子即可，因此我们设 f(i,S) 表示 [1,i] 中 A_i 的质因子集合为 S 的方案数，这样复杂度就降到了 O((2^{\omega(m)})^3 \log n)，由于 \omega(m)\sim (1+o(1)) \frac{\ln m}{\ln\ln m}（Hardy–Ramanujan Theorem, 1917。在 m\leq 10^{18} 的规模可以达到 15），这个做法还是不够好。

我们考察 m=30=2\times 3\times 5，那么此时 j=2\times 3 和 j'=3\times 5，f(i,j)=f(i,j') 是成立的，因为此时 2,3,5 是完全对称的。而对于 m=60=2^2\times 3\times 5，此时只有 3,5 等价，2 和它们不等价（指数不同），因此我们认为两个质因子 p_i,p_j 等价，当且仅当它们在 m 的指数 e_i,e_j 相等，那么我们可以定义根据此法对 \mathcal{P}[\{p_1,\cdots,p_{\omega(m)}\}] 划分为若干个等价类，等价类的规模 G(m) 如此定义：

一个较为简洁的上界估计是 2^{O(\sqrt{\ln m})}（证明需要一定的 analytic number theory 知识，从略），10^{18} 内 \leq 150 是可以预期的。

以上都是做法的大致轮廓，我们要将大致想法变成具体思路。首先对 m 质因子分解，由于 m 很大，O(\sqrt{m}) 的朴素算法不可行，可以使用 O(m^{1/4}) 的 Pollard-Rho + Miller-Rabin 算法分解质因数。然后统计得到 k 个等价类，等价类 i 指数为 e_i，质数数量 c_i。然后搜索得到 G(m) 个等价类。

尝试构建转移矩阵 F(S,T)，记 S_i,T_i 表示等价类 i 选择多少个质数。那么我们要求要选一个满足 T 的不与 S 互质，直接计算比较困难，先计算满足 T 的所有方案：

f_T = \prod_i \binom{c_i}{T_i} e_i^{T_i}

然后尝试减去与 S 互质的方案数，一个比较好的想法是用概率方法，乘上 1-\text{Pr}[\gcd(A_i,A_{i+1})=1] 这种，计算互质的概率也很简单，钦定 S 选择 T 剩余的即可：

\text{Pr}[\gcd(A_i,A_{i+1})=1]=\prod_i \dfrac{\binom{c_i-T_i}{S_i}}{\binom{c_i}{S_i}}

然后我们要统计答案，只需要计算 f' = f \cdot F^{n-1} 的各分量和即可。这样时间复杂度是 O(m^{1/4}+qG^3(m)\log n) 可能无法通过。

注意到我们要计算的是向量 * 矩阵的若干次幂的形式，并且询问较多（q\leq 200），因此可以预处理所有 F^{2^i}，然后询问只需要按照二进制位统计即可，这样时间复杂度就做到了 O(m^{1/4}+G^3(m)\log n+qG^2(m)\log n)，可以通过本题。

namespace solution {
    const int mod = 998244353;
    using mint = GF<int, mod>;
    using i64 = long long;
    using i128 = __int128;
    i64 n, m, q;
    mint fact[255], ifact[255];
    mint binom(int n, int m){ return (n < m || m < 0) ? 0 : fact[n] * ifact[m] * ifact[n - m]; }
    bool miller_rabin(i64 p){
        static constexpr int pris[] = {2, 3, 5, 7, 11, 13, 17, 19, 23, 29, 31, 37};
        static auto fastpow = [](i64 a, i64 b, i64 p){
            i64 res = 1; for(; b; b >>= 1, a = i128(a) * a % p) if(b & 1) res = i128(res) * a % p;
            return res;
        };
        i64 t, k; for(k = 0, t = p - 1; !(t & 1); k++, t >>= 1);
        for(int i : pris){
            if(p == i) return true;
            i64 a = fastpow(i, t, p), b = a;
            for(int j = 1; j <= k; j++){
                b = i128(a) * a % p;
                if(b == 1 && a != 1 && a != p - 1) return false;
                a = b;
            }
            if(a != 1) return false;
        }
        return true;
    }
    i64 pollard_rho(i64 n){
        if(!(n % 2)) return 2;
        if(!(n % 3)) return 3;
        static std::mt19937_64 rnd(std::random_device{}());
        std::uniform_int_distribution<i64> gen(1, n - 1);
        while(true){
            i64 c = gen(rnd), s = 0, t = 0, val = 1;
            auto f = [=](i64 x){ return (i128(x) * x + c) % n; };
            for(i64 i = 1; ; i <<= 1, s = t, val = 1){
                for (i64 j = 1; j <= i; j++){
                    t = f(t);
                    i64 diff = t > s ? t - s : s - t;
                    val = i128(val) * diff % n;
                    if (!(j & 127)) {
                        i64 d = std::__gcd(val, n);
                        if (d > 1 && d < n) return d;
                    }
                }
                i64 d = std::__gcd(val, n);
                if (d > 1 && d < n) return d;
                if(d == n) break;
            }
        }
    }
    std::map<i64, int> pri_fact;
    std::map<int, int> ed;
    std::basic_string<int> es;
    void dfs_pri_exp(i64 u){
        if(u == 1) return;
        if(pri_fact.count(u) || miller_rabin(u)) return pri_fact[u]++, void();
        i64 v = pollard_rho(u);
        dfs_pri_exp(v), dfs_pri_exp(u / v);
    }
    struct matrix {
        std::vector<std::vector<mint>> a;
        int n, m;
        matrix() = default;
        void resize(int _n, int _m){ n = _n, m = _m, a.clear(), a.resize(n, std::vector<mint>(m)); }
        matrix(int n, int m){ resize(n, m); }
        matrix(std::vector<std::vector<mint>> a) : a(a), n(a.size()), m(a[0].size()) {}
        std::vector<mint>& operator[](int i) { return a[i]; }
        const std::vector<mint>& operator[](int i) const { return a[i]; }
    };
    matrix operator*(const matrix& a, const matrix& b){
        matrix c(a.n, b.m);
        for(int i = 0; i < a.n; i++){
            for(int k = 0; k < a.m; k++){
                for(int j = 0; j < b.m; j++) c[i][j] += a[i][k] * b[k][j];
            }
        }
        return c;
    }
    matrix F, pw[65];
    std::basic_string<int> stat[255]; int k; mint f[255];
    void dfs_stat(size_t p, std::basic_string<int> s){
        if(p == es.size()) return stat[k++] = s, void();
        for(int i = 0; i <= ed[es[p]]; i++) dfs_stat(p + 1, s + i);
    }
    void solve(){
        cin >> m >> q, dfs_pri_exp(m); int allp = 0;
        for(auto [_, e] : pri_fact) ed[e]++, allp += e;
        for(auto [e, _] : ed) es += e;
        dfs_stat(0, {});
        fact[0] = 1; for(int i = 1; i <= allp; i++) fact[i] = fact[i - 1] * i;
        ifact[allp] = ~fact[allp]; for(int i = allp; i; i--) ifact[i - 1] = ifact[i] * i;
        F.resize(k, k);
        for(int S = 0; S < k; S++){
            f[S] = 1;
            for(int i = 0; i < es.size(); i++){
                int cnt = ed[es[i]], exp = es[i];
                f[S] *= binom(cnt, stat[S][i]) * (mint(exp) ^ stat[S][i]);
            }
        }
        for(int S = 1; S < k; S++){
            for(int T = 1; T < k; T++){
                mint coprime = 1;
                for(size_t i = 0; i < es.size(); i++){
                    int cnt = ed[es[i]];
                    coprime *= (binom(cnt - stat[T][i], stat[S][i]) / binom(cnt, stat[S][i]));
                }
                F[S][T] = f[T] * (mint(1) - coprime);
            }
        }
        pw[0] = F; for(int i = 1; i < 64; i++) pw[i] = pw[i - 1] * pw[i - 1];
        while(q--){
            cin >> n;
            matrix G(1, k);
            for(int i = 0; i < k; i++) G[0][i] = f[i];
            for(int i = 63; ~i; i--) if(((n - 1) >> i) & 1) G = G * pw[i];
            mint ans = 0; for(int i = 0; i < k; i++) ans += G[0][i];
            // There is a bug in the official data for this question, so we need to subtract 1 here. Please note that subtracting 1 is only to fit the result of the official data, not to subtract one to get the correct answer.
            if(n == 1) ans -= 1;
            std::cout << ans << '\n';
        }
    }
}

Student's Camp

有一个 (n+2) \times m 的网格图，其中第一层和最后一层的所有单元格始终存在。中间 n 层的每一层初始均有 m 个单元格。在 k 轮中，每一层左右两端的单元格各有 p 的概率被移除。求最终这 n+2 层单元格构成的图形是连通的概率，对 10^9+7 取模。

dp，设 f(i,l,r) 表示第 i 层剩余区间为 [l,r] 且前 i 层连通的概率。记 g(i) 表示 k 轮中某一层的某一端移除恰好 i 个单元格的概率，那么立即有：

\begin{aligned} g(i)&=\binom{k}{i}p^i(1-p)^{k-i}\\ f(i,l,r)&=g(l-1)g(m-r)\sum_{[l_0,r_0]\cap[l,r]\neq\varnothing} f(i-1,l_0,r_0) \end{aligned}

边界 f(0,1,m)=1，答案即为 \sum_{[l,r]} f(n+1,l,r)。时间复杂度 O(nm^4) 难以接受。

尝试优化 [l_0,r_0]\cap [l,r]\neq\varnothing 这一个条件，不过我们发现如果无交，那么只有 l_0>r 和 r_0<l 两种可能，并且这两种互斥，因此可以改为：

f(i,l,r)=g(l-1)g(m-r)\left(\sum_{[l_0,r_0]}f(i-1,l_0,r_0)-\sum_{[l_0,r_0],l_0>r} f(i-1,l_0,r_0)-\sum_{[l_0,r_0],r_0<l}f(i-1,l_0,r_0)\right)

那么我们实际上只需要处理这三个和式的信息就好了：

\begin{aligned} &F(i)=\sum_{[l_0,r_0]} f(i,l_0,r_0)\quad R(i,r)=\sum_{[l_0,r_0],l_0>r}f(i,l_0,r_0)\quad L(i,l)=\sum_{[l_0,r_0],r_0<l} f(i,l_0,r_0)\\ \implies& f(i,l,r)=g(l-1)g(m-r)(F(i-1)-R(i,r)-L(i,l)) \end{aligned}

根据对称性，很容易得到 L(i,l)=R(i,m-l+1)，因此只需要考虑 L，然后我们要削减状态数，由于 F,L 都只需要 \sum_l f(i,l,r)，因此我们可以转移这个：

\begin{aligned} h(i,r)&=\sum_{l=1}^{r} f(i,l,r)=\sum_{l=1}^{r} g(l-1)g(m-r)(F(i-1)-L(i,m-r+1)-L(i,l))\\ &=g(m-r)\left((F(i-1)-L(i,m-r+1))\sum_{l=1}^r g(l-1)-\sum_{l=1}^r g(l-1)L(i,l)\right) \end{aligned}

这样状态数就是 O(nm) 了，转移的话可以先预处理 g 前缀和，通过 L(i,l)=\sum_{i=1}^{l-1} h(i,l) 来处理出 L，通过 h 总和处理出 F。

然后处理出 g(l-1)L(i,l) 的前缀和，就可以 O(nm) 转移了。时间复杂度 O(nm+k)。

点亮 / The Mighty Spell

给定一个长度为 n 的颜色序列 c，其中 c_i \in [1, m]。序列中每个位置有 1/2 的概率被选中。

若存在某种颜色在所有被选中的位置中均未出现，则权值为 0。否则，权值为所有被选中位置组成的连续段的贡献之和。若一个连续段长度为 l，其贡献为 f(l)=2l^3 + 3l^2 + 3l + 3。

求所有选法的期望权值，答案乘上 2^n 后对 10^9 + 7 取模。m\leq 50。

利用期望的线性性，我们只需要研究每一个区间，钦定它成为连续段的贡献即可，设钦定的这个区间为 [l,r]，那么贡献由这几部分构成：

当然直接讲这三部分乘起来是不对的，因为我们还必须保证 [l-1,r],[l,r+1] 之类的区间不是连续段，就需要非常恶心的分讨，不过枚举区间复杂度已经高达 O(n^2)。

优化到更低复杂度前，我们必须规避这个恶心的分类讨论。一个高明的想法是，我们进行一种类似反演的操作，设：

f(r-l+1)=\sum_{[l_0,r_0]\subseteq [l,r]} g(r_0-l_0+1)\implies g(i)=f(i)-2f(i-1)+f(i-2)

由于我们实际上会将 [l,r] 的贡献，分摊到每一个子区间算一次，所以只需要将 f 换成 g，就可以去掉恶心分讨了，将式子写出来：

\left(\frac{1}{2}\right)^{r-l+1} g(r-l+1) \prod_{i=1}^m [i\not\in S[l,r]] \left(1-\frac{1}{2^{\sum_j [c_j=i]}}\right)

我们不妨枚举 r，后面的积式随着 l 的变化只会改变 O(m) 次（即当 S[l,r] 改变时），于是我们可以提前处理 (1/2)^i g(i) 的前缀和做到 O(nm)。

namespace solution {
    const int N = 2e5 + 5, M = 55, mod = 1e9 + 7;
    int n, m, a[N], cnt[M], bkt[M];
    using mint = GF<int, mod>;
    mint w[M], wi[M], F[N];
    mint f(mint x){ return  x * x * x * 2 + x * x * 3 + x * 3 + 3; }
    mint g(mint x){
        if(x == 1) return f(x);
        if(x == 2) return f(x) - f(x - 1) * 2;
        return f(x) - f(x - 1) * 2 + f(x - 2);
    }
    void solve(){
        cin >> n >> m; mint all = 1;
        for(int i = 1; i <= n; i++) cin >> a[i], cnt[a[i]]++;
        for(int i = 1; i <= m; i++){
            w[i] = mint(1) - ((~mint(2)) ^ cnt[i]), wi[i] = ~w[i];
            all *= w[i];
        }
        for(int i = 1; i <= n; i++) F[i] = F[i - 1] + ((~mint(2)) ^ i) * g(i);
        std::set<int, std::greater<>> s;
        mint ans = 0;
        for(int r = 1; r <= n; r++){
            if(bkt[a[r]]) s.erase(bkt[a[r]]);
            s.insert(bkt[a[r]] = r);
            mint tmp = all;
            for(auto it = s.begin(); it != s.end(); it++){
                int lr = *it, ll = std::next(it) == s.end() ? 1 : (*std::next(it) + 1);
                tmp *= wi[a[lr]], ans += (F[r - ll + 1] - F[r - lr]) * tmp;
            }
        }
        std::cout << (ans * (mint(2) ^ n)) << "\n";
    }
}

矩阵计数 / Lost Table

给定一个长度为 n 的序列 R 和一个长度为 m 的序列 C，计算非负整数构成的 n\times m 矩阵数量，满足行 i 最大值为 R_i，列 i 最大值为 C_i。答案对 10^9+7 取模。

交换 R,C 中的元素，对答案没有影响，因此不妨将 R,C 排序，然后记 A_{i,j}=\min(R_i,C_j) 表示每个位置的上界，那么容易发现 A 的分布构成了若干个 \Gamma 形。举个例子：这里以 C=\{1,1,2,3,4\},R=\{2,3,3,4\} 为例：

1 1 2 3 4
2 \color{red}1 \color{red}1 \color{yellow}2 \color{yellow}2 \color{yellow}2
3 \color{red}1 \color{red}1 \color{yellow}2 \color{green}3 \color{green}3
3 \color{red}1 \color{red}1 \color{yellow}2 \color{green}3 \color{green}3
4 \color{red}1 \color{red}1 \color{yellow}2 \color{green}3 \color{blue}4

考察每一个 \Gamma 形，设其 A 值为 v。记它的长（列数）为 N，宽（行数）为 M，横部分的宽为 a，纵部分的宽为 b，那么总单元格数为 S=Na+bM-ab。

我们要在每个单元格放 [0,v] 的数，但是要求 R_i=v 的 a 行和 C_i=v 的 b 列都至少有一个 v，求方案数。

至少有一个很难做（行列有重叠），但是一个也没有就很好做，那么不妨容斥，枚举 i,j 表示 i 行 j 列没有 v，其余位置随便填即可：

\sum_{i=0}^a\sum_{j=0}^b (-1)^{i+j}\binom{a}{i}\binom{b}{j}v^{Ni+Mj-ij}(v+1)^{S-(Ni+Mj-ij)}

直接做时间复杂度 \sum ab = O(n^2)，不过事实上可以优化，我们分离 i,j：

\begin{aligned} \implies&(v+1)^S\sum_{i=0}^a(-1)^i\binom{a}{i}v^{Ni}(v+1)^{-Ni}\sum_{j=0}^b(-1)^j\binom{b}{j}v^{j(M-i)}(v+1)^{-j(M-i)}\\ =&(v+1)^S\sum_{i=0}^a(-1)^i\binom{a}{i}\left(\frac{v}{v+1}\right)^{Ni}\sum_{j=0}^b \left(-\left(\frac{v}{v+1}\right)^{M-i}\right)^j\binom{b}{j}\\ =&(v+1)^S\sum_{i=0}^a(-1)^i\binom{a}{i}\left(\frac{v}{v+1}\right)^{Ni}\left(1-\left(\frac{v}{v+1}\right)^{M-i}\right)^b \end{aligned}

这样就可以做到 \sum a=O(n)（忽略快速幂、求逆等）了。

namespace solution {
    const int N = 2e5 + 5, mod = 1e9 + 7;
    using mint = GF<int, mod>;
    mint fact[N], ifact[N];
    int n, m, R[N], C[N], bkt[N << 1];
    mint binom(int n, int m){ return m < 0 || n < m ? 0 : fact[n] * ifact[m] * ifact[n - m]; }
    void solve(){
        cin >> n >> m; int _ = std::max(n, m);
        for(int i = 1; i <= n; i++) cin >> R[i], R[i]--, bkt[i] = R[i];
        for(int i = 1; i <= m; i++) cin >> C[i], C[i]--, bkt[i + n] = C[i];
        std::sort(R + 1, R + n + 1), std::sort(C + 1, C + m + 1);
        std::sort(bkt + 1, bkt + n + m + 1);
        int tot = std::unique(bkt + 1, bkt + n + m + 1) - bkt - 1;
        fact[0] = 1; for(int i = 1; i <= _; i++) fact[i] = fact[i - 1] * i;
        ifact[_] = ~fact[_]; for(int i = _; i; i--) ifact[i - 1] = ifact[i] * i;
        mint ans = 1;
        for(int i = 1, p = 1, q = 1; i <= tot; i++){
            int a = 0, b = 0, n_ = n - p + 1, m_ = m - q + 1;
            while(p <= n && R[p] == bkt[i]) b++, p++;
            while(q <= m && C[q] == bkt[i]) a++, q++;
            mint ret = 0, d = mint(bkt[i]) / mint(bkt[i] + 1);
            int S = (1ll * n_ * a + 1ll * b * m_ - 1ll * a * b) % (mod - 1);
            for(int i = 0; i <= a; i++){
                ret += binom(a, i) * (d ^ (1ll * n_ * i % (mod - 1))) * ((mint(1) - (d ^ (m_ - i))) ^ b) * (i & 1 ? mod -  1 : 1);
            }
            ans *= ret * (mint(bkt[i] + 1) ^ S);
        }
        std::cout << ans << '\n';
    }
}

棒棒糖

给定 N, M, x。对于所有 1 \le n \le N, 1 \le m \le M，求满足 \sum_{i=1}^n a_i = m 的自然数序列 a 的“前 x 大元素之和”的总和。

答案对 998244353 取模。

看不懂 dp，只能生成函数了，记答案为 f(n,m)，那么：

f(n,m)=\sum_a\sum_{v\geq 1}\min(x,\#\{a_i\geq v\})=\sum_{v\geq 1}\sum_{k=0}^{n}\min(x,k)F(n,m,k,v)

其中 F(n,m,k,v) 即长度为 n，总和为 m，恰好 k 个数 \geq v 的方案数。也就是：

F(n,m,k,v)=[z^m]\binom{n}{k}\left(\frac{z^v}{1-z}\right)^k\left(\frac{1-z^v}{1-z}\right)^{n-k}=[z^m]\binom{n}{k}z^{kv}\frac{(1-z^v)^{n-k}}{(1-z)^n}

顺势计算 \min(x,k)F(n,m,k,v) 固定 v 关于 m 的生成函数：

\begin{aligned} &\sum_{k=0}^n \min(x,k)\binom{n}{k}z^{kv}\frac{(1-z^v)^{n-k}}{(1-z)^n}=(1-z)^{-n}\sum_{k=0}^{n}\binom{n}{k}\min(x,k)(z^v)^k\sum_{i=0}^{n-k}\binom{n-k}{i}(-1)^i(z^v)^i\\ &=(1-z)^{-n}\sum_{S=0}^n(z^v)^S\sum_{k=0}^S \binom{n}{k}\min(x,k)\binom{n-k}{S-k}(-1)^{S-k}\\ &=(1-z)^{-n}\sum_{S=0}^n\binom{n}{S}z^{vS}\sum_{k=0}^S\binom{S}{k}\min(x,k)(-1)^{S-k} \end{aligned}

然后枚举 v：

\begin{aligned} &(1-z)^{-n}\sum_{v\geq 1}\sum_{S=0}^n\binom{n}{S}z^{vS}\sum_{k=0}^S\binom{S}{k}\min(x,k)(-1)^{S-k} \end{aligned}

时间复杂度 O(nm^2)。

namespace solution {
    const int N_ = 505, mod = 998244353;
    int N, M, x; using mint = GF<int, mod>;
    mint fact[N_ << 1], ifact[N_ << 1], f[N_], g[N_], h[N_], ans[N_];
    mint binom(int n, int m){ return n < m ? mint(0) : fact[n] * ifact[m] * ifact[n - m]; }
    void solve(){
        cin >> N >> M >> x; int _ = std::max(N, M) << 1;
        fact[0] = 1; for(int i = 1; i <= _; i++) fact[i] = fact[i - 1] * i;
        ifact[_] = ~fact[_]; for(int i = _; i; i--) ifact[i - 1] = ifact[i] * i;
        for(int S = 0; S <= N; S++){
            for(int k = 0; k <= S; k++) f[S] += binom(S, k) * std::min(x, k) * ((S - k) & 1 ? mod - 1 : 1);
        }
        for(int n = 1; n <= N; n++){
            for(int v = 1; v <= M; v++){
                for(int S = 0; S <= n && 1ll * S * v <= M; S++) g[v * S] += binom(n, S) * f[S];
            }
            for(int v = 0; v <= M; v++) h[v] = binom(v + n - 1, v);
            for(int a = 0; a <= M; a++) for(int b = 0; b <= a; b++) ans[a] += g[a - b] * h[b];
            for(int a = 1; a <= M; a++) std::cout << ans[a] << ' ';
            std::cout << '\n';
            for(int a = 0; a <= M; a++) ans[a] = g[a] = h[a] = 0;
        }
    }
}

异或图

给定一张 n 个点 m 条边的无向图和一个长度为 n 的数组 a_1, a_2, \cdots , a_n 以及一个整数 C，你需要求出有多少个长度为 n 的数组 b 满足：

  2. 对于每条边 (u, v)，b_u \ne b_v。

答案对 998244353 取模。

首先考虑 m=0 怎么做，对于 C\neq 0 的情况，我们在 b 序列末尾添加一个 C，计算 C=0 的答案，然后将末尾改为 C-1，再计算 C=0 的答案，两个答案作差就是最终答案（因为相当于钦定多填一个 C，剩余部分异或和就是 C 了）。所以我们只需要思考 C=0 的做法。这其实也是一道题，POI 2006 KRY-Crystals。

假如存在一个 a_i 足够大，以至于可以装下其他所有数的异或和，那么就可以简单用乘法原理计算其他位置的方案数。对于一般的情况，我们只需要拆位，即枚举 k 表示所有 b 的前 k 位与 a 均相同，我们希望存在一个 p，使得 a_p 在第 k 位是 1，但是 b_p 是 0，且 a,b 前 k 位都相同。那么此时对于低于 k 位的情况就可以随便选，相当于可以装下其他位置的异或和了。

枚举 k 后，设 \text{dp}(i,0/1,0/1) 表示考虑到前 i 个数，第 k 位异或和为 0/1，现在是否选出了目标的 p。转移的话分成如下两类：

最后统计 f(n,0,1) 即可（对于最后一位还要加上 f(n,0,0)，因为有没有 p 不重要），注意如果出现一位，所有数在这一位 1 的数量为奇数，那么就没有办法转移了（因为之后的转移，高位异或和都不再为 1），退出。

这样时间复杂度就可以做到 O(n\log V) 了，足够优秀而为我们所用。

现在把这张无向图的信息加上来。我们发现强制钦定两个点取不同的值极其困难，但是钦定取相同的值是容易的，因此可以 O(2^m) 枚举边集，钦定边集里的边都不满足两端点权值不同的限制，然后对每个边集 E 给一个 (-1)^{|E|} 的容斥系数即可，但是 2^m 很大难以接受。

对着边集容斥显然没有前途，不妨转为点集容斥，我们发现边集的限制本质上相当于给点集划分为若干个连通块，每个连通块要求填相同的数。那么不妨对连通块分开维护，考察连通块点集 S 的总容斥系数 f(S)=\sum_{G=(S,E) \text{ connected}} (-1)^{|E|}。

如何快速求 f(S)？这里有两个视角，最终会得到同样的结果：

因此可以 O(3^n) 或 O(n^2 2^n) 求出 f(S)，前者即足以通过。

最后我们要枚举点集的连通块划分方法，对于每一个连通块划分方法，大小为奇数的连通块相当于上限就是连通块中点的 a 的最小值，大小为偶数则可以随便填，简单乘上 (\min + 1) 作为贡献。

朴素枚举量为 \mathrm{BellB}(n) 难以接受，这一步可以 dp 设 F(S,T) 表示当前考虑完了点集 S 的元素，划分中所有大小为奇数的连通块中所有的 a 的最小值位置构成了集合 T。时间复杂度高达 O(4^n) 也难以接受。

将所有点按照 a 排序，然后枚举每个点 i 尝试作为连通块（可能大小为奇数，也可能是偶数）的最小元加进去，那么我们发现记录两个集合 S,T 是多余的，我们只需要记录 (S\cup T)\backslash ([1,i)\cap (S\backslash T)) 即可，因为 <i 的数已经考虑完他在不在 T 了，\geq i 的数如果在，那么一定只在 S 中。

转移的时候枚举集合 S，如果包含 i 就将它去掉（说明只在 S 中），否则枚举包含 i 的连通块集合 T 分大小的奇偶性进行转移即可。

最后会收集到一个一维数组 F(T)，实际上就是 F(U,T)，那么将 T 中的 a 取出来跑一遍那个 O(n\log V) 的 dp 即可。

总时间复杂度 O(n 3^n + 2^n n\log V) 可以通过本题。

namespace solution {
    using u64 = unsigned long long;
    const int mod = 998244353, N = (1 << 15) + 5;
    using mint = GF<int, mod>;
    u64 a[25], C;
    int n, m, e[25][25], id[25], rev[25];
    mint ind[N], f[N], g[25][N];
    mint crystal(std::vector<u64> a){
        static mint f[25][2][2];
        int n = a.size(); a.insert(a.begin(), 0ll);
        mint ans = 0;
        for(int k = 62; ~k; k--){
            memset(f, 0, sizeof(f)), f[0][0][0] = 1;
            u64 tag = 0; mint pw2 = mint(2) ^ k;
            for(int i = 1; i <= n; i++){
                u64 _ = ((a[i] >> k) & 1); mint c = mint::fit((a[i] % (1ull << k)) + 1);
                tag ^= a[i] >> (k + 1);
                for(int u : {0, 1}) for(int v : {0, 1}) f[i][u ^ _][v] += f[i - 1][u][v] * c;
                if(_) for(int u : {0, 1}) f[i][u][1] += (f[i - 1][u][1] * pw2) + f[i - 1][u][0];
            }
            if(tag) break;
            ans += f[n][0][1];
            if(!k) ans += f[n][0][0];
        }
        return ans;
    }
    mint solve(std::vector<u64> a){
        if(!C) return crystal(a);
        a.push_back(C); mint ans = crystal(a);
        a.back()--; ans -= crystal(a);
        return ans;
    }
    void solve(){
        cin >> n >> m >> C;
        for(int i = 1; i <= n; i++) cin >> a[i], id[i] = i;
        std::sort(id + 1, id + n + 1, [](int x, int y){ return a[x] < a[y]; });
        for(int i = 1; i <= n; i++) rev[id[i]] = i;
        std::sort(a + 1, a + n + 1);
        for(int i = 1, u, v; i <= m; i++){
            cin >> u >> v, u = rev[u], v = rev[v], e[u][v] = e[v][u] = 1;
        }
        for(int S = 0; S < (1 << n); S++){
            bool ok = true;
            for(int i = 1; i <= n; i++){
                if(!(S >> (i - 1) & 1)) continue;
                for(int j = i + 1; j <= n; j++){
                    if(!(S >> (j - 1) & 1)) continue;
                    if(e[i][j]){ ok = false; break; }
                }
                if(!ok) break;
            }
            ind[S] = ok;
        }
        for(int S = 0; S < (1 << n); S++){
            f[S] = ind[S]; int min = S & -S;
            for(int T = S & (S - 1); T; T = (T - 1) & S) if(T & min) f[S] -= ind[S ^ T] * f[T];
        }
        g[0][0] = 1;
        for(int i = 1; i <= n; i++){
            for(int S = 0; S < (1 << n); S++){
                if(S >> (i - 1) & 1){ g[i][S ^ (1 << (i - 1))] += g[i - 1][S]; continue; }
                int all = (S ^ ((1 << n) - 1)) & ((1 << n) - (1 << i));
                for(int T = all; ; T = (T - 1) & all){
                    mint c = g[i - 1][S] * f[T | (1 << (i - 1))];
                    if(__builtin_popcount(T) & 1) g[i][S | T] += c * (mint::fit(a[i]) + 1);
                    else g[i][S | T | (1 << (i - 1))] += c;
                    if(!T) break;
                }
            }
        }
        mint ans = 0;
        for(int S = 0; S < (1 << n); S++){
            std::vector<u64> vct;
            for(int i = 1; i <= n; i++) if(S >> (i - 1) & 1) vct.push_back(a[i]);
            ans += solve(vct) * g[n][S];
        }
        std::cout << ans << '\n';
    }
}

错排

我们从 n 个数中选出 m 个 >m 的数组成排列的前 m 个数（这一步方案数是 \binom{n-m}{m} m!）那么剩下的数组成的排列中 \leq m 的数相当于没有错位的限制，而其余的 n-2m 个 >m 的数有不能放到 p_i=i 的位置的限制，设 f(n,m) 表示长度为 n+m 的排列中，n 个数有错位限制，m 个数没有的方案数，那么可以写作 f(n-2m,m)。即答案就是 \binom{n-m}{m}m! f(n-2m, m)。

尝试求出 f 的递推公式。一种视角是考虑其中下一个没有限制的数，一种方法是直接放在对应编号处的 f(n,m-1)，另一种方法是放在其他位置，等价于钦定不能放在对应编号处，也就是改为了一个有限制的数，即 f(n+1,m-1)，得到：

f(n,m)=f(n,m-1)+f(n+1,m-1)

另一种视角是，我们考虑第一个位置怎么填，如果填一个没有限制的数，那么去掉第一个位置后，我们少了一个没有限制的数，但是原本限制位置为第一位的那个元素，可以随便放了，所以少了一个有限制的数，多了一个没有限制的数，总体来看就是少了一个有限制的数，即 m f(n-1,m)。如果填一个有限制，且限制位置部位第一位的元素，那么去掉这个元素后，限制位置为第一位的可以随便放了，相当于少了两个有限制的数，多了一个没有限制的数，即 (n-1)f(n-2,m+1)，得到：

f(n,m)=mf(n-1,m)+(n-1)f(n-2,m+1)=mf(n-1,m)+(n-1)(f(n-1,m)+f(n-2,m))

我们给出一个常数 C，尝试从 f(,m-C) 转移到 f(,m)。那么对着第一条递推式，我们可以枚举 f(n+1,m-1) 的转移次数得到：

f(n,m)=\sum_{k=0}^C \binom{C}{k}f(n+k,m-C)

这样只需要取 C=O(\sqrt{n}) 然后递推出 C,2C,3C,\cdots 这些列的 f，然后询问的时候 O(C) 用上面那个式子就可以求出答案。

时间复杂度 O((T+n)\sqrt{n})，假设 n,m 同阶。

namespace solution {
    const int N = 2e5, N_ = 2e5 + 5, B = 447, B_ = 455, mod = 998244353;
    using mint = GF<int, mod>;
    mint fact[N_], ifact[N_], f[B_][N_];
    mint binom(int n, int m){ return m < 0 || n < m ? 0 : fact[n] * ifact[m] * ifact[n - m]; }
    void solve(){
        fact[0] = 1; for(int i = 1; i <= N; i++) fact[i] = fact[i - 1] * i;
        ifact[N] = ~fact[N]; for(int i = N; i; i--) ifact[i - 1] = ifact[i] * i;
        for(int i = 0; i <= B; i++){
            f[i][0] = fact[i * B], f[i][1] = fact[i * B] * i * B;
            for(int j = 2; j <= N; j++) f[i][j] = f[i][j - 1] * i * B + (f[i][j - 1] + f[i][j - 2]) * (j - 1);
        }
        int T; cin >> T;
        while(T--){
            int n, m; cin >> n >> m; if(n < 2 * m){ std::cout << "0\n"; continue; }
            mint ans = 0; int k = m % B;
            for(int i = 0; i <= k; i++) ans += binom(k, i) * f[m / B][n - 2 * m + i];
            std::cout << ans * fact[m] * binom(n - m, m) << '\n';
        }
    }
}

匹配

给出一个 2n 个点（n 个左部点 n 个右部点）m 条边的二分图，每条边有黑白两种颜色，你需要选出一个黑色边数量为偶数的完美匹配，或指出无解。

先用 Dinic 等最大流算法随便选出一个完美匹配，如果不存在肯定无解，如果这组完美匹配恰好黑色边数量为偶数，那么直接输出这组解。

我们尝试对最大流的结果进行一次反悔来替换掉若干个匹配边，具体来说，我们在最大流的残量网络上找到一个环，那么环上肯定是被流过的边撤销流，没被流过的边添加流，dfs 找出这样的环就好了。

回响形态

给出一个长度为 n 的字符串 S，q 组询问，每组询问给出 k，你需要求出所有中心为 \frac{k}{2} 的子串的 Border 个数之和。

首先求出以 \frac{k}{2} 为根的最长子串 S’（假设长度为 m），那么问题就变成了计算按照 \frac{k}{2} 相等且位置对称的子串对数量。

我们尝试构造一个字符串 T=S’_1S’_mS’_2S’_{m-1}\cdots，那么我们断言，每一个长度为偶数，且起始位置为奇数的回文子串，都代表了一个这样的子串对，这是很自然的。

那么直接上 Manacher 即可。时间复杂度 O(nq) 考场没过真是遗憾。

Unpair Ampere

给出一个 n 个点 m 条边的 DAG，所有入度为 0 的点被分成两个集合 A,B，一个入度为 0 的点转换所属的集合代价为 a_u。其余点如果既能被 A 中的点到达，又能被 B 中的点到达，则代价为 a_u。求出最小代价和。

奇奇怪怪的最小代价，首选最小割建模。这道题建模也很直观，划分到集合为 A，我们认为是和 S 连通，划分到集合 B，我们认为是和 T 连通。

对于入度不为 0 的点，要刻画不能免费被两个集合的点连接的性质，可以建立一张反图来实现。

这样我们就得到了初步思路：对于每一个点 u 拆成两个点 u,u’ 分别表示正图和反图，然后：

然后就做完了。

音乐课

给定一个长度为 n 的序列 h，q 次操作，支持：

  • 1 p v：h_p\gets v。
  • 2 m k：构造一个长度为 n 的序列 a，使得 |a_i|\leq m, \sum a_i=k。最大化 \sum a_i h_i。求这个最大值。

考虑如果是 a_i\in[0,m] 怎么做，那么一个显然的方法是将 h_i 从小到大排序，然后从后往前，能放 m 就放 m，以此类推。

对于本题可以放负数，简单的策略不再适用，给式子变一下形 \sum a_ih_i=m\sum h_i-\sum (m-a_i)h_i，那么就有 m-a_i\in [0,2m],\sum m-a_i=nm-k，最小化 \sum h_i(m-a_i)，那么就转换成刚才的问题，将 h_i 较小的填上 m-a_i=2m 即可。

如果用数据结构维护，我们就需要支持查询全局和，以及前 k 小/大的和，用 fhq treap 维护即可。时间复杂度 O((n+q)\log n)。

Surprise me!

给定一个长度为 n 的序列 a 和一个 n 个点的树，计算：

\frac{1}{n(n-1)}\sum_{i=1}^n\sum_{j=1}^n \varphi(a_ia_j) \operatorname{dis}(i,j)

答案对 10^9+7 取模。

套路推式子题：

\begin{aligned} &\sum_{i=1}^n\sum_{j=1}^n \varphi(a_ia_j)\operatorname{dis}(i,j)\\ =&\sum_{i=1}^n\sum_{j=1}^n \frac{\varphi(a_i)\varphi(a_j) \gcd(a_i,a_j)}{\varphi(\gcd(a_i,a_j))} \operatorname{dis}(i,j)\\ =&\sum_{d=1}^n\frac{d}{\varphi(d)}\underbrace{\sum_{i=1}^n\sum_{j=1}^n[\gcd(a_i,a_j)=d]\varphi(a_i) \varphi(a_j) \operatorname{dis}(i,j)}_{f(d)}\\ \end{aligned}

因此可以设：

\begin{aligned} &F(d)=\sum_{k\mid d} f(k)=\sum_{i=1}^{n}\sum_{j=1}^{n}[d\mid \gcd(a_i,a_j)]\varphi(a_i)\varphi(a_j)\operatorname{dis}(i,j)\\ \implies&f(d)=\sum_{k\mid d}F(k)\mu\left(\frac{d}{k}\right) \end{aligned}

因此只需要计算 F(d)，我们将所有满足 d\mid a_i 的点找出来，建一个虚树，那么就要计算 \sum_i \sum_j \varphi(a_i) \varphi(a_j) \operatorname{dis}(i,j)：

\begin{aligned} &\sum_i\sum_j \varphi(a_i)\varphi(a_j)\operatorname{dis}(i,j)=\sum_i\sum_j \varphi(a_i)\varphi(a_j)(d_i+d_j-2d_{\text{lca}(i,j)})\\ =&\sum_i\sum_j\varphi(a_i)\varphi(a_j)d_i+\sum_i\sum_j \varphi(a_i)\varphi(a_j)d_j-2\sum_i\sum_j \varphi(a_i)\varphi(a_j)d_{\text{lca}(i,j)}\\ =&2\sum_i\varphi(a_i)d_i\sum_j\varphi(a_j)-2\sum_i\sum_j \varphi(a_i)\varphi(a_j)d_{\text{lca}(i,j)} \end{aligned}

前面一个很好求，后面那个双重和式，直接在虚树上做一个 dp，在 \text{lca}(i,j) 处统计点对 (i,j) 的贡献即可，时间复杂度 O(n\log^2 n)。

数字电路

给定一个 n+m 个点的树，其中 m 个点是叶子，每个叶子有一个颜色 A_u\in\{0,1\}，每个非叶子节点 u 有一个位置的权值 P_u\in [0,|\text{son}(u)|]。

对于一个非叶子节点 u，如果它的至少 P_u 个儿子颜色是 1，则 A_u=1，否则 A_u=0，记此时的答案为使得根节点颜色为 1 的分配各点权值方案数。

有 q 次修改，每次给出一个区间 [l,r]，你需要反转编号位于该区间的所有叶子的颜色，然后回答现在整棵树的答案。

设 f_u 表示点 u 为黑色的概率（因为概率乘上总方案数就是方案数），一个显然的状态转移方程是：

f_u=\sum_{S\subseteq \text{son}(u)} \frac{|S|}{|\text{son}(u)|}\prod_{v\in S} f_v\prod_{v\in \text{son}(u)\backslash S} (1-f_v)

但是直接做时间复杂度太高了，我们需要优化。注意到我们这个式子从另一个视角看，等同于每个儿子有一定概率被选择，计算选择数量除以 |\text{son}(u)| 的值，那么根据期望的线性性，我们很容易得到：

f_u=\frac{1}{|\text{son}(u)|} \sum_{v\in\text{son}(u)} f_v

边界是对于叶子 u，f_u=a_u，这样我们就可以方便的回答单次询问。如何回答多次询问？贡献的形式本质上相当于从叶子到根的链，因此 f_1 实际上就是：

f_1=\sum_{u\in\text{leaf}} a_u\prod_{v\in (1,u),v\neq u} \frac{1}{|\text{son}(v)|}

而后面那个乘积式可以树上 dfs 递推预处理出来，记为 g_u，那么问题就变成了反转 a_u 的一段区间，询问 \sum_u a_ug_u，使用任何一个你喜欢的数据结构即可。

时间复杂度 O(q\log m)。

白兔迷宫

在一个 n 个点 m 条边的有向图上进行随机游走。从起点 S 出发，每次等概率选择一条出边移动，到达终点 T 时停止。图中存在两类边：

  • 1 号边：经过后计数器加 1。
  • 0 号边：经过后计数器清零。

计数器初始值为 0。求游走停止时，计数器的期望和方差，答案对 998244353 取模。

根据随机变量的方差公式 \text{Var}(X)=E[X^2]-E^2[X]，我们只需要求出计数器的期望值和计数器平方的期望值即可。先来考虑一个（不算很经典）的问题：

维护一个初始值为 0 的计数器 x，不停地投掷一枚公平六面体骰子，如果投掷到：

求终止时计数器的期望值和计数器平方的期望值。

这道题有两种做法，都能得到答案。

最终可以求出 E[X]=\frac{3}{2},E[X^2]=3。

OK，我们可以通过这道题的方法来求出原题。先来尝试应用上一题的第二种做法，那么我们需要下面几个东西：

这样就可以通过 f(S)，求出 E[X]。列出 g,f 转移：

g(x)=\frac{1}{|\text{out}(x)|}\sum_{(x,y,1)\in G} g(y)\quad f(x)=\frac{1}{|\text{out}(x)|}\left(\sum_{(x,y,1)\in G} (f(y)+g(y))+\sum_{(x,y,0)\in G} f(y)\right)

特别的，f(T)=0,g(T)=1。

这样使用高斯消元即可解出全体 f,g。

然后考虑 E[X^2]，我们需要下面几个东西：

这样就可以列出转移：

G(x)=\frac{1}{|\text{out}(x)|}\sum_{(x,y,1)\in G} G(y)+g(y)\quad F(x)=\frac{1}{|\text{out}(x)|}\left(\sum_{(x,y,1)\in G} F(y)+2G(y)+g(y)+\sum_{(x,y,0)\in G} F(y)\right)

特别地，G(T)=F(T)=0，那么同样可以高斯消元求出来。时间复杂度 O(n^3)。

另外也可以使用上一题的第一种做法，可以发现从 x 出发，计数器值为 y，最终到达 T 的计数器值 / 值的平方期望是一次函数 / 二次函数，那么我们可以将各次项系数设出来，然后高斯消元求解，同样是 O(n^3) 的。

Rainbow Balls

给定一个包含 n 种颜色小球的可重集合，第 i 种颜色的小球数量为 a_i。令总球数为 S = \sum a_i。重复以下操作直至所有小球颜色相同：

  • 从中依次随机抽取两个不同的小球，设先后抽出的球颜色分别为 C_1 和 C_2。
  • 令 C_2\gets C_1，然后放回这两个球。

每次操作耗时 1 秒，求操作总时间的期望值。结果对 10^9 + 7 取模。

本题有基于 Martingale Theory 的做法但是严重超纲，不过 zak 给出了一个更加优秀的做法。

我们枚举最后是什么颜色，设 f_i 表示是当前这个颜色的有 i 个，停止事件期望。那么答案就是 \sum f_{a_i}，考虑如何设计这个 f 的递推：

f_i=(f_{i+1}+f_{i-1})p_i+f_i(1-2p_i)+v_i

其中 p_i 表示取出两个球，包含一个这个颜色的球的概率，v_i 表示此时最终会停留在这个颜色的概率（不能是 1 因为有可能并不算是这个颜色，就不能贡献）。

那么我们有：

\begin{aligned} &p_i=\frac{i}{S}\cdot \frac{S-i}{S-1} = \frac{i(S-i)}{S(S-1)}\\ &v_i=p_i(v_{i+1}+v_{i-1})+(1-2p_i)v_i\implies 2p_iv_i=p_i(v_{i+1}+v_{i-1})\implies v_i=\frac{v_{i+1}+v_{i-1}}{2} \end{aligned}

根据边界 v_0=0,v_S=1 可以解出 v_i=\frac{i}{S}。最后我们尝试求 f_i，为了方便我们先化简：

f_i=(f_{i+1}+f_{i-1})p_i+f_i(1-2p_i)+v_i\implies(f_{i+1}+f_{i-1}-2f_i)p_i+v=0\implies 2f_i-f_{i-1}-f_{i+1}=\frac{S-1}{S-i}

记 d_i=f_i-f_{i-1} 得到 d_{i+1}-d_i=-\frac{S-1}{S-i}。由于 f_0=0,f_S=0，因此 f_S-f_0=\sum_{i=1}^S f_i-f_{i-1}=\sum_{i=1}^S d_S=0，代入：

\begin{aligned} &Sd_1-\sum_{i=2}^S \sum_{t=1}^{i-1} \frac{S-1}{S-t}=Sd_1-(S-1)\sum_{t=1}^{S-1} \frac{1}{S-t} (S-t)=Sd_1-(S-1)^2=0\\ \implies& d_1=\frac{(S-1)^2}{S} \end{aligned}

因此 f_1=d_1=\frac{(S-1)^2}{S}，就可以直接递推出全体 f_i 了，时间复杂度 O(\max a_i)。

zbo

给定一棵包含 n 个节点的加权无向树，初始时 1 号点已存在城堡。现有 k 个依次给出的操作，每个操作会在指定节点 d_i 建造一座新城堡。

设当前已建造城堡的节点集合为 S。每次建造后，请输出 \sum_{x \in S} \sum_{y \in S} \text{dis}(x, y)。

注意到：

\sum_{x\in S}\sum_{y\in S}\mathrm{dis}(x,y)=\sum_{x\in S}\sum_{y\in S}d_x+d_y-2d_{\text{lca}(x,y)}=2(|S|-1)\sum_{x\in S}d_x-4\sum_{x\in S}\sum_{y\in S,y>x} d_{\mathrm{lca}(x,y)}

后面的和式应该如何维护？这也是一个经典问题，每次加上一个点 u，就查询 (1,u) 路径和，然后给 (1,u) 路径加 1，这样 u 和前面的点 v 的交集就是 d_{\mathrm{lca}(x,y)}。那么只需要重链剖分之后找一个数据结构维护 \sum a_ib_i 即可。

Nastia Plays with a Tree

Escaping on Beaveractor