QOJ 2603

这个博客出现完全是因为我嫌弃我 cnblogs 断更太久了。

首先发现 \(\mathcal{O}(\sum (n+k))\) 甚至 \(\mathcal{O}(\sum (n+k)\log n)\) 也是可以的，所以对于每一天求出 \(d\) 的期望值就可以了。

先考虑直接暴力求这个玩意。

首先 \(n=1\) 的情况答案就是 \((n+1)+(n+2)+\cdots +(n+k)\)。这边弄一个特殊情况是因为，后面的组合数会有 \(n-2\)。

对于天 \(i\)，我们枚举 \(d\) 的值算出概率再相加就可以了。那么设 \(d=j\)，则其他 \(n-1\) 个房间有 \(x=\frac{n+i-j}{2}\) 个细菌。枚举哪一个是众数（\(n\) 种），再把 \(x\) 个细菌分成 \(n-1\) 个正整数，就是 \(\binom{x-1}{n-2}\)。概率是什么？目前所有的情况有 \(cur=(n+1)\times (n+2)\times \cdots \times (n+i)\) 种，概率就是 \(\frac{n\times \binom{x-1}{n-2}}{cur}\times (i!)\)。要乘以阶乘是因为比如多出 \(2,1\) 个有多种实现顺序（如 \(\{1,1,2\}\),\(\{1,2,1\}\) 之类的）而我们对于每一个细菌都是看成一个不同的个体。那么这个 \(d=j\) 对答案的贡献就是 \(\frac{n\times \binom{x-1}{n-2}}{cur}\times (i!)\times j\)。发现 \(n\times \frac{1}{cur}\times (i!)\) 是可以提出来的，现在我们就是要快速求 \(\sum j\binom{x-1}{n-2}\)。

快速求 \(f_i=\sum j\binom{x-1}{n-2}\)：发现 \(n-2\) 这个是固定不动的，而 \(x-1\) 的范围是 \(0\sim r\)，\(r\) 是一个可以求出来的最大值（具体来说设 \(mn=(n+i)\bmod 2\)，\(r=\frac{n+i-mn}{2}-1\)）。当 \(i\rightarrow i+2\) 的时候，\(r\) 会增加 \(1\)，\(f_i\) 增加 \(mn\binom{r}{n-2}\)。而 \(0\sim r-1\) 的所有组合数都会相比 \(f_{i-2}\) 多出 \(2\) 次（\(j\rightarrow j+2\)），那么 \(f_{i}=f_{i-2}+\sum_{k=0}^{r-1} 2\binom{k}{n-2}+mn\binom{r}{n-2}\)。中间的一项可以用组合数前缀和预处理。

因此我们暴力求出 \(f_1,f_2\) 就可以线性递推出 \(f_3\sim f_k\) 了。因为要求逆元，所以时间复杂度多一个 \(\log\)。

#include <bits/stdc++.h>

using namespace std;

using ll = long long;

const int N = 2e6+6;
const ll mod = 998244353;
 
ll pw(ll x,ll y){
	ll res=1;
	while (y){
		if (y&1){
			res=res*x%mod;
		}
		x=x*x%mod;
		y>>=1;
	}
	return res;
}
 
ll fac[N],inv[N];
 
void init(){
	fac[0]=inv[0]=1;
	for (int i=1; i<N; i++){
		fac[i]=fac[i-1]*i%mod;
	}
	inv[N-1]=pw(fac[N-1],mod-2);
	for (int i=N-2; i; i--){
		inv[i]=inv[i+1]*(i+1)%mod;
	}
}
 
ll C(ll x,ll y){
	if (x==0 && x>=y) return 1;
	if (x<y || x<0 || y<0){
		return 0;
	}
	return fac[x]*inv[x-y]%mod*inv[y]%mod;
}

ll f[N],s[N];

int main(){
	ios::sync_with_stdio(false);
	cin.tie(0);

	int t;
	cin>>t;
	init();
	while (t--){
		int n,k;
		cin>>n>>k;
		if (n==1){
			ll ans=0;
			for (int i=1; i<=k; i++){
				ans=(ans+n+i)%mod;
			}
			cout<<ans<<"\n";
			continue;
		}
		s[0]=C(0,n-2);
		for (int i=1; i<=n+k; i++){
			s[i]=(s[i-1]+C(i,n-2))%mod;
		}
		ll ans=0,cur=1;
		for (int i=1; i<=2; i++){
			f[i]=0;
			for (int j=((n+i)&1)?1:2; j<=n+i; j+=2){
				int x=(n+i-j)/2;
				ll tmp=C(x-1,n-2);
				f[i]=(f[i]+tmp*j%mod)%mod;
			}
		}
		for (int i=3; i<=k; i++){
			int j=(((n+i)&1)?1:2);
			int r=(n+i-j)/2-1;
			f[i]=(f[i-2]+C(r,n-2)*j%mod)%mod;
			f[i]=(f[i]+s[r-1]*2%mod)%mod;
		}
		for (int i=1; i<=k; i++){
			cur=cur*(n+i-1)%mod;
			ll iv=pw(cur,mod-2);
			ll tmp=f[i]*iv%mod*fac[i]%mod*n%mod;
			ans=(ans+tmp)%mod;
		}
		cout<<ans<<"\n";
	}
	return 0;
}

暴力代码：

#include <bits/stdc++.h>

using namespace std;

using ll = long long;

const int N = 1e6+6;
const ll mod = 998244353;
 
ll pw(ll x,ll y){
	ll res=1;
	while (y){
		if (y&1){
			res=res*x%mod;
		}
		x=x*x%mod;
		y>>=1;
	}
	return res;
}
 
ll fac[N],inv[N];
 
void init(){
	fac[0]=inv[0]=1;
	for (int i=1; i<N; i++){
		fac[i]=fac[i-1]*i%mod;
	}
	inv[N-1]=pw(fac[N-1],mod-2);
	for (int i=N-2; i; i--){
		inv[i]=inv[i+1]*(i+1)%mod;
	}
}
 
ll C(ll x,ll y){
	if (x==0 && x>=y) return 1;
	if (x<y || x<0 || y<0){
		return 0;
	}
	return fac[x]*inv[x-y]%mod*inv[y]%mod;
}

int main(){
	ios::sync_with_stdio(false);
	cin.tie(0);

	int t;
	cin>>t;
	init();
	while (t--){
		int n,k;
		cin>>n>>k;
		if (n==1){
			ll ans=0;
			for (int i=1; i<=k; i++){
				ans=(ans+n+i)%mod;
			}
			cout<<ans<<"\n";
			continue;
		}
		ll ans=0,cur=1;
		for (int i=1; i<=k; i++){
			cur=cur*(n+i-1)%mod;
			ll iv=pw(cur,mod-2);
			for (int j=((n+i)&1)?1:2; j<=n+i; j+=2){
				int x=(n+i-j)/2;
				ll tmp=1ll*n*C(x-1,n-2)%mod;
				ans=(ans+tmp*iv%mod*j%mod*fac[i]%mod)%mod;
			}
		}
		cout<<ans<<"\n";
	}
	return 0;
}