首先这是一道dp题,对题意的把握和对状态的处理是解题关键。
题目给出的范围是n在1到1e11之间,由于在裂变过称中左儿子总是父亲节点的一个非平凡约数,容易看出裂变过程只与
素数幂有关,并且显然有素数不超过11个,幂指数不超过40,实际上可以用一个大小为11的数组来等价地表示状态,状态
与其内元素顺序无关,因此可以排序,压缩后的状态不超过3000个(准确地说是2957个,通过一个简单的dfs即可统计出此结果)。
以上解决了题目的规模问题。
这道题目我开始因为理解错题意wa了几次,不能通过统计儿子节点的期望高度的平均值再加1得到以父节点为根的数的期望高度,因为两颗子树
的高度在概率上对根树的影响不是相互独立的。
可以设dp(i, j)为以状态i为根深度为j的概率:
那么显然有dp(i, j) = sigma(dp(k, j - 1) * sigma(i / k, t) + sigma(k, t) * dp(i / k, j - 1) + dp(k, j - 1) * dp(i / k, j -1)) / (N - 2)
其中k | i, 且 k ≠ 1,k≠ i, t < j - 1。
以状态i为根的树的期望高度:exp(i) = sigma(j * dp(i, j)), 0 < j < 40.
因此状态转移可以通过枚举左儿子得到,整体复杂度O(3000 * 3000 * 40)。
#include <cstdio>
#include <cstring>
#include <algorithm>
#include <map>
#include <string>
#include <vector>
#include <set>
#include <cmath>
#include <ctime>
#include <cassert>
#pragma comment(linker, "/STACK:102400000,102400000")
#define lson (u << 1)
#define rson (u << 1 | 1)
#define cls(i, j) memset(i, j, sizeof i)
using namespace std;
typedef __int64 ll;
const double eps = 1e-;
const double pi = acos(-1.0);
const int maxn = 1e5 + ;
const int maxm = ;
const int inf = 0x3f3f3f3f;
const ll linf = 0x3fffffffffffffff;
const ll mod = 1e9 + ; int debug = ;
map<ll, int> mapi;
int buf[], k1;
ll S = (ll)1e11;
double dp[][]; int a[] = {, , , , , , , , , , };
bool cmp(int a, int b) { return a > b; }
int prime[maxn], k;
bool vis[]; ll Hash(int *t){
ll tem = , ans = ;
for(int i = ; i < ; i++){
ans += tem * t[i];
tem *= ;
}
return ans;
} int stack[], k2;
ll states[];
int ks;
int table[][]; void dfs(ll num, int limit, int next){
if(next){
memcpy(buf, stack, sizeof buf);
sort(buf, buf + , cmp);
ll hash_value = Hash(buf);
states[ks++] = hash_value;
memcpy(table[ks - ], buf, sizeof buf);
mapi[hash_value] = ks - ;
}
if(next > ) return;
for(int i = ; i <= limit; i++){
num *= a[next];
if(num > S) break;
stack[k2++] = i;
dfs(num, i, next + );
stack[--k2] = ;
}
} void shaffix(){
bool vis[ + ];
int mid = (int)5e5 + ;
cls(vis, );
for(int i = ; i < mid; i++){
if(vis[i]) continue;
prime[k++] = i;
for(ll j = (ll)i * ; j < mid; j += i) vis[j] = ;
}
} void cal(int id); void dfs1(int id, int next){
if(next > ){
int buf1[], buf2[];
memcpy(buf1, stack, sizeof stack);
for(int i = ; i < ; i++) buf2[i] = table[id][i] - buf1[i];
sort(buf1, buf1 + , cmp);
sort(buf2, buf2 + , cmp);
if(!buf1[] || !buf2[]) return;
ll hash_value1 = Hash(buf1), hash_value2 = Hash(buf2);
int id1 = mapi[hash_value1], id2 = mapi[hash_value2];
cal(id1), cal(id2);
double prefix_left = , prefix_right = ;
for(int i = ; i < ; i++){
dp[id][i + ] += dp[id1][i] * prefix_right +
prefix_left * dp[id2][i] + dp[id1][i] * dp[id2][i];
prefix_left += dp[id1][i], prefix_right += dp[id2][i];
}
return;
}
for(int i = ; i <= table[id][next]; i++){
stack[k2++] = i;
dfs1(id, next + );
stack[--k2] = ;
}
} void cal(int id){
if(vis[id]) return;
k2 = ;
dfs1(id, );
int sum = ;
for(int i = ; i < ; i++) sum *= + table[id][i];
sum -= ;
for(int i = ; i < ; i++) dp[id][i] /= sum;
vis[id] = ;
//printf("%d+\n", id);
} void init(){
shaffix();
cls(vis, );
vis[] = ;
mapi.clear();
//0...10 < 10^11
cls(stack, );
ks = k2 = ;
dfs(, , );
cls(dp, );
dp[][] = ;
for(int i = ; i < ks; i++) cal(i);
} double solve(ll num){
cls(buf, );
k1 = ;
int mid = (int)sqrt((double)num);
for(int i = ; i < k && prime[i] <= mid; i++){
if(num % prime[i]) continue;
while(num % prime[i] == ) ++buf[k1], num /= prime[i];
mid = (int)sqrt((double)num);
k1++;
}
if(num != ) buf[k1++] = ;
sort(buf, buf + , cmp);
ll hash_value = Hash(buf);
double ans = ;
int id = mapi[hash_value];
for(int i = ; i < ; i++) ans += dp[id][i] * i;
return ans;
} int main(){
//freopen("in.txt", "r", stdin);
//freopen("out.txt", "w", stdout);
int T, kase = ;
init();
scanf("%d", &T);
ll n;
while(T--){
scanf("%I64d", &n);
double ans = solve(n);
printf("Case #%d: %.6f\n", ++kase, ans);
}
return ;
}