跳转至

多线程计算

题意

主板上有 nm 列个处理器,每个处理器完成任务的时间都是 [0,1) 内均匀随机选择的一个实数。完成任务后,处理器会持续亮灯。

给定 h 组数对 (x_1,y_1),(x_2,y_2),\ldots,(x_h,y_j)。主板的一个亮灭状态被称为节能态,当且仅当 \exists i\in[1,h]\cup\mathbf{Z}恰好x_i 行和 y_i 列的处理器亮灯。

F_0=F_1=1,F_n=F_{n-1}+F_{n-2} (n>1)。在节能态触发时,若当时已有 k 个处理器亮灯,则每单位时间会节省 F_k 电量。求出节省的电量的期望值。

解析

由于 nm 个在 [0,1) 上独立分布的均匀随机变量的第 i 大的值的期望是 \dfrac{i}{nm+1},于是有至多 k 个处理器亮灯的期望时长为 \dfrac{k+1}{nm+1},容斥后可得有恰好 k 个处理器亮灯的期望时长为 \dfrac{1}{nm+1}

另一个角度的证明

x 时刻有恰好 k 个处理器亮灯的概率为 f(x)=\dbinom{n}{k}x^k(1-x)^{n-k}。不难得知有恰好 k 个处理器亮灯的期望时长就是 \int_0^1 f(x)\mathrm{d}x

f(x) 的一个原函数是 F(x)=\dbinom{n}{k}x^k\sum\limits_{i=0}^{n-k}\dfrac{\binom{n-k}{i}(-1)^ix^{i+1}}{k+i+1}。那么根据微积分基本定理

\begin{aligned} \int_0^1 f(x)\mathrm{d}x&=F(1)-F(0)\\ &=\binom{n}{k}\sum_{i=0}^{n-k}\dfrac{\binom{n-k}{i}(-1)^i}{k+i+1}\\ &=\sum_{i=0}^{n-k}(-1)^i\dfrac{n!}{k!i!(n-k-i)!(k+i+1)}\\ &=\sum_{i=0}^{n-k}(-1)^i\dfrac{n!}{(n+1)!}\cdot\dfrac{(n+1)!(k+i)!}{k!i!(n-k-i)!(k+i+1)!}\\ &=\dfrac{1}{n+1}\sum_{i=0}^{n-k}(-1)^i\binom{n+1}{k+i+1}\binom{k+i}{k}\\ &=\dfrac{1}{n+1}\sum_{i=0}^n(-1)^{i-k}\binom{n+1}{i+1}\binom{i}{k}\\ &=\dfrac{1}{n+1}\left((-1)^k\sum_i^n(-1)^i\binom{n+1}{i+1}\binom{i}{k}-(-1)^{-1-k}\binom{n+1}{0}\binom{-1}{k}\right) \end{aligned}

运用公式 \sum\limits_i(-1)^i\dbinom{l}{m+i}\dbinom{s+i}{n}=(-1)^{l+m}\dbinom{s-m}{n-l} 得到 \sum\limits_i^n(-1)^i\dbinom{n+1}{i+1}\dbinom{i}{k}=0,故原式等于 \dfrac{1}{n+1}

所用的公式可以从范德蒙德恒等式出发,结合上指标反转公式和对称定律证明。

那么每种恰好有 k 个处理器亮灯的状态的期望时长就是 l(k)=\dfrac{1}{(nm+1)\binom{nm}{k}}

f_{i,j} 表示有恰好 i 行和 j 列填满时节省电量的期望值,我们只需要求出 f_{x_i,y_i},就可以得到答案。

这个东西难以直接求出,考虑运用容斥原理。设 g_{i,j} 表示钦定 i 行和 j 列填满时节省电量的期望值,则

g_{i,j}=\binom{n}{i}\binom{m}{j}\sum_{k=0}^{nm-(im+jn-ij)}\binom{nm-(im+jn-ij)}{k}l(k+(im+jn-ij))

\pi(i)=\sum\limits_{j=0}^{nm-i}\dbinom{nm-i}{j}l(i+j),则 \pi 可以 O(nm\log nm) 卷积计算,g 可以通过 \pi 直接求出。

由于 g_{i,j}=\sum\limits_{x=i}^n\sum\limits_{y=j}^m\binom{x}{i}\binom{y}{j}f_{x,y},根据广义二项式反演公式,我们有

f_{i,j}=\sum_{x=i}^n\sum_{y=j}^m(-1)^{x+y-i-j}\binom{x}{i}\binom{y}{j}g_{x,y}

我们当然会用 NTT 将一维的二项式反演优化至 O(n\log n) 复杂度,但二维呢?我们可以分维反演。设 h_{i,j} 表示钦定 i 行,恰好 j 列填满时节省电量的期望值,则

h_{i,j}=\sum_{k=j}^m(-1)^{k-j}\binom{k}{j}g_{i,k}\\ f_{i,j}=\sum_{k=i}^n(-1)^{k-i}\binom{k}{i}h_{k,j}

于是我们以 O(nm\log nm) 的复杂度完成了本题。利用斐波那契数列的特殊性质 F_n=x_1^n+x_2^n 还可以做到线性复杂度,在此不提。

实现

#include <bits/stdc++.h>

template <typename Tp>
void read(Tp &res){
    static char ch; ch = std::getchar(), res = 0;
    while(!std::isdigit(ch)) ch = std::getchar();
    while(std::isdigit(ch)) res = res * 10 + ch - 48, ch = std::getchar();
}

typedef long long int ll;
const int maxn = 5e5 + 19, mod = 998244353;

int qpow(int a, int b){
    int res = 1;
    while(b){
        if(b & 1) res = (ll)res * a % mod;
        a = (ll)a * a % mod, b >>= 1;
    }
    return res;
}

namespace poly{
    int N, rev[maxn << 2], w[maxn << 2];
    int prepare(int n){
        N = 1;
        while(N < n) N <<= 1;
        for(int i = 0; i < N; ++i) rev[i] = (rev[i >> 1] >> 1) | (i & 1 ? N >> 1 : 0);
        return N;
    }
    void dft(int *f, int b = 1){
        for(int i = 0; i < N; ++i) if(i < rev[i]) std::swap(f[i], f[rev[i]]);
        for(int i = 2; i <= N; i <<= 1){
            w[0] = 1, w[1] = qpow(3, (mod - 1) / i); if(b == -1) w[1] = qpow(w[1], mod - 2);
            for(int j = 2; j < i / 2; ++j) w[j] = (ll)w[j - 1] * w[1] % mod;
            for(int j = 0; j < N; j += i){
                int *g = f + j, *h = f + j + i / 2;
                for(int k = 0; k < i / 2; ++k){
                    int p = g[k], q = (ll)h[k] * w[k] % mod;
                    g[k] = (p + q) % mod, h[k] = (p - q) % mod;
                }
            }
        }
    }
    void idft(int *f){
        dft(f, -1);
        for(int i = 0, inv = qpow(N, mod - 2); i < N; ++i) f[i] = (ll)f[i] * inv % mod;
    }
}

int n, m, h, x[maxn], y[maxn];
int pi[maxn << 2], tmp[maxn << 2];

int fact[maxn], ifact[maxn];
inline int binom(int n, int m){
    return (ll)fact[n] * ifact[n - m] % mod * ifact[m] % mod;
}
void init(int n){
    fact[0] = 1;
    for(int i = 1; i <= n; ++i) fact[i] = (ll)fact[i - 1] * i % mod;
    ifact[n] = qpow(fact[n], mod - 2);
    for(int i = n - 1; i >= 0; --i) ifact[i] = (ll)ifact[i + 1] * (i + 1) % mod;

    static int fib[maxn];
    fib[0] = fib[1] = 1;
    for(int i = 2; i <= n; ++i) fib[i] = (fib[i - 2] + fib[i - 1]) % mod;

    for(int i = 0; i <= n; ++i){
        pi[n - i] = (ll)qpow((ll)(n + 1) * binom(n, i) % mod * fact[n - i] % mod, mod - 2) * fib[i] % mod;
        tmp[i] = ifact[i];
    }
    poly::prepare(n + n + 1);
    poly::dft(pi), poly::dft(tmp);
    for(int i = 0; i < poly::N; ++i) pi[i] = (ll)pi[i] * tmp[i] % mod;
    poly::idft(pi);
    std::reverse(pi, pi + n + 1);
    for(int i = 0; i <= n; ++i) pi[i] = (ll)pi[i] * fact[n - i] % mod;
}

std::vector<std::vector<int> > f; 
int F[maxn << 2], G[maxn << 2];

int main(){
    read(n), read(m), read(h);
    for(int i = 1; i <= h; ++i)
        read(x[i]), read(y[i]);

    init(n * m);
    f.resize(n + 1, std::vector<int>(m + 1));

    for(int i = 0; i <= n; ++i)
        for(int j = 0; j <= m; ++j)
            f[i][j] = (ll)binom(n, i) * binom(m, j) % mod * pi[i * m + j * n - i * j] % mod;

    poly::prepare(m + m + 1);
    for(int i = 0; i <= n; ++i){
        for(int j = 0; j <= m; ++j){
            F[j] = j & 1 ? -ifact[j] : ifact[j];
            G[m - j] = (ll)f[i][j] * fact[j] % mod;
        }
        poly::dft(F), poly::dft(G);
        for(int j = 0; j < poly::N; ++j)
            F[j] = (ll)F[j] * G[j] % mod;
        poly::idft(F);
        for(int j = 0; j <= m; ++j)
            f[i][j] = (ll)F[m - j] * ifact[j] % mod;
        std::fill(F, F + poly::N, 0), std::fill(G, G + poly::N, 0);
    }

    poly::prepare(n + n + 1);
    for(int i = 0; i <= m; ++i){
        for(int j = 0; j <= n; ++j){
            F[j] = j & 1 ? -ifact[j] : ifact[j];
            G[n - j] = (ll)f[j][i] * fact[j] % mod;
        }
        poly::dft(F), poly::dft(G);
        for(int j = 0; j < poly::N; ++j)
            F[j] = (ll)F[j] * G[j] % mod;
        poly::idft(F);
        for(int j = 0; j <= n; ++j)
            f[j][i] = (ll)F[n - j] * ifact[j] % mod;
        std::fill(F, F + poly::N, 0), std::fill(G, G + poly::N, 0);
    }

    int ans = 0;
    for(int i = 1; i <= h; ++i)
        ans = (ans + f[x[i]][y[i]]) % mod;
    std::printf("%d\n", (ans + mod) % mod);
    return 0;
}

来源

UOJ Long Round #1 A 多线程计算

评论