HDU 6595 Everything Is Generated In Equal Probability (期望dp，线性推导)

Posted 2022-05-17 jiaaaaaaaqi

Everything Is Generated In Equal Probability

\[ Time Limit: 1000 ms\quad Memory Limit: 131072 kB \]

题意

给出一个 \(N\)，以相等的概率生成 \(n\) 且 \(n \in [1, N]\)，在以相等的概率生成长度为 \(n\) 的数组，最后将生成的数组扔到 \(CALCULATE\) 函数并返回一个数，问这个数的期望。

思路

先解释一下样例是怎么得来的。
令 \(dp[array]\) 表示数组 \(array\) 扔到函数里得到的期望，\(pair[array]\) 表示数组 \(array\) 中逆序对的数量。
则
\[ dp[array] = \frac1A_len(array)^len(array) \sum\left(dp[subsequence] + pair[subsequence]\right) \ ans = \frac1N\sum_n=1^N \frac1A_n^n\sum_len(array)=n\left(dp[array] + pair[array]\right) \]

\(N=2\)：

• \(\frac12\) 概率得到 \(n=1\)，只有一个 \(array=[1]\)，显然 \(dp[ [1] ] = 0\)
• \(\frac12\) 概率得到 \(n=2\)，此时有 \(\frac12\) 生成 \(array=[1,2]\)，\(\frac12\) 生成 \(array=[2,1]\)，则
  \[ \beginaligned &dp[[1,2]] = \frac14((dp[[1]]+0) + (dp[[2]]+0) + (dp[[1,2]]+0) + (dp[[\emptyset]]+0)) \\implies&dp[[1,2]] = 0\&dp[[2,1]] = \frac14((dp[[1]]+0) + (dp[[2]]+0) + (dp[[2,1]]+1) + (dp[[\emptyset]]+0)) \\implies&dp[[2,1]] = \frac13\&\emptyset 表示空集 \endaligned \]
  所以当 \(N=2\) 时的期望就是
  \[ ans = \frac12\left((dp[[1]]+0) + \frac12((dp[[1,2]]+0) + (dp[[2,1]]+1))\right)= \frac13 \]

\(N=3\) 同理，可以自行计算并算出每个序列的 \(dp\) 值。

计算 \(N=3\) 后，我们发现 \(dp[[1,2,3]] = \frac03， dp[[2,1,3]] = \frac13， dp[[2,3,1]]=\frac23，dp[[3,2,1]]=\frac33\)，在加上之前算出的 \(dp[[1,2]] = \frac03\)，\(dp[[2,1]] = \frac13\)，我们可以猜想 \(dp[array] = \fracpari[array]3\)，可以继续计算 \(n=4\) 的情况，同样满足猜想。

回到最初始的式子
\[ \beginaligned ans &= \frac1N\sum_n=1^N \frac1A_n^n\sum_len(array)=n\left(dp[array] + pair[array]\right) \ &= \frac1N\sum_n=1^N \frac1A_n^n \left(\frac4\sum_len(array)=n pair[array]3\right) \endaligned \]
令 \(f[i]\) 表示 \(\sum_len(array)=i pair[array]\)，只要计算出这个，最后的答案就可以 \(O\left(N\right)\) 得到。
对于一个长度为 \(n\) 的序列，我们假设把这个序列的最后一个数字拿掉，前面的 \(n-1\) 个数的 \(pair\) 数其实就可以看成 \(f[n-1]\) 的贡献，一共有 \(n\) 个数字可以拿掉，所以前 \(n-1\) 个数字这部分的总贡献就是\(n*f[n-1]\)。
现在把最后一个数字加进来，当加入的数字是 \(i\) 时，和其他 \(n-1\) 个数字会产生 \(\left(n-i\right)\) 对逆序对，剩余的 \(n-1\) 个数都在前面，可以随便排列，所以它的贡献就是 \(A_n-1^n-1\left(n-i\right)\)，则最后一个数字这部分的总贡献就是 \(\sum_i=1^n A_n-1^n-1\left(n-i\right) = A_n-1^n-1 \sum_i=1^n-1i\)。
现在就可以得到 \(f[i]\) 的线性递推式
\[ f[n] = n*f[n-1]+A_n-1^n-1\sum_i=1^n-1i \]
我们队最后把 \(f[i]\) 扔到 \(oeis\) 里面去....发现居然有 \(O\left(1\right)\) 公式 \(\fracn!n(n-1)4\)，打扰了....
最后只要把 \(f[i]\) 打表预处理出来，最后的答案就可以进一步化简
\[ ans = \frac1N \sum_n=1^N \frac1A_n^n\frac4f[n]3 \]
最后发现 \(ans\) 居然又有 \(O\left(1\right)\) 公式 \(\fracN^2-13\) ?? 打扰了...

/*************************************************************** 
    > File Name    : hdu6595.cpp
    > Author       : Jiaaaaaaaqi
    > Created Time : 2019年07月26日 星期五 09时17分58秒
 ***************************************************************/

#include <map>
#include <set>
#include <list>
#include <ctime>
#include <cmath>
#include <stack>
#include <queue>
#include <cfloat>
#include <string>
#include <vector>
#include <cstdio>
#include <bitset>
#include <cstdlib>
#include <cstring>
#include <iostream>
#include <algorithm>
#define  lowbit(x)  x & (-x)
#define  mes(a, b)  memset(a, b, sizeof a)
#define  fi         first
#define  se         second
#define  pii        pair<int, int>

typedef unsigned long long int ull;
typedef long long int ll;
const int    maxn = 1e5 + 10;
const int    maxm = 1e5 + 10;
const ll     mod  = 998244353;
const ll     INF  = 1e18 + 100;
const int    inf  = 0x3f3f3f3f;
const double pi   = acos(-1.0);
const double eps  = 1e-8;
using namespace std;

int n, m;
int cas, tol, T;

ll fac[maxn], facinv[maxn];
ll sum[maxn], inv[maxn], f[maxn];

ll fpow(ll a, ll b) 
    ll ans = 1;
    while(b) 
        if(b&1) ans = ans*a%mod;
        a = a*a%mod;
        b >>= 1;
    
    return ans;


void handle() 
    int mx = 3000;
    fac[1] = 1, inv[1] = 1, sum[1] = 1, f[1] = 0;
    for(int i=2; i<=mx; i++) 
        fac[i] = fac[i-1]*i%mod;
        sum[i] = (sum[i-1]+i)%mod;
        inv[i] = (mod-mod/i) * inv[mod%i]%mod;
        f[i] = f[i-1]*i%mod+fac[i-1]*sum[i-1]%mod;
        f[i] %= mod;
    
    facinv[mx] = fpow(fac[mx], mod-2);
    for(int i=mx-1; i>=1; i--) 
        facinv[i] = facinv[i+1]*(i+1)%mod;
    


int main() 
    // freopen("in", "r", stdin);
    handle();
    while(~scanf("%d", &n)) 
        ll ans = 0; 
        for(int i=1; i<=n; i++) 
            ll res = 4ll*f[i]*inv[3]%mod*facinv[i]%mod;
            ans = (ans+res)%mod;
        
        ans = ans*inv[n]%mod;
        printf("%lld\n", ans);
    
    return 0;