题意:

给出一个n*m(n<=10³,m<=10³)的黑白棋棋盘.

棋盘上有t(t<=10⁵)个位置坏掉了。

除了坏掉的位置之外，每个位置上有且只有1颗棋子，且初始时白面向上。

你可以进行k(k<=10⁹)次操作。

每次操作选择两颗棋子，将以它们为对角线的矩形内的所有棋子状态取反(不一定是主对角线).

问k次操作后，白面向上的棋子个数的期望值是多少。

分析:

由于期望具有线性性，我们可以分别计算每颗棋子的期望最后再求和。

假设某颗棋子在1次操作中被翻转的概率为p。

那么容易推出它k次操作后状态不变的概率为(1+(1-2*p)^k)/2.

现在问题就在于如何快速求出p.

由容斥得(x,y)被翻转的情况有:

|([1,x]×[1,y],[x,y]×[n,m])|+|([x,n]×[1,y],[1,x]×[y,m])|-|([x,x]×[1,y],[x,x]×[y,m])|-|([1,x]×[y,y],[x,n]×[y,y])|+1

用前缀和优化就好了。

#include <stdio.h>
#include <string.h>
#include <algorithm>
using namespace std;

#define g() getchar()
template<class Q>void Scan(Q&x){
	char c;while(c=g(),c<48||c>57);
	for(x=0;c>47&&c<58;c=g())x=10*x+c-48;
}

#define rep(i,a,b) for(int i=a;i<b;++i)

const int MAXN=1005;

int n,m,k,t;

bool a[MAXN][MAXN];

int sum[MAXN][MAXN];

inline int Sum(int x1,int y1,int x2,int y2){
	int res=sum[x2][y2]-sum[x2][y1-1];
	return res-(sum[x1-1][y2]-sum[x1-1][y1-1]);
}

inline void Init(){
	Scan(n),Scan(m),Scan(k),Scan(t);
	rep(i,1,n+1)rep(j,1,m+1)sum[i][j]=a[i][j]=1;
	rep(i,1,t+1){
		int x,y;
		Scan(x),Scan(y);
		sum[x][y]=a[x][y]=0;
	}
	rep(i,1,n+1)rep(j,1,m+1)sum[i][j]+=sum[i][j-1];
	rep(i,1,n+1)rep(j,1,m+1)sum[i][j]+=sum[i-1][j];
}

typedef long double ld;
typedef long long ll;

inline ld qpow(ld x,int k){
	ld ans=1;
	for(;k;x=x*x,k>>=1)
		if(k&1)ans=ans*x;
	return ans;
}

ld ans;

inline void Solve(){
	ll q=1ll*sum[n][m]*sum[n][m];
	rep(i,1,n+1)rep(j,1,m+1){
		if(!a[i][j])continue;
		ll p=1ll*Sum(1,1,i,j)*Sum(i,j,n,m)+1ll*Sum(1,j,i,m)*Sum(i,1,n,j);
		p-=1ll*Sum(i,1,i,j)*Sum(i,j,i,m)+1ll*Sum(1,j,i,j)*Sum(i,j,n,j);
		p=p<<1|1;
		ld pr=(ld)p/q;
		ans+=(qpow(1.-pr*2.,k)+1.)/2.;
	}
	printf("%.20f\n",(double)ans);
}

int main(){
	freopen("swap.in","r",stdin);
	freopen("swap.out","w",stdout);
	Init(); Solve(); return 0;
}
 

详见zhj的15年集训队论文

Portal

#include <math.h>
#include <stdio.h>
#include <string.h>
#include <algorithm>
using namespace std;

#define rep(i,a,b) for(int i=a;i<b;++i)
typedef long double ld; 

const int MAXN=1000; 

int n,k,cur;

ld dp[2][MAXN];

int main(){
	scanf("%d%d",&n,&k);
	for(int i=n;i;--i){
		cur^=1;
		int B=701;
		rep(j,1,B){
			dp[cur][j]=(dp[cur^1][j]+(1.+j)/2)*j;
			dp[cur][j]+=(dp[cur^1][j+1]+j);
			dp[cur][j]/=j+1;
			dp[cur][j]+=dp[cur^1][j]*(k-1);
			dp[cur][j]/=k;
		}
	}
	ld ans=dp[cur][1]*k;
	printf("%.20f\n",(double)ans);
	return 0;
}