【计算几何 06】计算几何的终极奥义

计算几何中由两个经典方法，也被称为最终奥义（其实是没办法的时候才能使用2333）——枚举和分治。最终奥义一般是在构成几何点数较少和其他算法无法正确解决的时候使用😀

枚举和计算几何

先引入一道经典例题：

[caioj 1211]统计正方形

题目描述

【题目描述】

给定平面上N个点,你需要计算以其中4个点为顶点的正方形的个数。注意这里的正方形边不一定需要和坐标轴平行。

【输入文件】

第一行一个数N,以下N个点的坐标。

【输出文件】

一个数正方形的个数。

【样例输入】

7
0 0
0 1
1 0
1 1
1 2
2 1
2 2

【样例输出】

3

【数据规模】
对于$20％$的数据,满足 $1\le N\le 20$；
对于$100％$的数据,满足 $1\le N\le 500,-50\le X[i],Y[i]\le 50$,点不会重叠。

这道题从算法角度看是使用不了叉积和旋转卡壳的，但注意到100%的数据范围：$N_{max}=500$,这样我们就可以考虑直接枚举了。

我们知道如果知道了正方形的任意两点的坐标，就可以推出其余两个点的坐标。为了方便起见，我们枚举正方形中靠下的两个点分别设为 $I(x_i,y_i)\mathrm{和}J(x_j,y_j)$ （规定其中 $J$的横坐标$x_j$ 大于$I$ 的横坐标$x_i$）。

然后我们令$X=x_j-x_i；Y=y_j-y_i$；
接着我们就可以推出剩下两点坐标：
$O(x_i-Y,y_i+x)，P(x_j-y,y_j+x)$;只需要用数组判断$O，P$两点是否存在即可，若存在则ans++，否则继续枚举下一个$I，J$。(如下图)

有一点需要注意的是，由于坐标可能为负数，为了方便后面的判断，需要让每个点的坐标加上一个比较大的正数让坐标全部变为正数，AC代码如下：

#include <bits/stdc++.h>
using namespace std;
typedef long long ll;
const int N = 1e3 + 50;
struct node {
    int x, y;
};
vector<node> a(N);
bool f[N][N];//标记点的存在
int n, x, y, ans;
int main() {
    // freopen("in.txt", "r", stdin);
    ios::sync_with_stdio(false), cin.tie(0), cout.tie(0);
    cin >> n;
    for (int i = 1; i <= n; ++i) {
        cin >> x >> y;
        a[i].x = x + 500, a[i].y = y + 500;  //避免负数的影响

        f[a[i].x][a[i].y] = true;
    }

    ans = 0;
    //开始枚举,先序定位 I，J点。由坐标推导正方形剩下的两个点。如果存在ans++
    for (int i = 1; i <= n; ++i)
        for (int j = 1; j <= n; ++j)
            if (i != j) {
                if (a[i].x < a[j].x && a[i].y <= a[j].y) {
                    int x = a[j].x - a[i].x;
                    int y = a[j].y - a[i].y;

                    if (f[a[i].x - y][a[i].y + x] && f[a[j].x - y][a[j].y + x])
                        ans++;
                }
            }

    cout << ans << endl;
}

分治与计算几何

在【计算几何 03】旋转卡壳算法的那篇博客中，我们使用了旋转卡壳求出了平面内最远点对，但现在重新来思考一下：如何求最近点对呢？

尽管我们能通过穷举$O(n^2)$ 出答案，但要想更优地解决这个问题，不妨先把问题简化，看看在一维上这个问题的求解方法。

一维内的最近点对

此时，$n$个点退化为 $x\mathrm{轴}$ 上的 $n$个实数 $x_1，x_2，...，x_n$ 。最接近点对即为这 $n$ 个实数中相差最小的两个实数。显然可以先将点排好序，然后用 $O(n)$ 线性扫描就可以了。但为了便于推广到二维的情形，尝试用分治法解决这个问题。

首先肯定要先对点进行排序，接着假设我们用 $m$点将 $S(x\mathrm{轴}\mathrm{点}\mathrm{序}\mathrm{列})$分为 $S_1$和 $S_2$ 两个集合，这样一来，对于所有的 $p(S_1\mathrm{中}\mathrm{的}\mathrm{点})\mathrm{和}q(S_2\mathrm{中}\mathrm{的}\mathrm{点})$ ，有 $p$ 小于 $q$。

根据分治思想：再对$S_1,S_2$重复操作，递归地在$S_1$和$S_2$ 上找出其最近点对$\{p1,p2\}\mathrm{和}\{q1,q2\}$，并设$d=min(|p1-p2|,|q1-q2|)$

由此易知，S中最近点对可能是$\{p1,p2\}$ ，可能是 $\{q1,q2\}$ ，也可能是某个 $\{q3,p3\}$ ，如上图所示。

如果此时最近点对是 $\{q3,p3\}$，即$|p3-q3|$小于 $d$ ，则 $p3$ 和$q3$ 两者与 $m$的距离都不超过$d$，且$p3$和$q3$一定分别出现在区间$(m-d,d]\mathrm{和}(d,m+d]$ 。此时，一维情形下的最近点对时间复杂度为$O(nlogn)$^[1]。

二维内的最近点对^[2]

现在可以对刚刚一维的做法进行推广，我们可以首先先对 $x$ 坐标进行排序，取出中点，然后对于左右两边的区域 $S_1，S_2$，重复上述过程，递归地求出左右两边分别的最近点对的距离 $d1，d2$，并且令$d=min(d1,d2)$ 。

@ 图片来自蒟蒻dalao

对于最多六个节点的证明

我们继续用回上一幅图进行证明；

@ 图片来自蒟蒻dalao

我们知道，左右两边分别的最近点对的距离$d1，d2$，并且令 $d=min(d1,d2)$，也就是右边这个$d\ast 2d$的矩形中，任意两点的之间的距离都要大于等于d。那么我们想在这个$d\ast 2d$的矩形中放尽量多的点，但彼此之间距离都要大于等于$d$，怎么放呢？

@ 图片来自蒟蒻dalao

聪明的你一定想到了，按照上图的红点位置放，是最优的策略了。但此时最多也就只能放下六个点。所以原结论是正确的。

再来道模板题

[caioj 1206]最近点对的距离

【题意】

给出n个点的坐标，求最近两点间的距离。

【输入格式】

第一行一个整数n(2 ≤ n ≤ 50000)。

下来n行，每行两个实数 x 和 y 表示点坐标。

【输出格式】

一行一个实数，表示最近两点间的距离（保留4位小数）。
【样例输入】

5
0 0
0 5
5 0
5 5
2 0

【样例输出】

2.0000

直接贴模板了：

#include <bits/stdc++.h>
using namespace std;
typedef long long ll;
const int N = 5e4 + 100;
const int inf = 21474836;
struct node {
    double x, y;
} p[N], f[N];
// vector<node> p(N), f(N);
double ans, Rans, Lans;
int n;

double dis(node x, node y) {
    return (sqrt((x.x - y.x) * (x.x - y.x) + (x.y - y.y) * (x.y - y.y)));
}
// 函数重构
// double min(double x, double y) {
//     if (x < y) return x;
//     return y;
// }
bool cmp(node x, node y) {
    if (x.x < y.x) return true;
    return false;
}

bool cmpp(node x, node y) {
    if (x.y < y.y) return true;
    return false;
}

double close(int left, int right) {
    if (left == right) return inf;
    if (left + 1 == right) return (dis(p[left], p[right]));
    int mid = (left + right) >> 1;
    //利用递归分别求出左右两边区域的最近距离
    double d1 = close(left, mid);
    double d2 = close(mid + 1, right);
    double d = min(d1, d2);
    int len = 0;

    for (int i = left; i <= right; ++i)
        if (fabs(p[i].x - p[mid].x) <= d) f[++len] = p[i];

    //选出mid左右d的条状区域的所有点并且按照y轴坐标排序
    sort(f + 1, f + len + 1, cmpp);

    for (int i = 1; i < len; i++) {
        for (int j = i + 1; j <= len; j++)
            if (f[j].y - f[i].y < d) {
                double d3 = dis(f[i], f[j]);
                if (d3 < d) d = d3;
            } else
                break;  //利用鹊巢原理进行优化
    }
    //判断是否有小于d的值存在
    return d;
}

int main() {
    // freopen("in.txt", "r", stdin);
    ios::sync_with_stdio(false), cin.tie(0), cout.tie(0);
    scanf("%d", &n);
    for (int i = 1; i <= n; i++) scanf("%lf%lf", &p[i].x, &p[i].y);
    sort(p + 1, p + n + 1, cmp);  //对于x轴坐标排序

    ans = close(1, n);  //计算答案

    printf("%0.4lf", ans);
}

结语

通过这篇博客相信你一定已经学会了计算几何的最终奥义！

同时如果本文中有什么错误的地方还请大神们指出，我会第一时间修改。

参考

平面内有N个点，如何快速求出距离最近的点对？

---

---

1. 关于时间复杂度的详细证明可以在OI wiki上查看click here ↩︎

2. 推广的思路来自蒟蒻dalao ↩︎