旋轉卡殼(1)--求凸包(點集)直徑 poj 2187


好早以前看的,現在再記下來吧,當做復習一遍。

那么,先提一下最基本最暴力的求凸包直徑的方法吧---枚舉。。。好吧。。很多問題都可以用 枚舉 這個“萬能”的方法來解決,過程很簡單方便是肯定的,不過在效率上就要差很遠了。  要求一個點集的直徑,即使先計算出這個點集的凸包,然后再枚舉凸包上的點對,這樣來求點集直徑的話依然會在凸包上點的數量達到O(n)級別是極大的降低它的效率,也浪費了凸包的優美性質。不過在數據量較小或者很適合時,何必要大費周折的用那些麻煩復雜的算法呢,枚舉依然是解決問題的很好的方法之一。

 

然后就是今天的旋轉卡殼算法了。(那個字念 qia 。。。搞了好久才發現讀都讀錯了。囧。。。)

旋轉卡殼可以用於求凸包的直徑、寬度,兩個不相交凸包間的最大距離和最小距離等。雖然算法的思想不難理解,但是實現起來真的很容易讓人“卡殼”。

其實簡單來說就是用一對平行線“卡”住凸包進行旋轉。

被一對卡殼正好卡住的對應點對稱為對踵點,對鍾點的具體定義不好說,不過從圖上還是比較好理解的。

可以證明對踵點的個數不超過3N/2個 也就是說對踵點的個數是O(N)

對踵點的個數也是我們下面解決問題時間復雜度的保證。

 

卡殼呢,具體來說有兩種情況:

1.

一種是這樣,兩個平行線正好卡着兩個點;

2.

 

一種是這樣,分別卡着一條邊和一個點。

而第二種情況在實現中比較容易處理,這里就只研究第二種情況。

 

在第二種情況中 我們可以看到 一個對踵點和對應邊之間的距離比其他點要大(借用某大神的圖··)

 

也就是一個對踵點和對應邊所形成的三角形是最大的 下面我們會據此得到對踵點的簡化求法。

下面給出一個官方的說明:

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])是為了處理凸包上兩條邊平行的特殊情況。

下面是道很基礎的旋轉卡殼求凸包直徑的題了:

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

用上面的知識很容易 解決,直接套兩個模板就可以了。 直接代碼:

 

 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 }

 

 

 

 

 

 

 

 

 

 

 

 

 


免責聲明!

本站轉載的文章為個人學習借鑒使用,本站對版權不負任何法律責任。如果侵犯了您的隱私權益,請聯系本站郵箱yoyou2525@163.com刪除。



 
粵ICP備18138465號   © 2018-2025 CODEPRJ.COM