c++ - 当存在 float 精度问题时，如何实际解决凸包？

Question

假设曲线 y = x^2 上有 100000 个点。你想找到这些点的凸包。所有的坐标都是 float 。

在我的 graham 扫描实现中，我对 float 进行操作的唯一地方是当我最初根据坐标对所有点进行排序，然后我有一个函数来确定三个点是左转还是右转。

积分:

struct point {
   double x; 
   double y;
};

排序比较器:

inline bool operator() (const point &p1, const point &p2) {
    return (p1.x < p2.x) || (p1.x == p2.x && p1.y > p2.y);
}

左转/右转:

inline int ccw(point *p1, point *p2, point *p3) {
  double left  = (p1->x - p3->x)*(p2->y - p3->y);
  double right = (p1->y - p3->y)*(p2->x - p3->x);
  double res = left - right;
  return res > 0;
}

我的程序说在 100 000 个点中只有 68894 个是凸包的一部分。但是因为它们在曲线上，所以它们都应该是凸包的一部分。

对于您的眼睛来说，这没有任何区别。见下图。红点是凸包的一部分。

但如果您看得足够近，并放大这些点，您会发现其中一些是蓝色的，因此它们不包含在凸包中。

现在我的初步假设是浮点错误导致了这个问题。

我想我可以使用一个对 float 具有任意精度的外部库，但我更感兴趣的是我们在 C++ 中拥有的简单数据类型。

我怎样才能提高准确性？我读过有关 epsilon 的内容，但是在这里使用 epsilon 有什么帮助呢？我仍然假设一些彼此接近的点是相同的，所以我不会得到接近 100% 的准确度。

解决这个问题的最佳方法是什么？

Answer 1

如果您确实使用 (x, x^2) 形式的点，那么所有点都应该在凸包上是正确的.但是，三个点可能共线。如果您要移动它们或做任何其他奇怪的事情，这就会消失。

如果您要选择 100000 点，我建议使用 [-50000,49999] 中的整数。你的ccw函数将计算 left和 right是绝对值小于 2.5e14 < 2^53 的整数，因此不会发生舍入。

无论输入如何，基于坐标的排序都将正常工作。

对于一般输入，以下 ccw谓词有问题:

inline int ccw(point *p1, point *p2, point *p3) {
  double left  = (p1->x - p3->x)*(p2->y - p3->y);
  double right = (p1->y - p3->y)*(p2->x - p3->x);
  double res = left - right;
  return res > 0;
}

减法和乘法都可以进行舍入。如果您的所有点都位于 H*W 边界框内，则 x 坐标差异将以 H*eps/2 附近的绝对误差计算，而 y 坐标差异将以 W*eps/附近的绝对误差计算2.因此，产品将以 H*W*eps/2 左右的绝对误差计算。如果fabs(left - right) < 3*H*W*eps/2 , 你需要评估 left和 right更确切地说。 eps这是 2^-52。

如果 double 我可能会推荐只使用 MPFR比较不会告诉你任何事情。但是，您可以在没有的情况下做到这一点。 Kahan 求和的技巧可以让您从差异中得到低位，而 2²⁷+1 技巧可以帮助您准确计算乘积。