Python实现最近点对求解

发表于 | 分类于 |

数学建模中遇到求最近点对的问题，按理说我们使用朴素法进行暴力求解答案也是可行的，但是由于数据过于庞大，我们最初使用朴素法在计算机上跑了6个小时都没有得到问题的解，最终只得作罢。后来我们利用分治法解决了最近点对问题，具体思路参看这里。

分治法

分治法解决最近点对问题的算法过程读者可以参见上文给出的参考博客。分治法的思想就是先将问题分解成一个个小的问题，最后将小问题的答案合并得到问题的解。这其中关键步骤是分解与合并，具体实现方法是使用递归，核心思想是空间换时间。

代码实现

上面给出的教程中是使用C++实现的，我这里给出Python实现。

由于我这里具体的问题是给出经纬度找出最近点对，所以我采用的距离为改进的球面余弦距离。

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
# -*- coding: utf-8 -*-
from collections import Counter
from math import sqrt,acos, sin, cos, radians
def nearest_dot(s):
    mid = int(len(s)/2)
    left = s[0:mid]
    right = s[mid:]
    mid_x = (left[-1][0]+right[0][0])/2.0
    if len(left) > 2:   lmin = nearest_dot(left)    #左侧部分最近点对
    else:   lmin = left
    if len(right) > 2:   rmin = nearest_dot(right)   #右侧部分最近点对
    else:   rmin = right
    if len(lmin) >1: dis_l = get_distance(lmin)
    else: dis_l = float("inf")
    if len(rmin) >1: dis_r = get_distance(rmin)
    else: dis_r = float("inf")
    d = min(dis_l, dis_r)   #最近点对距离
    mid_min=[]
    for i in left:
        if mid_x-i[0]<=d :   #如果左侧部分与中间线的距离<=d
            for j in right:
                if abs(i[0]-j[0])<=d and abs(i[1]-j[1])<=d:     #如果右侧部分点在i点的(d,2d)之间
                    if get_distance((i,j))<=d: mid_min.append([i,j])   #ij两点的间距若小于d则加入队列
    if mid_min:
        dic=[]
        for i in mid_min:
            dic.append({get_distance(i):i})
        dic.sort(key=lambda x: x.keys())
        return list(dic[0].values())[0]
    elif dis_l>dis_r:
        return rmin
    else:
        return lmin
        
# 求点对的距离
def get_distance(m):
    dx = m[0][1] - m[1][1]  # 经度差值
    dy = m[0][0] - m[1][0]  # 纬度差值
    # b = (lat1 + lat2) / 2.0  # 平均纬度
    # Lx = (a[3] * b*b*b  + a[2]* b*b  + a[1] * b + a[0]) * radians(dx) * 6367000.0  # 东西距离(单位：米)
    Lx = cos(dx) * radians(dx) * 6367000.0  # 东西距离(单位：米)
    Ly = 6367000.0 * radians(dy)  # 南北距离
    return sqrt(Lx * Lx + Ly * Ly)
def divide_conquer(s):
    s.sort()
    nearest_dots = nearest_dot(s)
    return nearest_dots
显示 Gitment 评论