【計算幾何 02】凸包問題(Convex Hull)


引言

首先介紹下什么是凸包?如下圖:

在一個二維坐標系中,有若干點雜亂排列着,將最外層的點連接起來構成的凸多邊型,它能包含給定的所有的點,這個多邊形就是凸包。

實際上可以理解為用一個橡皮筋包含住所有給定點的形態。

凸包用最小的周長圍住了給定的所有點。如果一個凹多邊形圍住了所有的點,它的周長一定不是最小,如下圖。根據三角不等式,凸多邊形在周長上一定是最優的。

凸包的求法

尋找凸包的算法有很多種,常用的求法有 Graham 掃描法和 Andrew 算法

Graham Scan 算法求凸包

Graham Scan 算法是一種十分簡單高效的二維凸包算法,能夠在 \(O(nlogn)\) 的時間內找到凸包。

Graham Scan 算法的做法是先確定一個起點(一般是最左邊的點和最右邊的點),然后一個個點掃過去,如果新加入的點和之前已經找到的點所構成的 "殼" 凸性沒有變化,就繼續掃,否則就把已經找到的最后一個點刪去,再比較凸性,直到凸性不發生變化。分別掃描上下兩個 "殼",合並在一起,凸包就找到了。這么說很抽象,我們看圖來解釋:

先找 "下殼",上下其實是一樣的。首先加入兩個點 A 和 B。

然后插入第三個點 C,並計算 \(\overrightarrow{AB}×\overrightarrow{BC}\) 的向量積,卻發現向量積系數小於(等於)0,也就是說 \(\overrightarrow{BC}\)\(\overrightarrow{AB}\) 的順時針方向上。

於是刪去 B 點。

按照這樣的方法依次掃描,找完 "下殼" 后,再找 "上殼"。

關於掃描的順序,有坐標序和極角序兩種,本文采用前者。坐標序是比較兩個點的 x 坐標,小的先被掃描(掃描上凸殼的時候反過來),如果兩個點 x 坐標相同,那么就比較 y 坐標,同樣的也是小的先被掃描(掃描上凸殼的時候也是反過來)。極角序使用 atan2 函數的返回值進行比較,讀者可以自己嘗試寫下。

下面貼下代碼:Graham Scan 算法

struct Point {
    double x, y;

    Point operator-(Point& p) {
        Point t;
        t.x = x - p.x;
        t.y = y - p.y;
        return t;
    }

    double cross(Point p)  // 向量叉積
    {
        return x * p.y - p.x * y;
    }
};

bool cmp(Point& p1, Point& p2) {
    if (p1.x != p2.x) return p1.x < p2.x;

    return p1.y < p2.y;
}

Point point[1005];  // 無序點
int convex[1005];   // 保存組成凸包的點的下標
int n;              // 坐標系的無序點的個數

int GetConvexHull() {
    sort(point, point + n, cmp);
    int temp;
    int total = 0;

    for (int i = 0; i < n; i++)  // 下凸包
    {
        while (total > 1 &&
               (point[convex[total - 1]] - point[convex[total - 2]])
                       .cross(point[i] - point[convex[total - 1]]) <= 0)
            total--;

        convex[total++] = i;
    }

    temp = total;

    for (int i = n - 2; i >= 0; i--)  // 上凸包
    {
        while (total > temp &&
               (point[convex[total - 1]] - point[convex[total - 2]])
                       .cross(point[i] - point[convex[total - 1]]) <= 0)
            total--;

        convex[total++] = i;
    }

    return total -
           1;  // 返回組成凸包的點的個數,實際上多了一個,就是起點,所以組成凸包的點個數是
               // total - 1
}

Andrew 算法求凸包

首先把所有點以橫坐標為第一關鍵字,縱坐標為第二關鍵字排序。

顯然排序后最小的元素和最大的元素一定在凸包上。而且因為是凸多邊形,我們如果從一個點出發逆時針走,軌跡總是“左拐”的,一旦出現右拐,就說明這一段不在凸包上。因此我們可以用一個單調棧來維護上下凸殼。

因為從左向右看,上下凸殼所旋轉的方向不同,為了讓單調棧起作用,我們首先 升序枚舉 求出下凸殼,然后 降序 求出上凸殼。

求凸殼時,一旦發現即將進棧的點( \(P\) )和棧頂的兩個點( \(S_1,S_2\) ,其中 \(S_1\) 為棧頂)行進的方向向右旋轉,即叉積小於 \(0\)\(\overrightarrow{S_2S_1}\times \overrightarrow{S_1P}<0\) ,則彈出棧頂,回到上一步,繼續檢測,直到 \(\overrightarrow{S_2S_1}\times \overrightarrow{S_1P}\ge 0\) 或者棧內僅剩一個元素為止。

通常情況下不需要保留位於凸包邊上的點,因此上面一段中 \(\overrightarrow{S_2S_1}\times \overrightarrow{S_1P}<0\) 這個條件中的“ \(<\) ”可以視情況改為 \(\le\) ,同時后面一個條件應改為 \(>\)

代碼實現
// stk[]是整型,存的是下標
// p[]存儲向量或點
tp = 0;                       //初始化棧
std::sort(p + 1, p + 1 + n);  //對點進行排序
stk[++tp] = 1;
//棧內添加第一個元素,且不更新used,使得1在最后封閉凸包時也對單調棧更新
for (int i = 2; i <= n; ++i) {
  while (tp >= 2  //下一行*被重載為叉積
         && (p[stk[tp]] - p[stk[tp - 1]]) * (p[i] - p[stk[tp]]) <= 0)
    used[stk[tp--]] = 0;
  used[i] = 1;  // used表示在凸殼上
  stk[++tp] = i;
}
int tmp = tp;  // tmp表示下凸殼大小
for (int i = n - 1; i > 0; --i)
  if (!used[i]) {
    //      ↓求上凸殼時不影響下凸殼
    while (tp > tmp && (p[stk[tp]] - p[stk[tp - 1]]) * (p[i] - p[stk[tp]]) <= 0)
      used[stk[tp--]] = 0;
    used[i] = 1;
    stk[++tp] = i;
  }
for (int i = 1; i <= tp; ++i)  //復制到新數組中去
  h[i] = p[stk[i]];
int ans = tp - 1;

根據上面的代碼,最后凸包上有 \(ans\) 個元素(額外存儲了 \(1\) 號點,因此 \(h\) 數組中有 \(ans+1\) 個元素),並且按逆時針方向排序。周長就是

\[\sum_{i=1}^{ans}\left|\overrightarrow{h_ih_{i+1}}\right| \]

參考


免責聲明!

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



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