题解 P1429 【平面最近点对（加强版）】
syksykCCC
·
·
题解
UPD：谢谢大家对之前题解的支持，这里我重新更换一下以前难看的图片，同时补补锅，回答一下疑点。
UPD2：根据新的题解规范修改了文章。另外大家喜欢的话不妨点点赞支持支持。
UPD3：应评论要求补充了时间复杂度分析。另外发现题解的图片中 \rm P_9 和 \rm P_{10} 全部标反了，但因为不影响算法理解，且原图源丢失，暂时请大家阅读时自行对换过来。
首先感谢：
- @plane 和他的题解；
- @frankchenfu 和他的题解。
下面我贴上 @frankchenfu 的代码，我将给出图示模拟，希望能帮助看不懂的同学理解
#include <cstdio>
#include <algorithm>
#include <cmath>
using namespace std;
const int maxn = 1000001;
const int INF = 2 << 20;
int n, temp[maxn];
struct Point
{
double x, y;
} S[maxn];
bool cmp(const Point &a, const Point &b)
{
if(a.x == b.x) return a.y < b.y;
else return a.x < b.x;
}
bool cmps(const int &a, const int &b) { return S[a].y < S[b].y; }
double min(double a, double b) { return a < b ? a : b; }
double dist(int i, int j)
{
double x = (S[i].x - S[j].x) * (S[i].x - S[j].x);
double y = (S[i].y - S[j].y) * (S[i].y - S[j].y);
return sqrt(x + y);
}
double merge(int left, int right)
{
double d = INF;
if(left == right) return d;
if(left + 1 == right) return dist(left, right);
int mid = left + right >> 1;
double d1 = merge(left, mid);
double d2 = merge(mid + 1, right);
d = min(d1, d2);
int i, j, k = 0;
for(i = left; i <= right; i++)
if(fabs(S[mid].x - S[i].x) < d) // 这里不太一样
temp[k++] = i;
sort(temp, temp + k, cmps);
for(i = 0; i < k; i++)
for(j = i + 1; j < k && S[temp[j]].y - S[temp[i]].y < d; j++)
{
double d3 = dist(temp[i], temp[j]);
if(d > d3) d = d3;
}
return d;
}
int main()
{
scanf("%d", &n);
for(int i = 0; i < n; i++) scanf("%lf%lf", &S[i].x, &S[i].y);
sort(S, S + n, cmp);
return !printf("%.4lf\n", merge(0, n - 1));
}
一个个来解释。
比如说我们拿到了这样一组样例：
\textbf{InputData:}
10
1 1
1 5
3 1
5 1
5 6
6 7
7 3
8 1
10 3
9 9
这组样例是排好序的，我们可以得到这样的图：
下面我们要开始 \mathrm{merge()}了。（注意我的图从 1 开始编号，略有不同）
```cpp
// 预处理
double d = INF;
if(left == right) return d;
if(left + 1 == right) return dist(left, right);
```
然后算出 $mid$ 并进行递归。
[](https://cdn.luogu.com.cn/upload/pic/61219.png )
也就是求出 $[\text P_1\sim \text P_5]$ 和 $[\text P_6\sim \text P_{10}]$ 中的最近点对。
显然 $d_1=2$（$\text P_1$ 与 $\text P_3$ 之间），$d_2=\sqrt{5}$（$\text P_7$ 和 $\text P_8$ 之间）。
那么 $d = \min(d_1,d_2)=2$。
```cpp
// 递归
int mid = left + right >> 1;
double d1 = merge(left, mid);
double d2 = merge(mid + 1, right);
d = min(d1, d2);
```
[](https://cdn.luogu.com.cn/upload/pic/61221.png)
**前方高能**
接下来就是这个算法的核心，如何求出跨越蓝线的最近点对了。
首先锁定点集 $\text{Temp}$。
```cpp
// 求出temp[]并排序
for(i = left; i <= right; i++)
if(fabs(S[mid].x - S[i].x) < d)
temp[k++] = i;
sort(temp, temp + k, cmps);
```
[](https://cdn.luogu.com.cn/upload/pic/61218.png)
由程序可知 $\text{Temp}$ 包含了离直线距离不超过 $d$ 的所有点。
为啥其它点不用考虑呢？
**因为其它点离该线的距离都大于等于 ${d}$ 了，就更不用说到另一侧的点的距离了，必定不如直接返回 ${d}$ 来的优**。
而按照纵坐标排序是为了后面的方便。
```cpp
// 求出最近点对（TLE）
for(i = 0; i < k; i++)
for(j = i + 1; j < k; j++)
{
double d3 = dist(temp[i], temp[j]);
if(d > d3) d = d3;
}
```
相信这段代码大家都能看懂，但是很可惜，开 `-O2` 才能过，这是为什么呢？
不难发现源代码是这样的：
```cpp
// 求出最近点对
for(i = 0; i < k; i++)
for(j = i + 1; j < k && S[temp[j]].y - S[temp[i]].y < d; j++)
{
double d3 = dist(temp[i], temp[j]);
if(d > d3) d = d3;
}
```
好像多了一个 `S[temp[j]].y - S[temp[i]].y < d`，它有什么用呢？
和上面同理，**如果两者的纵坐标之差都大于等于 ${d}$ 了，那再计算就没有丝毫意义了，因为如果加上水平距离的话，必定不如直接返回 ${d}$ 来的优**。
例如 $\text P_3$ 号点，它只需和 $\text P_4$ 匹配就行了，因为它们的纵坐标之差为 $0$，有可能会带来更好的答案，而它和 $\text P_7$ 匹配就没有任何意义，因为即便它们的纵坐标之差等于 $d$，但就算 $\text P_7$ 在 $\text P_3$ 的正上方，距离也**不会小于** $d$。
值得注意的是，这里比如 $\text P_3$ 和 $\text P_4$ 会被重新计算，但对时间复杂度没有影响。
[](https://cdn.luogu.com.cn/upload/pic/61220.png)
至此，我们求出了结果 $ans = \sqrt{2}$（$\text{P}_5$ 和 $\text P_6$ 之间）。
所以：
> $\textbf{OutputData:}
1.4142
这个算法的时间复杂度又是怎么样的呢？（不想看分析的可以跳过）
考虑对于左边的一个点 \text{P}_i，可能和它成为新的最近公共点对的点 \text{P}_j 有多少个。
// 求出最近点对
for(i = 0; i < k; i++)
for(j = i + 1; j < k && S[temp[j]].y - S[temp[i]].y < d; j++)
{
double d3 = dist(temp[i], temp[j]);
if(d > d3) d = d3;
}
不难发现 \text{P}_j 满足这些条件：
-
-
-
根据这三条，我们可以画出一个 d \times 2d 的矩形，合法的 \text{P}_j 一定是在这个矩形中的（不含边界）。
但是，我们又知道，如果两个点同在左侧，则距离 \ge d；如果两个点同在右侧，则距离 \ge d。
那么，无论是 mid 左边的正方形还是右边的正方形，每个里面都至多放 3 个点（包括 \text{P}_i）。
图中红色点即代表 \text{P}_j，展示了一种方案，其中每条虚线边的长度都 \ge d。
所以对于每个 \text{P}_i，检验至多其他 5 个点即可，所以求最近点对的部分时间复杂度是 O(n) 的！
（当然其实与 \text{P}_i 同侧的 \text{P}_j 是不可能对更新答案有用的，如果你愿意的话，检验至多 3 个点即可）
于是整个代码的瓶颈就在 sort 上了，时间复杂度为 T(n) = 2\cdot T(\frac n 2) + O(n \log n) = O(n \log^2 n)，当然如果使用归并排序替换 sort 的话，就可以做到严格的 O(n \log n)！
最后，如果大家发现了问题，也欢迎大家指正哦！