「SOL」开关(LOJ)

并不套路的概率生成函数

---

# 题面

\(n\) 个开关初始是关闭的，当某一些开关打开而另一些开关关闭时可以开门。

随机改变开关的状态——每次有 \(\frac{p_i}{\sum p}\) 的概率更改开关 \(i\) 的状态。求期望操作多少次第一次打开门。

数据规模：\(p_i\ge1,\sum p\le10^5\)。

---

# 解析

- P1. 一些函数的定义

可以拆解成每一个开关最后要变成指定的状态。

记 \(\sum p_i=P\)。由于对于同一个开关来说，操作是无序的，应用 EGF，设 \(F_i(x)\) 是「开关 \(i\) 操作 \(j\) 次的概率 \((\frac{p_i}{P})^j\)」的 EGF。

• 若 \(i\) 要求最后关闭，则操作偶数次：
  \[F_i(x)=\sum_{j=0}\frac{(\frac{p_i}{P})^{2j}x^{2j}}{(2j)!}=\frac{e^{\frac{p_i}{P}x}+e^{-\frac{p_i}{P}x}}{2} \]

• 同理，若要求最后打开，操作奇数次：
  \[F_i(x)=\frac{e^{\frac{p_i}{P}x}-e^{-\frac{p_i}{P}x}}{2} \]

从概率生成函数的角度考虑，要求期望，则设对应随机变量的概率——设 \(f_i\) 为「操作 \(i\) 次打开门」的概率，对应的 OGF 为 \(F(x)\)（普通概率生成函数才有 \(E(X)=F'(1)\)）。

先不管 OGF 怎么求，难道答案就是 \(F'(1)\) 了吗？注意到 \(f_i\) 只表示「操作 \(i\) 次打开了门」而并不是题目要求的「第一次打开门」。

设 \(g_i\) 表示操作了 \(i\) 次后，所有开关的状态恰好不变的概率，也即操作了每个开关偶数次，记其 OGF 为 \(G(x)\)。

最后设出 \(h_i\) 表示「操作 \(i\) 次后第一次打开门」的概率，其 OGF 为 \(H(x)\)，那么 \(H'(1)\) 才是我们真正要求的答案。根据实际意义可以推得

\[H(x)G(x)=F(x)\Rightarrow H(x)=\frac{F(x)}{G(x)} \]

但是我们现在将 \(F_i(x)\) 全部乘起来只能得到 \(f_i\) 的 EGF，记为 \(R(x)\)：

\[R(x)=\prod_{i=1}^nF_i(x) \]

同理，我们也只能快速得到 \(g_i\) 的 EGF 的表达式，记为 \(W(x)\)：

\[W(x)=\prod_{i=1}^n\frac{e^{\frac{p_i}{P}x}+e^{-\frac{p_i}{P}x}}{2} \]

考虑如何把 EGF 转成 OGF。

- P2. EGF to OGF

注意到 \(R(x)\) 和 \(W(x)\) 的形式是类似的。

不难发现 \(R(x)\) 可以表示成如下形式（\(W(x)\) 也可以，下面只通过 \(R(x)\) 举例）：

\[R(x)=\sum_{i}a_ie^{\frac{i}{P}x} \]

只需要把 \(F_i(x)\) 当成一对大小为 \(p_i\) 和 \(-p_i\) 的物品做一次 \(\mathcal{O}(nP)\) 的背包即可求出 \(a_i\)。

而知道了 \(a_i\) 我们就可以较为方便地求出其 OGF \(F(x)\)：

\[\begin{aligned} F(x)&=\sum_{k=0}^{+\infty}k!x^k[x^k]R(x)\\ &=\sum_{k=0}^{+\infty}k!x^k[x^k]\sum_{j}a_je^{\frac{j}{P}x}\\ &=\sum_{k=0}^{+\infty}k!x^k\sum_{j}a_j\frac{(\frac{j}{P})^k}{k!}\\ &=\sum_{j}a_j\sum_{k=0}^{+\infty}\left(\frac{jx}{P}\right)^k\\ &=\sum_{j}\frac{a_j}{1-\frac{jx}{P}} \end{aligned} \]

上式的 \(j\) 的取值根据背包可知应为 \([-P,P]\)。

- P3. 求导

根据之前的推导，我们已经可以求得 OGF \(F(x),G(x)\)，而最终目标

\[H(x)=\frac{F(x)}{G(x)} \]

对其求导并代入 \(x=1\) 可以得到答案：

\[\begin{aligned} &H'(x)=\frac{F'(x)G(x)-G'(x)F(x)}{G^2(x)}\\ \Rightarrow&H'(1)=\frac{F'(1)G(1)-G'(1)F(1)}{G^2(1)} \end{aligned} \]

注意 \(F(x),G(x)\) 都不涵盖所有可能情况，所以 \(F(1)\neq1,G(1)\neq1\)。

但是稍微细心一些，我们发现 \(F(x)\) 和 \(G(x)\) 在 \(1\) 处很有可能不收敛——在 P2 中求出的式子：

\[F(x)=\sum_{j}\frac{a_j}{1-\frac{jx}{P}} \]

若 \(a_P\neq 0\)，有 \(j=P\) 时 \(F(x)\) 在 \(1\) 处不收敛；\(G(x)\) 同理。

正无穷除以正无穷？洛必达

考虑 \(H(x)=\frac{F(x)}{G(x)}\) 分子分母同时乘以 \((1-x)\)，或者说 \(F(x)\) 和 \(G(x)\) 的表达式同时乘以 \((1-x)\)，可以消去不收敛的项。

又因为 \(H(x)\) 是一个形式幂级数，\(x\) 能取什么值并不重要，所以这样的操作是可以的。

乘上 \((1-x)\) 后，注意导函数在 \(x=1\) 处的取值会变成相反数。

---

# 源代码

/*Lucky_Glass*/
#include<cstdio>
#include<cstring>
#include<algorithm>
using namespace std;

inline int rin(int &r){
	int b=1,c=getchar();r=0;
	while(c<'0' || '9'<c) b=c=='-'?-1:b,c=getchar();
	while('0'<=c && c<='9') r=(r<<1)+(r<<3)+(c^'0'),c=getchar();
	return r*=b;
}

const int N=105,M=5e4+10,MOD=998244353;
#define con(type) const type &

inline int add(con(int)a,con(int)b){return a+b>=MOD?a+b-MOD:a+b;}
inline int sub(con(int)a,con(int)b){return a-b<0?a-b+MOD:a-b;}
inline int mul(con(int)a,con(int)b){return int(1ll*a*b%MOD);}
inline int ina_pow(con(int)a,con(int)b){return b?mul(ina_pow(mul(a,a),b>>1),(b&1)?a:1):1;}
inline int inv(con(int)k){return ina_pow(k,MOD-2);}

int n,sump;
int ars[N],arp[N],f[2][M<<1],g[2][M<<1];

int main(){
	const int INV2=ina_pow(2,MOD-2);
	rin(n);
	for(int i=1;i<=n;i++) rin(ars[i]);
	for(int i=1;i<=n;i++) sump+=rin(arp[i]);
	f[0][sump]=g[0][sump]=1;
	for(int i=1,now=0;i<=n;now+=arp[i],i++){
		int I=i&1,J=!I;
		fill(f[I]+(sump-now-arp[i]),f[I]+(sump+now+arp[i])+1,0);
		fill(g[I]+(sump-now-arp[i]),g[I]+(sump+now+arp[i])+1,0);
		int *fI=f[I]+sump,*fJ=f[J]+sump,*gI=g[I]+sump,*gJ=g[J]+sump;
		for(int j=-now;j<=now;j++){
			int val=mul(INV2,fJ[j]);
			fI[j+arp[i]]=add(fI[j+arp[i]],val);
			if(ars[i]) fI[j-arp[i]]=sub(fI[j-arp[i]],val);
			else fI[j-arp[i]]=add(fI[j-arp[i]],val);
			val=mul(INV2,gJ[j]);
			gI[j+arp[i]]=add(gI[j+arp[i]],val);
			gI[j-arp[i]]=add(gI[j-arp[i]],val);
		}
	}
	int *fi=f[n&1]+sump,*gi=g[n&1]+sump;
	int invs=ina_pow(sump,MOD-2),f0=fi[sump],g0=gi[sump],f1=0,g1=0;
	for(int i=-sump;i<=sump;i++){
		int ii=i<0? i+MOD:i,tmp=inv(sub(mul(ii,invs),1));
		f1=add(f1,mul(fi[i],tmp));
		g1=add(g1,mul(gi[i],tmp));
	}
	printf("%d\n",mul(sub(mul(f1,g0),mul(g1,f0)),inv(mul(g0,g0))));
	return 0;
}

---

THE END

Thanks for reading!

如若是狂风雪雨
亦或碧空如洗
此生百转千回阴晴只得由你
挥毫执墨笔 把你记作曲
也许 此中音律抚琴为君高歌一曲

——《琴弦上(vocaloid)》By 乐正绫/赤羽/星葵

> Link 琴弦上-Bilibili