Luogu P1742 minimum circle cover (computational geometry)

Posted by duckduckgoose on Sat, 07 Dec 2019 11:54:17 +0100

meaning of the title

Title Link

Sol

The violence method is \ (O(n^3) \) enumerate three points and check whether all points can be included

Consider a random algorithm. First, shuffle the sequence random.

Then we enumerate a point \ (i \) and maintain a current circle.

Enumerate another point \ (j \), if the point continues in the circle, otherwise replace the previous circle with the circle constructed by \ (i, j \).

Enumerate another point \ (k \), if the point continues in the circle, otherwise use \ (i, j, k \) to construct a new circle.

Such expected complexity is O(n)

At first, I thought that there was something wrong with the correctness of doing this, that is to say, it might find a bad solution. But it is obviously wrong, because if there is a better solution and the area is smaller than the current one, the solution should contain at least three points of the current one, which is contradictory.

#include<bits/stdc++.h>
using namespace std;
const int MAXN = 1e5 + 10;
int N;
double R;
struct Point {
    double x, y;
}p[MAXN], C;
double sqr(double x) {
    return x * x;
}
double dis(Point a, Point b) {
    return sqrt(sqr(a.x - b.x) + sqr(a.y - b.y));
}
void MakeC(Point p1, Point p2, Point p3) {
    double a = p2.x - p1.x,
           b = p2.y - p1.y,
           c = p3.x - p1.x,
           d = p3.y - p1.y,
           e = (sqr(p2.x) - sqr(p1.x) + sqr(p2.y) - sqr(p1.y)) / 2,
           f = (sqr(p3.x) - sqr(p1.x) + sqr(p3.y) - sqr(p1.y)) / 2;
    C.x = (e * d - b * f) / (a * d - b * c);
    C.y = (a * f - e * c) / (a * d - b * c);
    R = dis(C, p1);
}
int main() {
    cin >> N;
    for(int i = 1; i <= N; i++) scanf("%lf %lf", &p[i].x, &p[i].y);
    random_shuffle(p + 1, p + N + 1);
    for(int i = 1; i <= N; i++) {
        if(dis(p[i], C) < R) continue;
        C = p[i]; R = 0;
        for(int j = 1; j <= i - 1; j++) {
            if(dis(p[j], C) < R) continue;
            C.x = (p[i].x + p[j].x) / 2.0;
            C.y = (p[i].y + p[j].y) / 2.0;
            R = dis(C, p[j]);
            for(int k = 1; k <= j - 1; k++) {
                if(dis(p[k], C) < R) continue;
                MakeC(p[i], p[j], p[k]);
            }
        }
    }
    printf("%.10lf\n", R);
    printf("%.10lf %.10lf", C.x, C.y);
    return 0;
}

Topics: C++