最接近点对问题

一维最接近点问题

#!/usr/bin/env python
# -*- coding: utf-8 -*-
# @Time    : 2023/1/5 20:36
# @Author  : Lili
# @File    : 一维最接近点对问题.py
# @Description : 递归与分治，来探索二维中最接近点对问题的求解

from utils import lp_wrapper

# 统计运行时间的注解
@lp_wrapper()
def cpair1(S: list) -> (bool, float):
    """
    递归分治求解一维最接近点对问题

    :param S: 待求解点集合
    :param d: 最接近两点的距离
    :return: 成功与否
    """
    if len(S) < 3:
        d = float("inf")
        return False, d
    # mid: S中各点坐标的中位数位置
    mid = int(len(S)/2)
    S.sort()
    # print("中位数mid: ", mid)
    S1 = S[:mid]
    S2 = S[mid:]
    # print(S1, S2)
    d1 = cpair1(S1)[1]
    d2 = cpair1(S2)[1]
    # print("S1, d1", S1, d1)
    # print("S2, d2", S2, d2)
    return True, min(d1, d2, S[mid + 1] - S[mid])

二维最接近点问题

#!/usr/bin/env python
# -*- coding: utf-8 -*-
# @Time    : 2023/1/5 20:36
# @Author  : Lili
# @File    : 二维最接近点对问题.py, closest中对y的重构出错
# @Description :
from math import hypot
import sys
from random import randint

sys.setrecursionlimit(2400)

ID, X, Y = 0, 1, 2


def getx(point):
    return point[X]


def gety(point):
    return point[Y]


def distance(pointa, pointb):
    return hypot(pointa[X]-pointb[X], pointa[Y]-pointb[Y])


def closest(point_x, point_y, point_z, l, r):
    """
    计算最接近点对

    :param point_x: 原列表，按x坐标升序
    :param point_y: 按y坐标升序
    :param point_z: 同point_y用于存放分治后点集以及存放筛选后d矩形内的点集
    :param l: 起始位置
    :param r: 终止位置
    :return: (point_a, point_b, min_distance)
    """
    # 包含2点的情况
    if r-l == 1:
        return point_x[l], point_x[r], distance(point_x[l], point_x[r])
    elif r-l == 2:
        d1 = distance(point_x[l], point_x[l+1])
        d2 = distance(point_x[r], point_x[l+1])
        d3 = distance(point_x[l], point_x[r])
        if d1 <= d2 and d1 <= d3:
            return point_x[l], point_x[l+1], d1
        elif d2 <= d3:
            return point_x[r], point_x[l + 1], d2
        else:
            return point_x[l], point_x[r], d3
    # 分治法, 按照x坐标等分为左右两部分
    py_temp = point_y[l:r+1]
    pz_temp = point_z[l:r+1]
    mid = int((l+r)/2)
    f, g = l, mid + 1

    for py in point_y[l:r+1]:
        if py[ID] > mid:
            point_z[g] = py
            g += 1
        else:
            point_z[f] = py
            f += 1
    point_a, point_b, min_distance = closest(point_x, point_z, point_y, l, mid)
    ar, br, dr = closest(point_x, point_z, point_y, mid+1, r)
    if min_distance > dr:
        point_a, point_b = ar, br
        min_distance = dr
    # merge(z,y,l,m,r)
    point_y[l:r+1] = point_z[l:r+1]
    # point_y[l:r+1] = py_temp
    # point_z[l:r+1] = pz_temp
    k = l
    for i in range(l, r+1):
        if abs(point_y[mid][X]-point_y[i][X]) < min_distance:
            point_z[k] = point_y[i]
            k += 1
    for i in range(l, k):
        for j in range(i+1, k):
            if abs(point_z[j][Y]-point_z[i][Y]) >= min_distance:
                break
            dp = distance(point_z[i], point_z[j])
            if dp < min_distance:
                min_distance = dp
                point_a = point_x[point_z[i][ID]]
                point_b = point_x[point_z[j][ID]]
    return point_a, point_b, min_distance


def cpair2(point: [[]]):
    """
    递归与分治，求解二维最接近点问题算法

    :param point: [[id, x, y], ]点集列表，id: int
    :return: (bool, min_distance), 返回执行结果和最小距离值
    """
    min_dis = float("inf")
    if len(point) < 2:
        return False, None, None, min_dis
    point.sort(key=getx)
    point_x = point[:]
    i = 0
    for p in point_x:
        p[ID] = i
        i += 1
    point_y = point_x[:]
    point_y.sort(key=gety)
    point_z = point_y[:]
    point_a, point_b, min_dis = closest(point_x, point_y, point_z, 0, len(point)-1)
    del point_y
    del point_z
    return True, point_a, point_b, min_dis


# 测试
LEN = 20
BEGIN, END = 0, int(LEN*10)

points = []
for i in range(0, LEN):
    p = [i, randint(BEGIN, END), randint(BEGIN, END)]
    points.append(p)

print(points)
d, p1, p2 = float("inf"), 0, 0
for i in range(0, LEN):
    for j in range(i+1, LEN):
        dis = hypot(points[i][1]-points[j][1], points[i][2]-points[j][2])
        if dis < d:
            d = dis
            p1 = i
            p2 = j

print(d, points[p1], points[p2])
t = cpair2(points)[1:]
print(t[2], t[0], t[1])

# 错误项
# [[0, 69, 109], [1, 75, 128], [2, 87, 185], [3, 121, 107], [4, 152, 158], [5, 103, 42], [6, 145, 110], [7, 21, 122], [8, 176, 33], [9, 144, 121], [10, 194, 74], [11, 193, 22], [12, 110, 89], [13, 32, 61], [14, 135, 85], [15, 36, 5], [16, 138, 170], [17, 121, 141], [18, 55, 137], [19, 73, 15]]
# 11.045361017187261 [6, 145, 110] [9, 144, 121]
# 21.095023109728988 [9, 110, 89] [10, 121, 107]

# [[0, 59, 135], [1, 49, 96], [2, 16, 43], [3, 143, 131], [4, 192, 94], [5, 126, 82], [6, 66, 154], [7, 43, 154], [8, 71, 186], [9, 75, 173], [10, 161, 67], [11, 23, 54], [12, 6, 32], [13, 87, 159], [14, 124, 84], [15, 90, 128], [16, 36, 49], [17, 128, 156], [18, 70, 179], [19, 63, 75]]
# 2.8284271247461903 [5, 126, 82] [14, 124, 84]
# 7.0710678118654755 [9, 70, 179] [10, 71, 186]

# [[0, 171, 55], [1, 11, 133], [2, 6, 62], [3, 164, 193], [4, 127, 180], [5, 33, 152], [6, 188, 65], [7, 5, 78], [8, 161, 148], [9, 47, 16], [10, 68, 102], [11, 20, 62], [12, 78, 8], [13, 155, 182], [14, 1, 159], [15, 50, 119], [16, 187, 21], [17, 134, 80], [18, 174, 72], [19, 171, 105]]
# 14.0 [2, 6, 62] [11, 20, 62]
# 16.0312195418814 [2, 6, 62] [1, 5, 78]

# [[0, 34, 14], [1, 172, 198], [2, 150, 102], [3, 126, 15], [4, 66, 79], [5, 165, 97], [6, 38, 75], [7, 195, 158], [8, 195, 173], [9, 38, 91], [10, 146, 164], [11, 163, 84], [12, 19, 64], [13, 103, 139], [14, 116, 7], [15, 20, 3], [16, 110, 196], [17, 108, 130], [18, 154, 160], [19, 109, 6]]
# 7.0710678118654755 [14, 116, 7] [19, 109, 6]
# 10.295630140987 [7, 108, 130] [6, 103, 139]

二维算出来有问题，line70: # merge(z,y,l,m,r)这块重构point_y没看懂，《计算机算法设计与分析（王晓东）》p35