分治法找平面点集中最接近点对

分治法找平面中最近点对

写在前面

作为一个计算机系的菜鸟，这两年也没怎么好好写过代码，就从这篇博客开始好好记录一下平时作业中的代码吧，希望不会鸽的太狠。

正文

问题描述：利用分治法寻找平面上2048个点中距离最近的两个点，这里我们只找一对即可。

算法详解：分治法的思想就是将一个大问题不断分成小问题，一般情况我们分为两个小问题并且处理中间剩余情况所带来的新问题。对于这道题我们可以把大的点集分成两个小点集，剩余情况为距离最近的两个点分别位于两个小点集的情况。

​ 这里对于大点集的划分我们采用根据x轴坐标中位数划分的办法，即将点集中的点按照x坐标升序排序后前一半为第一个小点集，后一半为第二个小点集。这里我们就会得到
$$T(n)=2T(n/2)+?$$
​ 易知它的时间复杂度取决于剩余情况的算法，如果我们对两个点集中的点进行遍历搜索的话那时间复杂度将会达到$O(n^2)$，这显然是我们使用分治法时所不可接受的，经过分析我们可以知道，假设在小点集中的最小距离为$d$，大点集中x轴坐标中位数为$x_0$，则我们只需选取$|x_p-x_0|<d$的P点进行搜索即可在$O(n)$时间内完成对所有可能情况的搜索，这一段的证明我参考了其他大佬的博客，参考链接如下：

https://blog.csdn.net/liufeng_king/article/details/8484284

尚存在的问题：本题目要求将重复点去掉，但是我并没有想到太好的能在短时间复杂度内完成的算法，故暂时没有进行实现（我太菜了），希望有好心人指点一二。

附录：

#include <iostream>
#include <stdio.h>
#include <time.h>
#include <math.h>
#define file "C:/Users/12178/Desktop/point_2048.csv"
#define N 2048

using namespace std;

struct Point       //二维平面点
{
    int id;
    int x;
    int y;
};

float dis(Point a, Point b)     //计算距离
{
    float dx = a.x-b.x;
    float dy = a.y-b.y;
    return sqrt(dx*dx+dy*dy);
}

void Merge(Point a[], Point b[], int l, int m, int r)     //归并
{
    int i = l, j = m+1, k = l;
    while(i<=m && j<=r)
    {
        if(a[i].x < a[j].x)
        {
            b[k++] = a[i++];
        }
        else
        {
            b[k++] = a[j++];
        }        
    }

    if(i>m)
    {
        for(int h=j; h<=r; h++)
        {
            b[k++] = a[h];
        }
    }
    else
    {
        for(int h=i; h<=m; h++)
        {
            b[k++] = a[h];
        }
    }
}

void MergeSort(Point a[], Point b[], int l, int r)     //归并排序
{
    if(l<r)
    {
        int m = (l+r)/2;
        MergeSort(a,b,l,m);
        MergeSort(a,b,m+1,r);
        Merge(a,b,l,m,r);
        for(int i=l; i<=r; i++)
        {
            a[i] = b[i];
        }
    }
}

void Closest(Point p[], Point p1[], int l, int r, Point &a, Point &b, float &d)//分治算法
{
    if(r-l == 1)//只有两个点的情况
    {
        a = p[l];
        b = p[r];
        d = dis(a,b);
        return;
    }

    if(r-l == 2)//三个点的情况
    {
        float d1 = dis(p[l],p[l+1]);
        float d2 = dis(p[l+1],p[r]);
        float d3 = dis(p[l],p[r]);
        if(d1<=d2 && d1<=d3)
        {
            d = d1;
            a = p[l];
            b = p[l+1];
        }
        else if (d2<=d3)
        {
            d = d2;
            a = p[l+1];
            b = p[r];
        }
        else
        {
            d = d3;
            a = p[l];
            b = p[r];
        }
        return;  
    }

    int m = (l+r)/2;    //中点
    
    Closest(p,p1,l,m,a,b,d);    //左边
    Point a1,b1;
    float d1;
    Closest(p,p1,m+1,r,a1,b1,d1);   //右边
    if(d1<d)
    {
        a = a1;
        b = b1;
        d = d1;
    }

    //中间部分
    int h = l;
    for(int i=l; i<=r; i++)
    {
        if(abs(p[i].x-p[m].x) < d)
        {
            p1[h++] = p[i];
        }
    }

    //搜索中间的矩形区域
    for(int i=l; i<h; i++)
    {
        for(int j=i+1; j<h && abs(p1[j].y-p1[i].y)<d; j++)
        {
            float d2 = dis(p1[i], p1[j]);
            if(d2<d)
            {
                a = p1[i];
                b = p1[j];
                d = d2;
            }
        }
    }
}

bool DivideAndConquer(Point p[], int n, Point &a, Point &b, float &d)
{
    if(n<2) return false;

    Point *tmp = new Point[n];
    MergeSort(p, tmp, 0, n-1);
    Point *p1 = new Point[n];
    Closest(p, p1, 0, n-1, a, b, d);    //开始分治
    return true;
}

int main()
{
    FILE *fp = fopen(file, "r");
    
    int x=0,y=0,num=0;
    Point p[N];
    while(fscanf(fp, "%d,%d\n", &x, &y) != EOF)   
    {
        p[num].x = x;
        p[num].y = y;
        p[num].id = num;
        num++;
    }

    //clock_t begin = clock();
    Point a,b;
    float d;
    if(!DivideAndConquer(p,num,a,b,d))
    {
        cout << "Input Error!" << endl;
    }
    //clock_t end = clock();
    //double Runtime = (double)(end-begin)*1000/CLOCKS_PER_SEC;

    cout<<"最邻近点对为：("<<a.x<<","<<a.y<<")和("<<b.x<<","<<b.y<<") "<<endl;
    cout<<"最邻近距离为： "<<d<<endl;
    //cout<<"用时："<<Runtime<<"ms"<<endl;
    
    return 0;
}