图形处理单元上的高性能稀疏三角求解的制作方法

背景技术：

三角矩阵是一种类型的方阵，所述方阵在矩阵的主对角线上方或下方具有仅零元素。下三角矩阵在主对角线上方具有仅零元素，使得矩阵中的任何非零元素都在下三角中、在主对角线上或在主对角线下方。上三角矩阵在主对角线下方具有仅零元素，使得矩阵中的任何非零元素都在上三角中、在主对角线上或在主对角线上方。三角矩阵可用于在线性代数的领域中表示方程组。

稀疏三角矩阵是在填充的三角中具有大量零元素的三角矩阵；例如，稀疏的下三角矩阵在其下三角中具有一个或多个零值。稀疏三角求解(spts)是用于求解方程式ax＝y中的向量x的过程，其中a是具有n个行和n个列的稀疏三角矩阵，x是n个未知值的向量，且y是n个已知值的向量。如果矩阵a中唯一的非零值位于主对角线上和所述对角线的一侧，则可使用替换来求解向量x。求解向量条目x[n]依赖于已在下三角矩阵中进行前向替换的情况下求解所有先前向量条目(例如，x[0]-x[n-1])。然而，如果矩阵稀疏，则三角矩阵值中的一些也为零，并且可在并行处理器上并行地求解多行。

附图说明

在附图的图示中通过举例的方式而非限制的方式示出本公开。

图1示出根据实施方案的用于执行稀疏三角求解(spts)的并行计算系统的实施方案。

图2示出根据实施方案的计算装置。

图3示出根据实施方案的计算装置中的多个处理单元和存储器。

图4a示出根据实施方案的稀疏三角矩阵和向量的矩阵乘法。

图4b示出根据实施方案的spts的依赖图。

图4c示出根据实施方案的压缩稀疏行(csr)数据集和完成数组。

图5示出根据实施方案的用于标识向量中每个因子的完成标志的csr数据集中的元素。

图6示出根据实施方案的spts中的事件的时间线。

图7示出根据实施方案的spts中的事件的时间线。

图8是示出根据实施方案的用于在处理核心中执行spts的模块的框图。

图9a和图9b是示出根据实施方案的用于执行spts的过程的流程图。

具体实施方式

以下描述阐述了许多具体细节，诸如具体系统、部件、方法等的示例，以便提供对实施方案的良好理解。然而，本领域的技术人员将明白，可在没有这些具体细节的情况下实践至少一些实施方案。在其他情况下，未详细描述公知的部件或方法，或以简单的框图格式呈现公知的部件或方法，以避免对实施方案不必要的混淆。因此，阐述的具体细节仅仅是示例性的。特定实现方式可与这些示例性细节不同，并且仍被认为在实施方案的范围内。

稀疏三角求解(spts)试图求解方程式ax＝y中的向量x的未知值，其中a是稀疏三角矩阵，且y是已知值的向量。并行地求解稀疏三角矩阵中的行会导致一系列数据依赖性；向量x中的每个因子x[n]的求解取决于正在求解的在先因子x[0]-x[n-1]。求解可分成一系列层次，其中相同层次中的因子不直接或传递地彼此依赖，并因此可彼此并行地求解。快速执行并行spts的主要困难中的一个是找到任何特定输入矩阵的此数据依赖图，尤其是当输入矩阵包含成千上万个行和列时。因此，稀疏三角求解在高度并行的架构(诸如基于图形处理单元(gpu)的并行计算系统)上仍然表现不佳。基于在开始计算之前确定的依赖图，或通过在先前计算的结果可用时在并行工作器之间进行通信，并行计算系统可确定何时开始求解特定因子。

在并行计算系统上执行spts的一种方法包括：首先分析输入矩阵以确定可并行求解的行和因子，然后针对每个层次启动新内核，包括用于并行求解层次中的行中的每一个的线程。然而，将spts分成分析阶段和求解阶段会导致应用程序编程接口(api)更加繁琐，在所述api中用户在获得期望的求解之前另外调用分析(实现层次的细节)。此外，执行分析所需的时间可能会超过计算求解所需的时间。取决于所执行的分析，分析可能比求解阶段阶段花费多达数千倍的时间。如果矩阵不被重复使用，则执行分析所花费的时间可能不会被摊销。

在一些情况下，给定层次中的因子可能取决于先前层次中的已被求解的因子的子集；因此，应该能够进行因子求解。然而，当先前层次中尚未求解的其他因子阻止了先前层次的完成时，无法进行因子求解。因此，当将spts计算分成多个层次时，可能会丢失一些并行性。

在不需要单独分析阶段的一种方法中，可执行通过在稀疏三角矩阵上进行操作来动态地管理求解阶段期间的并行工作器的spts，所述稀疏三角矩阵根据压缩稀疏列(csc)格式进行存储。然而，许多应用程序以压缩稀疏行(csr)格式存储稀疏三角矩阵，并且将csr数据集转换为csc格式会消耗大量时间和存储器资源。

在一个实施方案中，并行计算系统可通过更新完成数组来对以csr格式存储的矩阵执行spts，所述完成数组指示向量x的因子何时被求解并且可用于后续计算。由计算系统启动的内核执行线程，所述线程用于使用矩阵的对应第n行中的元素来计算向量x中的每个因子x[n]。在每个线程中，执行自旋循环以重复监视完成数组中的完成标志，以确定因子x[n]所依赖于的在先因子(表示第n行的输入变量)是否已被求解。在在先因子已被求解之后，其用于计算因子x[n]的值的一部分。

为了减少由于对完成数组中的值进行自旋循环而导致的存储器争用，如果自旋循环的迭代次数或自旋循环所花费的时间量超过限度，则线程将启动依赖子内核，所述依赖子内核入队以在当前内核完成后开始。在其余线程已完成(即，已求解更多因子)后，依赖子内核将启动新线程，以对完成数组中的相同值恢复自旋循环。

因此，这种高性能spts机制允许并行计算系统对以csr格式存储的稀疏三角矩阵执行spts，而无需执行昂贵的转置操作以将其转换为诸如csc等另一种格式。因为不会将行及其对应因子分组到层次中，所以所述机制不会产生错误的依赖性，并因此能够在执行spts时找到更多动态并行性。高性能spts机制不需要单独的分析阶段；当在先因子已知或可通过求解前几行获得时，每个并行工作器开始求解一行以计算其对应因子的值。由于对完成数组中的多个完成标志进行自旋循环，允许自旋循环超时并通过启动子内核来稍后恢复的机制减少了存储器争用。在一个实施方案中，这种高性能spts机制在一些情况下比在单独阶段中执行分析和求解的spts机制执行起来快数千倍。

图1示出并行计算系统100的实施方案。计算系统100包括通过通信网络110彼此连接的多个计算装置101-103。计算装置101-103中的每一个具有处理和存储器存储能力。在一个实施方案中，计算系统100包含在单个物理外壳内，并且通信网络110是在外壳内连接计算装置101-103的总线或系统互连件。例如，计算装置101-103可包括在相同板上或在通过背板彼此连接的单独载体板上的诸如gpu、中央处理单元(cpu)、现场可编程门阵列(fpga)等处理单元。在一个实施方案中，计算系统100中的部件包含在单独的物理外壳中并且在地理上分布。例如，计算装置101-103可表示通过诸如互联网的广域网(wan)、局域网(lan)、无线网络或其他通信网络110彼此连接的个别服务器、个人计算机、移动装置等。在一个实施方案中，计算装置101-103表示相同类型或类似类型的装置；替代地，计算装置101-103是不同类型的装置。

图2示出实现高性能并行spts机制的计算装置101的实施方案。通常，计算装置101体现为许多不同类型的装置中的任一种，包括但不限于膝上型计算机或台式计算机、移动装置、服务器等。计算装置101包括通过总线201彼此通信的许多部件202-208。在计算装置101中，部件202-208中的每一个都能够直接通过总线201或者通过其他部件202-208中的一个或多个与任何其他部件202-208通信。计算装置101中的部件201-208包含在单个物理外壳内，诸如膝上型或台式计算机机架或移动电话外壳。在替代实施方案中，计算装置101的一些部件体现为外围装置，使得整个计算装置101不驻留在单个物理外壳内。

计算装置101还包括用于从用户接收信息或向用户提供信息的用户接口装置。具体地，计算装置101包括输入装置202，诸如键盘、鼠标、触摸屏或用于从用户接收信息的其他装置。计算装置101通过诸如监视器、发光二极管(led)显示器、液晶显示器或其他输出装置的显示器205向用户显示信息。

计算装置101另外包括用于通过有线或无线网络发射和接收数据的网络适配器207。计算装置101还包括一个或多个外围装置208。外围装置208可包括大容量存储装置、位置检测装置、传感器、输入装置、或由计算装置101使用的其他类型装置。

计算装置101包括一个或多个处理单元204，在多个处理单元204的情况下，它们能够并行操作。处理单元204被配置为接收并执行存储器子系统206中所存储的指令209。在一个实施方案中，处理单元204中的每一个包括驻留在共用集成电路衬底上的多个处理核心。存储器子系统206包括计算装置101使用的存储器装置，诸如随机存取存储器(ram)模块、只读存储器(rom)模块、硬盘和其他非暂时性计算机可读介质。

计算装置101的一些实施方案可包括比图2所示的实施方案更少或更多的部件。例如，在没有任何显示器205或输入装置202的情况下实现某些实施方案。其他实施方案具有多于一个特定部件；例如，计算装置101的实施方案可具有多个总线201、网络适配器207、存储器装置206等。

示出根据实施方案的包括计算装置101的所选择的部件的框图。图3示出处理单元204，所述处理单元204各自通过总线201连接至存储器206。尽管图3示出计算装置101中的一个，但计算系统100中的其他计算装置(例如，102-103)包括类似部件。

在一个实施方案中，处理单元204中的每一个是gpu、cpu、fpga或其他处理装置，并且位于与所述一组处理单元204中其他处理单元分离的集成电路管芯上。处理单元204中的每一个在单个集成电路管芯上包括一组处理核心。处理单元204(1)包括处理核心301-303，处理单元204(2)包括处理核心304-306，并且处理单元204(3)包括处理核心307-309。处理核心中的每一个被配置为执行计算机程序中的线程，如指令209所指引。处理核心301-309能够彼此独立地执行指令，并且因此能够在spts过程中执行并行线程，其中并行线程中的每一个在处理核心301-309中的一个中执行并计算向量x中的因子中的一个的值。

除了用于执行方程式ax＝y的高性能spts的指令集209，存储器206还以csr数据集321和输入数据322的形式存储稀疏三角矩阵a，所述稀疏三角矩阵a包括要求解的未知因子的向量x和已知乘积值的向量y。存储器206还存储完成数组323，所述完成数组323包括向量x中的因子中的每一个的完成标志，其中每个完成标志指示其对应因子是否已被求解。在一个实施方案中，存储器206中的信息存储在计算装置101中的单个存储器装置或子系统上。在替代实施方案中，信息跨相同计算装置101中或多个计算装置(例如，101-103)中的多个存储器装置分布。因此，用于更广泛的计算系统100的存储器系统可包括跨多个计算装置101-103分布的存储器装置。

图4a示出根据一个实施方案的稀疏三角矩阵410与向量x420的乘积，得到乘积向量y430(即，ax＝y)。矩阵a410用作spts的输入，并且是包括10个非零元素的下三角矩阵。矩阵a410的主对角线包括元素a、c、e、g和j。由于矩阵a410为下三角矩阵，因此可用前向替换将其求解。这意味着将使用先前求解的较高行的结果作为输入来求解一些行。因此，箭头411-415表示前向替换过程中的这些依赖性。例如，求解对应于相应的第1行、第2行和第3行的因子x[1]、x[2]和x[3]取决于与要求解的第0行相关联的因子x[0]。这些依赖性分别由箭头411、412和413指示。例如，依赖性箭头411指示计算来自第1行的项bx[0]的值取决于使用第0行中的项ax[0]来求解x[0]。依赖性箭头414和415指示通过第4行求解因子x[4]取决于分别求解对应于非零元素c和g的x[1]和x[3]因子。求解的x[1]和x[3]因子用于在求解因子x[4]时计算项hx[1]和ix[3]。尽管本文将高性能spts机制描述为被执行用于求解下三角矩阵，但所述机制可类似地用于执行向后替换以求解上三角矩阵。

图4b示出用于执行矩阵a410的spts的依赖图450。图中的节点中的每一个表示要求解的向量x中的因子中的一个。图4b中的依赖性411-415分别对应于图4a上的类似编号的依赖性411-415。因子x[0]-x[4]中的每一个的计算通过单独线程执行，其中单独线程中的每一个在单独处理核心(例如，处理核心301-309中的一个)中执行。

在依赖图450中，因子x[0]的值的计算没有依赖性，并且无需等待任何其他因子求解就首先加以计算。求解因子x[1]、x[2]和x[3]中的每一个分别通过依赖性411、412和413取决于x[0]；因此，在x[0]被求解并变得可用时开始这些计算。在用于求解x[1]、x[2]和x[3]的相应的单独线程中执行的计算至少部分地并行执行，但是可能需要花费不同的时间量来完成计算其相应因子。求解x[4]分别通过依赖性414和415取决于x[1]和x[3]。因此，当x[1]和x[3]都可用时，执行x[4]的完整计算。在一个实施方案中，多个线程中的每一个求解因子x[0]-x[4]中的一个。在替代性实施方案中，包括多个线程的并行工作组求解一个因子，单个线程求解多个因子，或多个线程求解向量x420中的多个因子。

图4c示出根据实施方案的表示稀疏三角矩阵a410的压缩稀疏行(csr)数据集321。csr数据集321包括三个数组：值数组(值[])、列数组(列[])和行指针数组(row_ptrs[])。值数组存储矩阵a410的非零元素。元素在矩阵中按从左到右(每行中的第0列至第4列)和从上到下(第0行至第4行)的顺序存储。针对值数组中的元素中的每一个，列数组标识矩阵的所述元素所在的列。列数组与值数组具有相同的条目数；列数组的每个元素都为值数组中具有相同数组索引的对应元素标识列。行指针数组标识哪些元素在矩阵的每一行中。具体地，行指针数组中的每个值都是在每一行的第一值处指向值数组和列数组的索引。行指针数组中的最终值比值数组或列数组中的最高索引大一。

根据用于对存储在csr格式中的矩阵执行spts的一种方法，用于求解x[0]的线程0通知依赖线程1、2和3(分别用于求解x[1]、x[2]和x[3])何时求解了x[0]。然而，当使用csr格式时，此类通知将需要遍历列数组以查找包含‘0’的所有条目，这指示所述行在第0列中具有非零值，并因此对x[0]具有数据依赖性。在列数组中找到‘0’条目后，通过执行行指针数组的搜索以确定来自列数组的‘0’条目位于行指针数组中的哪两个索引之间，线程0还将确定要唤醒哪个依赖线程。列数组的遍历和行指针数组的搜索在计算上是昂贵的，并且会导致高性能spts的无法实现的减慢。

图4c另外示出根据实施方案的完成数组323，所述完成数组323用于避免计算上昂贵的列数组的遍历和行指针数组的搜索。代替每个完成的线程唤醒其依赖线程，等待的依赖线程中的每一个检查完成数组323以确定其在先因子是否已被求解。完成数组323为向量x420中的因子中的每一个(以及因此矩阵a410中的行中的每一个)存储完成标志。向量x420中的每个因子对应于完成数组中的具有相同索引的标志(即，x[n]对应于完成[n])。在开始spts之前，完成数组中的所有标志都初始化为‘0’，从而指示因子都尚未被求解。每当线程完成为向量x420中的因子中的一个写入求解值时，线程也会断言完成数组323中的对应完成标志，以指示求解值可用于后续计算。在一个实施方案中，当完成标志的值为零时，将其断言，而当其值为非零值时，将其取消断言。

图5示出根据实施方案的使用csr数据集数组中的元素来确定针对被求解的向量x420的因子中的每一个将监视完成数组323中的哪些完成标志。在向量x420中，因子x[0]–x[4]中的每一个对应于矩阵a410的行中的一个，并因此转而对应于行指针数组中的行指针元素中的一个。因此，被执行用于求解因子中的一个的线程标识在行指针数组中与要求解的因子具有相同索引的行指针元素。行指针元素用作在列数组中定位列元素的索引，所述列元素标识矩阵a410的对应行中的第一非零值。所述线程将列元素用作在完成数组323中定位完成标志的索引。

在标识完成标志之后，线程执行重复检查完成标志直至断言完成标志为止的自旋循环。当断言完成标志时，线程通过在列数组中查找下一个列元素(例如，通过将索引递增1)来标识要监视的下一个完成标志。线程执行自旋循环以监视下一个完成标志。此过程重复，直至下一个列元素等于因子的索引并因此对应于矩阵a410的主对角线上的元素为止。受监视的完成标志中的每一个指示(所述因子的求解所依赖的)在先因子中的一个是否已被求解。因此，当到达主对角线时，行的所有在先因子都已被求解，并且线程能够求解其自有因子。

当对因子x[0]执行此过程时，值‘0’的对应列元素已表示主对角线上的元素。因此，线程已能够求解x[0]而无需检查任何完成标志。

图6示出根据实施方案的spts中的事件的时间线。如图6所示，向量x中的每个因子由并行线程0-4中的一个求解。在一个实施方案中，线程0-4在一个或多个处理单元(例如，一个或多个处理单元204)中执行，其中每个线程在处理核心中的一个(例如，内核301-309中的一个)中执行。在替代实施方案中，每个因子可由工作组的多个线程来求解，或工作组可求解多个因子。在图6的时间线中，时间从上到下进展。图解的右侧部分示出spts过程期间的不同时间的完成数组323中的完成标志的值。

在时间601处，初始化完成数组323，所有完成标志被设置为‘0’以指示尚未求解任何因子。线程0、1、2、3和4分别读取来自向量y430的乘积值y[0]、y[1]、y[2]、y[3]和y[4]以及来自矩阵主对角线的元素a、c、e、g和j。向量y430与主对角线元素的乘积值是已知的值，每个线程将使用这些值来根据向量x420求解其因子。

在时间602处，基于其如图5所示的确定，线程1-3(分别对应于因子x[1]-x[3])中的每一个执行自旋循环以监视完成数组的索引0(即，完成[0])处的完成标志的状态。根据如图5所示的确定，对应于因子x[4]的线程4执行自旋循环以监视完成[1]的状态。在此时间602期间，线程0没有依赖性(如先前参考图5讨论)，并继续求解x[0]。线程0将y[0]除以a，并将结果作为x[0]存储在向量x420中。在存储结果后，线程0通过将非零值写入完成[0]中来断言完成[0]。

在时间603处，完成[0]处的完成标志处于断言状态下。针对线程1-3中的每一个，列数组中的下一个列元素在主对角线上(参见图5)；因此，线程1-3的所有在先因子值已被求解。线程1-3退出其相应自旋循环并继续求解其相应因子。每个线程1-3读取x[0]的新计算值，连同x[0]在方程式ax＝y中相乘的矩阵元素(例如，b、d、f)。在时间604处，线程1-3分别基于所述行的x[0]、矩阵元素(b,d,f)和乘积值(y[1]-y[3])来求解因子x[1]-x[3]。线程1-3将x[1]-x[3]的所得计算值存储在向量x420中。

在存储计算值之后，每个线程1-3断言对应于其因子x[1]-x[3]的完成标志，以指示这些因子被求解。因此，通过在完成[1]、完成[2]和完成[3]处存储非零值，线程1、2和3分别断言这些位置。

在时间605处，完成[1]处于断言状态下；因此，线程4在完成[1]上停止自旋循环并读取x[1]的新计算值，连同x[1]在方程式ax＝y中相乘的矩阵元素h。在时间606处，线程4标识矩阵行中的包含非零元素的下一列。在列数组中，位置列[8](与第4行的初始位置列[7]相邻)指示第4行中的下一个非零元素是第3列中的非对角线元素。因此，线程4开始自旋循环以监视对应于第3列的完成[3]处的完成标志的状态。

在时间607处，先前求解x[3]，并且在时间603处通过线程3断言其完成标志。此外，第4行中的具有非零值(如由列[9]指定)的下一列在矩阵a410的主对角线上。因此，完成数组323指示已求解所有在先因子。作为响应，线程4退出自旋循环并开始读取将用于计算其因子x[4]的值的值x[3]和i。在时间608处，线程4基于在先因子x[1]和x[2]、矩阵元素h、i和j以及乘积值y[4]来计算因子x[4]的值。在已求解x[4]的情况下，线程4通过在完成[4]中存储非零值来断言完成[4]处的完成标志。

在一个实施方案中，线程0-4通过写入‘1’或其他固定数字作为非零值来在完成数组323中断言它们的相应完成标志。替代地，如图6所示，完成数组323的更新用于生成层次集合信息。断言线程不是简单地用‘1’值断言完成标志，而是通过递增在先因子的完成标志中的最高值来确定其求解因子的完成标志的值。然后，通过将完成数组323中的完成标志的所确定的值存储在对应于因子的位置处，线程断言因子的完成标志。

例如，针对这些线程1-3中的每一个，在先因子(即，x[0])的最高完成标志值为1。递增此值将产生新的完成标志值'2'。因此，线程1、2和3在时间604处分别针对其相应因子x[1]、x[2]和x[3]使用值‘2’来断言完成标志。当线程4在时间608处断言已求解因子x[4]的完成标志时，对应于在先因子x[1]和x[2]的完成标志中的最高完成标志的值为‘2’。因此，线程4使用递增值‘3’来断言x[4]的完成标志。

在spts的结尾(时间609)，断言了其所有完成标志元素的完成数组323指示层次集，所述层次集可任选地用于确定用于随后求解新向量x的因子的顺序。继续先前示例，包括元素[1,2,2,2,3]的完成数组指示首先在层次1中求解x[0]，随后在层次2中并行求解x[1]、x[2]和x[3]，最后在层次3中求解x[4]。完成数组323因此用于有效地生成可在后续计算中用于相同矩阵a410的层次集信息。

图7示出用于执行spts的实施方案的事件的时间线，其中执行自旋循环的线程可超时并启动子内核以减少由自旋循环引起的存储器争用。特别是在矩阵a具有大量行的情况下和/或长依赖性的链，用于重复检查大量对应因子的完成状态的自旋循环可能会由于用于读取完成数组323的重复存储器访问而导致大量存储器争用。存储器争用减慢试图求解其因子并在spts中取得进展的线程，从而进一步延长自旋循环所花费的时间。为了减少自旋循环的总数，如果迭代次数或自旋循环所花费的时间超过预定阈值，则执行自旋循环的每个线程都会终止自旋循环。在spts中已取得更多进展(即，求解更多因子)之后，自旋循环将在子内核的对应线程中恢复。

如图7所示，线程0-4中的每一个在第一内核711中执行，以求解向量x420中的因子x[0]–x[4]中的一个。每个线程0-4首先检查完成数组以确定其相应因子x[0]–x[4]已由前一线程求解。如果因子已被求解，则线程立即退出。在时间702处，完成数组323仅包括‘0’个元素，并因此指示向量x420中的所有因子尚未被求解。因此，线程1-3执行自旋循环以监视完成[0]以用于求解x[0]，而线程4执行自旋循环以监视完成[1]以用于求解x[1]。线程0求解x[0]并断言完成[0]处的相关联完成标志。

线程1-3中的每一个继续执行其自旋循环，直至断言由自旋循环监视的完成标志或自旋循环的迭代次数超过预定限度为止。在时间703处，完成[0]处的x[0]的完成标志处于断言状态下。因此，在其自旋循环迭代超过预定限度之前，线程1-3开始求解其相应因子x[1]-x[3]。然而，在求解x[1]之前，线程4对x[1]的超过预定迭代限度的完成标志执行自旋循环。线程4通过终止其自旋循环来放弃并检查指示子内核是否已被任何其他线程启动的全局“child_launched”变量。如果断言“child_launched”变量，则子线程已由另一个线程启动，并且线程4退出而不启动子内核。如果未断言“child_launched”变量，则先前未启动子内核，并且线程4启动第二内核712，所述第二内核712被入队以在第一内核711完成之后(即，当内核711中的所有其他线程完成时)启动。结合启动内核712，线程4断言“child_launched”变量，以防止第一内核711中的其他线程试图以此方式启动子内核。线程4然后退出。通过这种机制，执行太多迭代或花费太多时间自旋循环的线程暂时停止访问存储器系统，直至取得更多进展为止。

在时间703结束时，当所有并行线程0-4已完成时，第一内核711退出。在时间704处，启动第二内核712，并且分别为相同因子x[0]-x[4]再次调用并行线程0-4。结合子内核712的启动，将全局“child_launched”变量初始化为取消断言状态。线程0-4检查其相应因子x[0]-x[4]的完成标志。由于完成数组323指示因子x[0]-x[3]已被求解，因此线程0-3退出。

然而，完成[4]处的完成标志被取消断言，并指示因子x[4]尚未被求解。通过执行自旋循环来检查x[1](完成[1]处)和x[3](完成[3]处)的完成标志，线程4确定在先因子x[1]和x[3]是否已被求解。在时间705处，基于其在完成[1]处的完成标志，线程4确定x[1]已被求解。在时间706处，基于其在完成[3]处的完成标志，线程4确定x[3]已被求解。由于求解x[4]所取决于的所有在先因子都已被求解，因此线程4继续在时间707处求解x[4]。一旦所有因子x[0]-x[4]已被求解，将不再启动子内核并且完成spts。

图8示出根据实施方案的处理核心301中的用于对csr数据集321执行高性能spts的模块的框图。在一个实施方案中，处理核心301中的模块801-805是使用硬化电路模块实现的；然而，本发明不限于此；在替代实施方案中，使用可编程逻辑电路(例如，当使用fpga或其他可编程装置来实现处理核心301时)，软件模块或硬件、软件、可编程逻辑等的组合来实现模块801-805。在一个实施方案中，模块801-805执行在线程800中执行的用于求解因子x[n]的操作，所述因子x[n]表示向量x中的第n个元素。在一个实施方案中，线程800与如以参考图6或图7描述的线程0-4以类似方式起作用。

自旋循环模块801读取指针811(来自行指针数组和列数组)以确定完成数组323的要监视的适当完成标志，如先前参考图5描述，并执行自旋循环以重复检查来自完成数组323的完成标志810。因此，自旋循环模块801确定要求解的因子x[n]所取决于的任何在先因子的可用性。

在一个实施方案中，当自旋循环花费太多时间或已执行太多迭代时，线程800通过退出并启动子内核来减少存储器争用。计数器802对模块801执行的自旋循环的每次迭代进行计数。将迭代次数(或替代地，自旋循环所花费的时间)与预定限度803进行比较，并且当自旋循环的迭代数(或自旋循环所花费的时间)超过限度803时，子启动器804检查“child_launched”变量817。如果“child_launched”变量817指示子内核尚未由另一个线程启动，则子启动器804启动子内核712，所述子内核712在当前内核711中的所有线程完成后被入队以启动。线程800因此响应于超过限度803而停止自旋循环，然后在spts中已取得更多进展(即，已求解更多因子)之后，在子内核712中重新开始自旋循环。

求解器模块805响应于自旋循环模块801确定x[n]的所有在先因子已被求解而计算因子x[n]的值，如由它们在完成数组323中的相关联完成标志所指示。求解器805从csr数据集321在对应于x[n]的行中读取矩阵元素812，从对应于x[n]的乘积向量y430读取乘积值y[n]，并从向量x420读取已求解在先因子814。通过将在先因子814、乘积y[n]813和矩阵元素812带入如由方程式ax＝y定义的行的方程式然后代数求解因子x[n]，求解器805计算其因子x[n]的值。

求解器805将x[n]的求解值815存储在因子向量x420中，在所述因子向量x420中所述求解值815可用于求解其他线程中的因子。求解器805还针对x[n]确定完成标志816的非零值，并通过将所述值存储在完成数组323中来断言完成数组323中的完成标志。在一个实施方案中，求解器805通过递增在先因子x[n]的完成标志中的最大值来确定完成标志的值，以便计算因子x[n]的层次。替代地，求解器805将‘1’或另一个固定值用于完成标志。

图9a和图9b示出用于对以csr格式存储的稀疏三角矩阵执行高性能spts的过程900。过程900的操作由计算系统100的诸如存储器系统206等部件、处理核心301中的模块801-805等执行。

在框901处，存储器系统206将稀疏三角矩阵a410存储为csr数据集321。在csr数据集321内，值数组存储矩阵a410的元素，列数组针对存储在值数组中的元素标识中的每一个矩阵a410的列，并且行指针数组标识矩阵a410的每个行中的元素。在图4c中示出csr数据集321中的这些数组。

存储器系统206还存储包括向量x420中的因子x[0]–x[4]中的每一个的完成标志的完成数组323(完成[])。在框903处，将完成数组323中的完成标志中的每一个初始化为‘0’，从而指示因子都尚未被求解。在框905处，发起多个并行线程，其中每个线程计算向量x420中的因子中的一个。启动并行线程800中的一个以用于计算因子x[n]的值，其中n通常表示向量x420中的因子x[n]的索引。

在框907处，线程800通过读取完成[n]处的完成标志来确定因子x[n]是否已由先前过程求解。如果完成[n]为非零值(即，被断言)，则线程800在框933处退出，因为x[n]已被求解。如果完成[n]为零(即，无效)，则过程900在框909处继续。

x[n]的值是取决于向量x420中一组一个或多个在先因子中的每一个的数据；即，在先因子的值用于计算x[n]的值。因此，在框909处，自旋循环模块801基于csr数据集321中的行指针数组和列数组来标识要监视的下一个在先因子。与因子x[n]具有相同索引n的行指针用作标识列数组中位置的索引。列数组中的已标识位置标识矩阵a410的第n行中的对应于因子x[n]的非零元素的列。以图5为例，针对因子x[4]，row_ptrs[4]的值为‘7’，并且列[7]的值为‘1’。这指示在矩阵a410的第4行中，非零元素在第1列为中。如图4a所示，此非零元素为'h'。

在框910处，如果列数组中的已标识位置的索引不等于row_ptrs[n+1]，则尚未遍历第n行中的包含非零元素的所有列，并且尚未对所有在先因子执行自旋循环。继续第4行的先前示例，row_ptrs[4+1]为‘10’。对应于在先因子x[1]的已标识非零元素‘h’的列数组中的当前索引为‘7’。由于这些值不相等，因此尚未遍历第4行的所有列，并且过程900在框911处继续。

在框911处，自旋循环模块801确定已标识非零元素是否位于矩阵a410的主对角线上。在一个实施方案中，如果非零元素的列数等于因子的索引(即‘n’)，则所述非零元素在主对角线上。如果元素位于主对角线上，则过程900返回至框909，而不执行任何自旋循环。对应于主对角线上元素的因子是要由线程800求解的因子x[n]；因此，不执行自旋循环来监视其完成标志。相反，过程900在框910处继续至下一个在先因子。在框911处，如果元素位于非对角线(不在主对角线上)，则元素对应于在先因子且过程900在框913处继续。

在框913处，自旋循环模块801通过从存储器206中的完成数组323读取完成标志810来检查第一在先因子的完成标志。继续参考图5的先前示例，列数组中的已标识元素将在先因子的完成标志的索引存储在完成数组323中。因此，自旋循环模块801读取完成[1]处的完成标志，这指示在先因子x[1]是否已被求解。在框915处，如果完成标志未被断言，则在先因子尚未被求解，并且过程900在框917处继续。

在框917处，自旋循环模块801递增自旋循环计数器802，所述自旋循环计数器802对已执行的自旋循环的迭代次数进行计数。将所计数的迭代次数与限度803进行比较，并且在框919处，如果所计数的迭代次数未超过限度803，则过程900继续返回至框913。因此，过程900循环通过框913-919以执行自旋循环，以用于监视第一在先因子的完成标志，直至断言完成标志为止或直至超过限度803为止。

在自旋循环期间，如果完成标志被断言，则过程从框915进行至框909。在框909处，自旋循环模块801标识因子x[n]的求解所依赖于的下一个在先因子。列数组标识第n行中的非零元素；因此，所述行中的任何下一个非零元素将在与最近标识的非零列相邻的列数组中指示。因此，最近标识的列的索引递增一。继续参考图5的先前示例，在列[7]处指示第4行中的存储非零值的最近标识的列；因此，在列[7+1]或列[8]处指示第4行中的下一个非零元素。列数组中的此位置指示第4行在第3列处具有非零值。如图4a所示，此值为'i'。

第4行的第3列在列数组中的索引为'8'，这不等于row_ptrs[4+1]；因此，过程900从框910继续至911。第4行的第3列是矩阵a410中的非对角线元素。因此，从框911开始，过程900在框913处继续。在框913处，自旋循环模块801读取在列数组中的当前位置处指定的完成标志。继续先前示例，列[8]指定第3列。因此，执行自旋循环以用于监视完成[3]处的完成标志，所述完成标志指示x[3]是否被求解。因此，重复执行框909-919以执行自旋循环，以用于转而监视因子x[n]的多个在先因子(例如，x[1]和x[3])的完成标志。

在框910处，如果列数组中的已标识位置的索引等于row_ptrs[n+1]，则已遍历第n行中的包含非零元素的所有列。这意味着自旋循环过程已检测到所有在先因子的断言完成标志。继续先前示例，第4行在第1列、第3列和第4列中包括非零元素。最后一列‘4’在列数组中的索引为‘9’。针对第4行，row_ptrs[n+1]等于‘10’。因此，到列数组的索引递增至‘10’时，具有非零元素(索引从‘7’到‘9’)的所有列都已被自旋循环模块801遍历。即，自旋循环过程已确定所有在先因子x[1]和x[3]的完成标志已被断言。此时，已求解所有在先因子；因此，在框921处，求解器805基于求解的在先因子814、矩阵a410的第n行中的对应于x[n]的元素812以及对应于x[n]的乘积值y[n]813来计算因子x[n]的值。x[n]815的计算值被写入至存储器系统206中的向量x420。继续先前示例，其中n等于'4'，求解器805根据(y[4]-hx[1]-ix[3])/j计算x[4]，然后将结果写入至向量x420。

在框923处，求解器805计算非零值以断言x[n]的完成标志。求解器805使求解x[n]所依赖于的在先因子的完成标志中的最高值的完成标志递增。针对x[4]，在先因子x[1]和x[3]的完成标志分别为‘2’和‘2’，如图6在时间605处所示。由于最高值的完成标志为'2'，因此x[4]的完成标志为3。在框925处，求解器805在对应于完成数组323中的因子的位置处写入完成标志816。在一个实施方案中，完成标志在完成数组中与向量x420中的因子x[n]具有相同索引n(例如，对应于因子x[n]的完成[n])。完成标志指示x[n]被求解；因此，线程在框933处退出。

在框919处，如果自旋循环的由计数器802计数的迭代次数超过限度803，则过程900在框927处继续。在框927处，子启动器804检查child_started全局变量817，它指示当前内核的子内核先前是否已启动。如果child_started变量817被断言，线程在框933处退出而不启动子内核，并且因此终止自旋循环。在框927处，如果child_started变量817取消断言，则先前未启动子内核，并且过程900在框929处继续。在框929处，子启动器804启动子内核，并使子内核入队以在当前内核完成之后启动。因此，在当前执行的内核中的所有并行线程均已结束之后，将子内核入队以启动一组新的并行线程。

在框931处，子启动器804结合启动子内核及其并行线程并使其入队来断言子启动标志。这防止超过自旋循环的限度803的其他当前执行的线程启动另一个子内核并使其入队。如果子内核已入队，则当前内核中的由于自旋循环限度803而终止的任何线程将在相同子内核中的对应线程中继续。在框933处，线程退出。

在框935处，处理核心301等待当前内核中其余并行线程结束。并行线程中的每一个要么完成求解其因子，要么由于超过自旋循环限度803而终止。当所有线程已结束时，内核已完成且过程900在框937处继续。

在框937处，如果子内核已入队(即，在框929处)，则过程在框939处继续。在框939处，启动入队的子内核。子内核为向量x420中的因子x[n]中的每一个启动线程。通过框907和933终止用于求解先前线程中已求解的因子的线程(如完成数组指示)。如上所述求解其余未求解的因子。新线程执行自旋循环(即，框909-919)以确定其相应因子的所有在先因子何时可用，并在在先因子已被求解之后求解其因子。每个子内核都可通过其自旋循环限度803之外的线程中的一个，通过框919和927-933启动并入队其自己的子内核，直至向量x420中的所有因子都已求解并且整个spts都完整为止。

当求解向量x420中的所有因子时，所有线程都在框933处从框907或框925退出，并且没有子内核入队。因此，从框937开始，过程900在框941处继续，在所述框941处完成spts。

在一个实施方案中，框901-941表示用于执行矩阵a410的第一spts的操作，而框943-947表示用于执行矩阵a410的第二spts的操作，如图9b所示。针对第二spts，当求解方程式ax’＝y'中的第二向量x’时，使用通过执行第一spts生成的完成数组，其中a是来自第一spts的相同矩阵a410，x’是与向量x420不同的未知因子的第二向量，并且y'是与向量y不同的已知乘积值的第二向量430。

在框943处，处理单元(例如，处理单元204中的一个)确定向量x’中的要求解的因子的层次。向量x’中每个因子的层次由对应完成标志指示，所述对应完成标志与先前生成的完成数组中的因子具有相同索引。例如，当先前生成的完成数组包括完成标志[1,2,2,2,3]时，第1层次包括x’[0]，第2层次包括x’[1]、x’[2]和x’[3]，并且第3层次包括x’[4]。可并行求解相同层次中的因子。

在框945处，处理单元以对应于与所确定的层次的顺序计算向量x’中的因子x’[n]中的每一个的值，其中编号较低的层次的因子在编号较高的层次的因子之前被求解，并且相同层次中的因子被并行求解。基于其在先因子、矩阵a410的对应行中的元素以及对应乘积值求解y’[n]而求解每个因子x’[n]。如果根据层次依次求解因子x’[n]，则每个因子都将在其在先因子已被求解之后被求解。在所有因子被求解之后，在框947处完成第二spts。

过程900因此允许并行计算系统(例如，计算系统100)对以csr格式存储的矩阵执行spts，而无需将csr数据集昂贵地转换为其他格式(例如csc)，且无需执行任何昂贵的依赖性分析。因此，采用此方法的并行计算系统消耗更少功率，占用更少计算资源，并相比于用于实现相同结果的其他方法以更少的时间计算求解。

一种方法包括：将稀疏三角矩阵存储为压缩稀疏行(csr)数据集；并且针对第一向量中的多个因子中的每个因子计算所述因子的值。计算所述因子的所述值包括：针对所述因子标识所述第一向量中的一组一个或多个在先因子，其中所述因子的所述值取决于所述一个或多个在先因子中的每一个；响应于完成数组指示所述一个或多个在先因子值中的全部被求解，基于以下各项来计算所述因子的所述值：所述矩阵的行中的一个或多个元素；以及对应于所述行的乘积值；以及在所述完成数组中，针对所述因子断言指示所述因子被求解的第一完成标志。

所述方法还包括：将所述矩阵的元素存储在所述csr数据集的值数组中；在所述csr数据集的列数组中，针对存储在所述值数组中的所述元素中的每一个标识所述矩阵的列；以及在所述csr数据集的行指针数组中，标识所述矩阵的每一行中的元素。

所述方法还包括针对所述多个因子中的每个因子，发起用于计算所述因子的值的线程，其中所述线程是并行执行的多个线程中的一个。

所述方法还包括：针对所述多个线程中的每个线程：对在所述线程中执行的自旋循环的迭代次数进行计数，以用于监视所述完成数组中的第二完成标志，其中所述第二完成标志与所述在先因子中的一个相关联；以及响应于所述迭代次数超过限度，终止所述自旋循环；在所述多个并行处理线程中的所有线程已完成时将新线程入队以便执行；以及在所述新线程中监视所述完成标志。

所述方法还包括：检查子启动标志，其中在所述子启动标志被取消断言时执行所述新线程的所述入队；以及结合所述新线程的所述入队，断言所述子启动标志。

所述方法还包括：针对所述多个因子中的每个因子，执行第一自旋循环以监视所述完成数组中的第一完成标志，其中所述第一完成标志在所述csr数据集的列数组中的第一位置处指定，并且其中所述列数组中的所述第一位置由对应于所述因子的行指针指示；以及响应于确定所述第一完成标志被断言并且所述列数组中的第二位置对应于所述稀疏三角矩阵的非对角线元素，执行第二自旋循环以监视所述完成数组中的第二完成标志，其中所述第二完成标志在所述列数组的所述第二位置处指定。

在所述方法中，所述行指针在行指针数组中的位置对应于所述因子在所述向量中的位置，所述行指针是所述列数组中的所述第一位置的索引，并且所述列数组的所述第一位置存储所述完成数组中的所述完成标志的索引。

所述方法还包括：针对所述多个因子中的每个因子，存储所述在先因子中的每一个的完成标志；通过递增所述在先因子的所述完成标志中最高值的完成标志来确定所述因子的所述完成标志的值；以及通过将所述完成数组中的所述因子的所述完成标志的所确定值存储在对应于所述因子的位置处来断言所述因子的所述完成标志。

所述方法还包括：针对第二向量中多个因子中的每个因子，基于对应于所述完成数组中的所述因子的所述完成标志的值来确定所述因子的层次；以及根据所确定层次来依次计算所述第二向量中的所述因子中的每一个的值，其中并行计算相同层次上的所述多个因子中至少两个的值。

在所述方法中，所述csr数据集存储在存储器系统中，并且针对所述多个因子中的每个因子，在与所述存储器系统耦接的求解器电路中计算所述因子的所述值。

所述方法还包括：使用自旋循环电路来从所述完成数组读取所述第一完成标志，其中所述完成数组存储在所述存储器系统中；以及断言所述完成数组中的所述第一完成标志是由与所述自旋循环电路和所述存储器系统耦接的求解器电路执行的。

一种计算装置包括：存储器，所述存储器用于将稀疏三角矩阵存储为压缩稀疏行(csr)数据集；以及处理单元，所述处理单元与所述存储器耦接。所述处理单元用于：针对第一向量中的多个因子中的每个因子，通过执行以下各项来计算所述因子的值：针对所述因子标识所述第一向量中的一组一个或多个在先因子，其中所述因子的所述值取决于所述一个或多个在先因子中的每一个；响应于完成数组指示所述一个或多个在先因子值中的全部被求解，基于以下各项来计算所述因子的所述值：所述矩阵的行中的一个或多个元素；以及对应于所述行的乘积值；以及在所述完成数组中，针对所述因子断言指示所述因子被求解的完成标志。

在所述计算装置中，所述csr数据集还包括：值数组，所述值数组用于存储所述矩阵的元素；列数组，所述列数组用于针对存储在所述值数组中的所述元素中的每一个标识所述矩阵的列；以及行指针数组，所述行指针数组用于标识所述矩阵的每一行中的元素。

在所述计算装置中，所述处理单元包括多个处理核心，所述多个处理核心各自用于：执行多个并行线程中的一个线程，以用于计算所述向量中的所述因子中的一个的所述值。

在所述计算装置中，所述多个处理核心中的每一个进一步用于：通过在所述线程中执行自旋循环来监视所述完成数组中的完成标志；以及响应于所述自旋循环的迭代次数超过预定限度而：终止所述自旋循环；以及在所述多个并行处理线程中的全部已完成时，执行新线程以用于监视所述完成标志。

在所述计算装置中，所述处理单元包括多个处理核心，所述多个处理核心各自用于：执行第一自旋循环以监视所述完成数组中的第一完成标志，其中所述第一完成标志在所述csr数据集的列数组中的第一位置处指定，其中所述列数组中的所述第一位置由对应于第一向量中的所述多个因子中的一个的行指针指示；以及响应于确定所述第一完成标志被断言并且所述列数组中的第二位置对应于所述稀疏三角矩阵的非对角线元素，执行第二自旋循环以监视所述完成数组中的第二完成标志，其中所述第二完成标志在所述列数组的所述第二位置处指定。

所述计算装置还包括：所述存储器中的完成数组，其中所述完成数组用于存储所述多个因子中的每一个的完成标志；其中所述处理单元进一步用于针对所述多个因子中的每个因子，通过递增所述在先因子的所述完成标志中的最高值来确定所述因子的所述完成标志的值；以及通过将所述完成数组中的所述因子的所述完成标志的所确定值存储在对应于所述因子的位置处来断言所述因子的所述完成标志。

在所述计算装置中，所述处理单元进一步用于：针对第二向量中多个因子中的每个因子，基于对应于完成数组中的因子的完成标志的值来确定所述因子的层次；以及根据所确定层次来依次计算所述第二向量中所述因子中的每一个的值，其中并行计算相同层次上的所述多个因子中至少两个的值。

一种计算系统包括：存储器系统，所述存储器系统用于将稀疏三角矩阵存储为压缩稀疏行(csr)数据集；以及一组一个或多个处理单元，所述一组一个或多个处理单元与所述存储器系统耦接。所述一组处理单元中的每个处理单元用于：针对第一向量中的多个因子中的每个因子，通过执行以下各项来计算所述因子的值：针对所述因子标识所述第一向量中的一组一个或多个在先因子，其中所述因子的所述值取决于所述一个或多个在先因子中的每一个；响应于完成数组指示所述一个或多个在先因子值中的全部被求解，基于以下各项来计算所述因子的所述值：所述矩阵的行中的一个或多个元素；以及对应于所述行的乘积值；以及在所述完成数组中，针对所述因子断言指示所述因子被求解的完成标志。

在所述计算系统中，所述一组处理单元中的每个处理单元进一步用于：执行多个并行线程中的至少一个，其中所述多个并行线程中的每一个计算所述多个因子中的一个的值。

在所述计算系统中，所述一组处理单元中的每个处理单元还包括单个集成电路管芯上的多个处理核心，其中所述多个处理核心中的每一个用于执行所述多个并行线程中的一个。

在所述计算系统中，所述一组处理单元中的每个处理单元是与所述一组处理单元中的其他处理单元分离的集成电路管芯上的图形处理单元(gpu)。

如本文所用，术语“耦接至”可意味着直接或通过一个或多个中间部件间接地耦接。通过本文所述的各种总线提供的任何信号可与其他信号进行时分复用，并且可通过一根或多根公共总线提供。另外，电路部件或框之间的互连可显示为总线或单根信号线。总线中的每一个替代地是一根或多根单信号线，并且单信号线中的每一个替代地是总线。

某些实施方案可实现为计算机程序产品，其可包括存储在非暂时性计算机可读介质上的指令。这些指令可用于对通用或专用处理器进行编程以执行所描述的操作。计算机可读介质包括用于以机器(例如，计算机)可读的形式(例如，软件、处理应用)存储或传输信息的任何机制。非暂时性计算机可读存储介质可包括但不限于磁存储介质(例如，软盘)；光存储介质(例如，cd-rom)；磁光存储介质；只读存储器(rom)；随机存取存储器(ram)；可擦可编程存储器(例如，eprom和eeprom)；闪存或另一种类型的适合存储电子指令的介质。

另外，一些实施方案可在分布式计算环境中实践，计算机可读介质存储在其上和/或由多于一个计算机系统执行。另外，可在连接计算机系统的传输介质间拉动或推动计算机系统之间传递的信息。

通常，表示计算装置101的数据结构和/或其承载在计算机可读存储介质上的部分可以是数据库或其他数据结构，其可由程序读取并直接或间接地用于制造包括计算装置101的硬件。例如，数据结构可以是采用高级设计语言((hdl)，诸如verilog或vhdl)的硬件功能的行为级描述或寄存器传送级(rtl)描述。所述描述可由合成工具读取，所述合成工具可对描述进行合成以从合成库产生包括一系列门的网表。所述网表包括门集合，所述门集合还表示包括计算装置101的硬件的功能性。然后，可放置并且路由所述网表以产生描述要应用于掩模的几何形状的数据集。所述掩模然后可用于各种半导体制造步骤中，以产生对应于计算装置101的一个或多个半导体电路。替代地，计算机可读存储介质上的数据库可以是网表(带有或没有合成库)或数据集(根据需要)，或图形数据系统(gds)ii数据。

尽管以特定顺序显示和描述了本文方法的操作，但是可更改每种方法的操作顺序，以便可以相反的顺序执行某些操作或者以使得可至少部分地与其他操作同时执行某些操作。在另一个实施方案中，不同操作的指令或子操作可以是间歇和/或交替的方式。

在以上说明书中，已参考其具体示例性实施方案描述了实施方案。然而，将明显的是：在不脱离如在所附权利要求书中阐述的实施方案的更广范围的情况下，可对其做出各种修改和改变。因此，说明书和附图被认为是说明性的意义而不是限制性的意义。