poj 3608 Bridge Across Islands(旋转卡壳求凸包最短距离)

来源:互联网 发布:淘宝vr资源怎么搜索 编辑:程序博客网 时间:2024/04/27 13:59
Bridge Across Islands
Time Limit: 1000MS Memory Limit: 65536KTotal Submissions: 9843 Accepted: 2886 Special Judge

Description

Thousands of thousands years ago there was a small kingdom located in the middle of the Pacific Ocean. The territory of the kingdom consists two separated islands. Due to the impact of the ocean current, the shapes of both the islands became convex polygons. The king of the kingdom wanted to establish a bridge to connect the two islands. To minimize the cost, the king asked you, the bishop, to find the minimal distance between the boundaries of the two islands.

Input

The input consists of several test cases.
Each test case begins with two integers N, M. (3 ≤ N, M ≤ 10000)
Each of the next N lines contains a pair of coordinates, which describes the position of a vertex in one convex polygon.
Each of the next M lines contains a pair of coordinates, which describes the position of a vertex in the other convex polygon.
A line with N = M = 0 indicates the end of input.
The coordinates are within the range [-10000, 10000].

Output

For each test case output the minimal distance. An error within 0.001 is acceptable.

解题思路:

用旋转卡壳算法求凸包间最短距离:

两个凸多边形 PQ 之间的最小距离由多边形间的对踵点对确立。 存在凸多边形间的三种多边形间的对踵点对,因此就有三种可能存在的最小距离模式:

  1. “顶点-顶点”的情况
  2. “顶点-边”的情况
  3. “边-边”的情况

在P上找出最下方的点记作C,在Q上找最上方的点记作D,记C在逆时针方向的下一个顶点是A,同理定义B。那么我们得到了两个线段AC和DB,分别逆时针同时逐步旋转两线段,那么最短距离肯定在这两条线段上产生。

具体实现时,同时旋转体现为逐步选取逆时针方向上的下一个顶点作为C或D、A或B;两条线段的最短距离归结于四个顶点到另一线段的最短距离或直线的最短距离(平行),而点到线段的距离又归结于点到两个端点和点到直线(如果角度满足的话)的垂线距离,可见不平行的时候可以少考虑一种情况。

至于如何判断平行,代码中实际在用向量AC和DB的叉积判断两者是否平行。

代码:(套模板的)

#include<cstdio>#include<cmath>#include<algorithm>#define N 10005#define PI 3.14159265358979#define EPS 1e-7using namespace std;struct PO{    double x,y;}p1[N],p2[N],o;int n,m,top1,top2,stk1[N],stk2[N];inline int doublecmp(double x)//误差修正 {    if(x>EPS) return 1;    else if(x<-EPS) return -1;    return 0;}inline bool cmp(const PO &a,const PO &b)//点排序 {    if(doublecmp(a.x-b.x)==0) return a.y<b.y;    return a.x<b.x;}inline void read(){    for(int i=1;i<=n;i++) scanf("%lf%lf",&p1[i].x,&p1[i].y);    for(int i=1;i<=m;i++) scanf("%lf%lf",&p2[i].x,&p2[i].y);}inline PO operator +(PO a,PO b){    PO c;    c.x=a.x+b.x;    c.y=a.y+b.y;    return c;}inline PO operator -(PO a,PO b){    PO c;    c.x=a.x-b.x;    c.y=a.y-b.y;    return c;}inline double dot(PO &a,PO &b,PO &c)//求点积 {    return (b.x-a.x)*(c.x-a.x)+(b.y-a.y)*(c.y-a.y);}inline double cross(PO &a,PO &b,PO &c)//求叉积 {    return (b.x-a.x)*(c.y-a.y)-(b.y-a.y)*(c.x-a.x);}inline double getlen(PO &a)//求向量的模 {    return sqrt(a.x*a.x+a.y*a.y);}inline PO getty(PO b,PO a)//求投影向量 {    PO res=a;    double k=dot(o,a,b)/(getlen(a)*getlen(a));    res.x*=k; res.y*=k;    return res;}inline double getangle(PO &a,PO &b,PO &c,PO &d)//求两向量夹角{    PO t=c+(a-d);    return cross(a,b,t);}inline PO rotate(PO a,double hd)//点a逆时针转hd角度 {    PO c;    c.x=a.x*cos(hd)-a.y*sin(hd);    c.y=a.x*sin(hd)+a.y*cos(hd);    return c;}inline double getdis(PO &a,PO &b){    return sqrt((a.x-b.x)*(a.x-b.x)+(a.y-b.y)*(a.y-b.y));}inline double getdis_ps(PO &a,PO &b,PO &c){    PO t=c-b;    t=rotate(t,-PI*0.5);    PO s=c-a;    PO ty=getty(s,t);    s=a+ty;    if(doublecmp(dot(s,b,c))<=0) return getlen(ty);    else return min(getdis(a,b),getdis(a,c));}inline void graham(PO *p,int *stk,int &top,int gs)//求凸包 {    sort(p+1,p+1+gs,cmp);    top=-1;    stk[++top]=1; stk[++top]=2;    for(int i=3;i<=gs;i++)    {        while(top>=1&&doublecmp(cross(p[stk[top-1]],p[stk[top]],p[i]))<=0) top--;        stk[++top]=i;    }    int tmp=top;    for(int i=gs-1;i>=1;i--)    {        while(top>=tmp+1&&doublecmp(cross(p[stk[top-1]],p[stk[top]],p[i]))<=0) top--;        stk[++top]=i;    }}inline double rotating_calipers()//旋转卡壳 {    int s1=0,s2=0;    for(int i=1;i<top1;i++)        if(doublecmp(p1[stk1[i]].y-p1[stk1[s1]].y)<0||(doublecmp(p1[stk1[i]].y-p1[stk1[s1]].y)==0)&&doublecmp(p1[stk1[i]].x-p1[stk1[s1]].x)<0)//找到纵坐标最小的点,纵坐标相同取横坐标较小的; s1=i;    for(int i=1;i<top2;i++)        if(doublecmp(p2[stk2[i]].y-p2[stk2[s2]].y)>0||(doublecmp(p2[stk2[i]].y-p2[stk2[s2]].y)==0)&&doublecmp(p2[stk2[i]].x-p2[stk2[s2]].x)>0)//找到纵坐标最大的点,纵坐标相同取横坐标较小的;s2=i;    int t1=s1,t2=s2;    double ans=getdis(p1[stk1[s1]],p2[stk2[s2]]);    do    {        double af=getangle(p1[stk1[s1]],p1[stk1[(s1+1)%top1]],p2[stk2[s2]],p2[stk2[(s2+1)%top2]]);        if(doublecmp(af)==0)//卡壳到两条边        {            ans=min(ans,getdis_ps(p1[stk1[s1]],p2[stk2[s2]],p2[stk2[(s2+1)%top2]]));            ans=min(ans,getdis_ps(p1[stk1[(s1+1)%top1]],p2[stk2[s2]],p2[stk2[(s2+1)%top2]]));            ans=min(ans,getdis_ps(p2[stk2[s2]],p1[stk1[s1]],p1[stk1[(s1+1)%top1]]));            ans=min(ans,getdis_ps(p2[stk2[(s2+1)%top2]],p1[stk1[s1]],p1[stk1[(s1+1)%top1]]));            s1=(s1+1)%top1; s2=(s2+1)%top2;        }        else if(doublecmp(af)>0)//卡壳到第一个的边,第二个的点        {            ans=min(ans,getdis_ps(p2[stk2[s2]],p1[stk1[s1]],p1[stk1[(s1+1)%top1]]));            s1=(s1+1)%top1;        }        else//卡壳到第二个的边,第一个的点        {            ans=min(ans,getdis_ps(p1[stk1[s1]],p2[stk2[s2]],p2[stk2[(s2+1)%top2]]));            s2=(s2+1)%top2;        }    }while(t1!=s1||t2!=s2);    return ans;}inline void go(){    graham(p1,stk1,top1,n);    graham(p2,stk2,top2,m);    printf("%.5f\n",rotating_calipers());}int main(){    while(scanf("%d%d",&n,&m),n)read(),go();    return 0;}


0 0
原创粉丝点击