第二章 递归与分治

1.最接近点对问题（递归与分治）

【问题描述】

给定二维平面上n个点，找其中的一对点，使得在n个点组成的所有点对中，该点对间的距离最小。请使用递归与分治策略求解二维平面上的最接近点对问题，使得求解该问题的时间复杂度为。给出【问题分析】、【算法伪代码】，以及【算法时间复杂度分析】（包括计算时间T的递归方程和计算结果等）。

一、问题分析:

  首先,考虑一维情况下，此时空间S中的n个点可以看为x轴上的n个实数，最接近点对就变成了n个实数中最接近的2个数。

可以用暴力的方式，求没两个临近的点的距离，记录并比较找到最小值。但是，空间中的点如果进行暴力求解，其复杂度将大的惊人。因此，选择分治法来解决这一问题。

  将所给的平面上n个点的集合S分成两个子集S1和S2,每个子集中约有n/2个点，然后在每个子集中递归地求其最接近的点对，再分别对两部分递归求解最小距离然后比较取最小值。但是有一个问题，如果产生最小距离的两个点位于两个集合内呢？因此要考虑上两个集合边界上的点对距离，再跟集合内的最小距离比较就能得到最佳答案。

  抽象而言，即要递归地在S1和S2上找出其最接近点对{p1,p2}和{q1,q2}，并设d=min{|p1-p2|,|q1-q2|}，然后求解S1上最接近分界点的点到最接近分割点的S2上的点的距离，比较取最小值。

一维情况：

先将点的横坐标x1, x2,……,xn排好序，然后用一次线性扫描就可以找出最接近点对。这种方法的主要计算时间花在排序上，耗时为O(nlogn)。 

二维情况：

选取垂直线L：x=m来分割直线，其中m为S中个点的中位数，将空间S以直线x=m分割为S1和S2。

递归地在S1和S2上找出最小距离d1,d2设d=min（d1,d2），S中的最接近点对是d，或是某个点对{p,q}，其中p∈P1,q∈P2。

关于如何寻找p、q：

对于S1中任一点p，S2中最多只有n/2个点与它构成最接近点对的候选者，故需做n^2/4次计算和比较就能确定S的最接近点对。

二、算法伪代码:

1.    一维情况Cpair1：

 

2.    二维情况Cpair2：

bool Cpair2(Pointx X[], int n, Pointx& a, Pointx& b, float& d) {

    if (n < 2)

        return false;

    Pointx* tmpX = new Pointx[n];

    MergeSort(X, tmpX, 0, n - 1);

    Pointy* Y = new Pointy[n];

    for (int i = 0; i < n; i++) {

        Y[i].p = i;

        Y[i].x = X[i].x;

        Y[i].y = X[i].y;

    }

    Pointy* tmpY = new Pointy[n];

    MergeSort(Y, tmpY, 0, n - 1);

    Pointy* Z = new Pointy[n];

    closet(X, Y, Z, 0, n - 1, a, b, d);

    delete[]Y;

    delete[]Z;

    delete[]tmpX;

    delete[]tmpY;

    return true;

}

 

 

3.    任意两点距离公式

template<class Type>

float distance(const Type& u, const Type& v) {

    float dx = u.x - v.x;

    float dy = u.y - v.y;

    return sqrt(dx * dx + dy * dy);

}

4.    计算最接近点对

void closet(Pointx X[], Pointy Y[], Pointy Z[], int l, int r, Pointx& a, Pointx& b, float& d) {

//两点的情形

    if (r - l == 1) {

        a = X[l];

        b = X[r];

        d = distance(a, b);

        return;

    }

//三点的情形

    if (r - l == 2) {

        float d1 = distance(X[l], X[l + 1]);

        float d2 = distance(X[l + 1], X[r]);

        float d3 = distance(X[l], X[r]);

        if (d1 <= d2 && d1 <= d3) {

            a = X[l];

            b = X[l + 1];

            d = d1;

        }

        else if (d2 <= d1 && d2 <= d3) {

            a = X[l + 1];

            b = X[r];

            d = d2;

        }

        else {

            a = X[l];

            b = X[r];

            d = d3;

        }

        return;

    }

//多点的情形，用分治法

    int m = (l + r) / 2;

    int f = l, g = m + 1;

    for (int i = l; i <= r; i++) {

        if (Y[i].p > m)

            Z[g++] = Y[i];

        else

            Z[f++] = Y[i];

    }

    closet(X, Z, Y, l, m, a, b, d);               //S1

    if (l == 0 && m ==(N-1)/ 2) {

        cout << d << endl;

        if(a.x<b.x){

            cout << a.x << " " << a.y << endl;

            cout << b.x << " " << b.y << endl;

        }

        else{

            cout << b.x << " " << b.y << endl;

            cout << a.x << " " << a.y << endl;

        }

    }

    float dr;

    Pointx ar, br;

    closet(X, Z, Y, m + 1, r, ar, br, dr);    //S2

    if (r == N - 1 && m ==(N-1)/ 2) {

        cout << dr << endl;

        if(ar.x<br.x){

            cout << ar.x << " " << ar.y << endl;

            cout << br.x << " " << br.y << endl;

        }

        else{

            cout << br.x << " " << br.y << endl;

            cout << ar.x << " " << ar.y << endl;

        }

    }

    if (dr < d) {

        a = ar;

        b = br;

        d = dr;

    }

//重构数组Y

    Merge(Z, Y, l, m, r);

//d矩形条内的点置于Z中

    int k = l;

    for (int i = l; i <= r; i++) {

        if (fabs(X[m].x - Y[i].x) < d) {

            Z[k++] = Y[i];

        }

    }

//搜索Z[1:k-1]

    for (int i = l; i < k; i++) {

        for (int j = i + 1; j < k && Z[j].y - Z[i].y < d; j++) {

            float dp = distance(Z[i], Z[j]);

            if (dp < d) {

                d = dp;

                a = X[Z[i].p];

                b = X[Z[j].p];

            }

        }

    }

}

三、算法时间复杂度分析（计算时间T的递归方程、计算结果）:

  在分治算法中，当求解n个点的集合的最近点对时，对于上述三类情况中的前两者可由递归算得，而分析可得第三类情况的时间代价 ，合并后问题求解的总的时间按复杂度为：

  由此易知，T(n)=O(nlogn),预排序所需的计算时间为O(nlogn).因此，整个算法所需时间为O(nlogn)。在渐近的意义下，此算法已为最优算法。

2.最接近临点对问题

【问题描述】给定二维平面上n个点，找其中的一对点，使得在n个点组成的所有点对中，该点对间的距离最小。使用递归与分治策略求解二维平面上的最接近点对问题。假设所有点的集合为S，m为S中所有点的x坐标的中位数，垂直线x=m将集合S均匀分割为左右两个子集合S1和S2（当集合S中点的个数为奇数时，左集合S1中的点的个数比右集合S2中点的个数多1）。

【输入形式】在屏幕上输入点的个数，以及所有点的x和y坐标。

【输出形式】第一次分割时，将所有点集合S分割为左右两个子集合S1和S2，分别输出左右子集合S1和S2，以及所有点集合S的最接近点对的距离以及最接近点对。

【样例输入】

10

-15.4 -57.3

13.2 30.1

-87.5 93.2

47.6 -12.7

94.7 61.5

56.8 -57.1

27.8 43.5

-28.1 19.0

-96.2 47.5

55.5 -93.3

【样例输出】

42.8

-28.1 19.0

13.2 30.1

36.2

55.5 -93.3

56.8 -57.1

19.8

13.2 30.1

27.8 43.5

【样例说明】

输入：10个点，后续每行为每一点的x和y坐标。

输出：左右子集合S1和S2，以及所有点集合S的最接近点对的距离以及最接近点对。例如，前面三行中，S1的最接近点对的距离为42.8，最接近点对的x和y坐标分别为(-28.1,19.0)和(13.2,30.1)。输出最接近点对坐标时，先输出的点的x坐标小于后输出点的x坐标。中间三行和最后三行分别为子集合S2和集合S的最接近点对的距离以及最接近点对。

如下图所示，子集合S1点以蓝色表示，子集合S2以绿色表示。蓝色连线为子集合S1最接近点对间的线段；绿色连线为子集合S2最接近点对间的线段；紫色连线为集合S最接近点对间的线段。

代码实现：

#include<iostream>
#include<cmath>
using namespace std;
const int M = 50;
int N;
class Pointx {
public:
    int operator<=(Pointx a) const {
        return(x <= a.x);
    }
    int ID;
    float x, y;
};

class Pointy {
public:
    int operator<=(Pointy a) const {
        return(y <= a.y);
    }
    int p;
    float x, y;
};

template<class Type>
float distance(const Type& u, const Type& v) {
    float dx = u.x - v.x;
    float dy = u.y - v.y;
    return sqrt(dx * dx + dy * dy);
}

template<typename Type>                
void Copy(Type a[], Type b[], int left, int right) {
    for (int i = left; i <= right; i++) {
        a[i] = b[i];
    }
}

template<class Type>
void Merge(Type c[], Type d[], int l, int m, int r) {
    int i = l, j = m + 1, k = l;
    while ((i <= m) && (j <= r)) {
        if (c[i] <= c[j])
        {
            d[k++] = c[i++];
        }
        else
        {
            d[k++] = c[j++];
        }
    }
    if (i > m) {
        for (int q = j; q <= r; q++) {
            d[k++] = c[q];
        }
    }
    else {
        for (int q = i; q <= m; q++) {
            d[k++] = c[q];
        }
    }
}

template<class Type>
void MergeSort(Type a[], Type b[], int left, int right) {
    if (left < right) {
        int i = (left + right) / 2;
        MergeSort(a, b, left, i);
        MergeSort(a, b, i + 1, right);
        Merge(a, b, left, i, right);        //合并到b
        Copy(a, b, left, right);            //将b复制到a
    }
}

void closet(Pointx X[], Pointy Y[], Pointy Z[], int l, int r, Pointx& a, Pointx& b, float& d) {
    if (r - l == 1) {
        a = X[l];
        b = X[r];
        d = distance(a, b);
        return;
    }

    if (r - l == 2) {
        float d1 = distance(X[l], X[l + 1]);
        float d2 = distance(X[l + 1], X[r]);
        float d3 = distance(X[l], X[r]);
        if (d1 <= d2 && d1 <= d3) {
            a = X[l];
            b = X[l + 1];
            d = d1;
        }
        else if (d2 <= d1 && d2 <= d3) {
            a = X[l + 1];
            b = X[r];
            d = d2;
        }
        else {
            a = X[l];
            b = X[r];
            d = d3;
        }
        return;
    }

    int m = (l + r) / 2;
    int f = l, g = m + 1;
    for (int i = l; i <= r; i++) {
        if (Y[i].p > m)
            Z[g++] = Y[i];
        else
            Z[f++] = Y[i];
    }

    closet(X, Z, Y, l, m, a, b, d);                //S1
    if (l == 0 && m ==(N-1)/ 2) {
        cout << d << endl;
        if(a.x<b.x){
            cout << a.x << " " << a.y << endl;
            cout << b.x << " " << b.y << endl;
        }
        else{
            cout << b.x << " " << b.y << endl;
            cout << a.x << " " << a.y << endl;
        }
    }
    float dr;
    Pointx ar, br;
    closet(X, Z, Y, m + 1, r, ar, br, dr);        //S2
    if (r == N - 1 && m ==(N-1)/ 2) {
        cout << dr << endl;
        if(ar.x<br.x){
            cout << ar.x << " " << ar.y << endl;
            cout << br.x << " " << br.y << endl;
        }
        else{
            cout << br.x << " " << br.y << endl;
            cout << ar.x << " " << ar.y << endl;
        }
    }
    if (dr < d) {
        a = ar;
        b = br;
        d = dr;
    }

    Merge(Z, Y, l, m, r);

    int k = l;
    for (int i = l; i <= r; i++) {
        if (fabs(X[m].x - Y[i].x) < d) {
            Z[k++] = Y[i];
        }
    }
    for (int i = l; i < k; i++) {
        for (int j = i + 1; j < k && Z[j].y - Z[i].y < d; j++) {
            float dp = distance(Z[i], Z[j]);
            if (dp < d) {
                d = dp;
                a = X[Z[i].p];
                b = X[Z[j].p];
            }
        }
    }
}

bool Cpair2(Pointx X[], int n, Pointx& a, Pointx& b, float& d) {
    if (n < 2)
        return false;
    Pointx* tmpX = new Pointx[n];
    MergeSort(X, tmpX, 0, n - 1);

    Pointy* Y = new Pointy[n];
    for (int i = 0; i < n; i++) {
        Y[i].p = i;
        Y[i].x = X[i].x;
        Y[i].y = X[i].y;
    }

    Pointy* tmpY = new Pointy[n];
    MergeSort(Y, tmpY, 0, n - 1);

    Pointy* Z = new Pointy[n];
    closet(X, Y, Z, 0, n - 1, a, b, d);

    delete[]Y;
    delete[]Z;
    delete[]tmpX;
    delete[]tmpY;
    return true;

}

int main() {
    cin>>N;
    Pointx X[M];        ///
    for (int i = 0; i < N; i++) {
        X[i].ID = i;
        cin>>X[i].x>>X[i].y;
    }
    cout.setf(ios::fixed);
    cout.precision(1);
    Pointx a;
    Pointx b;
    float d=0;
    Cpair2(X, N, a, b, d);
    cout<<d<<endl;
    if(a.x<b.x){
        cout<<a.x<<" "<<a.y<<endl;
        cout<<b.x<<" "<<b.y<<endl;
    }
    else{
        cout<<b.x<<" "<<b.y<<endl;
        cout<<a.x<<" "<<a.y<<endl;
    }
    return 0;
}