【bzoj2510】【期望dp】【循环矩阵快速幂】弱题

桜のような恋でした

于 2019-09-13 22:36:40 发布

阅读量93

点赞数

分类专栏：期望dp

本文链接：https://blog.csdn.net/weixin_44111457/article/details/100812493

版权

期望dp 专栏收录该内容

8 篇文章 0 订阅

订阅专栏

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

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

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

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

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

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

【样例输入】
2 3 2
3 0

【样例输出】
1.667
1.333
$N \leq 1000, M \leq 100, 000, 000, K \leq 2, 147, 483, 647$

【思路】

转移十分显然，定义 $f [i] [k]$ 表示k次操作之后编号为i的球的期望个数:
$f [i] [k] = (1 / m) * f [i - 1] [k - 1] + ((m - 1) / m) * f [i] [k - 1]$
于是我们就可以矩阵快速幂转移啦。好像有什么不对…N<=1000。~~此时如果你有相当玄学的卡常技巧，你就可以去写代码了，相信你能把常数卡成1/1000~~ 。我们发现，转移矩阵是循环相同的。

循环矩阵例子:
$\begin{bmatrix} 1 & 2&3\\ 3& 1&2 \\2&3&1 \end{bmatrix}$
循环矩阵指的是上面这种矩阵。维护这样的矩阵很简单。我们发现，这样的矩阵自乘若干次以后仍然是循环矩阵，这意味着我们只需要维护一行即可。此时矩阵乘法类似于多项式乘法： $c [(i + j) m o d n] + = a [i] * b [j]$ 。很好记对吧？所以我们就可以在 $O(n^2logn)$ 内解决问题，甚至可以使用fft优化。(不过由于fft自带大常数，所以效率差不多)。
代码：

#include<bits/stdc++.h>
#define re register
using namespace std;
const int N=1e3+5;
inline int red(){
    int data=0;bool w=0; char ch=0;
    ch=getchar();
    while(ch!='-' && (ch<'0' || ch>'9')) ch=getchar();
    if(ch=='-') w=1,ch=getchar();
    while(ch>='0' && ch<='9') data=(data<<3)+(data<<1)+ch-'0',ch=getchar();
    return w?-data:data;
}
int n,m,c,k;
struct mat{
	double a[N];
	friend inline mat operator*(const mat&a,const mat&b){
		mat c;memset(c.a,0,sizeof(c.a));
		for(int re i=0;i<n;i++)
			for(int re j=0;j<n;j++)
				c.a[(i+j)%n]+=a.a[i]*b.a[j];
		return c;
	}
}ans,A;
void ksm(){while(k){if(k&1)ans=ans*A;A=A*A,k>>=1;}}
int main(){
	n=red();m=red();k=red();
	for(int re i=0;i^n;++i)ans.a[i]=red();
	A.a[0]=1.0*(m-1)/m,A.a[1]=1.0/m;ksm();
	for(int re i=0;i^n;++i)printf("%.3f\n",ans.a[i]);
}