刷题交流｜基于哈希表+二分答案的平面最近点对解法

3765

2022.01.09

发布于未知归属地

AtCoder ABC234Ex 是一道很有意思的题目。它要求找出平面上所有距离不超过 $K$ 的点对（题目保证这样的点对不超过 $4 \times 1 0^{5}$ 个）。

一种基于数据结构的解法是使用 K-D 树，但官方解答给出了一种更为优雅的解法，并给出了严格的复杂度证明。

简单来说，我们可以将平面划分成大小为 $K \times K$ 的网格，这样每个点必然落在一个网格中。容易证明，与一个点距离不超过 $K$ 的点必然落在该点所在网格及其周围八个网格的某一个之中。我们可以枚举这九个网格中的点，看它们是否能与当前点构成一个符合要求的点对。官方解答证明了，这一算法的时间复杂度为 $O (N + M)$ ，其中 $M$ 为符合要求的点对的总数。

由这一问题，很容易联想到经典的平面最近点对问题。我们有没有可能借鉴这题的方法，来解决平面最近点对问题呢？

答案是肯定的。我们可以使用二分答案的方法。对于某一个距离，如果存在距离不超过它的点对，就说明最近点对的距离小于它；否则说明最近点对的距离大于它。

假设当前需要检验的距离为 $K$ ，我们可以用上面的方法进行枚举。这时有两种情况：

符合要求的点对数不超过 $N$ ：此时我们可以直接使用上面的结论，得出检验的复杂度为 $O (N)$ 。
符合要求的点对数超过 $N$ ：注意到我们并不需要找出所有的点对，而只需要找到一个点对。下面考虑最坏的情况下我们要进行多少次枚举。在还没有找到符合要求的点对的情况下，每个格子里最多有四个点。这是因为我们可以将每个格子划分成四个小格，则任意小格中任取两个点距离必然不超过 $K$ ，根据抽屉原理可知我们最多在这四个小格中放入四个点。从而，对于任意一个点，我们最多进行 $9 * 4 - 1 = 35$ 次枚举。因此，总的时间复杂度也是 $O (N)$ 。

因此，总的时间复杂度就是 $O (N lo g C)$ ，其中 $C$ 取决于题目要求的精度。

与传统解法相比，这一方法的常数要大一些，但在理解和实现上有它的优势。

例题

SPOJ - CLOPPAIR

下面的解法利用了本题中最近点对唯一的条件。

参考代码（Rust）

use std::collections::HashMap;
use std::io;

fn read_line() -> String {
    let mut line = String::new();
    io::stdin().read_line(&mut line).unwrap();
    return line.trim().to_string();
}

fn read_ints() -> Vec<i64> {
    let line = read_line();
    return line
        .split(" ")
        .filter(|&s| s.len() > 0)
        .map(|s| s.parse::<i64>().unwrap())
        .collect();
}

fn main() {
    let n: usize = read_ints()[0] as usize;
    let mut points: Vec<(i64, i64)> = vec![];
    for _ in 0..n {
        let p = read_ints();
        points.push((p[0], p[1]));
    }

    let mut lo: f64 = 0.0;
    let mut hi: f64 = 4e6;
    let mut buckets: HashMap<(i64, i64), Vec<usize>> = HashMap::new();
    let mut ans: Vec<(usize, usize)> = vec![];
    let neighbors = [
        (0, 0),
        (-1, 0),
        (0, -1),
        (0, 1),
        (1, 0),
        (-1, -1),
        (-1, 1),
        (1, -1),
        (1, 1),
    ];

    while hi - lo >= 1e-8 {
        buckets.clear();
        ans.clear();
        let mid = (lo + hi) / 2.0;
        for i in 0..n {
            let x = (points[i].0 as f64 / mid).floor() as i64;
            let y = (points[i].1 as f64 / mid).floor() as i64;
            for &(dx, dy) in neighbors.iter() {
                let xp = x + dx;
                let yp = y + dy;

                if !buckets.contains_key(&(xp, yp)) {
                    continue;
                }

                for &j in buckets.get(&(xp, yp)).unwrap() {
                    if ((points[i].0 - points[j].0) * (points[i].0 - points[j].0)
                        + (points[i].1 - points[j].1) * (points[i].1 - points[j].1))
                        as f64
                        <= mid * mid
                    {
                        ans.push((j, i));
                        if ans.len() >= 2 {
                            break;
                        }
                    }
                }

                if ans.len() >= 2 {
                    break;
                }
            }

            if ans.len() >= 2 {
                break;
            }

            (*buckets.entry((x, y)).or_default()).push(i);
        }

        if ans.len() == 1 {
            let (i, j) = ans[0];
            println!(
                "{} {} {:.6}",
                i,
                j,
                (((points[i].0 - points[j].0) * (points[i].0 - points[j].0)
                    + (points[i].1 - points[j].1) * (points[i].1 - points[j].1))
                    as f64)
                    .sqrt()
            );
            std::process::exit(0);
        }

        if ans.is_empty() {
            lo = mid;
        } else {
            hi = mid;
        }
    }
}