Project Euler 816

发表于 更新于 分类于 阅读次数： 本文字数： 645 阅读时长 ≈ 1 分钟

Project Euler 816

题目

Shortest distance among points

We create an array of points \(P_n\) in a two dimensional plane using the following random number generator:

\(\begin{aligned} &s_0=290797\\ &s_{n+1}={s_n}^2 \bmod 50515093\\ &P_n=(s_{2n},s_{2n+1}) \end{aligned}\)

Let \(d(k)\) be the shortest distance of any two (distinct) points among \(P_0, \cdots, P_{k - 1}\).

E.g. \(d(14)=546446.466846479\)

Find \(d(2000000)\). Give your answer rounded to \(9\) places after the decimal point.

解决方案

经典问题平面最近点对问题。使用分治的思想可以解决。大概的做法是：

  • 分解：将一个平面上的所有点按照\(x\)坐标划分成相等的两部分；
  • 结局：然后递归求解这两部分的解；
  • 合并：假设这条划分接线为\(x=x_0\)，那么将\(x=x_0\)附近的点再按照\(y\)轴进行排序，再枚举这些点之间距离的最小值即可。最终得到了这一部分点的最小点对长度。

代码

本做法的时间复杂度为\(O(n\log^2 n)\)，虽然不是最优，但足以应对。

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
#include <bits/stdc++.h>
using namespace std;
typedef long long ll;
const int N=2e6;
ll INF = 1e17;
struct P {
    ll x, y;
    bool operator<(const P &p) const {
        return x < p.x || x == p.x && y < p.y;
    }
    ll dis2(P &p) const {
        ll dx = x - p.x, dy = y - p.y;
        return dx * dx + dy * dy;
    }
};
P p[N + 4];
int tp[N + 4];
ll dfs(int l, int r) {
    if (l == r)
        return INF;
    if (l + 1 == r)
        return p[l].dis2(p[r]);
    int mid = (l + r) >> 1;
    ll ans = min(dfs(l, mid - 1), dfs(mid, r));
    int k = 0;
    for (int i = l; i <= r; i++) {
        ll dx = p[mid].x - p[i].x;
        if (dx * dx <= ans)
            tp[++k] = i;
    }
    sort(tp + 1, tp + k + 1, [](const int &a,const int &b){return p[a].y < p[b].y; });
    for (int i = 1; i <= k; i++) {
        for (int j = i + 1; j <= k; j++) {
            ll dy = p[tp[j]].y - p[tp[i]].y;
            if(dy * dy > ans) break;
            ans = min(ans, p[tp[i]].dis2(p[tp[j]]));
        }
    }
    return ans;
}
void input(int N) {
    ll s = 290797;
    for (int i = 1; i <= N; i++) {
        p[i].x = s; s = s * s % 50515093;
        p[i].y = s; s = s * s % 50515093;
    }
}
int main()
{
    input(N);
    sort(p + 1, p + N + 1);
    double ans = sqrt(dfs(1, N));
    printf("%.9f\n", ans);
}