QOJ.ac
QOJ
ID | Problem | Submitter | Result | Time | Memory | Language | File size | Submit time | Judge time |
---|---|---|---|---|---|---|---|---|---|
#544976 | #7521. Find the Gap | asitshouldbe | WA | 3ms | 4652kb | C++17 | 40.0kb | 2024-09-02 21:23:34 | 2024-09-02 21:23:35 |
Judging History
answer
#include <bits/stdc++.h>
#define pi acos(-1.0)
// #define double long double
using namespace std;
const double eps = 1e-8, inf = 1e12;
const int N = 1e4 + 10;
int sgn(double x)
{
if (fabs(x) < eps) return 0;
if (x < 0) return -1;
else return 1;
}
int dsgn(double x, double y)
{
if (fabs(x - y) < eps) return 0;
if (x < y) return -1;
else return 1;
}
inline double sqr(double x) { return x * x; }
struct Point
{
double x, y;
Point() {}
Point(double _x, double _y) { x = _x;y = _y; }
void input() { cin >> x >> y; }
bool operator==(Point b) const { return sgn(x - b.x) == 0 && sgn(y - b.y) == 0; }
bool operator<(Point b) const { return sgn(x - b.x) == 0 ? sgn(y - b.y) < 0 : x < b.x; }
Point operator-(const Point &b) const { return Point(x - b.x, y - b.y); }
Point operator+(const Point &b) const { return Point(x + b.x, y + b.y); }
// 叉积
double operator^(const Point &b) const { return x * b.y - y * b.x; }
// 点积
double operator*(const Point &b) const { return x * b.x + y * b.y; }
// 数乘
Point operator*(const double &k) const {return Point(x * k, y * k);}
Point operator/(const double &k) const { return Point(x / k, y / k);}
// 返回长度
double len() { return hypot(x, y); }
// 返回长度的平方
double len2() { return x * x + y * y; }
// 返回两点的距离
double distance(Point p) { return hypot(x - p.x, y - p.y); }
int distance2(Point p) { return (x - p.x) * (x - p.x) + (y - p.y) * (y - p.y); }
//`计算pa 和 pb 的夹角 就是求这个点看a,b 所成的夹角`
//`测试 LightOJ1203`
double rad(Point a, Point b)
{
Point p = *this;
return fabs(atan2(fabs((a - p) ^ (b - p)), (a - p) * (b - p)));
}
//`绕着p点逆时针旋转angle`
Point rotate(Point p, double angle)
{
Point v = (*this) - p;
double c = cos(angle), s = sin(angle);
return Point(p.x + v.x * c - v.y * s, p.y + v.x * s + v.y * c);
}
//`化为长度为r的向量`
Point trunc(double r)
{
double l = len();
if (!sgn(l)) return *this;
r /= l;
return Point(x * r, y * r);
}
//`逆时针旋转90度`
Point rotleft() { return Point(-y, x); }
//`顺时针旋转90度`;
Point rotright() { return Point(y, -x); }
};
//`AB X AC`
double cross(Point A, Point B, Point C) { return (B - A) ^ (C - A); }
//`AB*AC`
double dot(Point A, Point B, Point C) { return (B - A) * (C - A); }
struct Line
{
Point s, e;
Line() {}
Line(Point _s, Point _e) { s = _s;e = _e; }
bool operator==(Line v) { return (s == v.s) && (e == v.e); }
//`根据一个点和倾斜角angle确定直线,0<=angle<pi`
Line(Point p, double angle)
{
s = p;
if (sgn(angle - pi / 2) == 0) e = (s + Point(0, 1));
else e = (s + Point(1, tan(angle)));
}
// ax+by+c=0
Line(double a, double b, double c)
{
if (sgn(a) == 0)
{
s = Point(0, -c / b);
e = Point(1, -c / b);
}
else if (sgn(b) == 0)
{
s = Point(-c / a, 0);
e = Point(-c / a, 1);
}
else
{
s = Point(0, -c / b);
e = Point(1, (-c - a) / b);
}
}
void input() { s.input();e.input(); }
void adjust() { if (e < s) swap(s, e); }
// 求线段长度
double length() { return s.distance(e); }
//`返回直线倾斜角 0<=angle<pi`
double angle()
{
double k = atan2(e.y - s.y, e.x - s.x);
// if(sgn(k) < 0)k += pi;
// if(sgn(k-pi) == 0)k -= pi;
return k;
}
bool operator<(Line &v) { return sgn(angle() - v.angle()) == 0 ? ((e - s) ^ (v.e - s)) < 0 : angle() < v.angle(); }
//`点和直线关系`
int relation(Point p)
{
int c = sgn((p - s) ^ (e - s));
if (c < 0) return 1; //`1 在左侧`
else if (c > 0) return 2; //`2 在右侧`
else return 3; //`3 在直线上`
}
//`点在线段上的判断`
bool pointonseg(Point p) { return sgn((p - s) ^ (e - s)) == 0 && sgn((p - s) * (p - e)) <= 0; }
//`两向量平行(对应直线平行或重合)`
bool parallel(Line v) { return sgn((e - s) ^ (v.e - v.s)) == 0; }
//`两向量正交`
bool orthogonal(Line v) { return sgn((e - s) * (v.e - v.s)) == 0; }
//`两直线关系`
int linecrossline(Line v)
{
if ((*this).parallel(v)) //`0 平行`
return v.relation(s) == 3; //`1 重合`
if ((*this).orthogonal(v)) return 3; //`3 正交`
else return 2; //`2 相交`
}
//`两线段相交判断`
int segcrossseg(Line v)
{
int d1 = sgn((e - s) ^ (v.s - s));
int d2 = sgn((e - s) ^ (v.e - s));
int d3 = sgn((v.e - v.s) ^ (s - v.s));
int d4 = sgn((v.e - v.s) ^ (e - v.s));
if ((d1 ^ d2) == -2 && (d3 ^ d4) == -2) return 2; //`2 规范相交`
return (d1 == 0 && sgn((v.s - s) * (v.s - e)) <= 0) || //`1 非规范相交`
(d2 == 0 && sgn((v.e - s) * (v.e - e)) <= 0) || //`0 不相交`
(d3 == 0 && sgn((s - v.s) * (s - v.e)) <= 0) ||
(d4 == 0 && sgn((e - v.s) * (e - v.e)) <= 0);
}
//`直线和线段相交判断`
int linecrossseg(Line v)
{
int d1 = sgn((e - s) ^ (v.s - s));
int d2 = sgn((e - s) ^ (v.e - s));
if ((d1 ^ d2) == -2) return 2; //`2 规范相交`
return (d1 == 0 || d2 == 0); //`1 非规范相交 0 不相交`
}
//`求两直线的交点 要保证两直线不平行或重合`
Point crosspoint(Line l)
{
Point u = e - s, v = l.e - l.s;
double t = (s - l.s) ^ v / (v ^ u);
return s + u * t;
}
//`返回点p在直线上的投影`
Point lineprog(Point p)
{
Point u = e - s, v = p - s;
return s + u * (u * v / u.len2());
}
// 点到直线的距离
double dispointtoline(Point p) { return fabs((p - s) ^ (e - s)) / length(); }
// 点到线段的距离
double dispointtoseg(Point p)
{
if (sgn((p - s) * (e - s)) < 0 || sgn((p - e) * (s - e)) < 0)
return min(p.distance(s), p.distance(e));
return dispointtoline(p);
}
//`返回线段到线段的距离 前提是两线段不相交,相交距离就是0了`
double dissegtoseg(Line v)
{
return min(min(dispointtoseg(v.s), dispointtoseg(v.e)), min(v.dispointtoseg(s), v.dispointtoseg(e)));
}
//`返回点p关于直线的对称点`
Point symmetrypoint(Point p)
{
Point pro = lineprog(p);
return pro * 2 - p;
}
};
struct circle
{
Point p; // 圆心
double r; // 半径
circle() {}
circle(Point _p, double _r) { p = _p;r = _r; }
circle(double x, double y, double _r) { p = Point(x, y);r = _r; }
//`三角形的外接圆`
//`需要Point的+ / rotate() 以及Line的crosspoint()`
//`利用两条边的中垂线得到圆心`
//`测试:UVA12304`
circle(Point a, Point b, Point c)
{
Line u = Line((a + b) / 2, ((a + b) / 2) + ((b - a).rotleft()));
Line v = Line((b + c) / 2, ((b + c) / 2) + ((c - b).rotleft()));
p = u.crosspoint(v);
r = p.distance(a);
}
//`三角形的内切圆`
//`参数bool t没有作用,只是为了和上面外接圆函数区别`
//`测试:UVA12304`
circle(Point a, Point b, Point c, bool t)
{
Line u, v;
double m = atan2(b.y - a.y, b.x - a.x), n = atan2(c.y - a.y, c.x - a.x);
u.s = a;
u.e = u.s + Point(cos((n + m) / 2), sin((n + m) / 2));
v.s = b;
m = atan2(a.y - b.y, a.x - b.x), n = atan2(c.y - b.y, c.x - b.x);
v.e = v.s + Point(cos((n + m) / 2), sin((n + m) / 2));
p = u.crosspoint(v);
r = Line(a, b).dispointtoseg(p);
}
// 输入
void input() { p.input();cin >> r; }
bool operator==(circle v) { return (p == v.p) && sgn(r - v.r) == 0; }
bool operator<(circle v) const { return ((p < v.p) || ((p == v.p) && sgn(r - v.r) < 0)); }
// 面积
double area() { return pi * r * r; }
// 周长
double circumference() { return 2 * pi * r; }
//`点和圆的关系`
int relation(Point b)
{
double dst = b.distance(p);
if (sgn(dst - r) < 0) return 2; //`2 圆内`
else if (sgn(dst - r) == 0) return 1; //`1 圆上`
return 0; //`0 圆外`
}
//`线段和圆的关系 比较的是圆心到线段的距离和半径的关系`
int relationseg(Line v)
{
double dst = v.dispointtoseg(p);
if (sgn(dst - r) < 0) return 2; //`2 相交`
else if (sgn(dst - r) == 0) return 1; //`1 相切`
return 0; //`0 相离`
}
//`直线和圆的关系 比较的是圆心到直线的距离和半径的关系`
int relationline(Line v)
{
double dst = v.dispointtoline(p);
if (sgn(dst - r) < 0) return 2; //`2 相交`
else if (sgn(dst - r) == 0) return 1; //`1 相切`
return 0; //`0 相离`
}
//`两圆的关系 需要Point的distance`
//`测试:UVA12304`
int relationcircle(circle v)
{
double d = p.distance(v.p);
if (sgn(d - r - v.r) > 0) return 5; //`5 相离`
if (sgn(d - r - v.r) == 0) return 4; //`4 外切`
double l = fabs(r - v.r);
if (sgn(d - r - v.r) < 0 && sgn(d - l) > 0) return 3; //`3 相交`
if (sgn(d - l) == 0) return 2; //`2 内切`
if (sgn(d - l) < 0) return 1; //`1 内含`
}
//`求两个圆的交点 需要relationcircle`
//`测试:UVA12304`
int pointcrosscircle(circle v, Point &p1, Point &p2)
{
int rel = relationcircle(v);
if (rel == 1 || rel == 5) return 0; //`0 无交点`
double d = p.distance(v.p);
double l = (d * d + r * r - v.r * v.r) / (2 * d);
double h = sqrt(r * r - l * l);
Point tmp = p + (v.p - p).trunc(l);
p1 = tmp + ((v.p - p).rotleft().trunc(h));
p2 = tmp + ((v.p - p).rotright().trunc(h));
if (rel == 2 || rel == 4) return 1; //`1 一个交点`
return 2; //`2 两个交点`
}
//`求直线和圆的交点,返回交点个数`
int pointcrossline(Line v, Point &p1, Point &p2)
{
if (!(*this).relationline(v)) return 0;
Point a = v.lineprog(p);
double d = v.dispointtoline(p);
d = sqrt(r * r - d * d);
if (sgn(d) == 0)
{
p1 = a;
p2 = a;
return 1;
}
p1 = a + (v.e - v.s).trunc(d);
p2 = a - (v.e - v.s).trunc(d);
return 2;
}
//`得到过a,b两点,半径为r1的两个圆`
int gercircle(Point a, Point b, double r1, circle &c1, circle &c2)
{
circle x(a, r1), y(b, r1);
int t = x.pointcrosscircle(y, c1.p, c2.p);
if (!t) return 0;
c1.r = c2.r = r;
return t; //`返回圆的个数`
}
//`得到与直线u相切,过点q,半径为r1的圆`
//`测试:UVA12304`
int getcircle(Line u, Point q, double r1, circle &c1, circle &c2)
{
double dis = u.dispointtoline(q);
if (sgn(dis - r1 * 2) > 0) return 0;
if (sgn(dis) == 0)
{
c1.p = q + ((u.e - u.s).rotleft().trunc(r1));
c2.p = q + ((u.e - u.s).rotright().trunc(r1));
c1.r = c2.r = r1;
return 2;
}
Line u1 = Line((u.s + (u.e - u.s).rotleft().trunc(r1)), (u.e + (u.e - u.s).rotleft().trunc(r1)));
Line u2 = Line((u.s + (u.e - u.s).rotright().trunc(r1)), (u.e + (u.e - u.s).rotright().trunc(r1)));
circle cc = circle(q, r1);
Point p1, p2;
if (!cc.pointcrossline(u1, p1, p2)) cc.pointcrossline(u2, p1, p2);
c1 = circle(p1, r1);
if (p1 == p2)
{
c2 = c1;
return 1;
}
c2 = circle(p2, r1);
return 2;
}
//`同时与直线u,v相切,半径为r1的圆`
//`测试:UVA12304`
int getcircle(Line u, Line v, double r1, circle &c1, circle &c2, circle &c3, circle &c4)
{
if (u.parallel(v)) return 0; // 两直线平行
Line u1 = Line(u.s + (u.e - u.s).rotleft().trunc(r1), u.e + (u.e - u.s).rotleft().trunc(r1));
Line u2 = Line(u.s + (u.e - u.s).rotright().trunc(r1), u.e + (u.e - u.s).rotright().trunc(r1));
Line v1 = Line(v.s + (v.e - v.s).rotleft().trunc(r1), v.e + (v.e - v.s).rotleft().trunc(r1));
Line v2 = Line(v.s + (v.e - v.s).rotright().trunc(r1), v.e + (v.e - v.s).rotright().trunc(r1));
c1.r = c2.r = c3.r = c4.r = r1;
c1.p = u1.crosspoint(v1);
c2.p = u1.crosspoint(v2);
c3.p = u2.crosspoint(v1);
c4.p = u2.crosspoint(v2);
return 4;
}
//`同时与不相交圆cx,cy相切,半径为r1的圆`
//`测试:UVA12304`
int getcircle(circle cx, circle cy, double r1, circle &c1, circle &c2)
{
circle x(cx.p, r1 + cx.r), y(cy.p, r1 + cy.r);
int t = x.pointcrosscircle(y, c1.p, c2.p);
if (!t) return 0;
c1.r = c2.r = r1;
return t; //`返回圆的个数`
}
//`过一点作圆的切线(先判断点和圆的关系)`
//`测试:UVA12304`
int tangentline(Point q, Line &u, Line &v)
{
int x = relation(q);
if (x == 2) return 0;
if (x == 1)
{
u = Line(q, q + (q - p).rotleft());
v = u;
return 1;
}
double d = p.distance(q);
double l = r * r / d;
double h = sqrt(r * r - l * l);
u = Line(q, p + ((q - p).trunc(l) + (q - p).rotleft().trunc(h)));
v = Line(q, p + ((q - p).trunc(l) + (q - p).rotright().trunc(h)));
return 2; //`返回切线的个数`
}
int tangentpoint(Point q, Point &u, Point &v)
{
int x = relation(q);
if (x == 2) return 0;
if (x == 1)
{
u = q;
v = u;
return 1;
}
double d = p.distance(q);
double l = r * r / d;
double h = sqrt(r * r - l * l);
u = p + ((q - p).trunc(l) + (q - p).rotleft().trunc(h));
v = p + ((q - p).trunc(l) + (q - p).rotright().trunc(h));
return 2; //`返回切点的个数`
}
//`求两圆相交的面积`
double areacircle(circle v)
{
int rel = relationcircle(v);
if (rel >= 4) return 0.0;
if (rel <= 2) return min(area(), v.area());
double d = p.distance(v.p);
double hf = (r + v.r + d) / 2.0;
double ss = 2 * sqrt(hf * (hf - r) * (hf - v.r) * (hf - d));
double a1 = acos((r * r + d * d - v.r * v.r) / (2.0 * r * d));
a1 = a1 * r * r;
double a2 = acos((v.r * v.r + d * d - r * r) / (2.0 * v.r * d));
a2 = a2 * v.r * v.r;
return a1 + a2 - ss;
}
//`求圆和三角形pab的相交面积`
//`测试:POJ3675 HDU3982 HDU2892`
double areatriangle(Point a, Point b)
{
if (sgn((p - a) ^ (p - b)) == 0) return 0.0;
Point q[5],p1, p2;
int len = 0;
q[len++] = a;
Line l(a, b);
if (pointcrossline(l, q[1], q[2]) == 2)
{
if (sgn((a - q[1]) * (b - q[1])) < 0) q[len++] = q[1];
if (sgn((a - q[2]) * (b - q[2])) < 0) q[len++] = q[2];
}
q[len++] = b;
if (len == 4 && sgn((q[0] - q[1]) * (q[2] - q[1])) > 0) swap(q[1], q[2]);
double res = 0;
for (int i = 0; i < len - 1; i++)
{
if (relation(q[i]) == 0 || relation(q[i + 1]) == 0)
{
double arg = p.rad(q[i], q[i + 1]);
res += r * r * arg / 2.0;
}
else res += fabs((q[i] - p) ^ (q[i + 1] - p)) / 2.0;
}
return res;
}
//`两圆公切线`
Point getpoint(double rad) { return Point(p.x + r * cos(rad), p.y + r * sin(rad)); }
int conmontangent(circle v, vector<Point> &p1, vector<Point> &p2)
{
bool flag = 0;
if (r < v.r) swap(*this, v), flag = 1;
double d = p.distance(v.p), rd = r - v.r, rs = r + v.r;
if (sgn(d - rd) < 0) return 0;
if (sgn(d) == 0) return -1;
double rad = Line(p, v.p).angle();
if (sgn(d - rd) == 0)
{
p1.push_back(getpoint(rad)), p2.push_back(getpoint(rad));
return 1; //`一条外公切线`
}
double rad1 = acos(rd / d);
p1.push_back(getpoint(rad + rad1)), p2.push_back(v.getpoint(rad + rad1));
p1.push_back(getpoint(rad - rad1)), p2.push_back(v.getpoint(rad - rad1));
if (sgn(d - rs) == 0)
{
p1.push_back(getpoint(rad)), p2.push_back(getpoint(rad));
if (flag) swap(p1, p2);
return 3; //`两条外公切线 一条内公切线`
}
else if (sgn(d - rs) > 0)
{
double rad2 = acos(rs / d);
p1.push_back(getpoint(rad + rad2)), p2.push_back(v.getpoint(rad + rad2 - pi));
p1.push_back(getpoint(rad - rad2)), p2.push_back(v.getpoint(rad - rad2 + pi));
if (flag) swap(p1, p2);
return 4; //`两条外公切线 两条内公切线`
}
else
{
if (flag) swap(p1, p2);
return 2; //`两条外公切线`
}
}
};
struct polygon
{
int n;
Point p[N];
Line l[N];
void input(int _n)
{
n = _n;
for (int i = 0; i < n; i++) p[i].input();
}
void add(Point q) { p[n++] = q; }
void getline()
{
for (int i = 0; i < n; i++)
l[i] = Line(p[i], p[(i + 1) % n]);
}
struct cmp
{
Point p;
cmp(const Point &p0) { p = p0; }
bool operator()(const Point &aa, const Point &bb)
{
Point a = aa, b = bb;
int d = sgn((a - p) ^ (b - p));
if (d == 0) return sgn(a.distance(p) - b.distance(p)) < 0;
return d > 0;
}
};
//`进行极角排序 mi为极点`
//`需要重载号好Point的 < 操作符(min函数要用) `
void norm(Point mi)
{
// Point mi = p[0];
// for(int i = 1;i < n;i++)mi = min(mi,p[i]);
sort(p, p + n, cmp(mi));
}
//`得到的凸包里面的点编号是0-n-1的`
//`注意如果有影响,要特判下所有点共点,或者共线的特殊情况`
//`测试 LightOJ1203 LightOJ1239`
void andrew(polygon &convex)
{
sort(p, p + n);
int &top = convex.n;
top = 0;
for (int i = 0; i < n; i++)
{
while (top >= 2 && sgn(cross(convex.p[top - 2], convex.p[top - 1], p[i])) <= 0) top--;
convex.p[top++] = p[i];
}
int temp = top;
for (int i = n - 2; i >= 0; i--)
{
while (top > temp && sgn(cross(convex.p[top - 2], convex.p[top - 1], p[i])) <= 0) top--;
convex.p[top++] = p[i];
}
top--;
}
//`判断是不是凸的`
bool isconvex()
{
bool s[5];
memset(s, false, sizeof(s));
for (int i = 0; i < n; i++)
{
int j = (i + 1) % n, k = (j + 1) % n;
s[sgn((p[j] - p[i]) ^ (p[k] - p[i])) + 1] = true;
if (s[0] && s[2]) return false;
}
return true;
}
//`判断点和任意多边形的关系`
int relationpoint(Point q)
{
for (int i = 0; i < n; i++)
if (p[i] == q) return 3; //` 3 点上`
for (int i = 0; i < n; i++)
if (l[i].pointonseg(q)) return 2; //` 2 边上`
int cnt = 0;
for (int i = 0; i < n; i++)
{
int j = (i + 1) % n;
int k = sgn((q - p[j]) ^ (p[i] - p[j]));
int u = sgn(p[i].y - q.y);
int v = sgn(p[j].y - q.y);
if (k > 0 && u < 0 && v >= 0) cnt++;
if (k < 0 && v < 0 && u >= 0) cnt--;
} //` 1 内部`
return cnt != 0; //` 0 外部`
}
//`直线u切割凸多边形左侧 注意直线方向`
//`测试:HDU3982`
void convexcut(Line u, polygon &po)
{
int &top = po.n; // 注意引用
top = 0;
for (int i = 0; i < n; i++)
{
int d1 = sgn((u.e - u.s) ^ (p[i] - u.s));
int d2 = sgn((u.e - u.s) ^ (p[(i + 1) % n] - u.s));
if (d1 >= 0) po.p[top++] = p[i];
if (d1 * d2 < 0) po.p[top++] = u.crosspoint(Line(p[i], p[(i + 1) % n]));
}
}
//`得到周长`
//`测试 LightOJ1239`
double getcircumference()
{
double sum = 0;
for (int i = 0; i < n; i++)
sum += p[i].distance(p[(i + 1) % n]);
return sum;
}
//`得到面积`
double getarea()
{
double sum = 0;
for (int i = 0; i < n; i++)
sum += (p[i] ^ p[(i + 1) % n]);
return fabs(sum) / 2;
}
//`得到方向`
bool getdir()
{
double sum = 0;
for (int i = 0; i < n; i++)
sum += (p[i] ^ p[(i + 1) % n]);
if (sgn(sum) > 0) return 1; //` 1 逆时针`
else return 0; //` 0 顺时针`
}
//`得到重心`
Point getbarycentre()
{
Point ret(0, 0);
double area = 0;
for (int i = 1; i < n - 1; i++)
{
double tmp = (p[i] - p[0]) ^ (p[i + 1] - p[0]);
if (sgn(tmp) == 0) continue;
area += tmp;
ret.x += (p[0].x + p[i].x + p[i + 1].x) / 3 * tmp;
ret.y += (p[0].y + p[i].y + p[i + 1].y) / 3 * tmp;
}
if (sgn(area)) ret = ret / area;
return ret;
}
//`多边形和圆交的面积`
//`测试:POJ3675 HDU3982 HDU2892`
double areacircle(circle c)
{
double ans = 0;
for (int i = 0; i < n; i++)
{
int j = (i + 1) % n;
if (sgn((p[j] - c.p) ^ (p[i] - c.p)) >= 0) ans += c.areatriangle(p[i], p[j]);
else ans -= c.areatriangle(p[i], p[j]);
}
return fabs(ans);
}
//`多边形和圆关系`
int relationcircle(circle c)
{
int x = 2; //` 2 圆完全在多边形内`
if (relationpoint(c.p) != 1) return 0; //` 0 圆心不在内部`
for (int i = 0; i < n; i++)
{
if (c.relationseg(l[i]) == 2) return 0; //` 0 其它`
if (c.relationseg(l[i]) == 1) x = 1; //` 1 圆在多边形里面,碰到了多边形边界`
}
return x;
}
//`旋转卡壳求凸包直径(最远点对)`
double rorating_calipers1()
{
double res = 0;
for (int i = 0, j = 1; i < n; i++)
{
while (dsgn(cross(p[i], p[i + 1], p[j]), cross(p[i], p[i + 1], p[j + 1])) < 0) j = (j + 1) % n;
res = max(res, max(p[i].distance(p[j]), p[i + 1].distance(p[j])));
}
return res;
}
//`旋转卡壳求最小矩形覆盖`
double rorating_calipers2(polygon &pt)
{
double res = 1e20;
for (int i = 0, a = 1, b = 1, c; i < n; i++)
{
while (dsgn(cross(p[i], p[i + 1], p[a]), cross(p[i], p[i + 1], p[a + 1])) < 0) a = (a + 1) % n;
while (dsgn(dot(p[i], p[i + 1], p[b]), dot(p[i], p[i + 1], p[b + 1])) < 0) b = (b + 1) % n;
if (!i) c = a;
while (dsgn(dot(p[i + 1], p[i], p[c]), dot(p[i + 1], p[i], p[c + 1])) < 0) c = (c + 1) % n;
double d = p[i].distance(p[i + 1]);
double H = cross(p[i], p[i + 1], p[a]) / d;
double R = dot(p[i], p[i + 1], p[b]) / d;
double L = dot(p[i + 1], p[i], p[c]) / d;
if (dsgn(res, (L + R - d) * H) > 0)
{
res = (L + R - d) * H;
pt.p[0] = p[i + 1] + (p[i] - p[i + 1]) * (L / d);
pt.p[1] = p[i] + (p[i + 1] - p[i]) * (R / d);
pt.p[2] = pt.p[1] + (p[i + 1] - p[i]).rotleft() * (H / d);
pt.p[3] = pt.p[0] + (p[i + 1] - p[i]).rotleft() * (H / d);
}
}
return res;
}
//`分治法求最近点对`
Point a[N];
double divide(int l, int r)
{
if (l == r) return 2e9;
if (l + 1 == r) return p[l].distance(p[r]);
int mid = l + r >> 1;
double d = min(divide(l, mid), divide(mid + 1, r));
int k = 0;
for (int i = l; i <= r; i++)
if (fabs(p[mid].x - p[i].x) < d) a[k++] = p[i];
// sort(a,a+k,[&](Point a,Point b)->bool {return a.y<b.y;});
for (int i = 0; i < k; i++)
for (int j = i + 1; j < k && a[j].y - a[i].y < d; j++)
d = min(d, a[j].distance(a[i]));
return d;
}
//`旋转卡壳求最大三角形面积`
double rotating_calipers3()
{
double res = 0;
for (int i = 0; i < n; i++)
{
int k = i + 1;
for (int j = i + 1; j < n; j++)
{
while (dsgn(cross(p[i], p[j], p[k]), cross(p[i], p[j], p[k + 1])) < 0) k = (k + 1) % n;
res = max(res, cross(p[i], p[j], p[k]));
}
}
return res / 2;
}
};
//`半平面交求凸多边形面积交`
double half_plane1(Line l[], int n)
{
double res = 0;
sort(l, l + n);
Line q[N];
Point p[N];
int h = 0, t = 0, k = 0;
q[t++] = l[0];
for (int i = 1; i < n; i++)
{
if (sgn(l[i].angle() - l[i - 1].angle()) == 0) continue;
while (h < t - 1 && l[i].relation(q[t - 1].crosspoint(q[t - 2])) == 2) t--;
while (h < t - 1 && l[i].relation(q[h].crosspoint(q[h + 1])) == 2) h++;
q[t++] = l[i];
}
while (h < t - 1 && l[h].relation(q[t - 1].crosspoint(q[t - 2])) == 2) t--;
q[t++] = q[h];
for (int i = h; i < t - 1; i++) p[k++] = q[i].crosspoint(q[i + 1]);
for (int i = 1; i < k - 1; i++) res += (p[i] - p[0]) ^ (p[i + 1] - p[0]);
return res / 2;
}
//`水平可见直线 从上向下看输出能看见哪些直线`
void half_plane2(Line l[], int n)
{
sort(l, l + n);
Line q[N];
int h = 0, t = 0, k = 0;
q[t++] = l[0];
for (int i = 1; i < n; i++)
{
if (sgn(l[i].angle() - l[i - 1].angle()) == 0) continue;
while (h < t - 1 && l[i].relation(q[t - 1].crosspoint(q[t - 2])) == 2) t--;
// while(h<t-1&&l[i].relation(q[h].crosspoint(q[h+1]))==2) h++;
q[t++] = l[i];
}
int ans[N];
for (int i = h; i < t; i++)
// for(auto j:q[i].id) ans[k++]=j;
sort(ans, ans + k);
cout << k << endl;
for (int i = 0; i < k; i++) cout << ans[i] << " ";
}
//`多边形内核`
bool half_plane3(Line l[], int n)
{
sort(l, l + n);
Line q[N];
int h = 0, t = 0;
q[t++] = l[0];
for (int i = 1; i < n; i++)
{
if (sgn(l[i].angle() - l[i - 1].angle()) == 0) continue;
while (h < t - 1 && l[i].relation(q[t - 1].crosspoint(q[t - 2])) == 2) t--;
while (h < t - 1 && l[i].relation(q[h].crosspoint(q[h + 1])) == 2) h++;
q[t++] = l[i];
}
while (h < t - 1 && l[h].relation(q[t - 1].crosspoint(q[t - 2])) == 2) t--;
return t - h >= 3;
}
//`最小圆覆盖`
circle increment(Point p[], int n)
{
random_shuffle(p, p + n);
circle ans;
ans.p = p[0], ans.r = 0;
for (int i = 1; i < n; i++)
if (ans.r < ans.p.distance(p[i]))
{
ans.p = p[i], ans.r = 0;
for (int j = 0; j < i; j++)
if (ans.r < ans.p.distance(p[j]))
{
ans.p = (p[i] + p[j]) / 2, ans.r = p[i].distance(p[j]) / 2;
for (int k = 0; k < j; k++)
if (ans.r < ans.p.distance(p[k]))
{
Point p1 = (p[i] + p[j]) / 2;
Point v1 = (p[i] - p[j]).rotright();
Point p2 = (p[i] + p[k]) / 2;
Point v2 = (p[i] - p[k]).rotright();
ans.p = Line(p1, p1 + v1).crosspoint(Line(p2, p2 + v2));
ans.r = ans.p.distance(p[i]);
}
}
}
return ans;
}
//`自适应辛普森积分`
inline double f(double x)
{ // 积分函数
return x;
}
double simpson(double l, double r)
{ // 辛普森公式
double mid = (l + r) / 2;
return (r - l) * (f(l) + f(r) + 4 * f(mid)) / 6;
}
double asr(double l, double r, double ans)
{ // 自适应
double mid = (l + r) / 2;
double a = simpson(l, mid), b = simpson(mid, r);
if (sgn(a + b - ans) == 0) return ans;
return asr(l, mid, a) + asr(mid, r, b);
}
struct Point3
{
double x, y, z;
//Point3(){}
Point3(double _x=0, double _y=0, double _z=0) { x = _x;y = _y;z = _z; }
void input() { cin>>x>>y>>z; }
bool operator==(const Point3 &b) const {return sgn(x-b.x)==0&&sgn(y-b.y)==0&&sgn(z-b.z)==0;}
bool operator<(const Point3 &b) const {return sgn(x-b.x)==0?(sgn(y-b.y)==0?sgn(z-b.z)<0:y<b.y):x<b.x;}
double len() { return sqrt(x * x + y * y + z * z); }
double len2() { return x * x + y * y + z * z; }
double distance(const Point3 &b) const {return sqrt((x-b.x)*(x-b.x)+(y-b.y)*(y-b.y)+(z-b.z)*(z-b.z));}
Point3 operator-(const Point3 &b) const { return Point3(x - b.x, y - b.y, z - b.z); }
Point3 operator+(const Point3 &b) const { return Point3(x + b.x, y + b.y, z + b.z); }
Point3 operator*(const double &k) const { return Point3(x * k, y * k, z * k); }
Point3 operator/(const double &k) const { return Point3(x / k, y / k, z / k); }
// 点乘
double operator*(const Point3 &b) const { return x * b.x + y * b.y + z * b.z; }
// 叉乘
Point3 operator^(const Point3 &b) const {return Point3(y*b.z-z*b.y,z*b.x-x*b.z,x*b.y-y*b.x);}
double rad(Point3 a, Point3 b)
{
Point3 p = (*this);
return acos(((a - p) * (b - p)) / (a.distance(p) * b.distance(p)));
}
//`化为长度为r的向量`
Point3 trunc(double r)
{
double l = len();
if (!sgn(l)) return *this;
r /= l;
return Point3(x * r, y * r, z * r);
}
};
struct Line3
{
Point3 s, e;
Line3() {}
Line3(Point3 _s, Point3 _e) { s = _s;e = _e; }
bool operator==(const Line3 v) { return (s == v.s) && (e == v.e); }
void input() { s.input();e.input(); }
double length() { return s.distance(e); }
// 点到直线距离
double dispointtoline(Point3 p) { return ((e - s) ^ (p - s)).len() / s.distance(e); }
// 点到线段距离
double dispointtoseg(Point3 p)
{
if (sgn((p - s) * (e - s)) < 0 || sgn((p - e) * (s - e)) < 0) return min(p.distance(s), e.distance(p));
return dispointtoline(p);
}
//`返回点p在直线上的投影`
Point3 lineprog(Point3 p) { return s + (((e - s) * ((e - s) * (p - s))) / ((e - s).len2())); }
//`p绕此向量逆时针arg角度`
Point3 rotate(Point3 p, double ang)
{
if (sgn(((s - p) ^ (e - p)).len()) == 0) return p;
Point3 f1 = (e - s) ^ (p - s);
Point3 f2 = (e - s) ^ (f1);
double len = ((s - p) ^ (e - p)).len() / s.distance(e);
f1 = f1.trunc(len);
f2 = f2.trunc(len);
Point3 h = p + f2;
Point3 pp = h + f1;
return h + ((p - h) * cos(ang)) + ((pp - h) * sin(ang));
}
//`点在直线上`
bool pointonseg(Point3 p) { return sgn(((s - p) ^ (e - p)).len()) == 0 && sgn((s - p) * (e - p)) == 0; }
};
struct Plane
{
Point3 a, b, c, o; //`平面上的三个点,以及法向量`
Plane() {}
Plane(Point3 _a, Point3 _b, Point3 _c)
{
a = _a;b = _b;c = _c;
o = pvec();
}
Point3 pvec() { return (b - a) ^ (c - a); }
//`ax+by+cz+d = 0`
Plane(double _a, double _b, double _c, double _d)
{
o = Point3(_a, _b, _c);
if (sgn(_a) != 0) a = Point3((-_d - _c - _b) / _a, 1, 1);
else if (sgn(_b) != 0) a = Point3(1, (-_d - _c - _a) / _b, 1);
else if (sgn(_c) != 0) a = Point3(1, 1, (-_d - _a - _b) / _c);
}
//`点在平面上的判断`
bool pointonplane(Point3 p) { return sgn((p - a) * o) == 0; }
//`两平面夹角`
double angleplane(Plane f) { return acos(o * f.o) / (o.len() * f.o.len()); }
//`平面和直线的交点,返回值是交点个数`
int crossline(Line3 u, Point3 &p)
{
double x = o * (u.e - a);
double y = o * (u.s - a);
double d = x - y;
if (sgn(d) == 0) return 0;
p = ((u.s * x) - (u.e * y)) / d;
return 1;
}
//`点到平面最近点(也就是投影)`
Point3 pointtoplane(Point3 p)
{
Line3 u = Line3(p, p + o);
crossline(u, p);
return p;
}
//`平面和平面的交线`
int crossplane(Plane f, Line3 &u)
{
Point3 oo = o ^ f.o;
Point3 v = o ^ oo;
double d = fabs(f.o * v);
if (sgn(d) == 0) return 0;
Point3 q = a + (v * (f.o * (f.a - a)) / d);
u = Line3(q, q + oo);
return 1;
}
};
struct CH3D
{
struct face
{
int a, b, c; // 表示凸包一个面上的三个点的编号
bool ok; // 表示该面是否属于最终的凸包上的面
};
int n;// 初始顶点数
Point3 P[N];
int num;// 凸包表面的三角形数
face F[8 * N];// 凸包表面的三角形
int g[N][N];
void input(int _n)
{
n = _n;
for (int i = 0; i < n; i++) P[i].input();
}
// 叉乘
Point3 cross(const Point3 &a, const Point3 &b, const Point3 &c) { return (b - a) ^ (c - a); }
//`三角形面积*2`
double area(Point3 a, Point3 b, Point3 c) { return ((b - a) ^ (c - a)).len(); }
//`四面体有向面积*6`
double volume(Point3 a, Point3 b, Point3 c, Point3 d) { return ((b - a) ^ (c - a)) * (d - a); }
//`正:点在面同向`
double dblcmp(Point3 &p, face &f)
{
Point3 p1 = P[f.b] - P[f.a];
Point3 p2 = P[f.c] - P[f.a];
Point3 p3 = p - P[f.a];
return (p1 ^ p2) * p3;
}
void deal(int p, int a, int b)
{
int f = g[a][b];
face add;
if (F[f].ok)
{
if (dblcmp(P[p], F[f]) > eps) dfs(p, f);
else
{
add.a = b;
add.b = a;
add.c = p;
add.ok = true;
g[p][b] = g[a][p] = g[b][a] = num;
F[num++] = add;
}
}
}
// 递归搜索所有应该从凸包内删除的面
void dfs(int p, int now)
{
F[now].ok = false;
deal(p, F[now].b, F[now].a);
deal(p, F[now].c, F[now].b);
deal(p, F[now].a, F[now].c);
}
bool same(int s, int t)
{
Point3 &a = P[F[s].a];
Point3 &b = P[F[s].b];
Point3 &c = P[F[s].c];
return fabs(volume(a, b, c, P[F[t].a])) < eps &&
fabs(volume(a, b, c, P[F[t].b])) < eps &&
fabs(volume(a, b, c, P[F[t].c])) < eps;
}
// 构建三维凸包
void create()
{
num = 0;
face add;
//***********************************
// 此段是为了保证前四个点不共面
bool flag = true;
for (int i = 1; i < n; i++)
{
if (!(P[0] == P[i]))
{
swap(P[1], P[i]);
flag = false;
break;
}
}
if (flag) return;
flag = true;
for (int i = 2; i < n; i++)
{
if (((P[1] - P[0]) ^ (P[i] - P[0])).len() > eps)
{
swap(P[2], P[i]);
flag = false;
break;
}
}
if (flag) return;
flag = true;
for (int i = 3; i < n; i++)
{
if (fabs(((P[1] - P[0]) ^ (P[2] - P[0])) * (P[i] - P[0])) > eps)
{
swap(P[3], P[i]);
flag = false;
break;
}
}
if (flag) return;
//**********************************
for (int i = 0; i < 4; i++)
{
add.a = (i + 1) % 4;
add.b = (i + 2) % 4;
add.c = (i + 3) % 4;
add.ok = true;
if (dblcmp(P[i], add) > 0) swap(add.b, add.c);
g[add.a][add.b] = g[add.b][add.c] = g[add.c][add.a] = num;
F[num++] = add;
}
for (int i = 4; i < n; i++)
for (int j = 0; j < num; j++)
if (F[j].ok && dblcmp(P[i], F[j]) > eps)
{
dfs(i, j);
break;
}
int tmp = num;
num = 0;
for (int i = 0; i < tmp; i++)
if (F[i].ok) F[num++] = F[i];
}
// 表面积
//`测试:HDU3528`
double area()
{
double res = 0;
if (n == 3)
{
Point3 p = cross(P[0], P[1], P[2]);
return p.len() / 2;
}
for (int i = 0; i < num; i++)
res += area(P[F[i].a], P[F[i].b], P[F[i].c]);
return res / 2.0;
}
double volume()
{
double res = 0;
Point3 tmp = Point3(0, 0, 0);
for (int i = 0; i < num; i++)
res += volume(tmp, P[F[i].a], P[F[i].b], P[F[i].c]);
return fabs(res / 6);
}
// 表面三角形个数
int triangle() { return num; }
// 表面多边形个数
//`测试:HDU3662`
int polygon()
{
int res = 0;
for (int i = 0; i < num; i++)
{
bool flag = true;
for (int j = 0; j < i; j++)
if (same(i, j))
{
flag = 0;
break;
}
res += flag;
}
return res;
}
// 重心
//`测试:HDU4273`
Point3 barycenter()
{
Point3 ans = Point3(0, 0, 0);
Point3 o = Point3(0, 0, 0);
double all = 0;
for (int i = 0; i < num; i++)
{
double vol = volume(o, P[F[i].a], P[F[i].b], P[F[i].c]);
ans = ans + (((o + P[F[i].a] + P[F[i].b] + P[F[i].c]) / 4.0) * vol);
all += vol;
}
ans = ans / all;
return ans;
}
// 点到面的距离
//`测试:HDU4273`
double ptoface(Point3 p, int i)
{
double tmp1 = fabs(volume(P[F[i].a], P[F[i].b], P[F[i].c], p));
double tmp2 = ((P[F[i].b] - P[F[i].a]) ^ (P[F[i].c] - P[F[i].a])).len();
return tmp1 / tmp2;
}
};
void solve()
{
int n;cin>>n;
Point3 p[N];
for(int i=0;i<n;i++) p[i].input();
Line3 l[N];int cnt=0;
for(int i=0;i<n;i++)
for(int j=0;j<i;j++) l[cnt++]=Line3(p[i],p[j]);
double res=1e20;
for(int i=0;i<cnt;i++)
for(int j=0;j<i;j++)
{
Point3 t=(l[i].e-l[i].s)^(l[j].e-l[j].s);
if(!t.len2()) continue;
double minv=1e20,maxv=-1e20;
for(int k=0;k<n;k++)
{
double dis=t*p[k];
minv=min(minv,dis);
maxv=max(maxv,dis);
}
double tt=(maxv-minv)/t.len();
res=min(res,tt);
}
cout<<fixed<<setprecision(15)<<res;
}
int main()
{
ios::sync_with_stdio(0), cin.tie(0), cout.tie(0);
int t = 1; //cin >> t;
while (t--) solve();
return 0;
}
Details
Tip: Click on the bar to expand more detailed information
Test #1:
score: 100
Accepted
time: 1ms
memory: 4500kb
input:
8 1 1 1 1 1 2 1 2 1 1 2 2 2 1 1 2 1 2 2 2 1 2 2 2
output:
1.000000000000000
result:
ok found '1.000000000', expected '1.000000000', error '0.000000000'
Test #2:
score: 0
Accepted
time: 1ms
memory: 4652kb
input:
5 1 1 1 1 2 1 1 1 2 1 2 2 2 1 1
output:
0.707106781186547
result:
ok found '0.707106781', expected '0.707106781', error '0.000000000'
Test #3:
score: -100
Wrong Answer
time: 3ms
memory: 4588kb
input:
50 973 1799 4431 1036 1888 4509 1099 1977 4587 1162 2066 4665 1225 2155 4743 1288 2244 4821 1351 2333 4899 1414 2422 4977 1540 2600 5133 1603 2689 5211 1666 2778 5289 1729 2867 5367 1792 2956 5445 1855 3045 5523 1918 3134 5601 1981 3223 5679 2044 3312 5757 2107 3401 5835 2170 3490 5913 2296 3668 606...
output:
100000000000000000000.000000000000000
result:
wrong answer 1st numbers differ - expected: '0.0000000', found: '100000000000000000000.0000000', error = '100000000000000000000.0000000'