P9972 [THUPC 2024 初赛] 勇闯末日塔题解

WaReTle · 2023-12-21 16:14:47 · 题解

本题解做法部分参考自官方题解，并对一些细节做了补充。

前置知识

P4001 [ICPC-Beijing 2006] 狼抓兔子（平面图最小割转对偶图最短路）。

P3249 [HNOI2016] 矿区（平面图转对偶图）。

三维向量的基本运算。

简要思路

考虑最大流转最小割。一个割对应着一个球面上的将 s 和 t 分开的环。建出题目中的平面图的对偶图，找到一个符合上述要求的最小环。将图分层，在与同一个点相邻的面之间连跨层的 0 权边可以解决删点的问题。

计算几何部分

基础操作

通过点乘可以算出两个向量夹角的余弦值。

通过叉乘可以求出两个向量确定的平面的法向量，并计算出两个向量夹角的正弦值的绝对值。

求两点的距离

设这两点到球心的向量的夹角为 \theta，则题目中的边的长度为 R\theta。

极角排序

平面图转对偶图需要对每个点出边进行极角排序。区别于平面上的极角排序，我们没有一个明确的向量可以作为基准（因为我们需要排序的向量不共面）。我使用了以下方法：

设当前正在排序出边的点为 u，它的相邻点为 v，球心到这些点的向量为 \vec{u},\vec{v}。我们求出 \vec{u}\times\vec{v}，这样就得到了若干平行于球在 u 点处的切面并与对应的 \vec{v} 垂直的向量.对这些向量进行排序即可。在这些向量中取一个向量作为基准向量，并计算其他向量与这个向量的夹角的正弦值与余弦值，最后通过 atan2 函数计算夹角的具体值。

细节：我们求正弦值时叉乘得到的是一个向量，难以判断正负性。观察到上面的向量平行于球在 u 点处的切面，因此它们叉乘得到的向量平行于 \vec{u}。因此我们可以认为叉乘得到的向量中与 \vec{u} 同向的向量为正，与 \vec{u} 反向的向量为负。

平面图转对偶图

极角排序后直接套用平面上的做法即可。

图论部分

下文中“面”指由原图的边围成的区域。

建图

将所有的面和所有的点作为点加入新图。如果不考虑删点，新图的边即为相邻面之间跨过原图边的连边，边权为跨过的原图边的边权。为了解决删点的问题，我们将新图分 L+1 层（编号 0\sim L），一层的面连向下一层相邻的点，每层的点连向本层相邻的面，边权均为 0。

最小环

我们要求求出的最小环将 s 和 t 分开。这可以通过任取一条 s 到 t 的路径设为关键路径，并钦定环必须跨过这条路径奇数次来解决。反映在图上，我们可以直接把新图上的点拆成奇偶点，所求环即为一个点在 0 层的偶点到这个点在任意层的奇点的最短路。

考虑奇偶点之后，跨过关键路径上的边的情况是容易处理的。对于从点跨过关键路径的情况，可以让这个点和关键路径一侧的面之间的边连接奇偶性相同的点，这个点和关键路径另一侧的面之间的边连接奇偶性不同的点。详见代码。

常数优化部分

最终答案的环必然经过与关键路径相邻的面。因此只需要以这些面为起点跑最短路。为了减小常数，最好使用 bfs 而不是 dfs 确定关键路径。

正常写法使用 double 即可，不需要 long double。

容易写错的地方

~~剧透警告~~

想清楚平面图转对偶图绕圈的过程。如果你是从一条边走到逆时针方向的下一条边，那么最后回到起点的时候走的应该是出发的边的上一条边。第一个样例输出答案的一半可能是这个原因。
如果你决定只以与关键路径相邻的面为起点跑最短路，记得跑仅和关键路径上的点相邻（不和边相邻）的面。否则如果你使用 bfs 确定关键路径可能会 wa 三个点（如果你用 dfs 确定关键路径，可能因为关键路径比较长误打误撞找到了正确答案）。
新图的点数大概有 5\times 10^4\sim 6\times 10^4 级别，数组不要开小。

代码实现

#include<bits/stdc++.h>
#define fi first
#define se second
using namespace std;
typedef long long ll;
typedef unsigned long long ull;
typedef pair<int,int> pii;
typedef pair<int,double> pid;
const int N=1005,M=3005;
const double pi=acos(-1);
int n,m,L,s,t,tot,vid[N][10][2],fid[N][10][2],mpe[N][N],fcnt;
int face[N][N];
double dis[50005],ef[N];
bool imp[N];
vector<int>fne[M],fnv[N];
bool vis[N][N],flag;
double R,K,ang[N][N];
struct vec{double x,y,z;}pt[N],base;
double operator*(vec a,vec b){return a.x*b.x+a.y*b.y+a.z*b.z;}
vec operator%(vec a,vec b){return {a.y*b.z-a.z*b.y,a.z*b.x-a.x*b.z,a.x*b.y-a.y*b.x};}
double operator!(vec a){return sqrt(a.x*a.x+a.y*a.y+a.z*a.z);}
struct edge{int u,v,ty;double w;}e[M];
double sgn(double x){return x>0?1:x<0?-1:0;}
double sin(vec a,vec b){return sgn(a%b*base)*!(a%b)/!a/!b;}
double cos(vec a,vec b){return a*b/!a/!b;}
vector<pii>og[N];
vector<pid>g[50005];
vector<int>pth;
bool inst[N];
double ans=1e9;
int q[N],hd,tl,pre[N];
void bfs(int u)
{
    q[hd=tl=1]=u;
    memset(pre,-1,sizeof(pre));
    pre[u]=0;
    while(hd<=tl)
    {
        u=q[hd++];
        for(auto i:og[u])
            if(!~pre[i.fi])
                pre[i.fi]=i.se,q[++tl]=i.fi;
    }
    u=t;
    while(pre[u])
        pth.push_back(pre[u]),u=e[pre[u]].u+e[pre[u]].v-u;
}
priority_queue<pair<double,int>,vector<pair<double,int>>,greater<pair<double,int>>>pq;
void ae(int u,int v,double w){g[u].emplace_back(v,w);}
void dij(int st)
{
    for(int i=1;i<=tot;++i)dis[i]=1e9;
    dis[st]=0,pq.emplace(dis[st],st);
    while(pq.size())
    {
        double d;int u;tie(d,u)=pq.top();pq.pop();
        if(dis[u]!=d)continue;
        if(dis[u]>ans)break;
        for(auto i:g[u])
        {
            int v;double w;tie(v,w)=i;
            if(dis[u]+w<dis[v])
                dis[v]=dis[u]+w,pq.emplace(dis[v],v);
        }
    }
    while(pq.size())pq.pop();
}
int main()
{
    scanf("%d%d%d%d%d%lf%lf",&n,&m,&L,&s,&t,&R,&K);
    for(int i=1;i<=n;++i)
    {
        double a,b;scanf("%lf%lf%lf",&a,&b,ef+i);
        a*=pi,b*=pi;
        pt[i]={R*sin(a)*cos(b),R*sin(a)*sin(b),R*cos(a)};
    }
    for(int i=1;i<=m;++i)
    {
        scanf("%d%d",&e[i].u,&e[i].v);
        double dist=R*acos(cos(pt[e[i].u],pt[e[i].v]));
        e[i].w=K*ef[e[i].u]*ef[e[i].v]/dist/dist;
        og[e[i].u].emplace_back(e[i].v,i);
        og[e[i].v].emplace_back(e[i].u,i);
        mpe[e[i].u][e[i].v]=mpe[e[i].v][e[i].u]=i;
    }
    for(int i=1;i<=n;++i)if(og[i].size())
    {
        static vec norm[N];
        for(auto j:og[i])
        {
            int v=j.fi;
            norm[v]=pt[i]%pt[v];
        }
        ang[i][og[i][0].fi]=0;
        base=pt[i];
        for(int j=1;j<og[i].size();++j)
            ang[i][og[i][j].fi]=
            atan2(sin(norm[og[i][0].fi],norm[og[i][j].fi]),
            cos(norm[og[i][0].fi],norm[og[i][j].fi]));
        sort(og[i].begin(),og[i].end(),
        [&](pii a,pii b){return ang[i][a.fi]<ang[i][b.fi];});
    }
    for(int i=1;i<=n;++i)if(og[i].size()>1)
        for(int j=0;j<og[i].size();++j)if(!vis[i][j]&&og[og[i][j].fi].size()>1)
        {
            vector<int>seqv,seqe;
            int u=i,pos=j;
            ++fcnt;
            while(1)
            {
                seqv.push_back(u);
                seqe.push_back(og[u][pos].se);
                vis[u][pos]=1,face[u][pos]=fcnt;
                int v=og[u][pos].fi;
                if(v==i)break;
                int poss=lower_bound(og[v].begin(),og[v].end(),make_pair(u,0),
                [&](pii a,pii b){return ang[v][a.fi]<ang[v][b.fi];})-og[v].begin();
                u=v,pos=(poss+1)%og[v].size();
            }
            for(int x:seqv)fnv[x].push_back(fcnt);
            for(int x:seqe)fne[x].push_back(fcnt);
        }
    bfs(s);
    for(int i:pth)e[i].ty=1,imp[e[i].u]=imp[e[i].v]=1;
    for(int i=1;i<=fcnt;++i)
        for(int j=0;j<=L;++j)
            for(int k=0;k<2;++k)
                fid[i][j][k]=++tot;
    for(int i=1;i<=n;++i)
        for(int j=0;j<=L;++j)
            for(int k=0;k<2;++k)
                vid[i][j][k]=++tot;
    for(int i=1;i<=n;++i)if(i!=s&&i!=t)
    {
        int rev=0;
        for(int j=0;j<og[i].size();++j)if(face[i][j])
        {
            for(int p=0;p<=L;++p)
                for(int q=0;q<2;++q)
                {
                    if(p<L)ae(fid[face[i][j]][p][q],vid[i][p+1][q^rev],0);
                    ae(vid[i][p][q],fid[face[i][j]][p][q^rev],0);
                }
            rev^=e[og[i][j].se].ty;
        }
    }
    for(int i=1;i<=m;++i)if(fne[i].size())
        for(int j=0;j<=L;++j)
            for(int k=0;k<2;++k)
                ae(fid[fne[i][0]][j][k],fid[fne[i][1]][j][k^e[i].ty],e[i].w),
                ae(fid[fne[i][1]][j][k],fid[fne[i][0]][j][k^e[i].ty],e[i].w);
    set<int>st;
    for(int i=1;i<=n;++i)if(imp[i])
        for(int j=0;j<fnv[i].size();++j)
        {
            int fr=fid[fnv[i][j]][0][0];
            if(st.count(fr))continue;
            st.insert(fr);
            dij(fr);
            for(int k=0;k<=L;++k)
                ans=min(ans,dis[fid[fnv[i][j]][k][1]]);
        }
    printf("%.12lf\n",ans);
    return 0;
}

P9972 [THUPC 2024 初赛] 勇闯末日塔 题解