「APIO2018」「LOJ #2586」選圓圈
阿新 • • 發佈：2020-08-19
Description
給定平面上的 \(n\) 個圓，用三個引數 \((x, y, R)\) 表示圓心座標和半徑。
每次選取最大的一個尚未被刪除的圓刪除，並同時刪除所有與其相切或相交的圓。
最後輸出每個圓分別是被那個圓所刪除的。
Hint
- \(1\le n\le 3\times 10^5\)
- \(0\le |x|, |y|, R \le 10^9\)
Solution 1
有一個非常簡單的 \(O(n^2)\) 暴力，由於每次都要掃一遍所有圓所以複雜度爆炸。
我們嘗試剪枝，縮小列舉的範圍。
若當前最大圓的半徑為 \(R\)，那麼我們做一個分塊操作：將整個平面 劃分為一個個方格，方格的邊長為 \(2R\)
對於當前這個圓，我們只要搜尋 其所在格子及其相鄰的 即可（兩個格子相鄰定義為所在行的距離不超過 1，且列距離也不超過 1，換言之，一個格子所有相鄰的就是周圍一圈 8 個）。顯然這樣是不會漏記的，因為所有圓的半徑都不超過方格大小，那麼一定不會出現兩個圓相交或相切，但卻不在相鄰兩個方格內。
但我們總不能每次都搜怎麼大的格子，最大圓的大小如果非常大而其他圓又很小的話這個剪枝沒有絲毫用處。因此我們引入一個 重構機制：設現在考慮到的圓的半徑為 \(R^\prime\)，原來方格大小為 \(L\)。若 \(R^\prime < L\)，那麼 重構整個方格，並以 \(2R^\prime\)
重構的複雜度會不會有問題？觀察到，當半徑不足方格大小的一半時才會重構，重構之後如果又來一次那又得一半。於是整個演算法不會有超過 \(O(\log R)\) 次重構，而一次重構的複雜度可以達到 \(O(n\log n)\)（排序時重構，搜尋時二分），所以總共最多也就兩隻 \(\log\)，不會有問題。
可以證明每個圓被檢查的次數為常數。但是我不會，這裡就不貼了，如果可以提供證明請私信，謝謝！
總複雜度為 \(O(n\log n\log R)\)。聽說雜湊表可以搞成一個 \(\log\)。
Code for Solution 1
/* * Author : _Wallace_ * Source : https://www.cnblogs.com/-Wallace-/ * Problem : APIO2018 LOJ #2586 選圓圈 */ #include <algorithm> #include <iostream> #include <vector> using namespace std; const int N = 3e5 + 5; const int inf = 1e9; int n; struct Triad { int x, y, r; inline void read() { cin >> x >> y >> r; x += inf, y += inf; } } circ[N]; int ans[N]; int cir_ord[N]; vector<Triad> f; vector<int> g[N]; int bsiz; inline long long sqr(int x) { return x * 1ll * x; } inline bool connect(Triad& a, Triad& b) { return sqr(a.x - b.x) + sqr(a.y - b.y) <= sqr(a.r + b.r); } typedef vector<Triad>::iterator vecIt; inline void init_blocks(int L) { if (bsiz && bsiz / 2 < L) return; bsiz = L; f.clear(); for (int i = 1; i <= n; i++) if (!ans[i]) f.push_back(Triad{circ[i].x / bsiz, circ[i].y / bsiz, i}); sort(f.begin(), f.end(), [](Triad& a, Triad& b) { return a.x != b.x ? a.x < b.x : a.y < b.y; }); int cnt = 0; for (vecIt i = f.begin(), j; i != f.end(); i = j, ++cnt) { for (j = i; j != f.end(); j++) if (i->x != j->x || i->y != j->y) break; f[cnt] = *i, g[cnt].clear(); for (vecIt k = i; k != j; k++) g[cnt].push_back(k->r); } f.resize(cnt); } void solve(int p) { if (ans[p]) return; init_blocks(circ[p].r * 2); int x = circ[p].x / bsiz; int y = circ[p].y / bsiz; for (int i = x - 1; i <= x + 1; i++) for (int j = y - 1; j <= y + 1; j++) { if (i < 0 || j < 0) continue; int c = lower_bound(f.begin(), f.end(), Triad{i, j, 0}, [](Triad& a, const Triad& b) { return a.x != b.x ? a.x < b.x : a.y < b.y; }) - f.begin(); if (c == int(f.size()) || f[c].x != i || f[c].y != j) continue; vector<int> buf; for (auto k : g[c]) if (connect(circ[k], circ[p])) ans[k] = p; else buf.push_back(k); g[c].swap(buf); } } signed main() { ios::sync_with_stdio(false); cin >> n; for (int i = 1; i <= n; i++) circ[i].read(), cir_ord[i] = i; sort(cir_ord + 1, cir_ord + 1 + n, [](int& a, int& b) { return circ[a].r != circ[b].r ? circ[a].r > circ[b].r : a < b; }); for (int i = 1; i <= n; i++) solve(cir_ord[i]); for (int i = 1; i <= n; i++) cout << ans[i] << ' '; return cout << endl, 0; }
Solution 2
既然是“二維平面”，那麼可以 KDT 暴力搞。
首先一個圓心為 \((x, y)\)，半徑為 \(R\) 的圓，我們可以粗略地認為是一個矩形區域 \((x-R, y-R)\sim (x+R, y+R)\)。
那麼考慮用 KDT 維護這些矩形。對於一個圓，搜尋可能與之有交集的區域即可。
但這樣複雜度是 \(O(\text{玄學})\) 而且可能被卡。
於是發揚人類智慧，將所有的點 隨機旋轉 一個角度。於是就可以過 LOJ 資料了。
不過Luogu 不需要旋轉 qwq。
Code for Solution 2
/*
* Author : _Wallace_
* Source : https://www.cnblogs.com/-Wallace-/
* Problem : APIO2018 LOJ #2586 選圓圈
*/
#include <algorithm>
#include <cmath>
#include <iostream>
using namespace std;
const double eps = 1e-3;
const int N = 3e5 + 5;
const int K = 2;
int n;
struct circle {
double p[K]; int r;
};
pair<circle, int> C[N];
int ans[N];
struct area {
double max[K], min[K];
};
inline area trans(circle a) {
return area{{a.p[0] + a.r, a.p[1] + a.r}, {a.p[0] - a.r, a.p[1] - a.r}};
}
struct node {
int lc, rc;
area range;
circle val;
int index, vaild;
} t[N];
int total = 0;
inline void pushup(int x) {
t[x].range = trans(t[x].val);
t[x].vaild = t[t[x].lc].vaild + t[t[x].rc].vaild + 1;
for (int i = 0; i < K; i++) {
if (t[x].lc) {
t[x].range.max[i] = max(t[x].range.max[i], t[t[x].lc].range.max[i]);
t[x].range.min[i] = min(t[x].range.min[i], t[t[x].lc].range.min[i]);
}
if (t[x].rc) {
t[x].range.max[i] = max(t[x].range.max[i], t[t[x].rc].range.max[i]);
t[x].range.min[i] = min(t[x].range.min[i], t[t[x].rc].range.min[i]);
}
}
}
int build(int l, int r, int d) {
if (l > r) return 0;
int mid = (l + r) >> 1, x = ++total;
nth_element(C + l, C + mid, C + r + 1, [&](pair<circle, int> a, pair<circle, int> b) {
return a.first.p[d] < b.first.p[d];
});
t[x].val = C[mid].first;
t[x].index = C[mid].second;
t[x].lc = build(l, mid - 1, d ^ 1);
t[x].rc = build(mid + 1, r, d ^ 1);
return pushup(x), x;
}
inline bool outside(circle& a, area& b) {
double x = a.p[0], y = a.p[1]; int r = a.r;
if (b.min[0] - (x + r) >= eps) return true;
if (b.min[1] - (y + r) >= eps) return true;
if ((x - r) - b.max[0] >= eps) return true;
if ((y - r) - b.max[1] >= eps) return true;
return false;
}
inline double sqr(double x) {
return x * x;
}
inline bool connect(circle& a, circle& b) {
return sqr(a.r + b.r) - (sqr(a.p[0] - b.p[0]) + sqr(a.p[1] - b.p[1])) >= -eps;
}
void select(int x, pair<circle, int>& c) {
if (!x || outside(c.first, t[x].range)) return;
if (!ans[t[x].index] && connect(c.first, t[x].val))
ans[t[x].index] = c.second;
if (t[t[x].lc].vaild) select(t[x].lc, c);
if (t[t[x].rc].vaild) select(t[x].rc, c);
t[x].vaild = t[t[x].lc].vaild + t[t[x].rc].vaild + (!ans[t[x].index] ? 1 : 0);
}
inline void rotate(double theta) {
double SIN = sin(theta), COS = cos(theta);
for (int i = 1; i <= n; i++) {
double x = C[i].first.p[0];
double y = C[i].first.p[1];
C[i].first.p[1] = x * SIN + y * COS;
C[i].first.p[0] = x * COS - y * SIN;
}
}
signed main() {
ios::sync_with_stdio(false);
cin >> n;
for (int i = 1; i <= n; i++) {
cin >> C[i].first.p[0] >> C[i].first.p[1] >> C[i].first.r;
C[i].second = i;
}
rotate(2.3);
int root = build(1, n, 0);
sort(C + 1, C + 1 + n, [](pair<circle, int>& a, pair<circle, int>& b) {
return a.first.r != b.first.r ? a.first.r > b.first.r : a.second < b.second;
});
for (int i = 1; i <= n; i++)
if (!ans[C[i].second]) select(root, C[i]);
for (int i = 1; i <= n; i++)
cout << ans[i] << ' ';
return cout << endl, 0;
}