0007算法笔记分治法最接近点对问题.docx

0007算法笔记分治法最接近点对问题.docx

《0007算法笔记分治法最接近点对问题.docx》由会员分享，可在线阅读，更多相关《0007算法笔记分治法最接近点对问题.docx（23页珍藏版）》请在冰豆网上搜索。

0007算法笔记分治法最接近点对问题

问题场景：

在应用中，常用诸如点、圆等简单的几何对象代表现实世界中的实体.在涉及这些几何对象的问题中，常需要了解其邻域中其他几何对象的信息.例如，在空中交通控制问题中，若将飞机作为空间中移动的一个点来看待,则具有最大碰撞危险的2架飞机，就是这个空间中最接近的一对点。

这类问题是计算几何学中研究的基本问题之一。

      问题描述：

给定平面上n个点，找其中的一对点,使得在n个点的所有点对中,该点对的距离最小。

严格地说，最接近点对可能多于1对。

为了简单起见，这里只限于找其中的一对.

     1、一维最接近点对问题

       算法思路：

      这个问题很容易理解，似乎也不难解决。

我们只要将每一点与其他n—1个点的距离算出，找出达到最小距离的两个点即可。

然而，这样做效率太低，需要O（n^2）的计算时间。

在问题的计算复杂性中我们可以看到,该问题的计算时间下界为Ω（nlogn）。

这个下界引导我们去找问题的一个θ（nlogn）算法.采用分治法思想，考虑将所给的n个点的集合S分成2个子集S1和S2,每个子集中约有n/2个点，然后在每个子集中递归地求其最接近的点对。

在这里，一个关键的问题是如何实现分治法中的合并步骤，即由S1和S2的最接近点对，如何求得原集合S中的最接近点对，因为S1和S2的最接近点对未必就是S的最接近点对.如果组成S的最接近点对的2个点都在S1中或都在S2中，则问题很容易解决.但是，如果这2个点分别在S1和S2中，则对于S1中任一点p，S2中最多只有n/2个点与它构成最接近点对的候选者，仍需做n^2/4次计算和比较才能确定S的最接近点对.因此,依此思路，合并步骤耗时为O（n^2）。

整个算法所需计算时间T（n）应满足:

　T（n）=2T（n/2）+O（n^2）。

它的解为T（n）=O（n^2），即与合并步骤的耗时同阶,这不比用穷举的方法好。

从解递归方程的套用公式法,我们看到问题出在合并步骤耗时太多.这启发我们把注意力放在合并步骤上。

     设S中的n个点为x轴上的n个实数x1，x2，。

。

xn.最接近点对即为这n个实数中相差最小的2个实数。

我们显然可以先将x1,x2，。

.,xn排好序，然后,用一次线性扫描就可以找出最接近点对.这种方法主要计算时间花在排序上,在排序算法已经证明,时间复杂度为O（nlogn）。

然而这种方法无法直接推广到二维的情形。

因此，对这种一维的简单情形，我们还是尝试用分治法来求解，并希望能推广到二维的情形.假设我们用x轴上某个点m将S划分为2个子集S1和S2，使得S1=｛x∈S｜x≤m｝；S2={x∈S|x>m｝.这样一来,对于所有p∈S1和q∈S2有p

如图所示。

    如果S的最接近点对是｛p3,q3}，即｜p3—q3｜〈d，则p3和q3两者与m的距离不超过d,即｜p3—m｜

同理,（m，m+d]中也至多包含S中的一个点。

由图可以看出,如果（m-d,m]中有S中的点，则此点就是S1中最大点。

同理，如果（m，m+d]中有S中的点,则此点就是S2中最小点.因此,我们用线性时间就能找到区间（m-d，m］和（m，m+d]中所有点，即p3和q3。

从而我们用线性时间就可以将S1的解和S2的解合并成为S的解。

也就是说，按这种分治策略，合并步可在O（n）时间内完成。

这样是否就可以得到一个有效的算法了呢？

还有一个问题需要认真考虑,即分割点m的选取，及S1和S2的划分.选取分割点m的一个基本要求是由此导出集合S的一个线性分割，即S=S1∪S2 ,S1∩S2=Φ,且S1=｛x|x≤m｝；S2={x｜x〉m｝。

容易看出,如果选取m=[max（S）+min（S）］/2，可以满足线性分割的要求.选取分割点后，再用O（n）时间即可将S划分成S1={x∈S｜x≤m｝和S2=｛x∈S｜x>m｝.然而，这样选取分割点m，有可能造成划分出的子集S1和S2的不平衡。

例如在最坏情况下,｜S1｜=1，|S2|=n—1，由此产生的分治法在最坏情况下所需的计算时间T（n）应满足递归方程:

  T（n）=T（n—1）+O（n）

   它的解是T（n）=O（n^2）。

这种效率降低的现象可以通过分治法中“平衡子问题”的方法加以解决。

即通过适当选择分割点m，使S1和S2中有大致相等个数的点。

自然地，我们会想到用S的n个点的坐标的中位数来作分割点.在选择算法中介绍的选取中位数的线性时间算法使我们可以在O（n）时间内确定一个平衡的分割点m.

     本程序确定平衡点采用m=[max（S）+min（S）］/2方法.如果需要利用中位数作分割点，看结合笔者博文《0005算法笔记——线性时间选择》改写。

   一维最接近临近点对问题程序清单如下：

［cpp］ viewplain copy

1.//2d10—1 一维最邻近点对问题  

2.#include "stdafx。

h”  

3.#include 〈ctime〉  

4.#include 〈iostream〉   

5.using namespace std;   

6.  

7.const int L=100；  

8.//点对结构体  

9.struct Pair  

10.｛  

11.    float d;//点对距离  

12.    float d1,d2；//点对坐标  

13.｝；  

14.float Random（）;  

15.int input（float s［]）；//构造S  

16.float Max（float s［]，int p,int q）;  

17.float Min（float s［］,int p，int q）；  

18.template 

19.void Swap（Type ＆x，Type &y）;  

20.template 

21.int Partition（Type s［]，Type x，int l，int r）;  

22.Pair Cpair（float s［]，int l，int r）;  

23.  

24.int main（）  

25.｛  

26.    srand（（unsigned）time（NULL））；  

27.    int m；  

28.    float s［L];  

29.    Pair d;  

30.    m=input（s）；  

31.    d=Cpair（s，0,m—1）；  

32.    cout〈〈endl〈〈”最近点对坐标为：

 （d1：

”〈

”〈

33.    cout<

 "〈〈d。

d<

34.    return 0;  

35.｝  

36.  

37.  

38.float Random（）  

39.｛  

40.    float result=rand（）%10000;  

41.     return result＊0。

01；  

42.}  

43.  

44.int input（float s［］）  

45.{  

46.    int length；  

47.    cout〈〈”输入点的数目:

 ”；  

48.    cin〉〉length；  

49.    cout<〈"点集在X轴上坐标为：

";  

50.    for（int i=0；i〈length；i++）  

51.    ｛  

52.        s[i]=Random（）；  

53.        cout<

54.    ｝  

55.      

56.    return length;  

57.｝  

58.  

59.  

60.float Max（float s［］,int l，int r）//返回s［]中的最大值  

61.｛  

62.    float s_max=s[l]；  

63.    for（int i=l+1;i〈=r；i++）  

64.        if（s_max〈s[i]）  

65.            s_max=s［i］；  

66.    return s_max;  

67.｝  

68.  

69.float Min（float s［],int l，int r）//返回s［］中的最小值  

70.{  

71.    float s_min=s[l］；  

72.    for（int i=l+1；i<=r;i++）   

73.        if（s_min〉s[i］）  

74.            s_min=s[i］；  

75.    return s_min;  

76.｝  

77.  

78.template 〈class Type〉  

79.void Swap（Type &x,Type ＆y）  

80.{  

81.    Type temp = x；  

82.    x = y；  

83.    y = temp;  

84.}  

85.  

86.template 〈class Type〉  

87.int Partition（Type s［]，Type x，int l，int r）  

88.｛  

89.    int i = l - 1，j = r + 1;  

90.  

91.    while（true）  

92.    ｛  

93.        while（s[++i］〈x ＆＆ i

94.        while（s［--j]>x）；  

95.        if（i>=j）  

96.        ｛  

97.            break；  

98.        }  

99.        Swap（s[i］，s[j］）；  

100.    ｝  

101.    return j；  

102.}  

103.  

104.//返回s[］中的具有最近距离的点对及其距离  

105.Pair Cpair（float s［］，int l，int r）  

106.｛  

107.    Pair min_d={99999,0，0｝;//最短距离  

108.  

109.    if（r—l〈1） return min_d；  

110.    float m1=Max（s，l,r）,m2=Min（s,l，r）；  

111.  

112.    float m=（m1+m2）/2；//找出点集中的中位数  

113.  

114.    //将点集中的各元素按与m的大小关系分组  

115.    int j = Partition（s,m,l,r）;  

116.  

117.    Pair d1=Cpair（s，l,j），d2=Cpair（s，j+1,r）；//递归  

118.    float p=Max（s,l，j），q=Min（s,j+1，r）;  

119.  

120.    //返回s[］中的具有最近距离的点对及其距离  

121.    if（d1。

d〈d2.d）  

122.    ｛  

123.        if（（q-p）〈d1。

d）  

124.        ｛  

125.            min_d。

d=（q-p）;  

126.            min_d。

d1=q;  

127.            min_d。

d2=p;  

128.            return min_d；  

129.        ｝  

130.        else return d1;  

131.    ｝  

132.    else  

133.    {  

134.        if（（q-p）〈d2。

d）  

135.        {  

136.            min_d。

d=（q—p）;  

137.            min_d.d1=q；  

138.            min_d。

d2=p;  

139.            return min_d；  

140.        ｝  

141.        else return d2;  

142.    ｝  

143.｝  

程序运行结果如下:

    该算法的分割步骤和合并步骤总共耗时O（n）。

因此，算法耗费的计算时间T（n）满足递归方程：

    解此递归方程可得T（n）=O（nlogn）.

     2、二维最接近点对问题

   将以上过程推广到二维最接近点对问题,设S中的点为平面上的点，它们都有2个坐标值x和y。

为了将平面上点集S线性分割为大小大致相等的2个子集S1和S2,我们选取一垂直线l：

x=m来作为分割直线。

其中m为S中各点x坐标的中位数。

由此将S分割为S1={p∈S|px≤m}和S2=｛p∈S｜px〉m｝.从而使S1和S2分别位于直线l的左侧和右侧，且S=S1∪S2.由于m是S中各点x坐标值的中位数,因此S1和S2中的点数大致相等。

递归地在S1和S2上解最接近点对问题，我们分别得到S1和S2中的最小距离d1和d2.现设d=min（d1，d2）.若S的最接近点对（p，q）之间的距离d（p，q）〈d则p和q必分属于S1和S2。

不妨设p∈S1,q∈S2.那么p和q距直线l的距离均小于d。

因此，我们若用P1和P2分别表示直线l的左边和右边的宽为d的2个垂直长条,则p∈S1，q∈S2，如图所示:

   距直线l的距离小于d的所有点

     在一维的情形,距分割点距离为d的2个区间（m—d,m］（m,m+d］中最多各有S中一个点。

因而这2点成为唯一的末检查过的最接近点对候选者。

二维的情形则要复杂些,此时，P1中所有点与P2中所有点构成的点对均为最接近点对的候选者。

在最坏情况下有n2/4对这样的候选者。

但是P1和P2中的点具有以下的稀疏性质，它使我们不必检查所有这n^2/4对候选者.考虑P1中任意一点p，它若与P2中的点q构成最接近点对的候选者，则必有d（p,q）〈d.满足这个条件的P2中的点有多少个呢?

容易看出这样的点一定落在一个d×2d的矩形R中,如下图所示：

包含点q的dX2d矩形R

   由d的意义可知P2中任何2个S中的点的距离都不小于d。

由此可以推出矩形R中最多只有6个S中的点。

事实上,我们可以将矩形R的长为2d的边3等分，将它的长为d的边2等分,由此导出6个（d/2）×（2d/3）的矩形。

如左图所示:

矩阵R中点的稀疏性

   若矩形R中有多于6个S中的点,则由鸽舍原理易知至少有一个δ×2δ的小矩形中有2个以上S中的点。

设u，v是这样2个点,它们位于同一小矩形中，则：

因此d（u,v）≤5d/6〈d.这与d的意义相矛盾.也就是说矩形R中最多只有6个S中的点。

图4（b）是矩形R中含有S中的6个点的极端情形。

由于这种稀疏性质，对于P1中任一点p，P2中最多只有6个点与它构成最接近点对的候选者。

因此，在分治法的合并步骤中，我们最多只需要检查6×n/2=3n对候选者，而不是n^2/4对候选者。

这是否就意味着我们可以在O（n）时间内完成分治法的合并步骤呢？

现在还不能作出这个结论,因为我们只知道对于P1中每个S1中的点p最多只需要检查P2中的6个点,但是我们并不确切地知道要检查哪6个点。

为了解决这个问题，我们可以将p和P2中所有S2的点投影到垂直线l上。

由于能与p点一起构成最接近点对候选者的S2中点一定在矩形R中,所以它们在直线l上的投影点距p在l上投影点的距离小于d。

由上面的分析可知,这种投影点最多只有6个。

因此,若将P1和P2中所有S的点按其y坐标排好序，则对P1中所有点p，对排好序的点列作一次扫描，就可以找出所有最接近点对的候选者,对P1中每一点最多只要检查P2中排好序的相继6个点.

程序清单如下:

[cpp］ viewplain copy

1.//2d10—2 二维最邻近点对问题  

2.＃include "stdafx。

h"  

3.＃include〈time.h〉  

4.＃include

5.#include

6.  

7.using namespace std;  

8.const int M=50；  

9.  

10.//用类PointX和PointY表示依x坐标和y坐标排好序的点  

11.class PointX {  

12.    public:

13.        int operator〈=（PointX a）const  

14.        ｛ return （x<=a。

x）; }  

15.        int ID； //点编号  

16.        float x，y； //点坐标   

17.｝；  

18.  

19.class PointY ｛   

20.    public：

21.        int operator〈=（PointY a）const  

22.        ｛ return（y〈=a.y）； ｝  

23.        int p； //同一点在数组x中的坐标   

24.        float x,y; //点坐标  

25.｝；  

26.  

27.float Random（）；  

28.template 

29.float dis（const Type&u，const Type&v）；   

30.  

31.bool Cpair2（PointX X[］， int n，PointX＆ a,PointX＆ b， float＆ d）;  

32.void closest（PointX X［］,PointY Y[]，PointY Z［］, int l, int r，PointX＆ a，PointX＆ b,float＆ d）；  

33.  

34.template 〈typename Type>   

35.void Copy（Type a［],Type b［］， int left，int right）；  

36.  

37.template 

38.void Merge（Type c[］，Type d[］，int l,int m，int r）；  

39.  

40.template 〈class Type〉  

41.void MergeSort（Type a［］，Type b［］,int left,int right）；  

42.  

43.int main（）  

44.｛   

45.    srand（（unsigned）time（NULL））;  

46.    int length；  

47.  

48.    cout〈<"请输入点对数:

”；  

49.    cin〉〉length;  

50.  

51.    PointX X［M］；  

52.    cout〈<”随机生成的二维点对为：

"〈

53.  

54.    for（int i=0；i

55.    ｛  

56.        X[i]。

ID=i；  

57.        X［i］.x=Random（）；  

58.        X［i].y=Random（）；  

59.        cout〈〈”（"〈〈X[i］。

x〈〈”，”〈〈X［i］.y<〈”） ";  

60.    }  

61.  

62.    PointX a；   

63.    PointX b;   

64.    float d；  

65.  

66.    Cpair2（X，length，a，b,d）;   

67.  

68.    cout〈〈endl；  

69.    cout〈〈”最邻近点对为：

（”<

x<〈",”〈

x〈<",”〈〈b.y〈〈"） ”<〈endl；  

70.    cout〈<”最邻近距离为：

 ”〈〈d〈〈endl；  

71.  

72.    return 0；  

73.｝  

74.  

75.float Random（）  

76.｛  

77.    float result=rand（）%10000；  

78.    return result＊0。

01；  

79.｝  

80.  

81.//平面上任意两点u和v之间的距离可计算如下  

82.template 〈class Type〉  

83.inline float dis（const Type& u,const Type＆ v）  

84.｛  

85.    float dx=u。

x—v.x；  

86.    float dy=u。

y—v。

y;   

87.    return sqrt（dx*dx+dy＊dy）；   

88.}  

89.  

90.bool Cpair2（PointX X[], int n，PointX＆ a,PointX& b，float& d）  

91.｛  

92.    if（n<2） return false;  

93.  

94.    PointX＊ tmpX = new PointX［n］;  

95.    MergeSort（X,tmpX,0,n-1）;  

96.  

97.    PointY＊ Y=new PointY［n];   

98.    for（int i=0；i〈n；i++） //将数组X中的点复制到数组Y中  

99.    {   

100.        Y［i］。

p=i;  

101.        Y［i].x=X［i］.x；  

102.        Y［i］。

y=X［i］.y；  

103.    ｝   

104.  

105.    PointY* tmpY = new PointY[n］;  

106.    MergeSort（Y,tmpY，0，n—1）;  

107.  

108.    PointY* Z=new PointY［n];  

109.    closest（X,Y，Z,0,n-1，a,b，d）;   

110.  

111.    delete []Y；   

112.    delete ［]Z;  

113.    delete ［]tmpX；  

114.    del