当存在浮动精度问题时，如何实际解决凸包？答案

【问题标题】：how can you practically solve the convex hull when there are floating precision issues?当存在浮动精度问题时，如何实际解决凸包？
【发布时间】：2014-11-25 04:58:29
【问题描述】：

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

在我的格雷厄姆扫描实现中，我对浮点数进行操作的唯一地方是，我最初按坐标对所有点进行排序，然后我有一个函数可以确定三个点是左转还是右转。

积分：

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 个是凸包的一部分。但是由于它们在曲线上，所以它们都应该是凸包的一部分。

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

但是如果你看的足够近，放大这些点，你会发现其中一些是蓝色的，所以它们不包含在凸包中。

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

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

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

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

【问题讨论】：

你试过long double吗？
您的凸包算法是否可能是正确的，但是在您最初评估 y = x^2 时发生了舍入？
首先，没有有限的精度可以让您表示实数。其次，您绘制的点是近似值。第三，真正的构造实数（无限精度）无法按照您想要比较它们的方式进行比较。第五，你为什么在乎？（你对如此重要的船体有什么目的）
mch，是的，但不幸的是我看不出有什么不同。 ajclinto，我没想到，大概也是这样。 Yakk，我的教授给了我们一个任务来实现不同的凸包算法，我已经实现了所有这些算法，但不知道如何处理浮点错误。
这不是一个通用的答案，但您可以使用两个整数来精确存储一个有理数，并在有理空间中进行所有数学运算。它不适用于sqrt 或sin，但只要自变量也是有理数，它就适用于任何具有有理系数的多项式。

标签： c++ floating-point precision floating-accuracy convex-hull

【解决方案1】：

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

如果您要选择 100000 点，我建议使用 [-50000,49999] 中的整数。您的 ccw 函数会将 left 和 right 计算为绝对值小于 2.5e14

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

对于一般输入，以下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 技巧可以帮助您准确计算乘积。

【讨论】：

【解决方案2】：

对于浮点数学，您经常需要引入“容差”的概念，有时表示为 epsilon。在您的情况下，您可以将 ccw() 函数设为三值：true/false/indeterminate。然后当你试图发现一个新点是否可以成为凸包的一部分时，你会问“它是 ccw=true 还是不确定的”，无论哪种方式你都接受了这个点。当斜率太接近直线而无法确定时，就会出现不确定的结果。

【讨论】：