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

【描述】
有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\le 1000,M\le 100,000,000,K\le 2,147,483,647$

【思路】

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

循环矩阵例子:
$$\left[\begin{array}{ccc}1& 2& 3\\ 3& 1& 2\\ 2& 3& 1\end{array}\right]$$
循环矩阵指的是上面这种矩阵。维护这样的矩阵很简单。我们发现，这样的矩阵自乘若干次以后仍然是循环矩阵，这意味着我们只需要维护一行即可。此时矩阵乘法类似于多项式乘法：$c[(i+j)modn]+=a[i]\ast b[j]$。很好记对吧？所以我们就可以在$O(n^{2}logn)$内解决问题，甚至可以使用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]);
}