【BZOJ 2510】 2510: 弱题 （矩阵乘法、循环矩阵的矩阵乘法）

2510: 弱题

Time Limit: 10 Sec  Memory Limit: 128 MB
Submit: 374  Solved: 196

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。

Source

2011福建集训

【分析】

　　综合题表里面目测的一道【稀少的】可做题。。

　　显然，矩阵乘法嘛。。$f[i][j]=f[i-1][j]*dfrac{m-1}{m}+f[i-1][j-1]*dfrac{1}{m}$

　　然后一个厉害的东西就是A,B是循环矩阵（就是矩阵的下一行是上一行循环右移一位得到的），那么A*B=C的C也是循环矩阵。

　　只需求第一行就可以求出整个C。

　　所以是$n^2 log$

#include<cstdio>
#include<cstdlib>
#include<cstring>
#include<iostream>
#include<algorithm>
using namespace std;
#define Maxn 1010

int a[Maxn];
int n,m,k;

struct node
{
    double w[Maxn][Maxn];
}t[5];

void mul(int x,int y,int z)
{
    for(int i=1;i<=n;i++)
     for(int j=1;j<=n;j++) t[2].w[i][j]=0;
    for(int k=1;k<=n;k++)
     for(int j=1;j<=n;j++)
      t[2].w[1][j]+=t[y].w[1][k]*t[z].w[k][j];
    for(int j=1;j<=n;j++) t[x].w[1][j]=t[2].w[1][j];
    for(int i=2;i<=n;i++)
    {
      for(int j=1;j<=n;j++)
        t[x].w[i][j]=t[x].w[i-1][j-1==0?n:j-1];
    }
}

void qpow(int b)
{
    for(int i=1;i<=n;i++) for(int j=1;j<=n;j++) t[1].w[i][j]=0;
    for(int i=1;i<=n;i++) t[1].w[i][i]=1.0;
    while(b)
    {
        if(b&1) mul(1,0,1);
        mul(0,0,0);
        b>>=1;
    }
}

int main()
{
    scanf("%d%d%d",&n,&m,&k);
    for(int i=1;i<=n;i++) scanf("%d",&a[i]);
    for(int i=1;i<=n;i++)
    {
        for(int j=1;j<=n;j++) t[0].w[i][j]=0;
        t[0].w[i][i]=1.0-1.0/m;
        t[0].w[i][i-1==0?n:i-1]=1.0/m;
    }
    qpow(k);
    for(int i=1;i<=n;i++)
    {
        double ans=0;
        for(int j=1;j<=n;j++)
        {
            ans+=a[j]*t[1].w[i][j];
        }
        printf("%.3lf
",ans);
    }
    return 0;
}

2017-04-10 20:47:14