凸包加kd树,求调.
查看原帖
凸包加kd树,求调.
817044
cjwdyzxfblzs楼主2023/6/12 21:38

这个kd树是拿的蒟蒻的k远点距离的那个题目的代码过来的,所以变量设置里面有一个k.但是不知道为什么样例只能输出1,得不到正确的答案,我感觉是精度问题,但是怎么改也改不出来欸.


#include <bits/stdc++.h>
using namespace std;
#define int long long
#define INF32_MAX 2147483647
namespace IO
{
    const int S = 1e5;
    char B[S], *H = B, *T = B;
    inline int gc() { return (H == T) && (T = (H = B) + fread(B, 1, S, stdin), H == T) ? EOF : *H++; }
    template <typename T = int>
    inline T read(T *p = nullptr)
    {
        static T *o = new T;
        !p && (p = o), *p = 0;
        int q = 1, c;
        while ((c = gc()) < '0' || c > '9')
            (c == '-') && (q = -1);
        *p = c ^ '0';
        while ((c = gc()) >= '0' && c <= '9')
            *p = *p * 10 + (c ^ '0');
        return *p *= q;
    }
    template <typename T, typename... Args>
    inline void read(T *p, Args... args) { read(p), read(args...); }
}
using IO::read;
const int N = 1e6;

/* Geometry objects */
struct Point { double x, y; };    /* point */
using Vec = Point;    /* vector */
struct Line { Point P; Vec v; };    /* line */
struct Seg { Point A, B; };    /* segment */
struct Circle { Point O; double r; };    /* circle */

/* constant */
const Point O = {0, 0};
const Line Ox = {O, {1, 0}}, Oy = {O, {0, 1}};
const double PI = acos(-1), eps = 1e-8;

/* Floating-point comparison */
// template<typename T> inline T abs(const T x) { return x >= 0 ? x : -x; } // |x|
bool eq(double a, double b) { return abs(a - b) < eps; }  // ==
bool gt(double a, double b) { return a - b > eps; }       // >
bool lt(double a, double b) { return a - b < -eps; }      // <
bool ge(double a, double b) { return a - b > -eps; }      // >= 
bool le(double a, double b) { return a - b < eps; }       // <=

/* Basic operations */
Vec r90a(Vec v) { return {-v.y, v.x}; } 
Vec r90c(Vec v) { return {v.y, -v.x}; }
Vec operator+(Vec u, Vec v) { return {u.x + v.x, u.y + v.y}; }
Vec operator-(Vec u, Vec v) { return {u.x - v.x, u.y - v.y}; }
Vec operator*(double k, Vec v) { return {k * v.x, k * v.y}; }
Vec operator/(Vec A, double p) { return {A.x / p, A.y / p}; }
double operator*(Vec u, Vec v) { return u.x * v.x + u.y * v.y; }
double operator^(Vec u, Vec v) { return u.x * v.y - u.y * v.x; } // cross ⇩
double Len(Vec v) { return sqrt(v.x * v.x + v.y * v.y); }
double slope(Vec v) { return v.y / v.x; }
double cross(Vec A, Vec B) { return A.x * B.y - A.y * B.x; } 

/* Vector correlation */
double cos_t(Vec u, Vec v) { return u * v / Len(u) / Len(v); }
Vec norm(Vec v) { return {v.x / Len(v), v.y / Len(v)}; }     // Unit vector
Vec pnorm(Vec v) { return (v.x < 0 ? -1 : 1) / Len(v) * v; }
Vec dirvec(Seg l) { return l.B - l.A; }                      // Segments'directions

/* Straight line correlation */
Line line(Point A, Point B) { return {A, B - A}; }           // 两点式直线
Line line(double k, double b) { return {{0, b}, {1, k}}; }   // 斜截式直线
Line line(Point P, double k) { return {P, {1, k}}; }         // 点截式直线
Line line(Seg l) { return {l.A, l.B - l.A}; }                // 线段所在的位置
double at_x(Line l, double x) { return l.P.y + (x - l.P.x) * l.v.y / l.v.x; }      // 给定直线的横坐标求纵坐标
double at_y(Line l, double y) { return l.P.x - (y + l.P.y) * l.v.x / l.v.y; }      // 给定直线的纵坐标求横坐标
Point pedal(Point P, Line l) { return l.P - (l.P - P) * l.v / (l.v * l.v) * l.v; }     // 求点到直线的垂足
Line perp(Line l, Point P) { return {P, r90c(l.v)}; }    // 过点做直线的垂线
Line bisec(Point P, Vec u, Vec v) { return {P, norm(u) + norm(v)}; }     // 角平分线

/* Line segment dependency */
Vec dvec(Seg l) { return l.B - l.A; } // Segment's direction
Point midp(Seg l) { return {(l.A.x + l.B.x) / 2, (l.A.y + l.B.y) / 2}; } // the middle of the segment
Line perp(Seg l) { return {midp(l), r90c(l.B - l.A)}; }     // The perpendicular line of the segment

/* The relationship between collection objects */
bool verti(Vec u, Vec v) { return eq(u * v, 0); }
bool paral(Vec u, Vec v) { return eq(u ^ v, 0); }
bool paral_x(Vec v) { return eq(v.y, 0); }
bool paral_y(Vec v) { return eq(v.x, 0); }
bool on(Point P, Line l) { return eq((P.x - l.P.x) * l.v.y, (P.y - l.P.y) * l.v.x); }
bool on(Point P, Seg l) { return eq(Len(P - l.A) + Len(P - l.B), Len(l.A - l.B)); }
bool operator==(Point A, Point B) { return eq(A.x, B.x) and eq(A.y, B.y); }
bool operator==(Line a, Line b) { return on(a.P, b) and on(a.P + a.v, b); }
bool operator==(Seg a,Seg b) { return (a.A == b.A and a.B == b.B) or (a.A == b.B and a.B == b.A); }
bool operator<(Point X, Point Y) { return lt(X.x, Y.x) or (eq(X.x, Y.x) and lt(X.y, Y.y)); }
bool tangency(Line l, Circle c) { return eq(abs((c.O ^ l.v) - (l.P ^ l.v)), c.r * Len(l.v)); }
bool tangency(Circle c1, Circle c2) { return eq(Len(c1.O - c2.O), c1.r + c2.r); }

/* Distance */
double dis(Point A, Point B) { return Len(A - B); }
double dis(Point P, Line l) { return abs((P ^ l.v) - (l.P ^ l.v)) / Len(l.v); }
double dis(Line a, Line b) { return abs((a.P ^ pnorm(a.v)) - (b.P ^ pnorm(b.v))); }

/* translate */
Line operator+(Line l, Vec v) { return {l.P + v, l.v}; }
Seg operator+(Seg l, Vec v) { return {l.A + v, l.B + v}; }

using Points = vector<Point>;
int n; Points P, H;
int qua(auto p) { return lt(p.y, 0) << 1 | lt(p.x, 0) ^ lt(p.y, 0); }
void __psort(Points &ps, Point c = O)   // 按照向量的叉乘排序
{
    sort(ps.begin(), ps.end(), [&](auto v1, auto v2) {
        return qua(v1 - c) < qua(v2 - c) || qua(v1 - c) == qua(v2 - c) && lt(cross(v1 - c, v2 - c), 0);
    });
}
double theta(Point p) { return p == O ? -1 / 0. : atan2(p.y, p.x); }    // 求极角
void psort(Points &ps, Point c = O)     // 极角排序
{
    sort(ps.begin(), ps.end(), [&](auto p1, auto p2) { return lt(theta(p1 - c), theta(p2 - c)); });
}
bool check(Point p, Point q, Point r) { return lt(0, (q - p) ^ (r - q)); }    // 检查是不是逆时针旋转的
Points chull(Points &ps)
{
    psort(ps, *min_element(ps.begin(), ps.end()));
    Points H{ps[0]};
    for (int i = 1; i < ps.size(); i ++ )
    {
        while (H.size() > 1 && !check(H[H.size() - 2], H.back(), ps[i]))
            H.pop_back();
        H.push_back(ps[i]);
    }
    return H;
}

/*=============================================================================================================================*/

#define pow(u) (u) * (u)
const int K = 2;
priority_queue<int, vector<int>, greater<int>> q;
int k = 2;
struct node
{
    int dim[K];
} s[N];
int lc[N], rc[N];
int L[N], R[N], D[N], U[N];
inline bool cmp1(node a, node b) { return a.dim[0] < b.dim[0]; }
inline bool cmp2(node a, node b) { return a.dim[1] < b.dim[1]; }
// inline int min(int a, int b) { return a < b ? a : b; }
// inline int max(int a, int b) { return a > b ? a : b; }
inline int dis(int a, int b){return max(pow(s[a].dim[0] - L[b]), pow(s[a].dim[0] - R[b])) + max(pow(s[a].dim[1] - D[b]), pow(s[a].dim[1] - U[b]));}
inline void update(int x)
{
    L[x] = R[x] = s[x].dim[0];
    D[x] = U[x] = s[x].dim[1];
    if (lc[x])
        L[x] = min(L[x], L[lc[x]]), R[x] = max(R[x], R[lc[x]]),
        D[x] = min(D[x], D[lc[x]]), U[x] = max(U[x], U[lc[x]]);
    if (rc[x])
        L[x] = min(L[x], L[rc[x]]), R[x] = max(R[x], R[rc[x]]),
        D[x] = min(D[x], D[rc[x]]), U[x] = max(U[x], U[rc[x]]);
}
int now = 1;
int build(int l, int r)
{
    if (l > r) return 0;
    if (l == r) { update(l);return l; }
    int mid = l + r >> 1;
    double avx = 0, avy = 0, vax = 0, vay = 0;
    for (int i = l; i <= r; i++) avx += s[i].dim[0], avy += s[i].dim[1];
    avx /= (r - l + 1), avy /= (r - l + 1);
    for (int i = l; i <= r; i++)
        vax += pow(avx - s[i].dim[0]),
        vay += pow(avy - s[i].dim[1]);
    if (vax > vay) nth_element(s + l, s + mid, s + r + 1, cmp1);
    else nth_element(s + l, s + mid, s + r + 1, cmp2);
    lc[mid] = build(l, mid - 1), rc[mid] = build(mid + 1, r);
    update(mid);
    return mid;
}
void query(int l, int r, int x)
{
    if (l > r) return;
    int mid = l + r >> 1;
    int t = pow(s[mid].dim[0] - s[x].dim[0]) + pow(s[mid].dim[1] - s[x].dim[1]);
    if (t > q.top()) q.pop(), q.push(t);
    double distl = dis(x, lc[mid]), distr = dis(x, rc[mid]);
    if (distl > q.top() and distr > q.top())
    {
        if (distl > distr)
        {
            query(l, mid - 1, x);
            if (distr > q.top()) query(mid + 1, r, x);
        }
        else
        {
            query(mid + 1, r, x);
            if (distl > q.top()) query(l, mid - 1, x);
        }
    }
    else
    {
        if (distl > q.top()) query(l, mid - 1, x);
        if (distr > q.top()) query(mid + 1, r, x);
    }
}

signed main()
{
    scanf("%d", &n);
    for (int i = 1; i <= n; i ++ )
    {
        double x, y;
        scanf("%lf%lf", &x, &y);
        P.push_back({x, y});
    }
    Points H = chull(P);

    for (int i = 1; i <= k; i ++ ) q.push(0);
    int cnt = 0;
    for (auto v : H) 
        s[++ cnt].dim[0] = v.x, 
        s[++ cnt].dim[1] = v.y;
    
    build(1, H.size());

    for (int i = 1; i <= H.size(); i ++ )
        query(1, H.size(), i);

    cout << ceil(q.top()) << endl;

    return 0;
}
2023/6/12 21:38
加载中...