ZZH與計數

題目

傳送門 to usOJ

題目描述

衆所周知，題目的難度與出題人的能力成正比。

現在「霸樹中學」的同學要集訓，要天天做題。由于 自願捐助款項太多 資金周轉，有時候會請來金牌出題人「标槍」，有時候請來廢鐵出題人「 O I D \tt OID OID」。

如果請來的是「标槍」，那麼會觸發效果「同學們的信心大幅下降了！」于是原本的信心值 v v v 的二進制（畢竟是資訊競賽）表示下為 1 1 1 的值，都會以 1 2 \frac{1}{2} 21​ 的機率變成 0 0 0 。

反之，如果請來「 O I D \tt OID OID」，會觸發效果「同學們 對出題人的鄙夷 的信心大幅增強了！」于是原來的信心值 v v v 的二進制表示下為 0 0 0 的位置，都會以 1 2 \frac{1}{2} 21​ 的機率變為 1 1 1 。

而請來的是誰，也是一個機率問題。有 p p p 的機率請來「标槍」，也有 1 − p 1-p 1−p 的機率請來「 O I D \tt OID OID」，并且進行了 m m m 天集訓（每天都重新請一個出題人，機率仍然為 p p p 與 1 − p 1-p 1−p 不變）。

告訴你最初每種信心值的數量，你能告訴我最後每種信心值的期望數量嗎？

資料範圍與提示

認為所有的二進制表示隻有 n n n 位，也就是說，任何數字在任意時刻均小于 2 n 2^{n} 2n 。

n ≤ 17 n\le 17 n≤17 但是天數 m ≤ 1 0 9 m\le 10^9 m≤109 。由于輸出要取模 998244353 998244353 998244353 ，是以 p = a b p=\frac{a}{b} p=ba​ 以及每種信心值的個數 v i v_i vi​ 滿足 max ⁡ ( a , b , v i ) < 998244353 \max(a,b,v_i)<998244353 max(a,b,vi​)<998244353 就不足為奇。

時限 3 s 3\text{s} 3s ，空間 1 GB 1\text{GB} 1GB 。

思路

是不是有人認為可以每一位分開考慮呢 😮 然後你就會發現不對，因為兩個二進制位之間有影響。

正确的姿勢是用 d p \tt dp dp 整體的處理一個數字。對于某個數字， f ( t , i , j ) f(t,i,j) f(t,i,j) 表示 t t t 天之後得到另一個數，并且原來為 0 0 0 的，現在有 i i i 個是 1 1 1 ；原來是 1 1 1 的，現在有 j j j 個還是 1 1 1 。 f f f 求的是 機率 。

v = ( 00000 ⋯ 00    11111 ⋯ 11 ) 2 v ′ = ( 111 ⋯ 1 ⏟ i 個 1 000    111 ⋯ 1 ⏟ j 個 1 000 ) 2 \begin{aligned} v&=(00000\cdots00\;11111\cdots11)_2\\ v'&=(\underbrace{111\cdots1}_{i個1}000\;\underbrace{111\cdots1}_{j個1}000)_2 \end{aligned} vv′​=(00000⋯0011111⋯11)2​=(i個1

111⋯1​​000j個1

111⋯1​​000)2​​

為了友善觀察，我在中間打了個空格，實際上沒有。而且 v , v ′ v,v' v,v′ 中的 0 , 1 0,1 0,1 不需要排列的那麼規整，不過為了觀察，我把它們放在一起。

然後轉移就比較簡單了。假設原來有 i 0 i_0 i0​ 個 0 0 0 ，有 j 0 = n − i 0 j_0=n-i_0 j0​=n−i0​ 個 1 1 1 。那麼

f ( t + 1 , i , j ) = ( 1 − p ) ∑ x = 0 i ∑ y = 0 j ( 1 2 ) n − x − y ( i x ) ( j y ) f ( t , x , y ) + p ∑ x = i i 0 ∑ y = j j 0 ( 1 2 ) x + y ( i 0 − i x − i ) ( j 0 − j y − j ) f ( t , x , y ) f(t+1,i,j)=\\ (1-p)\sum_{x=0}^{i}\sum_{y=0}^{j}\left({1\over 2}\right)^{n-x-y}{i\choose x}{j\choose y}f(t,x,y)\\ +p\sum_{x=i}^{i_0}\sum_{y=j}^{j_0}\left({1\over 2}\right)^{x+y}{i_0-i\choose x-i}{j_0-j\choose y-j}f(t,x,y) f(t+1,i,j)=(1−p)x=0∑i​y=0∑j​(21​)n−x−y(xi​)(yj​)f(t,x,y)+px=i∑i0​​y=j∑j0​​(21​)x+y(x−ii0​−i​)(y−jj0​−j​)f(t,x,y)

（易錯點：為什麼沒有一些你們認為的組合數？因為這裡計算的是到達某類數字中 任意一個 的機率，而不是所有這一類數的機率的總和。為什麼有一些沒想到的組合數？因為兩類數字之間可能有多條邊。）

枚舉一個 i 0 i_0 i0​ ，然後做 O ( n 6 log ⁡ m ) \mathcal O(n^6\log m) O(n6logm) 的矩陣加速，複雜度 O ( n 7 log ⁡ m ) \mathcal O(n^7\log m) O(n7logm) 。蛤？感覺不行？但是 i 0 × j 0 ≤ n 4 i_0\times j_0\le\frac{n}{4} i0​×j0​≤4n​ ，于是有了 1 64 \frac{1}{64} 641​ 的超級牛逼的常數！

然而我們隻求得了機率。 O ( 2 n × 2 n ) \mathcal O(2^n\times 2^n) O(2n×2n) 的計算貢獻是不可接受的。是以我們還得 d p \tt dp dp 求。（一個值得思考的問題，為啥這裡 d p \tt dp dp 可以做到更快？因為這裡的貢獻計算是有過程的，是 0 0 0 與 1 1 1 轉化的數量決定。）

用 g ( t , S , i , j ) g(t,S,i,j) g(t,S,i,j) 表示，考慮了前 t t t 個二進制位之後，得到一個過程中的數字 S S S ，但是還有 i i i 個 0 0 0 要變成 1 1 1 ，有 j j j 個 1 1 1 要變成 0 0 0 。（關于這個，可以了解為在 D A G \tt DAG DAG 上推貢獻。）

初始時 g ( 0 , S , i , j 0 − j ) = f ( m , i , j ) × v S g(0,S,i,j_0-j)=f(m,i,j)\times v_S g(0,S,i,j0​−j)=f(m,i,j)×vS​

（如果你仔細思考會明白，我還要指出 j 0 = b i t c o u n t ( S ) j_0=bitcount(S) j0​=bitcount(S) 才算嚴謹。）

複雜度 O ( n 3 × 2 n ) \mathcal O(n^3\times 2^n) O(n3×2n) 。然而同樣也有 1 4 \frac{1}{4} 41​ 的常數（實際上可能更小）。

答案當然是 g ( n , 0 , 0 , 0 ) , … , g ( n , 2 n − 1 , 0 , 0 ) g(n,0,0,0),\dots,g(n,2^n-1,0,0) g(n,0,0,0),…,g(n,2n−1,0,0) ，推到底了。然後你發現 g g g 開不下，要把第一位滾動掉。然後你發現非常容易滾動，因為 S S S 隻會與 S ⊕ 2 t S\oplus 2^t S⊕2t 發生關聯。

代碼

卡常技巧：計算 g g g 的複雜度極高！你的循環是怎樣枚舉的， g g g 的變量就要怎麼開。

是以我開成了 g ( i , j , S ) g(i,j,S) g(i,j,S) 而不是 g ( S , i , j ) g(S,i,j) g(S,i,j) 。

#include <cstdio>
#include <iostream>
#include <vector>
#include <cstring>
using namespace std;
typedef long long int_;
const int __M__ = 20000000;
char inBuf[__M__], *iS = inBuf, *iT = inBuf;
inline char getChar(){
	if(iS == iT){
		iT = fread(inBuf,1,__M__,stdin)
			+ (iS = inBuf);
	}
	return *(iS ++);
}
inline int readint(){
	int a = 0; char c = getChar(), f = 1;
	for(; c<'0'||c>'9'; c=getChar())
		if(c == '-') f = -f;
	for(; '0'<=c&&c<='9'; c=getChar())
		a = (a<<3)+(a<<1)+(c^48);
	return a*f;
}
inline void writeint(int x){
	if(x > 9) writeint(x/10);
	putchar((x-x/10*10)^48);
}

const int Mod = 998244353;
inline int qkpow(int_ b,int q){
	int ans = 1; b %= Mod;
	for(; q; q>>=1,b=b*b%Mod)
		if(q&1) ans = ans*b%Mod;
	return ans;
}

const int MaxN = 17;
int N; // 矩陣的大小
struct Matrix{
	int a[90][90];
	void clear(){
		for(int i=0; i<N; ++i)
		for(int j=0; j<N; ++j)
			a[i][j] = 0;
	}
	Matrix operator * (const Matrix &b) const {
		Matrix c; c.clear();
		for(int i=0; i<N; ++i)
		for(int j=0; j<N; ++j)
		if(a[i][j] != 0)
		for(int k=0; k<N; ++k)
			c.a[i][k] = (c.a[i][k]+1ll
				*a[i][j]*b.a[j][k])%Mod;
		return c; // 複雜度還蠻高的
	}
};
Matrix qkpow(Matrix b,int q){
	Matrix ans = b; -- q;
	for(; q; q>>=1,b=b*b)
		if(q&1) ans = ans*b;
	return ans;
}

const int inv2 = (Mod+1)>>1;
int n, m, p, c[MaxN+1][MaxN+1];
int f[MaxN+1][MaxN+1][MaxN+1];
int powTwo[MaxN+1]; // (1/2)^x
void solveF(){
	powTwo[0] = 1;
	for(int i=1; i<=n; ++i)
		powTwo[i] = 1ll*powTwo[i-1]*inv2%Mod;
	for(int i=0; i<=n; ++i)
	for(int j=c[i][0]=1; j<=i; ++j)
		c[i][j] = c[i-1][j-1]+c[i-1][j];
	Matrix S; // S 為轉移矩陣
	for(int i0=0; i0<=n; ++i0){
		int j0 = n-i0+1; // 開區間
		N = (i0+1)*j0; S.clear();
		for(int i=0; i<=i0; ++i)
		for(int j=0; j<j0; ++j){
			int zxy = i*j0+j; // 卡常
			for(int x=0; x<=i; ++x)
			for(int y=0; y<=j; ++y)
				S.a[x*j0+y][zxy] =
					1ll*powTwo[n-x-y]*
					c[i][x]%Mod*c[j][y]%Mod
					*(Mod+1-p)%Mod;
			for(int x=i; x<=i0; ++x)
			for(int y=j; y<j0; ++y)
				S.a[x*j0+y][zxy] =
					1ll*powTwo[x+y]*
					c[i0-i][x-i]%Mod*
					c[j0-1-j][y-j]%Mod
					*p%Mod;
			S.a[zxy][zxy] += (Mod+1ll-p)
				*powTwo[n-i-j]%Mod; // 不變
			S.a[zxy][zxy] %= Mod;
		}
		if(m) S = qkpow(S,m);
		else{
			S.clear(); // 機關矩陣
			for(int i=0; i<N; ++i)
				S.a[i][i] = 1;
		}
		for(int i=0; i<=i0; ++i)
		for(int j=0; j<j0; ++j)
			f[i0][i][j] = S.a[j0-1][i*j0+j];
	}
}

int g[MaxN+1][MaxN+1][1<<MaxN];
int cnt[1<<MaxN]; // bitcount
int num[1<<MaxN]; // amount of number
void solveG(){
	for(int S=0; S<(1<<n); ++S){
		if(S) cnt[S] = cnt[S-(S&-S)]+1;
		for(int i=0; i<=n-cnt[S]; ++i)
		for(int j=0; j<=cnt[S]; ++j)
			g[i][cnt[S]-j][S] = 1ll*num[S]
				*f[n-cnt[S]][i][j]%Mod;
	}
	for(int t=0; t<n; ++t)
	for(int i=0; i<n-t; ++i)
	for(int j=0; i+j<n-t; ++j){
		int *g_ = g[i][j], *gi = g[i+1][j], *gj = g[i][j+1];
		for(int S=0; S<(1<<n); ++S)
			if(!(S>>t&1))
				g_[S] = (g_[S]+gj[S^(1<<t)])%Mod;
			else // if(S>>t&1)
				g_[S] = (g_[S]+gi[S^(1<<t)])%Mod;
	}
}

int main(){
	n = readint(), m = readint();
	p = readint(); p = 1ll*p*
		qkpow(readint(),Mod-2)%Mod;
	solveF(); // 預處理 F
	for(int i=0; i<(1<<n); ++i)
		num[i] = readint();
	solveG();
	writeint(g[0][0][0]); // 去掉行末空格
	for(int i=1; i<(1<<n); ++i){
		putchar(' ');
		writeint(g[0][0][i]);
	}
	putchar('\n');
	return 0;
}