如何找到矩形内到点集的最大距离最小的点

阿华AIGC实验室

2026-5-11

带矩形约束的最小最大距离点求解方案

这个问题本质是带矩形约束的最小最大距离优化问题——我们要在指定矩形区域内找一点 ( P(x,y) )，使得 ( P ) 到点集所有点的最远欧氏距离尽可能小。用数学语言描述就是：
[
\min_{P \in \text{矩形}} \left( \max_{p_i \in pts} | P - p_i |_2 \right)
]

核心思路

无约束情况下，这个最优点是点集的**最小包围圆（Minimal Enclosing Circle, MEC）**的圆心——因为最小包围圆的圆心是到所有点的最大距离最小的点。但现在多了矩形约束，所以要分两种情况处理：

如果MEC圆心在矩形内部，那它直接就是我们要的最优解；
如果MEC圆心在矩形外部，我们需要把最优解限制在矩形的边界上，在四条边（线段）上分别求解最优点，再取全局最小值。

具体实现步骤

1. 计算点集的最小包围圆（MEC）

最小包围圆可以用Welzl算法高效求解，这是一个递归随机算法，核心逻辑是：

递归处理点集，逐步构建包含所有点的最小圆；
圆最多由3个点确定（不共线的三点），如果点集为空或已找到3个边界点，直接构造对应圆；
随机选取点，判断是否在当前圆内，不在则加入边界点列表重新计算。

2. 判断MEC圆心是否在矩形内

假设矩形范围是 ( [x_{min}, x_{max}] \times [y_{min}, y_{max}] )，判断圆心 ( (cx, cy) ) 是否满足：
[
x_{min} \leq cx \leq x_{max} \quad \text{且} \quad y_{min} \leq cy \leq y_{max}
]
如果满足，这个圆心就是最优解。

3. 当MEC圆心在矩形外时，在矩形边界找最优点

此时最优解必然在矩形的四条边上（因为矩形是凸集，约束下的最优解会落在边界上）。我们需要对每条边单独求解：

对于一条线段（比如水平底边 ( y = y_{min}, x \in [x_{min}, x_{max}] )），问题转化为单变量优化：找 ( x ) 使得 ( \max_{p_i} \sqrt{(x - p_{i,x})^2 + (y_{min} - p_{i,y})^2} ) 最小。

由于这个目标函数是凸函数（多个二次函数的上包络），可以用三分法高效找到最小值点：

不断将当前区间三等分，比较两个三等分点的目标值，舍弃不可能包含最小值的区间；
重复迭代直到区间长度小于精度阈值（比如 ( 1e-6 )），或者迭代固定次数（比如100次）保证收敛。

同样的方法适用于矩形的顶边、左边、右边。

4. 筛选全局最优解

收集所有候选点：

MEC圆心（如果在矩形内）；
四条边上用三分法找到的最优点；
矩形的四个顶点（极端情况下最优解可能在顶点，比如所有给定点都在矩形外的同一侧）。

计算每个候选点的Dist值，选择Dist值最小的点作为最终结果。

代码示例（伪代码）

def minimal_enclosing_circle(points):
    # 实现Welzl算法，返回(圆心x, 圆心y, 半径)
    # 处理边界情况：空点集、单点、两点、三点
    pass

def is_point_in_rect(cx, cy, x_min, x_max, y_min, y_max):
    return x_min <= cx <= x_max and y_min <= cy <= y_max

def optimize_on_segment(seg_x_start, seg_x_end, seg_y_start, seg_y_end, pts):
    # 用三分法优化线段上的点，返回最优(x,y)
    def max_dist_sq(x, y):
        max_dsq = 0.0
        for (px, py) in pts:
            dsq = (x - px)**2 + (y - py)**2
            if dsq > max_dsq:
                max_dsq = dsq
        return max_dsq

    # 处理水平线段（y固定）
    if seg_y_start == seg_y_end:
        left, right = seg_x_start, seg_x_end
        for _ in range(100):
            if right - left < 1e-8:
                break
            m1 = left + (right - left)/3
            m2 = right - (right - left)/3
            d1 = max_dist_sq(m1, seg_y_start)
            d2 = max_dist_sq(m2, seg_y_start)
            if d1 < d2:
                right = m2
            else:
                left = m1
        best_x = (left + right)/2
        return (best_x, seg_y_start)
    # 处理垂直线段（x固定）
    else:
        top, bottom = seg_y_start, seg_y_end
        for _ in range(100):
            if bottom - top < 1e-8:
                break
            m1 = top + (bottom - top)/3
            m2 = bottom - (bottom - top)/3
            d1 = max_dist_sq(seg_x_start, m1)
            d2 = max_dist_sq(seg_x_start, m2)
            if d1 < d2:
                bottom = m2
            else:
                top = m1
        best_y = (top + bottom)/2
        return (seg_x_start, best_y)

def find_optimal_point(x_min, x_max, y_min, y_max, pts):
    cx, cy, r = minimal_enclosing_circle(pts)
    candidates = []
    
    # 加入MEC圆心（如果在矩形内）
    if is_point_in_rect(cx, cy, x_min, x_max, y_min, y_max):
        candidates.append((cx, cy))
    
    # 加入四条边的最优点
    candidates.append(optimize_on_segment(x_min, x_max, y_min, y_min, pts))
    candidates.append(optimize_on_segment(x_min, x_max, y_max, y_max, pts))
    candidates.append(optimize_on_segment(x_min, x_min, y_min, y_max, pts))
    candidates.append(optimize_on_segment(x_max, x_max, y_min, y_max, pts))
    
    # 加入四个顶点
    candidates.extend([
        (x_min, y_min), (x_min, y_max),
        (x_max, y_min), (x_max, y_max)
    ])
    
    # 找到Dist最小的点
    min_dist = float('inf')
    best_p = None
    for p in candidates:
        current_max = 0.0
        for pt in pts:
            dx = pt[0] - p[0]
            dy = pt[1] - p[1]
            dist = (dx**2 + dy**2)**0.5
            if dist > current_max:
                current_max = dist
        if current_max < min_dist:
            min_dist = current_max
            best_p = p
    
    return best_p, min_dist