CF1153F Serval and Bonus Problem

官方的解法是\(O(n ^ 2)\)的，這里給出一個\(O(n \log n)\)的做法。

首先對於長度為\(l\)的線段，顯然它的答案就是長度為\(1\)的線段的答案\(\times l\)，這樣做只是為了方便計算。

考慮對於數軸上區間\([0,1]\)內任意一個點\(x\)，它被一條隨機線段覆蓋的概率是多少：線段的兩個端點都在它左邊的概率是\(x ^ 2\)、都在它右邊的概率是\((1 - x) ^ 2\)，那么它被覆蓋的概率即為\(p(x) = 1 - x^2 - (1 - x) ^ 2 = 2 x (1 - x)\)。

那么他被\(\ge k\)條線段覆蓋的概率為\(f(x) = \sum \limits _{i = k} ^ n \binom{n}{i} p(x) ^ i (1 - p(x)) ^ {n - i}\)。

根據定積分的定義就得到區間\([0,1]\)內被\(\ge k\)條線段覆蓋的長度期望為\(\int _{0} ^{1} f(x) \mathrm{d} x\)。

現在我們的問題就是怎么算這個東西了：

\[\begin{align*} \int _ 0 ^ 1 f(x) \mathrm{d} x & = \int _ 0 ^ 1 \sum \limits _{i = k} ^ n \binom {n} {i} (2x (1 - x)) ^ i (1 - 2x (1 - x)) ^ {n - i} \mathrm{d} x \\ & = \sum \limits _{i = k} ^ n \binom {n} {i} \int _ 0 ^ 1 (2x (1 - x)) ^ i \sum \limits _{j = 0} ^ {n - i} \binom {n - i} {j} (-2x(1 - x)) ^ j \mathrm{d} x \\ & = \sum \limits _{i = k} ^n \binom {n} {i} \sum \limits _{j = 0} ^ {n - i} \binom {n - i} {j} (-1) ^ j 2 ^ {(i + j)} \int _ 0 ^ 1 x ^ {i + j} (1 - x) ^ {i + j} \mathrm{d} x \end{align*} \]

推到這里，就可以把積分去掉了，這是一個Beta Function的形式：記結論吧

\[B(x, y) = \int _ 0 ^ 1 t ^ {x - 1} (1 - t) ^ {y - 1} \mathrm{d} t = \frac {(x - 1) ! (y - 1) !} {(x + y - 1) !} \]

代入得：

\[\begin{align*} \int _ 0 ^ 1 f(x) \mathrm{d} x & = \sum \limits _{i = k} ^n \binom {n} {i} \sum \limits _{j = 0} ^ {n - i} \binom {n - i} {j} (-1) ^ j 2 ^ {(i + j)} \frac {((i + j) ! ) ^ 2} {(2(i + j) + 1) !} \\ & = n ! \sum \limits _{i = k} ^ n \sum \limits _ {j = 0} ^ {n - i} \frac {1} {i !} \frac {(-1) ^ j} {j !} \frac {2 ^ {i + j} ((i + j) !) ^ 2} {(2(i+ j) + 1) ! (n - (i + j)) !} \end{align*} \]

令\(t = i + j\)，\(f[i] = \frac {1} {i}\)，\(g[j] = \frac {(-1) ^ j} {j !}\)，\(h[t] = \frac{2^t (t !) ^ 2} {(2 t + 1) ! (n - t) !}\)，考慮枚舉\(t\)：

\[\begin{align*} \int _ 0 ^ 1 f(x) \mathrm{d} x & = n ! \sum \limits _{t = k} ^ n h[t] \sum \limits _ {i = k} ^ t f[i] g[t - i] \end{align*} \]

這后面就是一個非常顯然的卷積式子了，直接FFT即可。

答案記得\(\times l\)。

//written by newbiechd
#include <cstdio>
#include <cctype>
#include <algorithm>
#define BUF 1000000
using namespace std;
const int N = 100003, yyb = 998244353, Gg = 3, Gi = 332748118;
char buf[BUF], *p1, *p2;
inline char gc() { return p1 == p2 && (p2 = (p1 = buf) + fread(buf, 1, BUF, stdin), p1 == p2) ? EOF : *p1++; }
inline int rd() {
    register int f = 0;
    register char c;
    while (!isdigit(c = gc())) {}
    do
        f = f * 10 + (c ^ 48);
    while (isdigit(c = gc()));
    return f;
}
int rev[N], G[2][N];
inline int power(int x, int y) {
    register int o = 1;
    for (; y; y >>= 1, x = 1ll * x * x % yyb)
        if (y & 1)
            o = 1ll * o * x % yyb;
    return o;
}
inline void ntt(int *f, int len, int opt) {
    register int i, j, k, x, y, p, q;
    for (i = 1; i < len; ++i)
        if (i < rev[i])
            swap(f[i], f[rev[i]]);
    for (i = 1; i < len; i <<= 1) {
        p = G[opt][i];
        for (j = 0; j < len; j += i << 1)
            for (k = 0, q = 1; k < i; ++k, q = 1ll * p * q % yyb)
                x = f[j | k], y = 1ll * q * f[i | j | k] % yyb, f[j | k] = (x + y) % yyb, f[i | j | k] = (x - y) % yyb;
    }
}
int fac[N], invFac[N], f[N], g[N], h[N];
int main() {
#ifndef ONLINE_JUDGE
    freopen("a.in", "r", stdin);
    freopen("a.out", "w", stdout);
#endif
    int n = rd(), m = n << 1 | 1, k = rd(), l = rd(), i, len, tmp, ans = 0;
    fac[0] = 1;
    for (i = 1; i <= m; ++i)
        fac[i] = 1ll * fac[i - 1] * i % yyb;
    invFac[m] = power(fac[m], yyb - 2);
    for (i = m; i; --i)
        invFac[i - 1] = 1ll * invFac[i] * i % yyb;
    for (i = k; i <= n; ++i)
        f[i] = invFac[i];
    for (i = 0; i <= n; ++i)
        g[i] = i & 1 ? yyb - invFac[i] : invFac[i];
    for (len = 1; len <= m; len <<= 1) {}
    for (i = 1; i < len; ++i)
        rev[i] = (rev[i >> 1] >> 1) | (i & 1 ? len >> 1 : 0);
    for (i = 1; i < len; i <<= 1)
        G[0][i] = power(Gg, (yyb - 1) / (i << 1)), G[1][i] = power(Gi, (yyb - 1) / (i << 1));
    ntt(f, len, 0), ntt(g, len, 0);
    for (i = 0; i < len; ++i)
        f[i] = 1ll * f[i] * g[i] % yyb;
    ntt(f, len, 1), tmp = power(len, yyb - 2);
    for (i = 1; i <= n; ++i)
        f[i] = 1ll * f[i] * tmp % yyb;
    for (i = 1, tmp = 2; i <= n; tmp = (tmp << 1) % yyb, ++i)
        h[i] = 1ll * fac[i] * fac[i] % yyb * tmp % yyb * invFac[i << 1 | 1] % yyb * invFac[n - i] % yyb;
    for (i = k; i <= n; ++i)
        ans = (1ll * f[i] * h[i] % yyb + ans) % yyb;
    ans = 1ll * ans * fac[n] % yyb * l % yyb, printf("%d\n", (ans + yyb) % yyb);
    return 0;
}