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

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] = \frac{1}{A_{len(array)}^{len(array)}} \sum\left(dp[subsequence] + pair[subsequence]\right) \\ ans = \frac{1}{N}\sum_{n=1}^{N} \frac{1}{A_n^n}\sum_{len(array)=n}\left(dp[array] + pair[array]\right) \]

\(N=2\)：

• \(\frac{1}{2}\) 概率得到 \(n=1\)，只有一個 \(array=[1]\)，顯然 \(dp[ [1] ] = 0\)
• \(\frac{1}{2}\) 概率得到 \(n=2\)，此時有 \(\frac{1}{2}\) 生成 \(array=[1,2]\)，\(\frac{1}{2}\) 生成 \(array=[2,1]\)，則

\[\begin{aligned} &dp[[1,2]] = \frac{1}{4}((dp[[1]]+0) + (dp[[2]]+0) + (dp[[1,2]]+0) + (dp[[\emptyset]]+0)) \\ \implies&dp[[1,2]] = 0\\ &dp[[2,1]] = \frac{1}{4}((dp[[1]]+0) + (dp[[2]]+0) + (dp[[2,1]]+1) + (dp[[\emptyset]]+0)) \\ \implies&dp[[2,1]] = \frac{1}{3}\\ &\emptyset 表示空集 \end{aligned} \]

所以當 \(N=2\) 時的期望就是

\[ ans = \frac{1}{2}\left((dp[[1]]+0) + \frac{1}{2}((dp[[1,2]]+0) + (dp[[2,1]]+1))\right)= \frac{1}{3} \]

\(N=3\) 同理，可以自行計算並算出每個序列的 \(dp\) 值。

計算 \(N=3\) 后，我們發現 \(dp[[1,2,3]] = \frac{0}{3}， dp[[2,1,3]] = \frac{1}{3}， dp[[2,3,1]]=\frac{2}{3}，dp[[3,2,1]]=\frac{3}{3}\)，在加上之前算出的 \(dp[[1,2]] = \frac{0}{3}\)，\(dp[[2,1]] = \frac{1}{3}\)，我們可以猜想 \(dp[array] = \frac{pari[array]}{3}\)，可以繼續計算 \(n=4\) 的情況，同樣滿足猜想。

回到最初始的式子

\[\begin{aligned} ans &= \frac{1}{N}\sum_{n=1}^{N} \frac{1}{A_n^n}\sum_{len(array)=n}\left(dp[array] + pair[array]\right) \\ &= \frac{1}{N}\sum_{n=1}^{N} \frac{1}{A_n^n} \left(\frac{4\sum_{len(array)=n} pair[array]}{3}\right) \end{aligned} \]

令 \(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-1}i\)。
現在就可以得到 \(f[i]\) 的線性遞推式

\[f[n] = n*f[n-1]+A_{n-1}^{n-1}\sum_{i=1}^{n-1}i \]

我們隊最后把 \(f[i]\) 扔到 \(oeis\) 里面去....發現居然有 \(O\left(1\right)\) 公式 \(\frac{n!n(n-1)}{4}\)，打擾了....
最后只要把 \(f[i]\) 打表預處理出來，最后的答案就可以進一步化簡

\[ ans = \frac{1}{N} \sum_{n=1}^{N} \frac{1}{A_n^n}\frac{4f[n]}{3} \]

最后發現網上說 \(ans\) 居然又有 \(O\left(1\right)\) 公式 \(\frac{N^2-1}{9}\) 😫 打擾了...

/*************************************************************** 
    > 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;
}