天天看點

程式設計之美之尋找最近點對

轉載于  http://sxnuwhui.blog.163.com/blog/static/137068373201264104915935/?COLLCC=3097019025&COLLCC=3097019024&

在二維平面上的n個點中,如何快速的找出最近的一對點,就是最近點對問題。

    一種簡單的想法是暴力枚舉每兩個點,記錄最小距離,顯然,時間複雜度為O(n^2)。

    在這裡介紹一種時間複雜度為O(nlognlogn)的算法。其實,這裡用到了分治的 思想。将所給平面上n個點的集合S分成兩個子集S1和S2,每個子集中約有n/2個點。然後在每個子集中遞歸地求最接近的點對。在這裡,一個關鍵的問題是 如何實作分治法中的合并步驟,即由S1和S2的最接近點對,如何求得原集合S中的最接近點對。如果這兩個點分别在S1和S2中,問題就變得複雜了。

    為了使問題變得簡單,首先考慮一維的情形。此時,S中的n個點退化為x軸上的n個實數x1,x2,...,xn。最接近點對即為這n個實數中相差最小的兩個實數。顯然可以先将點排好序,然後線性掃描就可以了。但我們為了便于推廣到二維的情形,嘗試用分治法解決這個問題。

    假設我們用m點将S分為S1和S2兩個集合,這樣一來,對于所有的p(S1中的點)和q(S2中的點),有p<q。

    遞歸地在S1和S2上找出其最接近點對{p1,p2}和{q1,q2},并設

d = min{ |p1-p2| , |q1-q2| }

    由此易知,S中最接近點對或者是{p1,p2},或者是{q1,q2},或者是某個{q3,p3},如下圖所示。

程式設計之美之尋找最近點對

    如果最接近點對是{q3,p3},即|p3-q3|<d,則p3和q3兩者與m的距離都不超過d,且在區間(m-d,d]和(d,m+d]各有且僅有一個點。這樣,就可以線上性時間内實作合并。

    此時,一維情形下的最近點對時間複雜度為O(nlogn)。

    在二維情形下,類似的,利用分治法,但是難點在于如何實作線性的合并?

程式設計之美之尋找最近點對

    由上圖可見,形成的寬為2d的帶狀區間,最多可能有n個點,合并時間最壞情況下為n^2,。但是,P1和P2中的點具有以下稀疏的性質,對于P1中的任意一點,P2中的點必定落在一個d X 2d的矩形中,且最多隻需檢查六個點(鴿巢原理)。

    這樣,先将帶狀區間的點按y坐标排序,然後線性掃描,這樣合并的時間複雜度為O(nlogn),幾乎為線性了。

算法描述:已知集合S中有n個點,分治法的思想就是将S進行拆分,分為2部分求最近點對。算法每次選擇一條垂線L,将S拆分左右兩部分為SL和SR,L一般取點集S中所有點的中間點的x坐标來劃分,這樣可以保證SL和SR中的點數目各為n/2,

(否則以其他方式劃分S,有可能導緻SL和SR中點數目一個為1,一個為n-1,不利于算法效率,要盡量保持樹的平衡性)

依次找出這兩部分中的最小點對距離:δL和δR,記SL和SR中最小點對距離δ = min(δL,δR),如圖1:

程式設計之美之尋找最近點對

     以L為中心,δ為半徑劃分一個長帶,最小點對還有可能存在于SL和SR的交界處,如下圖2左圖中的虛線帶,p點和q點分别位于SL和SR的虛線範圍内,在這個範圍内,p點和q點之間的距離才會小于δ,最小點對計算才有意義。

程式設計之美之尋找最近點對

Figure 2

      對于SL虛框範圍内的p點,在SR虛框中與p點距離小于δ的頂多隻有六個點,就是圖二右圖中的2個正方形的6的頂點。這個可以反推證明,如果右邊這2個正方形内有7個點與p點距離小于δ,例如q點,則q點與下面正方形的四個頂點距離小于δ,則和δ為SL和SR中的最小點對距離相沖突。是以對于SL虛框中的p點,不需求出p點和右邊虛線框内所有點距離,隻需計算SR中與p點y坐标距離最近的6個點,就可以求出最近點對,節省了比較次數。

(否則的話,最壞情形下,在SR虛框中有可能會有n/2個點,對于SL虛框中的p點,每次要比較n/2次,浪費了算法的效率)

HDOJ1007:http://acm.hdu.edu.cn/showproblem.php?pid=1007

代碼如下:

#include <iostream>
#include <cstdio>
#include <cstring>
#include <cmath>
#include <algorithm>
using namespace std;
const double INF = 1e20;
const int N = 100005;

struct Point
{
    double x;
    double y;
}point[N];
int n;
int tmpt[N];

bool cmpxy(const Point& a, const Point& b)
{
    if(a.x != b.x)
        return a.x < b.x;
    return a.y < b.y;
}

bool cmpy(const int& a, const int& b)
{
    return point[a].y < point[b].y;
}

double min(double a, double b)
{
    return a < b ? a : b;
}

double dis(int i, int j)
{
    return sqrt((point[i].x-point[j].x)*(point[i].x-point[j].x)
                + (point[i].y-point[j].y)*(point[i].y-point[j].y));
}

double Closest_Pair(int left, int right)
{
    double d = INF;
    if(left==right)
        return d;
    if(left + 1 == right)
        return dis(left, right);
    int mid = (left+right)>>1;
    double d1 = Closest_Pair(left,mid);
    double d2 = Closest_Pair(mid+1,right);
    d = min(d1,d2);
    int i,j,k=0;
    //分離出寬度為d的區間
    for(i = left; i <= right; i++)
    {
        if(fabs(point[mid].x-point[i].x) <= d)
            tmpt[k++] = i;
    }
    sort(tmpt,tmpt+k,cmpy);
    //線性掃描
    for(i = 0; i < k; i++)
    {
        for(j = i+1; j < k && point[tmpt[j]].y-point[tmpt[i]].y<d; j++)//  了解!!
        {
            double d3 = dis(tmpt[i],tmpt[j]);
            if(d > d3)
                d = d3;
        }
    }
    return d;
}


int main()
{
    while(true)
    {
        scanf("%d",&n);
        if(n==0)
            break;
        for(int i = 0; i < n; i++)
            scanf("%lf %lf",&point[i].x,&point[i].y);
        sort(point,point+n,cmpxy);
        printf("%.2lf\n",Closest_Pair(0,n-1)/2);
    }
    return 0;
}

           

繼續閱讀