并不套路的概率生成函数
# 题面
$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$ 为「操作 $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)$ 才是我们真正要求的答案。根据实际意义可以推得
但是我们现在将 $F_i(x)$ 全部乘起来只能得到 $f_i$ 的 EGF,记为 $R(x)$:
同理,我们也只能快速得到 $g_i$ 的 EGF 的表达式,记为 $W(x)$:
考虑如何把 EGF 转成 OGF。
- P2. EGF to OGF
注意到 $R(x)$ 和 $W(x)$ 的形式是类似的。
不难发现 $R(x)$ 可以表示成如下形式($W(x)$ 也可以,下面只通过 $R(x)$ 举例):
只需要把 $F_i(x)$ 当成一对大小为 $p_i$ 和 $-p_i$ 的物品做一次 $\mathcal{O}(nP)$ 的背包即可求出 $a_i$。
而知道了 $a_i$ 我们就可以较为方便地求出其 OGF $F(x)$:
上式的 $j$ 的取值根据背包可知应为 $[-P,P]$。
- P3. 求导
根据之前的推导,我们已经可以求得 OGF $F(x),G(x)$,而最终目标
对其求导并代入 $x=1$ 可以得到答案:
注意 $F(x),G(x)$ 都不涵盖所有可能情况,所以 $F(1)\neq1,G(1)\neq1$。
但是稍微细心一些,我们发现 $F(x)$ 和 $G(x)$ 在 $1$ 处很有可能不收敛——在 P2 中求出的式子:
若 $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$ 处的取值会变成相反数。
# 源代码
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56
| #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