地球上大约有7.5x10 ^ 18的沙粒。这个故事是关于在数以千亿计的计算中找到相差约1e-18的方程式的变化。这就是7粒沙粒,与我们对整个地球的预期不同。
花了几天时间生成了数十亿字节的调试日志和GDB断点后,我终于在编译器中发现了一个非常特殊的错误。我认为这是一个有趣的故事。
早在2008年,我就开始开发名为空间人口模型(SPM)的科学建模平台。该软件旨在对海洋环境进行建模或模拟,以近似鱼类种群的健康状况。 SPM的输出用于创建科学报告,并提交给各国政府以设定商业捕鱼配额。
由于这些报告的政治性质,结果的可重复性和完整性至关重要。来自世界各地的科学家将使用该软件在不同的操作系统和编译器之间重新运行模型,以验证结果。
从2008年开始,就决定SPM将不再是多线程的。这是因为:没有C ++标准线程,普通计算机具有2个内核,与自动分化库不兼容)
快进到2020/21,在这种情况下,多核系统司空见惯。我最近为自己购买了AMD 5950X 16c / 32t CPU,并渴望将其应用于某些建模工作。 SPM仍被用作世界领先的空间建模平台,并且一直在不断更新。我一直对将并发引入SPM来发明科学建模的新方法感兴趣,这些方法可以改进基础数学。作为概念验证,我想开始对内部梯度计算进行处理。 SPM中的模型是用户定义的,并且可能很大(输入文本超过10万行),因此我们无法通过自动微分来解析确定梯度函数。我们使用迭代方法调整模型参数来计算梯度。这些调整是独立的,因此可以并行化。
在主线程(IEEE-754)上进行计算时,保持浮点运算的顺序完全相同。
SPM等科学建模平台的开发需要围绕可重复性的大量思考和工作。给定相同的输入文件,任何用户都必须能够重新运行模型并获得完全相同的答案。这似乎非常明显,但事实并非如此……
IEEE-754是浮点运算的标准。计算机中的浮点数是一个近似值。每当您看到浮点值时,都不太可能完全准确。这意味着浮点算法是非关联的。这样,A + B!= B +A。模型中操作发生的顺序将影响输出。对于小型程序或游戏而言,这无关紧要;对于操作数以千亿计的科学模型,这种不断增加的误差是必然的。
我们产生的每个线程必须以完全相同的顺序运行相同的操作。当线程将结果返回给主线程时,它们必须以确保将来所有操作均以完全相同的顺序进行的方式进行操作。这与模型中的线程数或参数数无关。
旁注:在科学模型中,我们使用双精度数字(double)而不是单精度(float)。单精度没有足够的精度来处理所需的计算数量。在模型的最后,浮点数学近似值所增加的误差非常大。这就限制了我们只能在GPU上使用CPU,因为当您不得不考虑编写GPU特定代码的复杂性时,GPU的双精度性能并不比现代CPU更好。
准备好任务列表后,我开始着手重构代码。一切都进行了修改,以使用中央模型类作为系统的父级。这使我可以生成与线程一样多的模型类。一切都在编译和运行而没有崩溃。是时候检查输出了。
略有差异,无需担心,因为我可能已经通过对它们进行变换来更改了某些方程的执行顺序。我正在使用运行时间少于5秒的小型模型进行测试,因此遍历每个实例化的类以检查代码不会花费太长时间。第二天下来,没有发现明显的东西。是时候通过GDB和消毒程序运行应用程序以查找问题了……那里什么也没有。
一个简单的模型将具有募集/繁殖,衰老和死亡。我开始减少和简化模型中的流程,以查看是否可以识别任何关键流程或参数。我删除了随机数生成并注释掉了大量复杂的数学运算法则,但是还是没有运气。
经过三天的尝试,发现问题后,我添加了大量的调试日志,并在GDB中逐步进行了查找,以寻找可以发现结果变化的地方。日志文件超过20MB,并且该模型仅在几次迭代后才会显示问题……为什么不立即呢?
我稍后在出现问题的模型中调用函数吗?还是该问题开始得较早,但以我正在打印的精度更高。我以15的精度打印调试输出。是时候将精度提高到20了。这基本上是您可以从精度上翻倍打印出来的最大精度。我再进行几次模型运行,结果发现结果的变化已经移到模型的更早位置……几乎立即。
该模型将加载配置文件,然后在内存中构造所有用户定义的对象。作为构建这些对象的一部分,将运行一堆计算以构建缓存。为了节省线程之间的数据移动,每个线程将重复相同的过程。每个新线程都会获取已加载的配置文件,创建所有对象并运行初始计算。 GDB显示,即使这些初始计算也有不同的结果。
我的工作文件夹中到处都是名为“ fuck”,“ fuck1”和“ fuckN”的文件,每个文件都具有用于小型模型运行的超过20MB的调试输出。在文件中运行了Diff,以查看发生了什么更改以及更改的数量有多小...但是还是没有运气..几乎所有的错误都立即发生了,但是我不知道为什么。
SPM使用GCC和MinGW64编译器。这使我们可以为Windows和Linux使用通用的Makefile创建几乎相同的代码。我一直喜欢TDM-GCC,因为它具有最新版本。在这种情况下,从2020年3月开始使用TDM-GCC 9.2.0。使用GCC 10.2.1在OpenSuSe Tumbleweed Linux上运行清理程序。这些都是GCC编译器套件的非常现代的版本。我从来没有将代码错误归咎于工具。
在这一点上,我正在努力寻找这个错误的原因。一些测试向我表明:
我的新二进制文件每次运行都会产生相同的结果,但与原始结果不同
我的新二进制文件将产生相同的结果,而不管线程数..从1到100
最后的观察使我感到好奇。当我仅使用一个线程运行新代码时,为什么我的新代码会产生不同的结果?是什么使其与原始二进制文件不同。仅执行一个线程将消除任何竞争条件。
在原始二进制文件中,我们不创建任何线程。我们加载模型,构建并运行它。在新代码中,即使只有一个线程,我们也将创建该线程,加载模型,构建并运行它。我在线程创建过程中正在做什么会导致此问题? GB日志文件和GDB的时间。
逐步执行每项操作以验证输出,将我带入所谓的选择性。这是一段静态代码,需要输入,调用标准的简单数学段并返回输出。这几行代码在一个线程中返回了不同的结果,而没有返回。是时候分析以下代码了:
double CLogisticSelectivity :: calculateResult(int Age){double dRet = 0.0;双倍dTemp =(dA50-年龄)/ dAto95;如果(dTemp> 5.0)dRet = 0.0;否则(dTemp< -5.0)dRet = dAlpha;否则dRet = dAlpha /(1.0 + pow(19.0,dTemp));返回dRet;}
如果我给代码Age = 1然后运行它,我得到以下结果:
是什么赋予了?这里没有什么应该在线程中有所不同。我没有传入值1的整数。我编辑了代码,使dA50和dAto95局部变量消除了任何潜在的竞争条件或线程异常,但在输出中仍然存在相同的错误。一切都在同一函数中初始化,唯一的外部值是Integer。奇怪的…
查看代码,唯一可能导致错误的是操作员和pow()调用。这些是C ++标准产品的一部分,因此它们中的任何一个极不可能出现错误。接下来是什么?让GodBolt尝试其他一些编译器,看看我是否得到了奇怪的结果。 MSVC ++…不,Clang…不,Linux上的GCC-Trunk…不。我只在TDM-GCC中得到这种奇怪的行为。是时候尝试使用其他MingW64编译器了。 Nuwen ...是的,Mingw64 ..是的,TDM-GCC ...是的。
当我创建一个新线程并在该线程中运行浮点运算时,得到的答案略有不同。
我去SourceForge记录一个错误。在等待一个月并且没有任何回应之后,我跳上IRC来问同样的问题。在24小时内,我的问题发表了评论,表明已解决。任何预先构建的发行版都没有意识到这一点。
新线程在构造时未调用_fpreset(),这“重新初始化了浮点数学包”。这意味着在线程内部时,浮点算术中使用的默认值会略有不同。
由于没有任何MinGW64的下游版本具有此修复程序,因此我不得不求助于自己在线程内调用_fpreset()。肮脏但有效的修复方法。
全部n全部。我花了大约5天的时间在我的代码中追踪这个错误。我生成了千兆字节的日志文件,并且不得不下降到地球上7.5粒沙子的精度。导致缺少关键函数调用的编译器原来是问题的原因。很多时候,在尝试寻找根本原因时,我发现自己在质疑自己编写代码,诊断错误并保持理智的能力。我很高兴找到答案,并有前进的道路。
使用AMD 5950X和32个线程测试代码时,求解模型的时间从23.5小时缩短到90分钟。这是一项重大的改进,仅此一项工作就为开发科学报告时的模型迭代提供了重要的新机会。
由安全架构师Zaita撰写 软件架构师-在LinkedIn上查看他