OpenBLAS项目与矩阵乘法优化

本篇文章7869字，读完约20分钟

说到矩阵计算，学过高等数学的人可能听说过，但如果他们不是这个领域的研究者，他们可能只会停留在“听说”的层面。在矩阵计算领域，开源项目openblas影响很大。除了被ibm和华为这样的大公司使用之外，它还吸引了全世界研究机构和开发者的注意。

最近，雷锋。(公开号码:雷锋。人工智能研究学会荣幸地邀请到了彭锋科技公司的创始人、开源项目的创始人和主要维护者张宪义，他将向我们介绍开源项目和矩阵乘法的优化。

嘉宾介绍了中国科学院张宪义博士、麻省理工学院博士后、开源项目的创始人和主要维护者、perfxlab彭锋技术公司的创始人。2016年获中国计算机联合会科技进步奖二等奖。

张宪义师从张云泉教授，获得中国科学院高性能计算博士学位。在他的博士学位期间，基于got扁球体的原始基础，他创建了openblas，一个开源矩阵计算库，并领导团队不断地修复和维护它。目前，他已经成为一个在矩阵计算领域有着巨大影响力的开源项目。

他想出了在麻省理工学院创业的主意。回到中国后，他创立了perfxlab彭锋技术公司进行“深度学习”，为计算机视觉和语音识别等公司提供集成的性能优化解决方案。雷锋[新知照]频道此前采访并报道了彭锋科技。

露天爆破工程课程内容介绍

矩阵乘法优化算法

逐步调优实现

以下是公开课的完整视频，共64分钟:

下面是公开课内容的正文和ppt安排。

雷锋的朋友们，你们好。我是张宪义。今天，我主要介绍我们的开源矩阵计算库openblas和矩阵乘法的优化。

首先，什么是爆破？

Blas是基本线性代数子程序的缩写，主要用于基本矩阵计算或向量计算。它分为三个层次:

Blas 1级，主要完成矢量之间的点运算或乘加运算，并计算相应的元素；

在blas 2级别，它主要生成矩阵和向量，如ppt的蓝色部分所示，矩阵a*向量x，并得到向量Y..此外，可能存在对称矩阵变形；

在blas级别3，它主要是矩阵和矩阵的计算，最典型的是矩阵*b矩阵，得到c矩阵。m*n的C矩阵由矩阵的宽度和高度得到。

为什么blas是一个非常重要的库或接口，是因为它是许多科学计算的核心之一。每年，当我们列出超级计算机的排名表时，我们都要做linpack测试，这个测试的很多部分都要做BLAS级矩阵和矩阵计算。此外，还有许多科学和工程模拟，它们在转换后成为矩阵运算。如果你能很好地优化矩阵，这将对整个应用程序的改进非常有帮助。

布拉斯和深度学习的关系，在近年来深度学习的普及之后，比如阿列克谢网，这一点大家都很清楚，如果我们做具体的性能划分，ppt上的这幅图来自于某篇论文。在削减之后，我们可以看到每个部分花费了多少时间，并且发现它的大部分时间花费在conv层上，而大量时间花费在fc层上。目前，体基层的常用实现是生成矩阵，即矩阵与矩阵的乘法，即BLAS级。全连通层一般转化为矩阵和向量的乘积，并完成blas运算。

也就是说，基于矩阵学习的深度学习，90%或更多的时间是由blas操作的。当然，随着新算法的出现，卷积层有一种特殊的3×3卷积核算法，或者可以通过fft类算法来实现，但是一般来说，扩展矩阵也是非常广泛的。

目前，blas只是一个定义好的实现接口，但是有很多具体的实现。从商业角度来看，有许多商业版本，如英特尔、amd、英伟达、ibm等。，它们基本上是经过优化以匹配自己的硬件。

至于开放源代码，有以下几种类型:gotoblas，这是众所周知的，与openblas有着深刻的关系，但它的开发在2010年被终止。有时间为每个人分析背后的原因。主要的开发者Goto离开了ut奥斯汀的研究小组，进入公司，然后停止了开发。地图集是由一所美国学校制作的，openblas是基于Goto Blast的，blis是Goto Blast离开后基于Goto Blast的扩展项目，目前还处于相对早期的阶段，成熟度还很低。

露天爆破已经有好几年的历史了。它开始进入2011年初。最初的原因是got opals放弃了。我们重新划分了社区分布，并继续开发和维护它。它的维护不是一个简单的错误修复过程。如果你想获得更好的性能，你需要继续跟踪硬件。例如，如果您创建了一个新的硬件体系结构或适应了一个更广泛的硬件体系结构，您必须优化它以获得更好的性能。Openblas被认为是目前世界上最好的开源矩阵计算库，并于去年获得中国计算机联合会科技进步奖二等奖。同时，它也进入了许多主流的linux安装包。例如，ubuntu有我们的openblas包，几乎所有你能想到的linux发行版都已经进入了，但这不是我们主动要做的。还有一套openhpc，这也是最近高性能计算的来源。

目前，openblas的开发是为了支持几乎所有主流的cpu处理器，并实现更好的优化性能。从操作系统来看，它基本上是由普通主流操作系统支持的。总的来说，从自适应处理器范围和支持的操作系统来看，它是开源库中实现最广泛的。

因此，露天爆破有许多用户。例如，有开源项目julia language、gnu octave等等；深入学习的人都熟悉Mxnet和caffe，它们可以选择openblas作为cpu端的计算后端；Ibm和arm也在他们的产品中使用OpenBlast。nvidia将OpenBlast列为与cpu进行比较测试的基准。还有一些以色列制造编译器的初创公司，以及一些国内互联网公司，如搜狗。前一段时间，我和搜狗公司的人谈过。我们的图书馆已经在网上稳定运行了一年多，所以工程质量还可以。

不久前，ibm制作了一个power ai软件框架，因为深度学习非常流行。我们可以看到在顶层有一些开源框架，我们的openblas在底层计算中。当然是为了携带他的服务器。

简而言之，blas库的性能越高越好。在英特尔的sandy bridge平台上，与mkl的性能相比，它基本上是相互重叠的，并且达到了相当的性能。

此图显示了龙芯上的一个结果，相对完整，整个blas是多线程的，性能经过充分测试，我们是性能更高的，提高了一到两倍。这是因为我们优化了龙芯3a，所以我们取得了非常好的效果。

刚才我们主要介绍了露天爆破的性能和效果。我们在github上举办的。欢迎对矩阵乘法或优化感兴趣的学生加入并贡献代码。我们公司愿意出一笔钱来支持这个项目向前发展。接下来，我将从一些技术性的干货开始，主要谈论每个人都对优化感兴趣的部分。我参考了这些矩阵乘法的教程，以及utaustinflame集团制作的教程。我基本上挖掘出了他的内容，并一步一步地带领每个人。如果我们从最简单的矩阵乘法到高性能的矩阵乘法，可能需要几个步骤。我们是怎么得到它的？或者它为什么被优化，以及每一步可以获得多少性能增益。这样，我们希望能为其他一些优化方案提供一些帮助。

让我们首先看看基本的实现。

我认为只要我学过线性代数等等，这种矩阵乘法就是一个很简单的问题。如果它被转换成C代码，它就是一个三重循环。我在这张图片中列出了一个三重循环。矩阵乘法的代码已经存在。它的函数是矩阵a*矩阵B，矩阵B加到矩阵C，其中C是结果矩阵。找到第一行，对应于这个位置的结果C，取出第一行中的每个元素，乘以第二列，相应地乘以，然后相加得到这个结果。但是，如果在当前处理器上运行，这种实现的性能将比优化的blas库差很多倍，甚至10倍。

为什么？我们做了最后的性能测试

这个数字也是从教程中截取的，它代表了一个性能图表，越高越好。这里的测试平台是英特尔酷睿i5，它只测试单线程，但忽略多线程。该初始实施可能为1千兆位/秒。随着规模变大，矩阵的性能正在下降。为什么？因为在实现的过程中，没有考虑缓存的原因，当矩阵小的时候，速度可以更快，而当矩阵大的时候，它就会掉下来，所以在图中有一个下降的过程。

这是一个非常原始和基本的实现。

再往前走一步，我们怎样才能稍微优化一下呢？我们需要调整矩阵的乘法顺序。我们在这里做了一个小块，在一个函数中单独提到了P，以点乘的形式写了它，每次做一个1*4的结果，然后把它分开出来成为一个函数。在P的这一步，我们需要稍微改变计算顺序，把I放在里面，J放在外面。我们为什么要改变这个背景？事实上，我们假设矩阵首先存储在列中，列项的值连续存储，并且行与行之间有间隔，这对于复制更有利。在这个实现之后，它对整体性能没有帮助，所以为什么为了以后的优化而以这种形式编写它，留下一个特定的空.

当矩阵很小时，在缓存中，性能基本不变，但是当它超过缓存时，性能稍微好一点，从刚才很低的值开始，也达到了接近1 gflop/s，主要原因是a(0，p)在一定程度上被重用，节省了一些缓存。另一方面，就循环利用而言，它相当于比这部分有一定的循环扩张，因此其效率有所提高。

该块的重用只从内存中读取一次，因此超出缓存的情况得到了一定程度的改善。

在此基础上，我们需要看看有什么更好的方法来做优化。我们的基准是如何在adddot1*4的基准上做到这一点。我们首先要考虑的是我们是否可以使用寄存器变量来代替操作内存。我可以申请一堆寄存器变量，比如c 00，01，它们在C语言中是双寄存器的，还有矩阵A的一部分，它也使用寄存器变量。

其余的操作是一些寄存器的变量。当然，B还是和以前一样，最后写回C。它完成的过程基本上和以前的做法一样，只是我们引入了寄存器变量，这样可以把更多的数据保存在寄存器中，而不是缓存中，这样可以减轻缓存的压力，提高一些性能。

可以看出，红线表示我们的优化性能，达到2个触发器/秒，蓝线表示之前的性能。这里的优化是使用寄存器来减少缓存消耗，并实现大约50%的性能提升。完成这一步后，我们如何对其进行优化？

我们刚刚用寄存器替换了A和C的部分，所以当B模仿这个部分时，我们有可能做一些优化吗？在前面的实现中，需要计算零件B的每个坐标计算，并且零件B访问的每个位置都需要计算一次，这引入了额外的开销，这是可以避免的。

我们使用指针的方式，在开始时指出相应的指针位置，只要每次计算时指针都连续移动，而不是每次读取一个位置再重新计算，这样速度会更快。

让我们看看。经过这种小规模的优化后，b指标的消耗降低了，4g的性能由2g提高到2g。这种改进对于小矩阵是有效的，但是如果所有的大矩阵都在高速缓存范围之外，它就不会有效。因此，假设矩阵都在缓存中，该块的性能也得到了极大的提高。

当你完成这一部分时，你会想，我已经用寄存器代替了A矩阵，B矩阵已经通过索引得到了改进。接下来我们还能做什么？

事实上，展开是一种常见的方式。

在最里面的循环中，它可以扩展到4次吗？这样做时，我们可以降低整个循环的成本，并使其更好地流动。

这一部分也将提高性能。这个数字比较了与初始阶段相比，在中间阶段取得了多大的改进。通过使用寄存器变量和指针，红线的性能在某个底部循环展开后得到了实现，与蓝线相比有了明显的提高，但这还没有完成，它只是一个基础。但是在1*4的级别，没有什么可以挖掘的了，所以我们需要在更高的级别找到一些方法。

在1*4处，a的值有一些重用，但是如果块被放大一点，例如，当它变成4*4时，也就是说，在每次计算中它是一个正方形而不是一列。对于整个优化，您可以重用您的访问内存，并充分发挥您的计算能力。

当我们变成这样一个4*4的小正方形时，加点函数也应该写成这个模式。当然，这部分也应该使用刚刚完成的1*4方法。以前有1个值，但现在有4个值。关于寄存器变量，在C语言中有16个变量，它们都是寄存器变量，并且都是通过指针优化的。

然而，如果我们这样做，整体性能的提高将是有限的。因为这只是一个初步的结果，我们可以看到对于小矩阵，在缓存的范围内没有改进。但是，在超过缓存后，大规模矩阵的性能会有一定的提高。

在优化结果为4*4的基础上，我们可以做进一步的优化和开发。我们为什么要切换到4*4优化是因为我们当前的cpu处理器基本上想获得高性能，必须使用矢量化指令，无论是旧的sse2、avx还是avx 2.0等等。对于cpu优化，如果我们想获得高性能，我们必须使用单指令多数据的矢量化指令。

经过4×4块后，在这种情况下对矢量指令更有利。这里，这部分循环的内容用向量指令重写。当然，我对这部分指令没有任何嵌入式汇编或纯汇编操作，所以我直接以英特尔内蕴的形式编写它。您可以看到，这部分被写成一个128位的sse，它是一个用于进行双精度浮点双精度运算的矩阵，数据是双精度的。从a中加载这两个值。最后两次加载被加载到这个向量寄存器中，并且B部分每次都加载由load复制的指令。剩下的一个基本上是用向量的内部变量来操作的，它看起来类似于前一个，所以它在编译时变成了矢量化指令。这两行被视为前一部分的值，而这一部分被视为后一部分的值

使用这个矢量寄存器指令后，他的性能明显提高，从刚才的4g可以达到6g以上，这是一个整体的变化。缓存中的矢量加速效果非常明显，基本上是翻倍的。

下一步要解决的是这个缓存的问题。问题是没有大块，超过缓存大小后性能会下降。为了解决这个问题，我们需要在更高的层次上做阻塞。

要转换成代码，在这一层做k分割，在下层做m分割。kc和mc都是参数。这些参数都是可调整的，应该根据l2cache的大小进行调整，然后每次对一个小块C进行矩阵乘法，相当于一个子问题。这个子问题的实现与刚才的4*4基本相同。

这部分阻塞的性能如图所示。蓝线表示没有阻塞的性能，红线表示完成后的性能。当问题的规模在缓存中时，它的性能提升基本上很小，但是当大规模矩阵被阻塞时，可以看出它得到了非常明显的提升，快了2倍。

对于大型矩阵，有必要优化其他问题，以充分利用缓存，使子问题更小，并提高其数据局部性。接下来，当我们进行阻塞时，如果只有代码级别改变，它基本上就完成了。

之后，引入一些操作以进一步优化。当我们分析程序的性能瓶颈时，A和B的内存访问都比较慢，而且矩阵中的很多内存访问都是不连续的，所以内存访问的性能比较差。一方面，缓存不能使用，另一方面，它也对tlb有影响。当然，C部分也有一些影响。矩阵c通常很大，所以没有办法打包，只有a和b。打包意味着，我有一部分连续记忆空，m*k，对应于以前的mc和kc。在这个记忆中，空，在每次计算之前，我从原始矩阵中读取一个需要的矩阵，并将它存储在一个连续的块中。packedmatrix a函数的具体实现非常简单。基本上，它从相应的位置被取出，并且当它被放置在连续的存储器地址中时结束。

你为什么要做这一步？这一步的意义在于，在传递包之后，下一步在adddot4*4中读取的元素不是从矩阵A中读取的，而是从包之后的缓冲区位置中读取的。一个优点是A矩阵已经被预热并放入中央处理器高速缓存；第二个优点是，您可以看到，当我存储时，这种连续存储也是连续读取，这非常高效，缓存效率也非常高。此外，传递包的步骤对于性能的提高非常明显。

此图是上一步的操作。打包后，大部分性能改善非常明显，超过50%，达到9gflop/s，除了最小矩阵大小没有影响，或引入额外开销，效果更差。

对于矩阵乘法的实现，填充矩阵是一种非常重要的优化方法。将来，每个人都会想，如果a已经做了包装，b是否也能做包装，同样的理由也是可能的，那就是把它复制到一个连续的空房间里。b部分的打包操作类似于a部分的打包操作，a部分从原始矩阵中读取其数据，并将其放入连续的空房间中，因此其内存访问是连续的。当然，在这一部分，由于B访问内存是一个流访问内存，这一部分的改进会稍微小一些，只比A高10%左右

对于矩阵乘法的实现，填充矩阵是一种非常重要的优化方法。将来，每个人都会想，如果a已经做了包装，b是否也能做包装，同样的理由也是可能的，那就是把它复制到一个连续的空房间里。b部分的打包操作类似于a部分的打包操作，a部分从原始矩阵中读取其数据，并将其放入连续的空房间中，因此其内存访问是连续的。当然，在这一部分，由于B访问内存是一个流访问内存，这一部分的改进会稍微小一些，只比A高10%左右

当您达到这一点时，与初始三重循环的性能改进相比，您的矩阵乘法的性能明显提高。如果你想做得更好，内核不仅可以写内部指令，还可以写嵌入式汇编，重新排列流水线，更多地使用硬件资源，这可能会增加10%。当然，这一部分对于实现blas是很重要的，将会详细介绍。

让我们整体回顾一下矩阵乘法的算法。我将把算法的这一部分放在最后，并从头开始一步一步地实现它，以便每个人在最后都能清楚地看到它。矩阵a *矩阵b得到矩阵c，它对应于最外面的环。当每走一步，它实际上都是阻塞的。刚才提到了阻塞的原因，这是为了匹配硬件结构，因为内存、缓存等。是分层的，越来越小。阻塞实际上提高了每层缓存的利用率。

让我们今天在这里分享它，谢谢。

问答问题1:什么是访问优化？

张宪义:内存访问优化解决了处理器的数据读取性能问题。经过计算，优化相对容易，但优化内存访问非常困难。密集矩阵乘法的数据比较规则，读取数据的顺序也比较规则，这就更容易优化。然而，我们已经对稀疏矩阵做了大量的优化，例如稀疏矩阵乘以向量的优化，这对于内存访问来说更加困难，因为没有办法预测下一次访问的位置，这使得优化变得困难。

问题2:2:openblas和其他矩阵库之间是什么关系？

张宪义:openblas和其他blas实现实际上完成了接口。blas只是接口的定义，它可以通过多种方式实现。我们认为我们在开源中的openblas实现非常好。如果它是标准的blas，它有一个参考实现，但它是一个非常简单的fortran实现，而且它的性能非常差，所以我们比他们快得多。Mkl是英特尔公司制造的炸弹，我们和他们是平等的。我们还没有完全测试过的Engien声称做得很好，但是他们所做的可能对x86平台有一些影响，但是我怀疑其他平台的效果。但是，我没有做具体的对比测试，所以我没有多少发言权。

问题3:从开始到掌握需要多长时间？

张宪义:如果我指导你，你可以在几个月内开始做一些事情。欢迎大家。

问题4:与高通相比，该库的表现如何？

张宪义:老实说，高通公司的图书馆还没有经过测试。我认为它自称更快，因为在32位arm上，我们的openblas没有进行矢量化优化，但是高通做了，所以它可能比我们更快，但是我们公司的内部perfblas是优化的。

问题5:分区的目的是什么？

张宪义:这是为了优化访问和内存。分区后，访问和内存集中在一定的区域，可以提高数据的局部性，从而提高缓存利用率和性能。

问题6:如果硬件很弱，你能玩神经网络吗？

张宪义:我们给出的数据之一是，在arm cortexa57平台上，我们在1.7ghz运行带有4个内核的alexnet模型，图片是150毫秒，相对较快。另一方面，我们仍然在做一些其他的模型，精简和优化，然后配合底层库的优化。在一个小模型下，运行一个小图片的印象只需要50毫秒。因此，在arm处理器上，仍有可能实现神经网络推理的实时定位。

问题7:内部版本和开源版本有很大区别吗？

张宪义:内部版本为深度学习做了一些区分，它的高性能可能会翻倍。优化的一部分是矩阵的大小，它不同于刚才提到的正方形矩阵。深度学习的矩阵大多是中小型的，某个维度会相对较小，所以使用的优化策略或阻塞策略会有所不同。有时直接做更好，不要阻塞或包装。

雷锋原创文章。严禁擅自转载。详情请参考转载说明。

来源：搜狐微门户