FT03:苦熬Bentley-Ottmann算法的那些日子

发生了什么

为什么会讲这么一个内容,本心是上个月忙完了甲方的项目,这个月在做公司内部的算法库。苦熬苦业一个月其实还是蛮失败的,姑且记录一下自己都干了些什么。

起点其实是出于某些目的,手头的算法需要用到基于扫描线的线段求交,而这个算法的核心结构需要一个二叉搜索树来暂存线段。众所都周知,算法论文是不会替你考虑如何实现的。那些完美的复杂度分析几乎全部建立在“我认为该操作的复杂度为O(1)”这种其实并不符合现实的前提下。而令人感动的是,Bentley-Ottmann的这篇论文已经非常贴心的处理掉了和扫描线平行,以及多点段共交点这种非常规情况。但在处理退化情况的时候,常规的二叉搜索树,比如std::set就变的完全力不从心。简单来讲就是这样,现在说一下算法具体细节。

算法基于扫描线,这意味着在平面内存在一个扫描线,每扫到一个特定的位置,在扫描线的缓存区就会添加或是删除一些线段,并且将最近更新的线段和它左右两边的其他线段进行判断,确认有没有相交点。如果有,则将计算出的相交点也作为一个“特定的位置”留待后续扫描,直到所有待扫描位置全部处理完。在每次处理完线段时,扫描线上的线段要么增加——此时用增加的线段和两侧的线段进行判定,要么减少——此时用被删除的位置的两侧的线段进行判断,或是交换位置——因为两个线段有了新的相邻线段,因此用新的相邻线段进行相交判定。

为了处理多线共交点的退化情况,该算法通过将每个事件点绑定所有与之相关的线段进行进一步操作。具体来讲,每扫描到一个新的事件点,以该位置作为终点和交点的线段会被从扫描线上去掉,然后以该位置作为交点和起点的线段会被加入扫描线上。可以看出,作为交点的线段被去掉后重新插入,由于在插入时是按线段的延伸方向排序的(可以简单的通过延伸一段很小的距离来比较位置),因此作为交点的两个线段会被事实上交换位置从而获得新的相邻线段。

难点正在这里:如何通过事件点索引线段?强行搜索不是不行,但显然会因为在事件点上重复储存线段信息大量浪费空间。另一方面,即便可以在渐进意义上满足处理事件点的O(\log n)复杂度,但实际上缺引入了不大不小的常系数。如果扫描线的二叉搜索树允许直接操作树的节点,就可以让事件点在预处理阶段就和特定的节点进行绑定,进而在删除和添加过程中直接以O(1)的代价找到对应节点开始后续操作,最终省去处理事件点时对每一个关联的线段进行搜索的工作。

AVL树这部分内容可以说是纯抄了,蓝本是韦易笑的这个对linux红黑树和AVL树的项目做benchmark的项目。不得不说,AVL树真的比红黑树又好理解,又好实现,性能还有来有回。其实从这个角度来看,计算机行业的老古董其实也是又多又臭的。话题调转回来,做完这个之后肯定第一件事是做性能测试,而对比对象自然是std::set。为了公平竞争,自然就要封装一个类似的壳。其实这方面的内容真没啥好说的,主要就一个内容:迭代器。

迭代器

首先,迭代器是一个总的概念,STL的各种容器所用的迭代器并不是按照面向对象继承自同一个迭代器基类的类型,而是符合同一个规格的一系列类型,这个规格称之为iterator requirements,是C++标准中一系列具名要求中的一个,而针对迭代器的要求并不是只有这一个基础版本。由于迭代器本身根据所工作的容器的不同,是被分类为前向迭代器、双向迭代器、随机访问迭代器等一系列类别。STL用这些标签完成了TMP中非常经典的tag dispatch模式。在静态期通过tag选择特定的函数实现以达成优化要求。而这个tag本身也是requirement的一部分,与之相关的还有容器的内容类型、指针类型和引用类型。而对于内部实现,迭代器其实就是一个指针:在线性容器中就是数组中的某个位置的指针,在树和链表中就是一个节点的指针。通过重载解引用操作符来访问实际数据的位置。这里值得一提的是,range-based iteration固定使用的是begin()而不是cbegin()。这意味着哪怕你使用的是

for (const auto& v : vec)

该语句调用的也是begin()而不是cbegin()。这也就意味着begin()虽然按语义返回的是一个可修改的迭代器,但依然必须有一个const版本的重载。而cbegin()的实现其实就是对const重载版本的begin()的调用。

既然想到哪说到哪,这里额外提一句关于类函数const重载的问题。保持重用是一个基本原则,在const重载的场景下有一个比较典型的应用。首先将this指针转换成const版本,然后调用const版本的重载成员函数,接着将调用的结果用非const化,这是为数不多的建议的用例。好处是可以只保持在一处进行实现,其余的地方仅进行调用。在修改内容时可以保证逻辑统一,在修改中出现不易发现的错误。以一个向量类来讲

struct Vector {
	double data[Dimension];
   const double& operator[](size_t i) const {
       return data[i];
   }

   double& operator[](size_t i) {
       return const_cast<double&>(const_cast<const Vector&>(*this)[i]);
   }
};

这个技术往常是非常常用的,而在const_iterator和iterator中,由于语义上表示const和非const,我也曾一度想到过这个方法。但最终想明白的是,由于两个iterator使用的是同一个节点指针,因此其实根本不需要先const化再去除,而可以直接大胆的直接使用const_cast。提这几句姑且算是标记一下const_cast的正确使用方式吧。

后续

在有了一个类似STL的封装之后,我就进行月初那些诡异的benchmark结果。看来这个实现的效果相比于std::set有着非常显著,甚至可以达到平均15%以上的效率提升。虽然我依然怀疑内部有些逻辑错误或者其他问题就是了。在有了这个底层框架之后,算法的实现就比较直接了。抛开优先队列如何实现不谈,相比于STL的优先队列缺少的唯一接口就是搜索。事实上AVL树作为底层也可以实现一个优先队列,无非是有些杀鸡用牛刀了而已。

关于线段求交的整个流程,算法首先会在初始阶段首先将输入的一系列线段对象转化为包含线段的节点指针,而这些指针是放在std::vector中的。这一方面可以保证内存连续从而对缓存友好,另一方面在内存安全性上也是有保证的。在初始化事件队列时,遍历每一个端点,并在现存的队列中查询该位置是否已经被保存,以保证多个端点只要位置一样就可以被储存为同一个事件点。在处理事件点时,由于事件点已经绑定了与之相关的线段的节点,因此可以在线段插入扫描线后直接进行前后搜寻,以没有常系数的O(\log n)的时间找到相邻线段从而开始进行求交测试。

算法的内容基本到此为止了,在测试了几组基本情况、几组和扫描线平行线段以及几组共点线段的用例之后,大致认为该算法没什么问题了。但很可惜的是,该算法在两条线段共线重合处与第三条线段相交时,没有办法正确汇报所有相交的所有线段,甚至因为我一开始只实现了set而没有实现multiset,对于共线情况,靠后加入扫描线的线段甚至会被二叉搜索树拒绝加入。即便是允许加入,由于multiset的实现依然是二叉搜索树,相同key的不同节点依然处于逻辑上连续的几个节点上,从而在逻辑排序上,共线的N个线段中的N-1个必然被边缘的那一个线段阻挡,没有办法和相交的那条线段进行求交测试。而在交点处线段交换位置之后,一方面,第三条线段的位置会被交换到整个重合线段的另一边,依然丧失了和其他重合线段的求交机会;另一方面,由于处理事件点的流程是先汇报当前位置是否存在交点,再寻找该事件点之后的位置是否存在交点,因此这一流程本身就否定了在事件点出汇报本事件点上新发现的交点的可能性。

结果

最终,该方案还是落选了。这么一个不太成熟的set重制方案可以说是非常的失败了。在最终实现上,由于必须输出交点,因此共线的情况都被输出成共线范围的首尾两个交点。这部分的分支是独立于论文思路外单独做的魔改流程,说的更明确一点就是共线情况下并不会在事件队列里添加所谓的“相交事件点”而是直接汇报。值得说的是该实现在没有共线的情况下的确可以正确的输出所有想要的结果,而在存在共线的情况下,至少输出的交点位置是一个不少的,损失精确性的是哪些线段产生了相交这个信息。

这个项目本身是我上个月交接完了手头的工作,而预计工作两个月的项目提前完成了,所以手贱自己给自己找的事,因此无论如何这个月得给圆上。如果继续在Bentley-Ottmann这棵树上吊死属于是纯纯跟自己找不自在了,再加上上班偶尔摸个鱼,其实进度并不算快。所以最后只能找一个相对naive但是工程上实打实的提升性能的方案暂且丢上去用。只能说别再手贱了吧。

后记

对于这个笔记,我觉得我还是太刻意了。终归这是个写着玩的内容,好歹能让我以后忘了某些东西的时候有个参照,明白大概该去哪里复习。严格的分段落,希望能够每个段落一个主题对自己还是要求太高了。加上这次的内容其实是从公司的项目里折腾出来,领导让我“准备一下技术分享”的内容,所以还是更随着性子一点比较好。毕竟别搞成念幻灯片是才能说是演讲入门了。

另一方面,我决定还是把TMP从标题里拿掉。从实践的角度来讲,不能二进制化就是模板的原罪,公司内部要迭代内部代码库的版本或是为客户升级新版本,最省力最容易的方式就是保持接口不变的情况下只更新作为二进制的库文件。而不能二进制的模板就意味着对于底层代码库的每一层调用可能都需要重新编译重新发布。我开始可以理解那么多人期待module的理由了。然而现实是,C++20的module就是完全不可用的一坨,而另一方面,基层大量开发人员仍然在用C++11,甚至14都不愿意主动开(VS默认开14我觉得是一个很好的事情),这就意味着更新到更现代的C++真的就不止10年那么长。当然,想明白这个问题之后其实也就明白了至少在公司范围内,可能是很难有一组人和我一起做一个至少形式上类似CGAL那样的泛用算法库了。

这周某天午休时突然想到,身边这群做应用层的人都在用C++11。相比于C++98的老古董固然是先进很多了,但这些人最终也会工作十几年,也会在这一层形成进步的壁垒,这是无可否认的现实,我不能指望同事,至少说希望每一个同事和我一样想要“懂得更多”。我很高兴我处在一个允许我用C++20的公司,这至少允许我在可能的范围内优化一些东西。但我所在的终究还是传统行业,即便是传统行业里面已经算是很前卫的,但终归越传统的行业也会越保守。某种程度上一直都想去想去学Rust了,我还是比较欣赏那个社群的进步性的,无奈有些人的无视现实的进步性又太过空中楼阁。这是成长吗?我不好说。

下一篇应该是在两个月之后。公司项目要求我做一个关于no fit polygon的算法,大概最终也会形成这么一篇流水账式的笔记吧。

1. Bentley-Ottmann线段求交算法 doi:10.1109/TC.1979.1675432

2. "Computational Geometry 2nd edition" Mark de Berg "Chapter 2 Line segment intersection"

3. 改进参考书中的事件点到线段的索引设计 https://www.youtube.com/watch?v=I9EsN2DTnN8

4. AVL树的实现 https://github.com/skywind3000/avlmini  作者 https://www.zhihu.com/people/skywind3000

QR Code
微信扫一扫,欢迎咨询~

联系我们
武汉格发信息技术有限公司
湖北省武汉市经开区科技园西路6号103孵化器
电话:155-2731-8020 座机:027-59821821
邮件:tanzw@gofarlic.com
Copyright © 2023 Gofarsoft Co.,Ltd. 保留所有权利
遇到许可问题?该如何解决!?
评估许可证实际采购量? 
不清楚软件许可证使用数据? 
收到软件厂商律师函!?  
想要少购买点许可证,节省费用? 
收到软件厂商侵权通告!?  
有正版license,但许可证不够用,需要新购? 
联系方式 155-2731-8020
预留信息,一起解决您的问题
* 姓名:
* 手机:

* 公司名称:

姓名不为空

手机不正确

公司不为空