Uva11762 Race to 1 数学期望

/*
    题目描述：给定一个整数n，每次随机选择一个小于n的素数y，如果n能整除y就耗费一单位时间将n变为n/y，否则耗费
            一单位时间，n不发生变化，问平均花费多少时间能把n变为1

    思路：设dp[x]表示将x变为1的时间数学期望，p(x)表示小于等于x的素数个数，则有
                             dp[x] = Σ 1/p(x) * (1 + dp[x / y])  + Σ 1/p(x) * (1 + dp[x])
                       =>    dp[x] = (Σdp[x/y] + p(x))/(p(x) - num)        其中num为小于x的素数中不能整除x的素数个数
         可以先线性筛素数，得到所有1e6内的素数，对于一个n，二分得到p（n），然后通过记忆化搜索的方法计算
*/
#pragma warning(disable:4786)
#pragma comment(linker, "/STACK:102400000,102400000")
#include<iostream>
#include<cstdio>
#include<cstring>
#include<algorithm>
#include<stack>
#include<queue>
#include<map>
#include<set>
#include<vector>
#include<cmath>
#include<string>
#include<sstream>
#define LL long long
#define FOR(i,f_start,f_end) for(int i=f_start;i<=f_end;++i)
#define mem(a,x) memset(a,x,sizeof(a))
#define lson l,m,x<<1
#define rson m+1,r,x<<1|1
using namespace std;
const int INF = 0x3f3f3f3f;
const int mod = 1e9 + 7;
const double PI = acos(-1.0);
const double eps=1e-6;
const int maxn = 1e6 + 5;
const int maxv = 1e5 ;
double dp[maxn];
int isprime[maxn],prime[maxv];
int pnum ;
void find_prime(int n)
{
    pnum = 0;
    mem(prime, 0);
    mem(isprime, 1);
    for(int i = 2 ; i < n ; i++){
        if(isprime[i])      prime[pnum++] = i;
        for(int j = 0; j<pnum && i * prime[j] < n ; j++){
            isprime[i * prime[j] ] = 0;
            if(i % prime[j] == 0)       break;
        }
    }
}
int bsearch(int left , int right , int val)
{
    int mid , ret = -1 ;
    while(left <= right){
        mid = left + (right - left ) / 2;
        if( prime[mid] <= val){
            ret = mid ;
            left = mid + 1;
        }
        else{
            right = mid - 1;
        }
    }
    return ret;
}
double dfs(int x)
{
    if(x == 1)
        return dp[1] = 0 ;
    if( fabs(dp[x]) > eps)
        return dp[x];
    int n = bsearch(0 , pnum - 1, x ) , cnt = 0 ;       //cnt为小于x的素数中x不能整除的个数
    double sum = 0 ;                                            //sum为y是小于x的素数且x%y==0时dfs(x / y)的和
    for(int i = 0 ; i <= n ; i++){
        if(x % prime[i])        ++cnt;
        else        sum += dfs( x / prime[i] );
    }
    double p = sum + n + 1 , q = n + 1 - cnt ;
    double ret = 1.0 * p / q ;
    return dp[x] = ret;
}
int main()
{
    find_prime(1e6 + 2);
    mem(dp , 0);
    int T , kase = 0 , n;
    scanf("%d",&T);
    while(T--){
        scanf("%d",&n);
        double ans = dfs(n);
        printf("Case %d: %.12lf\n",++kase , ans);
    }
    return 0;
}