【BZOJ2510】弱题 期望DP+循环矩阵乘法

【BZOJ2510】弱题

Description

有M个球，一开始每个球均有一个初始标号，标号范围为1～N且为整数，标号为i的球有a_i个，并保证Σa_i = M。

每次操作等概率取出一个球（即取出每个球的概率均为1/M），若这个球标号为k（k < N），则将它重新标号为k + 1；若这个球标号为N，则将其重标号为1。（取出球后并不将其丢弃）

现在你需要求出，经过K次这样的操作后，每个标号的球的期望个数。

Input

第1行包含三个正整数N，M，K，表示了标号与球的个数以及操作次数。

第2行包含N个非负整数a_i，表示初始标号为i的球有a_i个。

Output

应包含N行，第i行为标号为i的球的期望个数，四舍五入保留3位小数。

Sample Input

2 3 2
3 0

Sample Output

1.667
1.333

HINT

【样例说明】
第1次操作后，由于标号为2球个数为0，所以必然是一个标号为1的球变为标号为2的球。所以有2个标号为1的球，有1个标号为2的球。
第2次操作后，有1/3的概率标号为2的球变为标号为1的球（此时标号为1的球有3个），有2/3的概率标号为1的球变为标号为2的球（此时标号为1的球有1个），所以标号为1的球的期望个数为1/3*3+2/3*1 = 5/3。同理可求出标号为2的球期望个数为4/3。
【数据规模与约定】
对于10%的数据，N ≤ 5, M ≤ 5, K ≤ 10；
对于20%的数据，N ≤ 20, M ≤ 50, K ≤ 20；
对于30%的数据，N ≤ 100, M ≤ 100, K ≤ 100；
对于40%的数据，M ≤ 1000, K ≤ 1000；
对于100%的数据，N ≤ 1000, M ≤ 100,000,000, K ≤ 2,147,483,647。

题解：一开始too naive，以为同样用一个期望DP的黑科技就能过（k=min(k,5000)）~

发现正解又是矩阵乘法，但是矩乘不是n^3的吗？本题有特殊性质。

我们的DP方程长这样：f[i][j]=f[i的上一个][j-1]/m+f[i][j-1]*(m-1)/m

所以我们的转移矩阵的每一行都是循环相同的，将转移矩阵自乘若干次后，每一行仍然是循环相同的，所以我们的矩阵实际上只需要维护一行，那么转移一次的代价自然就是O(n^2)的。

那么具体实现呢？其实非常简单，直接c[i+j]+=a[i]*b[j]。

#include <cstdio>
#include <cstring>
#include <iostream>
using namespace std;
int n,m,k;
struct M
{
	double v[1010];
	M (){memset(v,0,sizeof(v));}
	double& operator [] (int x)	{return v[x];}
	M operator * (M a)
	{
		M ret;
		for(int i=0;i<n;i++)	for(int j=0;j<n;j++)	ret[(i+j)%n]+=v[i]*a[j];
		return ret;
	}
};
M ans,x;
void pm(int y)
{
	while(y)
	{
		if(y&1)	ans=ans*x;
		x=x*x,y>>=1;
	}
}
int main()
{
	scanf("%d%d%d",&n,&m,&k);
	int i;
	for(i=0;i<n;i++)	scanf("%lf",&ans[i]);
	x[0]=(double)(m-1)/m,x[1]=(double)1/m;
	pm(k);
	for(i=0;i<n;i++)	printf("%.3lf
",ans[i]);
	return 0;
}