模板 旋转卡壳 凸包

模板 旋转卡壳 凸包

 

 

好早之前看的，如今再记下来吧，当作复习一遍。

那么，先提一下最基本最暴力的求凸包直径的方法吧---枚举。。。好吧。。不少问题均可以用 枚举 这个“万能”的方法来解决，过程很简单方即是确定的，不过在效率上就要差很远了。  要求一个点集的直径，即便先计算出这个点集的凸包，而后再枚举凸包上的点对，这样来求点集直径的话依然会在凸包上点的数量达到O（n）级别是极大的下降它的效率，也浪费了凸包的优美性质。不过在数据量较小或者很适合时，何须要大费周折的用那些麻烦复杂的算法呢，枚举依然是解决问题的很好的方法之一。

 

而后就是今天的旋转卡壳算法了。

旋转卡壳能够用于求凸包的直径、宽度，两个不相交凸包间的最大距离和最小距离等。虽然算法的思想不难理解，可是实现起来真的很容易让人“卡壳”。

其实简单来讲就是用一对平行线“卡”住凸包进行旋转。

被一对卡壳正好卡住的对应点对称为对踵点，对锺点的具体定义很差说，不过从图上仍是比较好理解的。

能够证实对踵点的个数不超过3N/2个 也就是说对踵点的个数是O(N)的

对踵点的个数也是咱们下面解决问题时间复杂度的保证。

 

卡壳呢，具体来讲有两种状况：

1.

一种是这样，两个平行线正好卡着两个点；

2.

 

一种是这样，分别卡着一条边和一个点。

而第二种状况在实现中比较容易处理，这里就只研究第二种状况。

 

在第二种状况中 咱们能够看到 一个对踵点和对应边之间的距离比其余点要大

 

也就是一个对踵点和对应边所造成的三角形是最大的 下面咱们会据此获得对踵点的简化求法。

下面给出一个官方的说明：

http://cgm.cs.mcgill.ca/~orm/rotcal.frame.html

Compute the polygon's extreme points in the y direction. Call them ymin and ymax. Construct two horizontal lines of support through ymin and ymax. Since this is already an anti-podal pair, compute the distance, and keep as maximum. Rotate the lines until one is flush with an edge of the polygon. A new anti-podal pair is determined. Compute the new distance, compare to old maximum, and update if necessary. Repeat steps 3 and 4 until the anti-podal pair considered is (ymin,ymax) again. Output the pair(s) determining the maximum as the diameter pair(s). 

要是真的按这个实现起来就麻烦到吐了。。

根据上面的第二种状况，咱们能够获得下面的方法：

若是qa，qb是凸包上最远两点，必然能够分别过qa，qb画出一对平行线。经过旋转这对平行线，咱们可让它和凸包上的一条边重合,如图中蓝色直线，能够注意到，qa是凸包上离p和qb所在直线最远的点。因而咱们的思路就是枚举凸包上的全部边，对每一条边找出凸包上离该边最远的顶点，计算这个顶点到该边两个端点的距离，并记录最大的值。

直观上这是一个O(n2)的算法，和直接枚举任意两个顶点同样了。

然而咱们能够发现 凸包上的点依次与对应边产生的距离成单峰函数（以下图：）

 

这个性质就很重要啦。

 

根据这个凸包的特性，咱们注意到逆时针枚举边的时候，最远点的变化也是逆时针的，这样就能够不用从头计算最远点，而能够紧接着上一次的最远点继续计算。因而咱们获得了O(n)的算法。这就是所谓的“旋转”吧！

利用旋转卡壳，咱们能够在O（n）的时间内获得凸包的对锺点中的长度最长的点对。

 

又因为最远点对必然属于对踵点对集合 ，那么咱们先求出凸包 而后求出对踵点对集合 而后选出距离最大的便可。

 

 

那么具体的代码就很容易实现了，利用叉积，代码只有这么几行的长度：

 1 double rotating_calipers(Point *ch,int n) 
 2 {
 3     int q=1;
 4     double ans=0;
 5     ch[n]=ch[0];
 6      for(int p=0;p<n;p++)
 7 {
 8       while(cross(ch[p+1],ch[q+1],ch[p])>cross(ch[p+1],ch[q],ch[p]))
 9             q=(q+1)%n;
10     ans=max(ans,max(dist(ch[p],ch[q]),dist(ch[p+1],ch[q+1])));
11      }
12 return ans;
13 }

其中cross()是计算叉积，由于凸包上距离一条边最远的点和这条边的两个端点构成的三角形面积是最大的。之因此既要更新(ch[p],ch[q])又要更新(ch[p+1],ch[q+1])是为了处理凸包上两条边平行的特殊状况。

下面是道很基础的旋转卡壳求凸包直径的题了：

http://poj.org/problem?id=2187

http://acm.csu.edu.cn/OnlineJudge/problem.php?id=1393

用上面的知识很容易 解决，直接套两个模板就能够了。

直接代码：

 1 #include<iostream>
 2 #include<cstdio>
 3 #include<cstring>
 4 #include<algorithm>
 5 #include<cmath>
 6 #define N 50005
 7 using namespace std;
 8 struct point
 9 {
10        int x,y;
11 }p[N];
12 int top,n,s[N];
13 int cross(point a,point b,point c)//叉积
14 {
15     return (b.x - a.x) * (c.y - a.y) - (c.x - a.x) * (b.y - a.y);
16 }
17 
18 int dis(point a,point b)//距离
19 {
20        return (a.x - b.x) * (a.x - b.x) + (a.y - b.y) * (a.y - b.y);
21 }
22 
23 bool cmp(point a,point b)
24 {
25      int ans = cross(p[0],a,b);
26      if( ans > 0 || (ans == 0 && dis(p[0],a) > dis(p[0],b) )) return true;
27      return false;
28 }
29 
30 void graham()//凸包模板。。
31 {
32      s[0] = 0;
33      s[1] = 1;
34      top = 1;
35      for(int i = 2;i != n ;i ++)
36      {
37              while(top && cross(p[s[top - 1]],p[s[top]],p[i] ) < 0) top--;
38              s[++top] = i;
39      }
40      top ++;
41 }
42 
43 void RC()//旋转卡壳
44 {
45      int q,p1,pp,qq,rr,r,ans = 0;
46      q = 1;
47      ans = dis(p[s[0]],p[s[1]]);
48      for(int i = 0; i != top ;i ++)
49      {
50            while(abs(cross(p[s[(i+1)%top]],p[s[i%top]],p[s[(q+1)%top]])) > abs(cross(p[s[(i+1)%top]],p[s[i%top]],p[s[q%top]]))) q = (q + 1)%top;
51              ans = max(ans , max(dis(p[s[(i+1)%top]],p[s[q]]),dis(p[s[i%top]],p[s[q]])));
52      }
53      printf("%d\n",ans);
54 }
55      
56 int main()
57 {
58     while(~scanf("%d",&n))
59     {
60                           for(int i =0 ;i != n;i ++) scanf("%d%d",&p[i].x,&p[i].y);
61                           int u = 0;
62                           for(int i = 1;i != n;i ++)//寻找基点
63                           {
64                                   if(p[u].y > p[i].y || (p[u].y == p[i].y && p[u].x > p[i].x)) u = i;
65                           }
66                           if(u)
67                           {
68                                swap(p[u].y,p[0].y);
69                                swap(p[u].x,p[0].x);
70                           }
71                           sort(p + 1,p + n,cmp);
72                           graham();
73                           RC();
74     }
75 return 0;
76 }

View Code

 

下一页：