[ACM] poj 2096 Collecting Bugs （機率DP，期望）

Collecting Bugs

Time Limit: 10000MS	Memory Limit: 64000K
Total Submissions: 2026	Accepted: 971
Case Time Limit: 2000MS	Special Judge

Description

Ivan is fond of collecting. Unlike other people who collect post stamps, coins or other material stuff, he collects software bugs. When Ivan gets a new program, he classifies all possible bugs into n categories. Each day he discovers exactly one bug in the program and adds information about it and its category into a spreadsheet. When he finds bugs in all bug categories, he calls the program disgusting, publishes this spreadsheet on his home page, and forgets completely about the program.

Two companies, Macrosoft and Microhard are in tight competition. Microhard wants to decrease sales of one Macrosoft program. They hire Ivan to prove that the program in question is disgusting. However, Ivan has a complicated problem. This new program has s subcomponents, and finding bugs of all types in each subcomponent would take too long before the target could be reached. So Ivan and Microhard agreed to use a simpler criteria --- Ivan should find at least one bug in each subsystem and at least one bug of each category.

Macrosoft knows about these plans and it wants to estimate the time that is required for Ivan to call its program disgusting. It's important because the company releases a new version soon, so it can correct its plans and release it quicker. Nobody would be interested in Ivan's opinion about the reliability of the obsolete version.

A bug found in the program can be of any category with equal probability. Similarly, the bug can be found in any given subsystem with equal probability. Any particular bug cannot belong to two different categories or happen simultaneously in two different subsystems. The number of bugs in the program is almost infinite, so the probability of finding a new bug of some category in some subsystem does not reduce after finding any number of bugs of that category in that subsystem.

Find an average time (in days of Ivan's work) required to name the program disgusting.

Input

Input file contains two integer numbers, n and s (0 < n, s <= 1 000).

Output

Output the expectation of the Ivan's working days needed to call the program disgusting, accurate to 4 digits after the decimal point.

Sample Input

1 2      

Sample Output

3.0000      

Source

Northeastern Europe 2004, Northern Subregion

解題思路：

一個軟體有 s 個子系統，存在 n 種 bug。某人一天能找到一個 bug。問，在這個軟體中找齊 n 種 bug，并且每個子系統中至少包含一個 bug 的時間的期望值（機關：天）。注意：bug 是無限多的，每個 bug 屬于任何一種 bug 的機率都是 1/n；出現在每個系統是等可能的，為 1/s。

令 ex[i][j] 表示已經找到了 i 種 bug，且 j 個子系統至少包含一個 bug，距離完成目标需要的時間的期望。

這樣說有點不好了解。換一種方法。有n個bug1有s個bug2，ex[i][j]表示為已經發現了i個bug1，j個bug2,離目标發現s個bug1,n個bug2還需要的期望天數。這樣所求就是ex[0][0]，而ex[n][s]則為0,因為已經全部發現，不需要額外的天數了。到下一天，可能産生四種情況，ex[i][j]，沒有發現任何bug，ex[i+1][j]，發現一個bug1, ex[i][j+1]，發現一個bug2，ex[i+1][j+1]，發現兩個bug，那麼

ex[i][j]= 1+p1*ex[i][j] +p2*ex[i+1][j]+p3*ex[i][j+1]+p4*ex[i+1][j+1]   ，需要加1，因為是下一天

（後來右想了想，為什麼由ex[i][j]可以推出四種情況，然後就可以說ex[i][j]=後面的式子呢，其實這裡ex[i][j] “ 實際上”，這裡說的實際，是在實際工作中，你要發現兩個bug，肯定首先得發現一個bug吧，是以ex[i][j]其實是由後面那四個式子推出來的，因為我們這裡ex[i][j]表示的意義就是距離完成目标還需要的時間期望。)

p1= (i/n)*(j/s)=(i*j)/(n*s)

p2=((n-i)/n)*(j/s)=(n-i)*j/(n*s)

p3=(i/n)*((s-j)/s)=i*(s-j)/(n*s)

p4=((n-i)/n)*((s-j)/s)=(n-i)*(s-j)/(n*s)

這樣整理出來就是

ex[i][j]=1+ (i*j)/(n*s)*ex[i][j]+(n-i)*j/(n*s)*ex[i+1][j]+i*(s-j)/(n*s)*ex[i][j+1]+(n-i)*(s-j)/(n*s)*ex[i+1][j+1]

右邊(i*j)/(n*s)*ex[i][j]移到左邊合并同類項。得到：

ex[i][j]=(  1+(n-i)*j/(n*s)*ex[i+1][j]+i*(s-j)/(n*s)*ex[i][j+1]+(n-i)*(s-j)/(n*s)*ex[i+1][j+1]  )     /   (1-(i*j)/(n*s))

參考資料：

http://kicd.blog.163.com/blog/static/126961911200910168335852/

代碼：

#include <iostream>
#include <string.h>
#include <iomanip>
using namespace std;
const int maxn=1005;
double ex[maxn][maxn];

int main()
{
    memset(ex,0,sizeof(ex));
    int n,s;
    cin>>n>>s;
    double m=n*s*1.0;
    for(int i=n;i>=0;i--)
        for(int j=s;j>=0;j--)
    {
        if(i==n&&j==s)
        continue;
        ex[i][j]=(1+((n-i)*(s-j)/m)*ex[i+1][j+1]+((n-i)*j/m)*ex[i+1][j]+(i*(s-j)/m)*ex[i][j+1])/(1-i*j/m);
    }
    cout<<setiosflags(ios::fixed)<<setprecision(4)<<ex[0][0];
    return 0;
}