Loj #3045. 「ZJOI2019」開關

題目描述

九條可憐是一個貪玩的女孩子。

這天，她和她的好朋友法海哥哥去玩密室逃脫。在他們面前的是 \(n\) 個開關，開始每個開關都是關閉的狀態。要通過這關，必須要讓開關達到指定的狀態。目標狀態由一個長度為 \(n\) 的 \(01\) 數組 \(s\) 給出，\(s_i = 0\) 表示第 \(i\) 個開關在最后需要是關着的，\(s_i = 1\) 表示第 \(i\) 個開關在最后需要被打開。

然而作為闖關者，可憐和法海並不知道 \(s\)。因此他們決定采用一個比較穩妥的方法：瞎按。他們根據開關的外形、位置，通過一些玄學的方法給每一個開關賦予了一個值 \(p_i(p_i > 0)\)。每一次，他們會以正比於 \(p_i\) 的概率（第 \(i\) 個開關被選中的概率是 \(\frac{p_i}{\sum_{j=1}^n p_j}\)）選擇並按下一個開關。開關被按下后，狀態會被反轉，即開變關，關變開。注意,每一輪的選擇都是完全獨立的。

在按開關的過程中，一旦當前開關的狀態達到了 \(s\)，那么可憐和法海面前的門就會打開，他們會馬上停止按開關的過程並前往下一關。作為一名歐皇，可憐在按了\(\sum_{i=1}^n s_i\) 次后，就打開了大門。為了感受一下自己的運氣是多么的好，可憐想要讓你幫她計算一下，用這種隨機按的方式，期望需要按多少次開關才能通過這一關。

輸入格式

第一行輸入一個整數 \(n\)，表示開關的數量。

第二行輸入 \(n\) 個整數 \(s_i(s_i \in \{0, 1\})\)，表示開關的目標狀態。第三行同樣輸入 \(n\) 個整數 \(p_i\)，表示每個開關的權值。

輸出格式

輸出一行一個整數，表示答案對 \(998244353\) 取模后的值。即如果答案的最簡分數表示為 \(\frac{x}{y}\ (x \ge 0, y \ge 1, \gcd(x, y) = 1)\)，你需要輸出 \(x \times y^{−1} \bmod 998244353\)。

數據范圍與提示

對於 \(100\%\) 的測試數據，保證 \(n\ge 1,\sum_{i=1}^n p_i\le 5\times 10^4,p_i\ge 1\)。

\(\\\)

完全不會啊，還是\(\text{Orz}\)別人的題解好了。https://www.cnblogs.com/zhoushuyu/p/10687696.html

首先設\(P=\sum_{i=1}^np_i\)

為了方便，下面的\(p_i=\frac{p_i}{P}\)。

先求出走了\(n\)步到達結束狀態的概率的\(EGF\)：

考慮當\(s_i\)為\(1\)的時候，結束狀態時應該按了\(i\)奇數次，否則按了偶數次。

我們也知道：

\[\frac{e^x-e^{-x}}{2}=\sum_{i\ is\ odd}\frac{x^i}{i!}\\ \frac{e^x+e^{-x}}{2}=\sum_{i\ is\ even}\frac{x^i}{i!}\\ \]

所以概率的\(EGF\)為：

\[F(x)=\prod_{i=1}^n\frac{e^{p_ix}+(-1)^{s_i}e^{-p_ix}}{2} \]

但是有可能在走完\(n\)步之前就經過了終止狀態。全過程可以認為是先走到了終止狀態，再繞了若干圈回到終止狀態。

考慮走了\(n\)步回到出發點的概率的\(EGF\)：

\[G(x)=\prod_{i=1}^n\frac{e^{p_ix}+e^{-p_ix}}{2} \]

設它們的\(OGF\)分別為\(h,f,g\)，則：

\[h(x)*g(x)=f(x)\\ h(x)=\frac{f(x)}{g(x)} \]

下面考慮求\(OGF\)。

我們將\(F\)換一種寫法：

\[F(x)=\sum_{i=-P}^Pa_ie^{\frac{i}{P}x}\\ \]

則：

\[f(x)=\sum_{i=-P}^Pa_i(\sum_{j=0}^{\infty}(\frac{i}{P})^jx^j)\\ =\sum_{i=-P}^Pa_i\frac{1}{1-\frac{i}{P}x} \]

答案為\(h'(1)\)（帶入一下就知道了）：

\[h'(x)=\frac{f'(x)*g(x)-f(x)*g'(x)}{g^2(x)} \]

但是帶入\(x=1\)時\(f(x)\)和\(g(x)\)不收斂，於是考慮\(f(x)\)與\(g(x)\)同時乘上\(\prod_{i}(1-\frac{i}{P}x)\)。此時

\[f(x)=\sum_ia_i\prod_{j\neq i}(1-\frac{j}{P}x) \]

容易得到：

\[f(1)=a_P\prod_{j\neq P}(1-\frac{j}{P}) \]

接下來我們要求\(f'(1)\)。首先我們要利用下面的公式：

\[(\prod_i (1+a_ix))'=\sum_ia_i\prod_{j\neq i}(1+a_jx) \]

證明的話可以將左邊部分暴力展開，求導，就會發現等於右邊的部分。

對\(f(x)\)求導的話可以分為兩部分：

1.首先是\(i\neq P\)的部分：

\[A(x)=\sum_{i\neq P}a_i\sum_{k\neq i}(-\frac{k}{P})\prod_{j\neq i,j\neq k}(1-\frac{j}{P}x) \]

因為我最后要的是\(f'(1)\)啊，所以將\(x=1\)帶入就會發現，當\(k\neq P\)的時候后面的\(\prod_{j\neq i,j\neq k}(1-\frac{j}{P}x)=0\)。所以：

\[A(1)=-\sum_{i\neq P}a_i\prod_{j\neq i,j\neq P}(1-\frac{j}{P}) \]

2.然后是\(i=P\)的部分：

\[B(1)=a_p\sum_{i\neq P}(-\frac{i}{P}) \prod_{j\neq P,j\neq i}(1-\frac{j}{P}) \]

\(\\\)

於是

\[f'(1)=A(1)+B(1)\\ =-\prod_{i\neq P}(1-\frac{i}{P})(\sum_{j\neq P}\frac{a_j}{1-\frac{j}{P}}+a_P\sum_{j\neq P}\frac{\frac{j}{P}}{1-\frac{j}{P}}) \]

帶入\(h'(1)=\frac{f'(1)g(1)-f(1)g'(1)}{g^2(1)}\)就可以了。
設\(f(x)=\sum a_i \prod_{j\neq i}(1-\frac{j}{P}x)，g(x)=\sum b_i \prod_{j\neq i}(1-\frac{j}{P}x)\):

又因為\(a_P=b_P=\frac{1}{2^n}\)，所以

\[Ans=2^n\sum_{i\neq P}\frac{b_i-a_i}{1-\frac{i}{P}} \]

下面考慮計算\(a_i,b_i\)。

以\(G(x)\)為例：

\[\prod_i\frac{e^{\frac{p_i}{P}x}+e^{-\frac{p_i}{P}x}}{2}=\sum_{i=-P}^P b_ie^{\frac{i}{P}x} \]

可以發現這就是個背包。

代碼：

#include<bits/stdc++.h>
#define ll long long
#define N 50005

using namespace std;
inline int Get() {int x=0,f=1;char ch=getchar();while(ch<'0'||ch>'9') {if(ch=='-') f=-1;ch=getchar();}while('0'<=ch&&ch<='9') {x=(x<<1)+(x<<3)+ch-'0';ch=getchar();}return x*f;}

const ll mod=998244353;
ll ksm(ll t,ll x) {
	ll ans=1;
	for(;x;x>>=1,t=t*t%mod)
		if(x&1) ans=ans*t%mod;
	return ans;
}

int n;
ll sum,invs;
int s[N],p[N];
ll f[N<<1],g[N<<1];
ll tem[N<<1];
const ll inv2=mod+1>>1;

int main() {
	n=Get();
	for(int i=1;i<=n;i++) s[i]=Get();
	for(int i=1;i<=n;i++) p[i]=Get();
	f[N]=g[N]=1;
	for(int i=1;i<=n;i++) {
		sum+=p[i];
		ll flag=s[i]?mod-inv2:inv2;
		memset(tem,0,sizeof(tem));
		for(int j=sum;j>=p[i]-sum;j--) tem[N+j]=f[N+j-p[i]]*inv2%mod;
		for(int j=-sum;j<=sum-p[i];j++) (tem[N+j]+=f[N+j+p[i]]*flag)%=mod;
		memcpy(f,tem,sizeof(f));
		memset(tem,0,sizeof(tem));
		for(int j=sum;j>=p[i]-sum;j--) tem[N+j]=g[N+j-p[i]]*inv2%mod;
		for(int j=-sum;j<=sum-p[i];j++) (tem[N+j]+=g[N+j+p[i]]*inv2)%=mod;
		memcpy(g,tem,sizeof(g));
	}
	ll ans=0;
	for(int i=-sum;i<sum;i++) (ans+=sum*ksm(sum-i,mod-2)%mod*(mod-f[N+i]+g[N+i]))%=mod;
	cout<<ans*ksm(2,n)%mod;
	return 0;
}