QOJ.ac

QOJ

ID题目提交者结果用时内存语言文件大小提交时间测评时间
#672018#7521. Find the GapForever_Young#WA 423ms3964kbC++2020.4kb2024-10-24 15:22:082024-10-24 15:22:10

Judging History

你现在查看的是最新测评结果

  • [2024-10-24 15:22:10]
  • 评测
  • 测评结果:WA
  • 用时:423ms
  • 内存:3964kb
  • [2024-10-24 15:22:08]
  • 提交

answer

#include <bits/stdc++.h>
using namespace std;
const double EPS = 1e-8;
struct Point3D {
  double x, y, z;
};
struct Line3D {
  Point3D a, b;
};
struct Plane {
  Point3D a, b, c;
};
struct PlaneF {
  // ax + by + cz + d = 0
  double a, b, c, d;
};
inline bool zero(double x) { return (x > 0 ? x : -x) < EPS; }
// 平方
inline double sqr(double d) { return d * d; }
// 计算cross product U x V
inline Point3D xmult(const Point3D& u, const Point3D& v) {
  Point3D ret;
  ret.x = u.y * v.z - v.y * u.z;
  ret.y = u.z * v.x - u.x * v.z;
  ret.z = u.x * v.y - u.y * v.x;
  return ret;
}
// 计算dot product U . V
inline double dmult(const Point3D& u, const Point3D& v) {
  return u.x * v.x + u.y * v.y + u.z * v.z;
}
// 矢量差 U - V
inline Point3D subt(const Point3D& u, const Point3D& v) {
  Point3D ret;
  ret.x = u.x - v.x;
  ret.y = u.y - v.y;
  ret.z = u.z - v.z;
  return ret;
}
// 取平面法向量
inline Point3D pvec(const Plane& s) {
  return xmult(subt(s.a, s.b), subt(s.b, s.c));
}
inline Point3D pvec(const Point3D& s1, const Point3D& s2, const Point3D& s3) {
  return xmult(subt(s1, s2), subt(s2, s3));
}
inline Point3D pvec(const PlaneF& p) {
  Point3D ret;
  ret.x = p.a;
  ret.y = p.b;
  ret.z = p.c;
  return ret;
}
// 两点距离
inline double dis(const Point3D& p1, const Point3D& p2) {
  return sqrt((p1.x - p2.x) * (p1.x - p2.x) + (p1.y - p2.y) * (p1.y - p2.y) +
              (p1.z - p2.z) * (p1.z - p2.z));
}
// 向量大小
inline double vlen(const Point3D& p) {
  return sqrt(p.x * p.x + p.y * p.y + p.z * p.z);
}
// 向量大小的平方
inline double sqrlen(const Point3D& p) {
  return (p.x * p.x + p.y * p.y + p.z * p.z);
}
// 判三点共线
bool dotsInline(const Point3D& p1, const Point3D& p2, const Point3D& p3) {
  return sqrlen(xmult(subt(p1, p2), subt(p2, p3))) < EPS;
}
// 判四点共面
bool dotsOnplane(const Point3D& a, const Point3D& b, const Point3D& c,
                 const Point3D& d) {
  return zero(dmult(pvec(a, b, c), subt(d, a)));
}
// 判点是否在线段上, 包括端点和共线
bool dotOnlineIn(const Point3D& p, const Line3D& l) {
  return zero(sqrlen(xmult(subt(p, l.a), subt(p, l.b)))) &&
         (l.a.x - p.x) * (l.b.x - p.x) < EPS &&
         (l.a.y - p.y) * (l.b.y - p.y) < EPS &&
         (l.a.z - p.z) * (l.b.z - p.z) < EPS;
}
bool dotOnlineIn(const Point3D& p, const Point3D& l1, const Point3D& l2) {
  return zero(sqrlen(xmult(subt(p, l1), subt(p, l2)))) &&
         (l1.x - p.x) * (l2.x - p.x) < EPS &&
         (l1.y - p.y) * (l2.y - p.y) < EPS && (l1.z - p.z) * (l2.z - p.z) < EPS;
}

// 判点是否在线段上, 不包括端点
bool dotOnlineEx(const Point3D& p, const Line3D& l) {
  return dotOnlineIn(p, l) &&
         (!zero(p.x - l.a.x) || !zero(p.y - l.a.y) || !zero(p.z - l.a.z)) &&
         (!zero(p.x - l.b.x) || !zero(p.y - l.b.y) || !zero(p.z - l.b.z));
}
bool dotOnlineEx(const Point3D& p, const Point3D& l1, const Point3D& l2) {
  return dotOnlineIn(p, l1, l2) &&
         (!zero(p.x - l1.x) || !zero(p.y - l1.y) || !zero(p.z - l1.z)) &&
         (!zero(p.x - l2.x) || !zero(p.y - l2.y) || !zero(p.z - l2.z));
}
// 判点是否在空间三角形上, 包括边界, 三点共线无意义
bool dotInplaneIn(const Point3D& p, const Plane& s) {
  return zero(vlen(xmult(subt(s.a, s.b), subt(s.a, s.c))) -
              vlen(xmult(subt(p, s.a), subt(p, s.b))) -
              vlen(xmult(subt(p, s.b), subt(p, s.c))) -
              vlen(xmult(subt(p, s.c), subt(p, s.a))));
}
bool dotInplaneIn(const Point3D& p, const Point3D& s1, const Point3D& s2,
                  const Point3D& s3) {
  return zero(vlen(xmult(subt(s1, s2), subt(s1, s3))) -
              vlen(xmult(subt(p, s1), subt(p, s2))) -
              vlen(xmult(subt(p, s2), subt(p, s3))) -
              vlen(xmult(subt(p, s3), subt(p, s1))));
}
// 判点是否在空间三角形上, 不包括边界, 三点共线无意义
bool dotInplaneEx(const Point3D& p, const Plane& s) {
  return dotInplaneIn(p, s) &&
         sqrlen(xmult(subt(p, s.a), subt(p, s.b))) > EPS &&
         sqrlen(xmult(subt(p, s.b), subt(p, s.c))) > EPS &&
         sqrlen(xmult(subt(p, s.c), subt(p, s.a))) > EPS;
}
bool dotInplaneEx(const Point3D& p, const Point3D& s1, const Point3D& s2,
                  const Point3D& s3) {
  return dotInplaneIn(p, s1, s2, s3) &&
         sqrlen(xmult(subt(p, s1), subt(p, s2))) > EPS &&
         sqrlen(xmult(subt(p, s2), subt(p, s3))) > EPS &&
         sqrlen(xmult(subt(p, s3), subt(p, s1))) > EPS;
}
// 判两点在线段同侧, 点在线段上返回0, 不共面无意义
bool sameSide(const Point3D& p1, const Point3D& p2, const Line3D& l) {
  return dmult(xmult(subt(l.a, l.b), subt(p1, l.b)),
               xmult(subt(l.a, l.b), subt(p2, l.b))) > EPS;
}
bool sameSide(const Point3D& p1, const Point3D& p2, const Point3D& l1,
              const Point3D& l2) {
  return dmult(xmult(subt(l1, l2), subt(p1, l2)),
               xmult(subt(l1, l2), subt(p2, l2))) > EPS;
}
// 判两点在线段异侧, 点在线段上返回0, 不共面无意义
bool oppositeSide(const Point3D& p1, const Point3D& p2, const Line3D& l) {
  return dmult(xmult(subt(l.a, l.b), subt(p1, l.b)),
               xmult(subt(l.a, l.b), subt(p2, l.b))) < -EPS;
}
bool oppositeSide(const Point3D& p1, const Point3D& p2, const Point3D& l1,
                  const Point3D& l2) {
  return dmult(xmult(subt(l1, l2), subt(p1, l2)),
               xmult(subt(l1, l2), subt(p2, l2))) < -EPS;
}
// 判两点在平面同侧, 点在平面上返回0
bool sameSide(const Point3D& p1, const Point3D& p2, const Plane& s) {
  return dmult(pvec(s), subt(p1, s.a)) * dmult(pvec(s), subt(p2, s.a)) > EPS;
}
bool sameSide(const Point3D& p1, const Point3D& p2, const Point3D& s1,
              const Point3D& s2, const Point3D& s3) {
  return dmult(pvec(s1, s2, s3), subt(p1, s1)) *
             dmult(pvec(s1, s2, s3), subt(p2, s1)) >
         EPS;
}
bool sameSide(const Point3D& p1, const Point3D& p2, const PlaneF& s) {
  return (s.a * p1.x + s.b * p1.y + s.c * p1.z + s.d) *
             (s.a * p2.x + s.b * p2.y + s.c * p2.z + s.d) >
         EPS;
}
// 判两点在平面异侧, 点在平面上返回0
bool oppositeSide(const Point3D& p1, const Point3D& p2, const Plane& s) {
  return dmult(pvec(s), subt(p1, s.a)) * dmult(pvec(s), subt(p2, s.a)) < -EPS;
}
bool oppositeSide(const Point3D& p1, const Point3D& p2, const Point3D& s1,
                  const Point3D& s2, const Point3D& s3) {
  return dmult(pvec(s1, s2, s3), subt(p1, s1)) *
             dmult(pvec(s1, s2, s3), subt(p2, s1)) <
         -EPS;
}
bool oppositeSide(const Point3D& p1, const Point3D& p2, const PlaneF& s) {
  return (s.a * p1.x + s.b * p1.y + s.c * p1.z + s.d) *
             (s.a * p2.x + s.b * p2.y + s.c * p2.z + s.d) <
         -EPS;
}
// 判两直线平行
bool parallel(const Line3D& u, const Line3D& v) {
  return sqrlen(xmult(subt(u.a, u.b), subt(v.a, v.b))) < EPS;
}
bool parallel(const Point3D& u1, const Point3D& u2, const Point3D& v1,
              const Point3D& v2) {
  return sqrlen(xmult(subt(u1, u2), subt(v1, v2))) < EPS;
}
// 判两平面平行
bool parallel(const Plane& u, const Plane& v) {
  return sqrlen(xmult(pvec(u), pvec(v))) < EPS;
}
bool parallel(const Point3D& u1, const Point3D& u2, const Point3D& u3,
              const Point3D& v1, const Point3D& v2, const Point3D& v3) {
  return sqrlen(xmult(pvec(u1, u2, u3), pvec(v1, v2, v3))) < EPS;
}
bool parallel(const PlaneF& u, const PlaneF& v) {
  return sqrlen(xmult(pvec(u), pvec(v))) < EPS;
}
// 判直线与平面平行
bool parallel(const Line3D& l, const Plane& s) {
  return zero(dmult(subt(l.a, l.b), pvec(s)));
}
bool parallel(const Point3D& l1, const Point3D& l2, const Point3D& s1,
              const Point3D& s2, const Point3D& s3) {
  return zero(dmult(subt(l1, l2), pvec(s1, s2, s3)));
}
bool parallel(const Line3D& l, const PlaneF& s) {
  return zero(dmult(subt(l.a, l.b), pvec(s)));
}
// 判两直线垂直
bool perpendicular(const Line3D& u, const Line3D& v) {
  return zero(dmult(subt(u.a, u.b), subt(v.a, v.b)));
}
bool perpendicular(const Point3D& u1, const Point3D& u2, const Point3D& v1,
                   const Point3D& v2) {
  return zero(dmult(subt(u1, u2), subt(v1, v2)));
}
// 判两平面垂直
bool perpendicular(const Plane& u, const Plane& v) {
  return zero(dmult(pvec(u), pvec(v)));
}
bool perpendicular(const Point3D& u1, const Point3D& u2, const Point3D& u3,
                   const Point3D& v1, const Point3D& v2, const Point3D& v3) {
  return zero(dmult(pvec(u1, u2, u3), pvec(v1, v2, v3)));
}
bool perpendicular(const PlaneF& u, const PlaneF& v) {
  return zero(dmult(pvec(u), pvec(v)));
}
// 判直线与平面垂直
bool perpendicular(const Line3D& l, const Plane& s) {
  return sqrlen(xmult(subt(l.a, l.b), pvec(s))) < EPS;
}
bool perpendicular(const Point3D& l1, const Point3D& l2, const Point3D& s1,
                   const Point3D& s2, const Point3D& s3) {
  return sqrlen(xmult(subt(l1, l2), pvec(s1, s2, s3))) < EPS;
}
bool perpendicular(const Line3D& l, const PlaneF& s) {
  return sqrlen(xmult(subt(l.a, l.b), pvec(s))) < EPS;
}
// 判两线段相交, 包括端点和部分重合
bool intersectIn(const Line3D& u, const Line3D& v) {
  if (!dotsOnplane(u.a, u.b, v.a, v.b)) {
    return 0;
  } else if (!dotsInline(u.a, u.b, v.a) || !dotsInline(u.a, u.b, v.b)) {
    return !sameSide(u.a, u.b, v) && !sameSide(v.a, v.b, u);
  } else {
    return dotOnlineIn(u.a, v) || dotOnlineIn(u.b, v) || dotOnlineIn(v.a, u) ||
           dotOnlineIn(v.b, u);
  }
}
bool intersectIn(const Point3D& u1, const Point3D& u2, const Point3D& v1,
                 const Point3D& v2) {
  if (!dotsOnplane(u1, u2, v1, v2)) {
    return 0;
  } else if (!dotsInline(u1, u2, v1) || !dotsInline(u1, u2, v2)) {
    return !sameSide(u1, u2, v1, v2) && !sameSide(v1, v2, u1, u2);
  } else {
    return dotOnlineIn(u1, v1, v2) || dotOnlineIn(u2, v1, v2) ||
           dotOnlineIn(v1, u1, u2) || dotOnlineIn(v2, u1, u2);
  }
}
// 判两线段相交, 不包括端点和部分重合
bool intersectEx(const Line3D& u, const Line3D& v) {
  return dotsOnplane(u.a, u.b, v.a, v.b) && oppositeSide(u.a, u.b, v) &&
         oppositeSide(v.a, v.b, u);
}
bool intersectEx(const Point3D& u1, const Point3D& u2, const Point3D& v1,
                 const Point3D& v2) {
  return dotsOnplane(u1, u2, v1, v2) && oppositeSide(u1, u2, v1, v2) &&
         oppositeSide(v1, v2, u1, u2);
}
// 判线段与空间三角形相交, 包括交于边界和(部分)包含
bool intersectIn(const Line3D& l, const Plane& s) {
  return !sameSide(l.a, l.b, s) && !sameSide(s.a, s.b, l.a, l.b, s.c) &&
         !sameSide(s.b, s.c, l.a, l.b, s.a) &&
         !sameSide(s.c, s.a, l.a, l.b, s.b);
}
bool intersectIn(const Point3D& l1, const Point3D& l2, const Point3D& s1,
                 const Point3D& s2, const Point3D& s3) {
  return !sameSide(l1, l2, s1, s2, s3) && !sameSide(s1, s2, l1, l2, s3) &&
         !sameSide(s2, s3, l1, l2, s1) && !sameSide(s3, s1, l1, l2, s2);
}
// 判线段与空间三角形相交, 不包括交于边界和(部分)包含
bool intersectEx(const Line3D& l, const Plane& s) {
  return oppositeSide(l.a, l.b, s) && oppositeSide(s.a, s.b, l.a, l.b, s.c) &&
         oppositeSide(s.b, s.c, l.a, l.b, s.a) &&
         oppositeSide(s.c, s.a, l.a, l.b, s.b);
}
bool intersectEx(const Point3D& l1, const Point3D& l2, const Point3D& s1,
                 const Point3D& s2, const Point3D& s3) {
  return oppositeSide(l1, l2, s1, s2, s3) && oppositeSide(s1, s2, l1, l2, s3) &&
         oppositeSide(s2, s3, l1, l2, s1) && oppositeSide(s3, s1, l1, l2, s2);
}
// 计算两直线交点, 注意事先判断直线是否共面和平行!
// 线段交点请另外判线段相交(同时还是要判断是否平行!)
#include <algorithm>
using namespace std;
Point3D intersection(Point3D u1, Point3D u2, Point3D v1, Point3D v2) {
  double dxu = u2.x - u1.x;
  double dyu = u2.y - u1.y;
  double dzu = u2.z - u1.z;
  double dxv = v2.x - v1.x;
  double dyv = v2.y - v1.y;
  double dzv = v2.z - v1.z;
  double t;
  if (!zero(dxu * dyv - dyu * dxv)) {
    t = (dyv * (v1.x - u1.x) + dxv * (u1.y - v1.y)) / (dxu * dyv - dyu * dxv);
  } else if (!zero(dxu * dzv - dzu * dxv)) {
    t = (dzv * (v1.x - u1.x) + dxv * (u1.z - v1.z)) / (dxu * dzv - dzu * dxv);
  } else {
    t = (dzv * (v1.y - u1.y) + dyv * (u1.z - v1.z)) / (dyu * dzv - dzu * dyv);
  }
  Point3D ret;
  ret.x = u1.x + dxu * t;
  ret.y = u1.y + dyu * t;
  ret.z = u1.z + dzu * t;
  return ret;
}
// 计算直线与平面交点, 注意事先判断是否平行, 并保证三点不共线!
// 线段和空间三角形交点请另外判断
Point3D intersection(const Line3D& l, const Plane& s) {
  Point3D ret = pvec(s);
  double t = (ret.x * (s.a.x - l.a.x) + ret.y * (s.a.y - l.a.y) +
              ret.z * (s.a.z - l.a.z)) /
             (ret.x * (l.b.x - l.a.x) + ret.y * (l.b.y - l.a.y) +
              ret.z * (l.b.z - l.a.z));
  ret.x = l.a.x + (l.b.x - l.a.x) * t;
  ret.y = l.a.y + (l.b.y - l.a.y) * t;
  ret.z = l.a.z + (l.b.z - l.a.z) * t;
  return ret;
}
Point3D intersection(const Point3D& l1, const Point3D& l2, const Point3D& s1,
                     const Point3D& s2, const Point3D& s3) {
  Point3D ret = pvec(s1, s2, s3);
  double t =
      (ret.x * (s1.x - l1.x) + ret.y * (s1.y - l1.y) + ret.z * (s1.z - l1.z)) /
      (ret.x * (l2.x - l1.x) + ret.y * (l2.y - l1.y) + ret.z * (l2.z - l1.z));
  ret.x = l1.x + (l2.x - l1.x) * t;
  ret.y = l1.y + (l2.y - l1.y) * t;
  ret.z = l1.z + (l2.z - l1.z) * t;
  return ret;
}
Point3D intersection(const Line3D& l, const PlaneF& s) {
  Point3D ret = subt(l.b, l.a);
  double t = -(dmult(pvec(s), l.a) + s.d) / (dmult(pvec(s), ret));
  ret.x = ret.x * t + l.a.x;
  ret.y = ret.y * t + l.a.y;
  ret.z = ret.z * t + l.a.z;
  return ret;
}
// 计算两平面交线, 注意事先判断是否平行, 并保证三点不共线!
Line3D intersection(const Plane& u, const Plane& v) {
  Line3D ret;
  ret.a = parallel(v.a, v.b, u.a, u.b, u.c)
              ? intersection(v.b, v.c, u.a, u.b, u.c)
              : intersection(v.a, v.b, u.a, u.b, u.c);
  ret.b = parallel(v.c, v.a, u.a, u.b, u.c)
              ? intersection(v.b, v.c, u.a, u.b, u.c)
              : intersection(v.c, v.a, u.a, u.b, u.c);
  return ret;
}
Line3D intersection(const Point3D& u1, const Point3D& u2, const Point3D& u3,
                    const Point3D& v1, const Point3D& v2, const Point3D& v3) {
  Line3D ret;
  ret.a = parallel(v1, v2, u1, u2, u3) ? intersection(v2, v3, u1, u2, u3)
                                       : intersection(v1, v2, u1, u2, u3);
  ret.b = parallel(v3, v1, u1, u2, u3) ? intersection(v2, v3, u1, u2, u3)
                                       : intersection(v3, v1, u1, u2, u3);
  return ret;
}
// 点到直线距离
double disptoline(const Point3D& p, const Line3D& l) {
  return vlen(xmult(subt(p, l.a), subt(l.b, l.a))) / dis(l.a, l.b);
}
double disptoline(const Point3D& p, const Point3D& l1, const Point3D& l2) {
  return vlen(xmult(subt(p, l1), subt(l2, l1))) / dis(l1, l2);
}
// 点到直线最近点
Point3D ptoline(const Point3D& p, const Line3D& l) {
  Point3D ab = subt(l.b, l.a);
  double t = -dmult(subt(p, l.a), ab) / sqrlen(ab);
  ab.x *= t;
  ab.y *= t;
  ab.z *= t;
  return subt(l.a, ab);
}
// 点到平面距离
double disptoplane(const Point3D& p, const Plane& s) {
  return fabs(dmult(pvec(s), subt(p, s.a))) / vlen(pvec(s));
}
double disptoplane(const Point3D& p, const Point3D& s1, const Point3D& s2,
                   const Point3D& s3) {
  return fabs(dmult(pvec(s1, s2, s3), subt(p, s1))) / vlen(pvec(s1, s2, s3));
}
double disptoplane(const Point3D& p, const PlaneF& s) {
  return fabs((dmult(pvec(s), p) + s.d) / vlen(pvec(s)));
}
// 点到平面最近点
Point3D ptoplane(const Point3D& p, const PlaneF& s) {
  Line3D l;
  l.a = p;
  l.b = pvec(s);
  l.b.x += p.x;
  l.b.y += p.y;
  l.b.z += p.z;
  return intersection(l, s);
}
// 直线到直线距离
double dislinetoline(const Line3D& u, const Line3D& v) {
  Point3D n = xmult(subt(u.a, u.b), subt(v.a, v.b));
  return fabs(dmult(subt(u.a, v.a), n)) / vlen(n);
}
double dislinetoline(const Point3D& u1, const Point3D& u2, const Point3D& v1,
                     const Point3D& v2) {
  Point3D n = xmult(subt(u1, u2), subt(v1, v2));
  return fabs(dmult(subt(u1, v1), n)) / vlen(n);
}
// 直线到直线的最近点对
// p1在u上,p2在v上,p1到p2是uv之间的最近距离
// 注意,保证两直线不平行
void linetoline(const Line3D& u, const Line3D& v, Point3D& p1, Point3D& p2) {
  Point3D ab = subt(u.b, u.a), cd = subt(v.b, v.a), ac = subt(v.a, u.a);
  double r = (dmult(ab, cd) * dmult(ac, ab) - sqrlen(ab) * dmult(ac, cd)) /
             (sqrlen(ab) * sqrlen(cd) - sqr(dmult(ab, cd)));
  p2.x = v.a.x + r * cd.x;
  p2.y = v.a.y + r * cd.y;
  p2.z = v.a.z + r * cd.z;
  p1 = ptoline(p2, u);
}
// 两直线夹角cos值
double angleCos(const Line3D& u, const Line3D& v) {
  return dmult(subt(u.a, u.b), subt(v.a, v.b)) / vlen(subt(u.a, u.b)) /
         vlen(subt(v.a, v.b));
}
double angleCos(const Point3D& u1, const Point3D& u2, const Point3D& v1,
                const Point3D& v2) {
  return dmult(subt(u1, u2), subt(v1, v2)) / vlen(subt(u1, u2)) /
         vlen(subt(v1, v2));
}
// 两平面夹角cos值
double angleCos(const Plane& u, const Plane& v) {
  return dmult(pvec(u), pvec(v)) / vlen(pvec(u)) / vlen(pvec(v));
}
double angleCos(const Point3D& u1, const Point3D& u2, const Point3D& u3,
                const Point3D& v1, const Point3D& v2, const Point3D& v3) {
  return dmult(pvec(u1, u2, u3), pvec(v1, v2, v3)) / vlen(pvec(u1, u2, u3)) /
         vlen(pvec(v1, v2, v3));
}
double angleCos(const PlaneF& u, const PlaneF& v) {
  return dmult(pvec(u), pvec(v)) / (vlen(pvec(u)) * vlen(pvec(v)));
}
// 直线平面夹角sin值
double angleSin(const Line3D& l, const Plane& s) {
  return dmult(subt(l.a, l.b), pvec(s)) / vlen(subt(l.a, l.b)) / vlen(pvec(s));
}
double angleSin(const Point3D& l1, const Point3D& l2, const Point3D& s1,
                const Point3D& s2, const Point3D& s3) {
  return dmult(subt(l1, l2), pvec(s1, s2, s3)) / vlen(subt(l1, l2)) /
         vlen(pvec(s1, s2, s3));
}
double angleSin(Line3D l, const PlaneF& s) {
  return dmult(subt(l.a, l.b), pvec(s)) /
         (vlen(subt(l.a, l.b)) * vlen(pvec(s)));
}
// 平面方程形式转化 Plane -> PlaneF
PlaneF planeToPlaneF(const Plane& p) {
  PlaneF ret;
  Point3D m = xmult(subt(p.b, p.a), subt(p.c, p.a));
  ret.a = m.x;
  ret.b = m.y;
  ret.c = m.z;
  ret.d = -m.x * p.a.x - m.y * p.a.y - m.z * p.a.z;
  return ret;
}

int n;
Point3D p[51];

double eq(double a, double b) { return fabs(a - b) < 1e-7; }

bool cmp(const Point3D& a, const Point3D& b) {
  if (eq(a.x, b.x)) {
    if (eq(a.y, b.y)) {
      return a.z < b.z;
    }
    return a.y < b.y;
  }
  return a.x < b.x;
}

double gao(const Point3D& vec) {
  Line3D line = {Point3D({0, 0, 0}), vec};
  Point3D maxp {-1e9, -1e9, -1e9}, minp {1e9, 1e9, 1e9};
  for (int i = 1; i <= n; i++) {
    auto pp = ptoline(p[i], line);
    if (cmp(pp, minp)) {
      minp = pp;
    }
    if (cmp(maxp, pp)) {
      maxp = pp;
    }
  }
  /*
  printf("line = (%.5f,%.5f,%.5f) (%.5f,%.5f,%.5f)\n", line.a.x, line.a.y,
         line.a.z, line.b.x, line.b.y, line.b.z);
  printf("min=(%.5f,%.5f,%.5f), max=(%.5f,%.5f,%.5f)\n", minp.x, minp.y, minp.z,
         maxp.x, maxp.y, maxp.z);
  printf("d = %.5f\n", dis(maxp, minp));
  */
  return dis(maxp, minp);
}

int main() {
  scanf("%d", &n);
  for (int i = 1; i <= n; i++) {
    scanf("%lf%lf%lf", &p[i].x, &p[i].y, &p[i].z);
  }
  double ans = 1e9;
  for (int i = 1; i <= n; i++) {
    for (int j = i + 1; j <= n; j++) {
      for (int k = j + 1; k <= n; k++) {
        Plane pl;
        pl.a = p[i];
        pl.b = p[j];
        pl.c = p[k];
        Point3D vec = pvec(pl);
        ans = min(ans, gao(vec));
      }
    }
  }
  for (int i = 1; i <= n; i++) {
    for (int j = i + 1; j <= n; j++) {
      for (int k = 1; k <= n; k++) {
        if (k == i || k == j) continue;
        for (int l = k + 1; l <= n; l++) {
          if (l == i || l == j) continue;
          Point3D delta = subt(p[l], p[k]);
          Point3D p3 = {p[i].x + delta.x, p[i].y + delta.y, p[i].z + delta.z};
          if (zero(dis(p3, p[j]))) continue;
          Plane pl = {p[i], p[j], p3};
          Point3D vec = pvec(pl);
          ans = min(ans, gao(vec));
        }
      }
    }
  }
  printf("%.12f\n", ans);
}

详细

Test #1:

score: 100
Accepted
time: 1ms
memory: 3892kb

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.000000000000

result:

ok found '1.000000000', expected '1.000000000', error '0.000000000'

Test #2:

score: 0
Accepted
time: 0ms
memory: 3940kb

input:

5
1 1 1
1 2 1
1 1 2
1 2 2
2 1 1

output:

0.707106781187

result:

ok found '0.707106781', expected '0.707106781', error '0.000000000'

Test #3:

score: -100
Wrong Answer
time: 423ms
memory: 3964kb

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:

1000000000.000000000000

result:

wrong answer 1st numbers differ - expected: '0.0000000', found: '1000000000.0000000', error = '1000000000.0000000'