几个月前,我决定在曲线区域上实现集合运算。我有一本关于计算几何的规范教科书,它描述了由直线组成的多边形的方法,似乎其他项目已经将这些技术扩展到参数曲线。我想可能需要几个星期。
不幸的是,计算几何领域体现了一个根本性的矛盾。在几何学中,三角形的角度加起来正好等于π弧度,如果u与v成顺时针方向,那么v必然与u成逆时针方向。另一方面,计算机使用浮点表示法来嘲弄这些简单的欧几里得真理。
学术文献在很大程度上忽略了这一点。算法被证明是几何合理的,健壮的实现留给读者作为练习。这类似于散列冲突导致不可避免的数据丢失的世界,学术界的反应是隐含地假设存在完美的散列函数。如果一篇论文是建立在不切实际的假设上的,我们就不能根据它自己的条件来评价它;我们必须从经验上理解,当这些假设被打破时,它发挥的作用有多好。
考虑到这一点,我们现在可以看看现有的多边形裁剪算法,多边形裁剪是多边形集操作的术语。每种技术都是一个共同主题的变体:
关于这个主题的十几篇论文只有第三步不同。由于我们包含线段的决定通常在交点处反转,因此它们描述了各种方法,用于使用我们对大约一个线段的决定来通知我们关于相邻线段的决定。
我的教科书描述了一种使用双连通边列表的方法,这是一种通用的几何数据结构。我认为这意味着它可以被重用来解决其他问题,所以我开始了我的实现。
我在第一周就完成了实现,但它并不可靠。DCEL是可以增量更新的链接循环的集合。在执行集合操作时,我们递增地将原始循环集一分为二,然后确定哪些应该包含,哪些不应该包含。尽管我尽了最大的努力,我还是不断地寻找新的形状,使相邻的面纠缠在一起,创造出一条同时位于预期结果内外的莫比乌斯条带。
慢慢地,我意识到问题不在于数据结构,而是每篇论文都粉饰的第一步:找到所有的交叉点。DCEL假设顶点上的边具有总顺序:前一条边直接顺时针方向,下一条边直接逆时针方向。如果我们错过了交叉点,我们可能会得出结论,两条曲线相对于另一条曲线都是顺时针方向的,导致一切都分崩离析。
我开始寻找更好的方法来寻找交叉点,希望如果我找到一种足够准确的方法,我在DCEL上的工作可以挽救下来。不幸的是,我在文献中找到并在野外实现的方法并不比我一直使用的更好。我的数据结构需要精确的输入,没有内部矛盾,我无法交付。
在这一点上,我开始怀疑我是不是错过了什么根本的东西。我想,也许如果我仔细分析一下其他更成熟的库是如何处理我的病态形状的,我就可以倒着工作,看看我哪里出了问题。但是,当我开始将这些形状引入像Pap.js这样成熟的项目时,我发现它们同样失败了。
为了找到我的病理输入,我一直在使用基于属性的测试。给定形状的随机组合,我将执行区域中的点测试,并将其与通过单独查询每个形状并根据操作合并结果而生成的参考结果进行比较。大多数输入都运行良好,但是在几十万次输入之后,它不可避免地会发现一些故障。
事实证明,其他项目并不那么急于找到自己的失败模式。有些只有几个基于示例的测试,另一些则有一个包含数千个输入的静态套件,用于验证他们的更改。如果说我漏掉了什么,那就是似乎没有人预料到这些操作会特别强劲。
通过一个简单而令人抓狂的事实来理解浮点算术是最好的:a+(b-a)不一定等于b。它可能相等,也可能只相差一点点,其中“小”是相对于两个数字中较大的一个。这意味着当我们比较两个浮点数时,我们不能进行精确的比较,我们必须询问它们的差异是否小于某个ε值。
这个ε代表了我们已经听天由命的数字不确定性的程度。关于如何将这种不确定性降到最低,民间有大量的智慧,但事实仍然是,每个算术操作都会注入一点不确定性,而且随着每一次连续的操作,不确定性会累积增长。在处理小数字时,这种不确定性可能会使值本身相形见绌。
精确数学意义上的相交发生在曲线之间的距离正好为零的地方:
但在实际意义上,它是在曲线之间的距离足够接近于零的地方:
这至少有一个交集,但我们也可以很容易地在该重叠范围内返回三个或十个。这种不确定性令已发表的技术深恶痛绝,这些技术依赖于计算这些交叉点的数量来确定我们是在另一个形状的内部还是外部。单个虚假交叉点可能会导致整个结果消失。如果具有相似曲率的两个对象相互移动,结果将倾向于闪烁。
这些技术可能适用于需要较小ε的直线,但它们完全不适合参数曲线相对不精确的情况。
几个星期过去了,我的测试发现的错误从感觉上的随机变成了感觉上的恶意。浮点算术可能是确定性的,但当我看着屏幕,等待测试失败时,我想象着FPU中的一个恶魔在它们流过时轻推它们,试图将它们移出自我一致性的阈值。
有一天,我意识到这正是纠错码背后的叙述;我们假设我们的数据在途中被随机更改,我们希望将其恢复到一致的状态。第一次我好像不能把它做对,为什么不事后再修呢?
表面上看,应该只有三个交叉点,一个在左边,两个在右边。但由于上述原因,当曲线彼此收敛时,任何相交例程都可能在左侧找到多个交点:
左边的部分足够小,因此不够精确,以至于我们对这个问题的空间直觉只会误导我们。因此,最好将其视为图表:
现在,我们必须决定要包括哪些边,以及要排除哪些边。由于我们正在尝试查找并集,因此我们希望保留其他形状之外的所有线段:
这是一个结构良好的结果;只有一个循环,一旦我们移除该循环,就没有剩余的边了。但我们不能指望得到这样的幸运;左边的边缘可能已经屈服于浮点错误,并认为它们都在另一个里面:
这不是一个结构良好的结果;没有循环,只有一堆剩余的边。要使这一点保持一致,我们需要要么关闭循环,要么删除剩余部分。这些更改的成本由我们添加或删除的边的聚合长度来衡量。
最小更改集等效于悬挂顶点之间的最短路径。找到路径后,我们将其经过的每条边的包含内容反转。在本例中,它通过我们最初排除的一条边,因此我们将该边添加回来,并返回我们刚刚创建的循环。
在本例中,我们有一个完整的循环,但是一旦我们提取了它,就会有一个剩余的边:
在这里,我们再次搜索所有剩余的边,以寻找两个悬挂顶点之间的最短路径。因为它穿过我们剩下的边缘,所以我们把它移走。
每个浮点不一致都会作为这两种情况中的一种或它们的某种组合浮出水面。通过搜索悬挂边之间的最短路径,我们可以找到产生一致结果的尽可能小的编辑。当然,一致的结果不一定是正确的,但是浮点错误往往聚集在最小的边,这使得这是一个合理的启发式方法。更重要的是,它毫无问题地经受住了数以千万计的生成性测试用例。
我不确定这是不是一种新颖的方法,但至少它代表了对开源技术水平的有意义的改进。直观地说,在广泛的几何和数值算法中,这可能是避免ε地狱的一种手段。如果有人知道这方面的现有技术,我会很感兴趣的。
我在Artifex库上的工作正在进行中,但我希望它能证明它对其他人有用,并期待在不久的将来分享我自己的项目。
感谢Alex Engelberg、Elana Hashman、Angus Fletcher、Reid McKenzie和Zack Maril对本文早期草稿的反馈。