[WerKeyTom_FTD的模拟赛]记忆的轮廓

题目背景

四次死亡轮回后，昴终于到达了贤者之塔，当代贤者夏乌拉一见到昴就上前抱住了昴“师傅！你终于回来了！你有着和师傅一样的魔女的余香，肯定是师傅”。
众所周知，大贤者是嫉妒魔女沙提拉的老公，400年前与神龙、剑圣一起封印魔女因子暴走的莎缇拉。在魔女茶会的时候，莎缇拉也表示过对昴浓浓的爱意，昴便是被莎缇拉召唤来异世界的。
而贤者之塔中的资料与试炼，似乎都指向同一种可能性……记忆的轮廓，逐渐显形……

题目描述

通往贤者之塔的路上，有许多的危机。
我们可以把这个地形看做是一颗树，根节点编号为1，目标节点编号为n，其中1-n的简单路径上，编号依次递增，在[1,n]中，一共有n个节点。
我们把编号在[1,n]的叫做正确节点，[n+1,m]的叫做错误节点。一个叶子，如果是正确节点则为正确叶子，否则称为错误叶子。
莎缇拉要帮助昴到达贤者之塔，因此现在面临着存档位置设定的问题。为了让昴成长为英雄，因此一共只有p次存档的机会，其中1和n必须存档。被莎缇拉设置为要存档的节点称为存档位置。
当然不能让昴陷入死循环，所以存档只能在正确节点上进行，而且同一个节点不能存多次档。因为通往贤者之塔的路上有影响的瘴气，因此莎缇拉假设昴每次位于树上一个节点时，都会等概率选择一个儿子走下去。每当走到一个错误叶子时，再走一步就会读档。
具体的，每次昴到达一个新的存档位置，存档点便会更新为这个位置（假如现在的存档点是i，现在走到了一个存档位置j>i，那么存档点便会更新为j）。读档的意思就是回到当前存档点。
初始昴位于1，当昴走到正确叶子n时，便结束了路程。莎缇拉想知道，最优情况下，昴结束路程的期望步数是多少？

输入格式

第一行一个正整数T表示数据组数。
接下来每组数据，首先读入三个正整数n,m,p。
接下来m-n行，描述树上所有的非正确边（正确边即连接两个正确节点的边），用两个正整数j，k表示j与k之间有一条连边，j和k可以均为错误节点，也可以一个为正确节点另一个为错误节点。数据保证j是k的父亲。

输出格式

T行每行一个实数表示每组数据的答案。请保留四位小数。

样例输入

1
3 7 2
1 4
2 5
3 6
3 7

样例输出

9.000

数据范围及约定

50%，n=p
70%，50<=p<=n<=500
100%，50<=p<=n<=700，m<=1500，T<=5
数据保证每个除了n的正确节点均有至少2个儿子，至多3个儿子。

题目来源

原创

题目背景相关

题目名选自《Re：从零开始的异世界生活》第六章标题“记忆的轮廓”，本章讲述了水门都市战后为了解决魔女教遗留的难题，人工精灵狐引导大家去向贤者之塔中的贤者请教，揭开秘密的故事。目前更新到24节。

50%算法

n=p时显然每个正确节点都是存档位置是最优的。
那么问题在于如何计算期望。
首先设d[i]表示i的儿子数。设g[i]表示对于一个错误节点i，期望走多少步会读档。那么g[i]=1+1/d[i]*sigma{g[j]}其中j是i的儿子。对于每个正确节点i预处理s[i]表示i的错误儿子的g值和，那么s[i]=sigma{g[j]}，j是i的错误儿子。
设f[i]表示正确节点i走到n的期望步数，显然f[n]=0，我们倒着递推。
f[i]=1+1/d[i]*f[i+1]+1/d[i]*sigma{g[j]+f[i]}[j是i的错误儿子]
移项得f[i]=d[i]+f[i+1]+s[i]
复杂度线性。

70%算法

我们设dp，f[i,j]表示当前存档点为i，还剩j次存档机会。
首先我们需要预处理一个a[i,j]，表示存档点为i，从i开始走到正确节点j的期望步数（中间不能存档）。
显然有边界条件a[i,i]=0。对于i<j，可以列出递推式：
a[i,j]=a[i,j-1]+1+1/d[j-1]*0+1/d[j-1]*sigma{g[k]+a[i,j]}[k是j-1的错误儿子]
移项得a[i,j]=a[i,j-1]*d[j-1]+d[j-1]+s[j-1]
可以用n^2的时间预处理a，然后做dp就很好转移了。
枚举下一次的存档点k，那么f[i,j]可以由f[k,j-1]+a[i,k]转移而来。
复杂度O(n^2p)

100%算法

我们首先可以从小到大枚举存档次数，现在只需要考虑怎么优化这个dp（显然可以看出dp是把一个过程重复做了p次）。
观察两个可转移状态j和l，其中j<l，j和l哪一个更优呢？
我们来观察a[i-1,j]-a[i,j]的值，归纳可得a[i-1,j]-a[i,j]=
这说明了什么呢？无论i为多少，i每次左移，对于j<l，a[i,j]的增量会严格小于a[i,k]的增量。
假设定义函数Sj(i)=f[j]+a[i,j]，那么Sj和Sl这两个函数都是单调函数，而且至多有一个交点。
因此便可以考虑使用单调队列优化，默认从队首开始决策越来越劣，同时保证队列中两两元素的函数交点在单调递减，每次在队尾加入元素，然后取最优决策时检测队首是否比第二个优否则踢出队首，这样决策的选择均摊O(1)。
然后考虑如何求两个函数的交点，可以考虑使用二分。
复杂度O(np log n)

不过，有没有发现，70%和100%好像都没考虑一个问题。观察a数组，可以看到它是恐怖的增长的，我们最终答案会不会爆炸？
我们来估计答案的上界。考虑一种可行方案，每n/p个正确节点就设立一次存档位置，那么答案最大是多少呢？考虑最坏情况，观察a的转移，应该每变换一次存档点，大约需要3^(n/p)s[i]+3(n/p-1)*s[i+1]+3^(n/p-2)*s[i+2]+……
因为最多m个节点，s的上限是1500（实际上也远远达不到），把所有s都视为这个上限，提取公因数，计算一下那个等比数列求和，由于p是有下界的，因此n/p有上界14，发现最后也就是个12位数的样子，那么我们估计出答案最大也不会超过这个，可以放心做了。而至于a会爆炸的问题，double是可以存很多位的，而且太大的a肯定不可能被用上。

那么其实，针对答案不会特别大，a的增长又很恐怖，我们还可以思考对70%的算法优化。那就是设定一个常数step，每次转移最多从距当前step步远的位置转移过来。step取40多基本不会有问题了，因为a的下界已经是2^40了，而答案的上界远远没有达到，经过精确计算还可以再把step调小一点。
复杂度O(40np)

标程

正解一

#include<cstdio>#include<algorithm>#define fo(i,a,b) for(i=a;i<=b;i++)#define fd(i,a,b) for(i=a;i>=b;i--)using namespace std;typedef double db;const int maxn=2000+10,maxm=2000+10;const db inf=100000000000000;int d[maxn];db a[maxn][maxn],f[maxn][2],g[maxm],s[maxn];int dl[maxn],xy[maxn],h[maxm],go[maxm],next[maxm];bool bz[maxm];int i,j,k,l,r,t,n,m,p,tot,top,head,tail,ca;db ans;int read(){    int x=0;    char ch=getchar();    while (ch<'0'||ch>'9') ch=getchar();    while (ch>='0'&&ch<='9'){        x=x*10+ch-'0';        ch=getchar();    }    return x;}void add(int x,int y){    d[x]++;    go[++tot]=y;    next[tot]=h[x];    h[x]=tot;}void dfs(int x){    bz[x]=1;    int t=h[x];    g[x]=1;    while (t){        dfs(go[t]);        g[x]+=(db)1/d[x]*g[go[t]];        t=next[t];    }}int getxy(int j,int k,int x){    int l=0,r=j,mid;    while (l<r){        mid=(l+r+1)/2;        if (f[j][x]+a[mid][j]<=f[k][x]+a[mid][k]) l=mid;else r=mid-1;    }    return l;}void solve(int x){    int i,j,k;    fo(i,1,n) f[i][x]=inf;    dl[head=tail=1]=n;    fd(i,n-1,1){        while (head<tail&&xy[head+1]>=i) head++;        k=dl[head];        f[i][x]=f[k][1-x]+a[i][k];        while (head<tail&&getxy(i,dl[tail],1-x)>=xy[tail]) tail--;        dl[++tail]=i;        xy[tail]=getxy(i,dl[tail-1],1-x);    }}int main(){    freopen("memory.in","r",stdin);freopen("memory.out","w",stdout);    ca=read();    while (ca--){        n=read();m=read();p=read();        fo(i,1,m) h[i]=d[i]=bz[i]=0;        tot=0;        fo(i,1,m-n){            j=read();k=read();            add(j,k);        }        fo(i,1,n-1) d[i]++;        fo(i,n+1,m)            if (!bz[i]) dfs(i);        fo(i,1,n){            t=h[i];            s[i]=0;            while (t){                s[i]+=g[go[t]];                t=next[t];            }        }        fo(i,1,n){            a[i][i]=0;            fo(j,i+1,n)                a[i][j]=(db)d[j-1]*a[i][j-1]+s[j-1]+d[j-1];        }        fo(i,1,n-1) f[i][0]=inf;        f[n][0]=0;        fo(i,2,p) solve(1-i%2);        ans=f[1][1-p%2];        printf("%.4lf\n",ans);    }}

正解二

#include<cstdio>#include<algorithm>#define fo(i,a,b) for(i=a;i<=b;i++)#define fd(i,a,b) for(i=a;i>=b;i--)using namespace std;typedef double db;const int maxn=2000+10,maxm=100000+10;const db inf=1000000000000000;int d[maxn];db a[maxn][maxn],f[maxn][maxn],g[maxm],s[maxn];int dl[maxn],xy[maxn],h[maxm],go[maxm],next[maxm];bool bz[maxm];int i,j,k,l,r,t,n,m,p,ca,tot,top,head,tail;db ans;int read(){    int x=0;    char ch=getchar();    while (ch<'0'||ch>'9') ch=getchar();    while (ch>='0'&&ch<='9'){        x=x*10+ch-'0';        ch=getchar();    }    return x;}void add(int x,int y){    d[x]++;    go[++tot]=y;    next[tot]=h[x];    h[x]=tot;}void dfs(int x){    bz[x]=1;    int t=h[x];    g[x]=1;    while (t){        dfs(go[t]);        g[x]+=(db)1/d[x]*g[go[t]];        t=next[t];    }}int main(){    freopen("memory.in","r",stdin);freopen("memory.out","w",stdout);    ca=read();    while (ca--){    n=read();m=read();p=read();    fo(i,1,m) h[i]=d[i]=bz[i]=0;    tot=0;    fo(i,1,m-n){        j=read();k=read();        add(j,k);    }    fo(i,1,n-1) d[i]++;    fo(i,n+1,m)        if (!bz[i]) dfs(i);    fo(i,1,n){        t=h[i];        s[i]=0;        while (t){            s[i]+=g[go[t]];            t=next[t];        }    }    fo(i,1,n){        a[i][i]=0;        fo(j,i+1,n)            a[i][j]=(db)d[j-1]*a[i][j-1]+s[j-1]+d[j-1];    }    fo(i,1,n)        fo(j,1,p)            f[i][j]=inf;    f[n][1]=0;    fo(j,2,p)        fo(i,1,n)            fo(k,i+1,n){                if (k-i>40) break;                if (f[k][j-1]+a[i][k]<f[i][j])                    f[i][j]=f[k][j-1]+a[i][k];            }    ans=f[1][p];    printf("%.4lf\n",ans);    }}

题目情况

50分：3人
0分：3人
平均分：25

想说的话

原创的压轴题，fanvree好强，可惜考试时没对拍。