[UVA][概率DP] 11762 - Race to 1
B |
Race to 1 Input: Standard Input Output: Standard Output |
Dilu have learned a new thing about integers, which is - any positive integer greater than 1 can be divided by at least one prime number less than or equal to that number. So, he is now playing with this property. He selects a number N. And he calls this D.
In each turn he randomly chooses a prime number less than or equal to D. If D is divisible by the prime number then he divides D by the prime number to obtain new D. Otherwise he keeps the old D. He repeats this procedure until D becomes 1. What is the expected number of moves required for N to become 1.
[We say that an integer is said to be prime if its divisible by exactly two different integers. So, 1 is not a prime, by definition. List of first few primes are 2, 3, 5, 7, 11, …]
Input
Input will start with an integer T (T <= 1000), which indicates the number of test cases. Each of the next T lines will contain one integer N (1 <= N <= 1000000).
Output
For each test case output a single line giving the case number followed by the expected number of turn required. Errors up to 1e-6 will be accepted.
Sample Input Output for Sample Input
3 1 3 13
|
Case 1: 0.0000000000 Case 2: 2.0000000000 Case 3: 6.0000000000 |
Problemsetter: Md. Arifuzzaman Arif
Special Thanks: Sohel Hafiz
給一個數字 N,從 [2, N] 隨機挑一個質數 p,如果 N%p = 0,N = N/p。求嘗試次數的期望值。
題目解法:
很久沒寫概率 DP 了,由於有一次可能嘗試一次不變。
dp[n] = p0*dp[n] + 1 + sigma(p1*dp[n/prime_factor]);
移項之後便可以得到遞迴公式。
其中 p0 表示沒選到質因數的機率,p1 表示每個質因數選到的機率。
#include <stdio.h>
#include <string.h>
#include <vector>
#define maxL (1000000>>5)+1
#define GET(x) (mark[x>>5]>>(x&31)&1)
#define SET(x) (mark[x>>5] |= 1<<(x&31))
using namespace std;
int mark[maxL];
vector<int> pfactor[1000005];
int psum[1000005] = {};
double dp[1000005];
char used[1000005] = {};
void sieve() {
register int i, j, k;
SET(1);
int n = 1000000;
for(i = 2; i <= n; i++) {
psum[i] = psum[i-1];
if(!GET(i)) {
psum[i]++;
pfactor[i].push_back(i);
for(j = i+i; j <= n; j += i) {
SET(j);
pfactor[j].push_back(i);
}
}
}
}//d[n]=p0*d[n]+1+p1*d[n/d1]+p1*d[n/d2]...
double dfs(int n) {
if(n == 1) return 0;
if(used[n]) return dp[n];
used[n] = 1;
double p0, p1, e = 0;
p0 = 1-(double)pfactor[n].size()/psum[n];
p1 = 1.0/psum[n];
for(int i = pfactor[n].size()-1; i >= 0; i--)
e += dfs(n/pfactor[n][i])*p1;
e = (e+1)/(1-p0);
dp[n] = e;
return e;
}
int main() {
sieve();
int testcase, n, cases = 0;
scanf("%d", &testcase);
while(testcase--) {
scanf("%d", &n);
printf("Case %d: %.10lf\n", ++cases, dfs(n));
}
return 0;
}