电路仿真时电路稀疏矩阵的基于消去图的并行分解方法，属于电子设计自动化(EDA)领域，其特征在于，从电路矩阵的符号分析结果中提取出消去图，用来表示LU分解过程中的数据依赖性，根据消去图的不同结构使用两种不同的并行方法：基于群的并行和流水线的并行，从而降低LU分解的运算时间、加快电路仿真的速度，在一系列测试电路矩阵上的测试结果表明本发明专利技术的LU分解时间在并行线程数为1～8个时比LU分解软件KLU快1.66～7.72倍。
【技术实现步骤摘要】
本专利技术涉及一种针对电路仿真时电路稀疏矩阵的基于消去图的并行LU分解方法，属于电子设计自动化(EDA)
技术介绍
在科学计算中，求解线性方程组Ax = b(Ax = b的形式如图1所示，A为nXn矩阵，b和χ都是η维列向量，A和b已知，χ未知待求)具有举足轻重的地位。Ax = b的求解方法主要可分为两大类迭代法和直接法。迭代法如雅可比迭代法、高斯——赛德尔迭代法、超松驰迭代法等，都是先设定一组解的初值，然后通过一个迭代公式进行迭代，使得解逐渐向真实的解靠近，当误差小于给定值时，迭代收敛。但是迭代法通常对矩阵A的要求较高，在一般的计算问题中可能难以满足。不同于迭代法，直接法通过固定的求解公式直接对 A进行有限步骤的计算，从而求得所需要的解。相比于迭代法，直接法的适用性更强，稳定性更好。因此直接法在各类科学计算中用得更为普遍。传统的直接法如高斯消元法、LU分解法(也称为三角分解法)、乔里斯基分解法等。在一般的问题中，矩阵A通常是非常稀疏的。比如对于电路仿真问题，通常电路矩阵中每行的非零元个数仅为5个左右，而且与矩阵规模无关。针对稀疏矩阵特别设计的求解算法，可以大大降低运算复杂度而减少运算时间。在各类直接解法中，运用最多的是LU分解法。LU分解法指的是将一个nXn的方矩阵A分解成一个下三角矩阵L和一个上三角矩阵U的乘积(其中U矩阵的对角线元素都为1)，即A = LU,其中L和U也都是nXn的矩阵。从而求解线性方程组Ax = b的问题转化成求解两个三角方程Ly = b和te = y (y是η维列向量)。图2显示了 LU分解的基本形式，其中L和U未写数值的地方都是0。稀疏矩阵的LU分解过程包括预处理和数值分解两部分。其中预处理是指使用一定的算法对矩阵进行行列交换以达到在数值分解过程中减少运算量的目的，有些预处理方法中还引入了选主元(选主元指的是将绝对值大的元素通过矩阵行列交换操作，交换到对角线的位置上)的步骤，以保证分解过程中的数值稳定性。不管采用哪种预处理方法，预处理之后的矩阵只是对原矩阵进行了一些行、列的交换，而没有其他的变化。本专利技术主要针对 LU分解算法的后者，即数值分解部分，如无特别说明，后文提及的LU分解皆指数值分解部分。本专利技术并不局限于某种特定的预处理方法，对于所有预处理方法都能适用。目前针对稀疏矩阵的LU分解算法大致可分为向左看算法(Left-looking Algorithm)和向右看算法(Right-looking Algorithm)两大类，其中向左看算法由于对稀疏矩阵存储结构有着很好的适应性，在LU分解软件中得到的广泛的应用。目前加州大学伯克利分校(University of California, Berkeley) Sherry Li 开发的 SuperLU 和弗罗里达大学(University of Florida) Tim Davis开发的KLU等软件的LU分解部分都是以向左看算法为基础，针对大规模稀疏矩阵特点进行了优化。在向左看LU分解算法中，从第1列到第η列依次分解A的每一列，即每次计算一个列向量。对每一列的分解可以概括为三个步骤符号分析、数值分解和数值分配，如图3所示。以正在分解A的第k列为例说明。符号分析指从矩阵A的第k列的非零元结构(非零元结构指的是一个集合，这个集合包括第k列上所有非零元素所在的行的行号)计算出LU分解完成后第k列上的非零元结构。具体的符号分析所采用的方法，参见文献 J. R. Gilbert and Ε. Ng, Predicting structure in nonsymmetric sparse matrix factorizations, Graph Theory and Sparse Matrix Computation， Springer Verlag， 1993 (该方法可称为“Gilbert符号分析方法”)。数值分解步骤是根据这一列的符号分析结果进行数值计算，从而获得第k列的所有非零位置上的数值结果。数值分析所采用的算法， 参见文献 J. R. Gilbert and Τ. Peierls, Sparse partial pivoting in time proportional to arithmetic operations, SIAM J. Sci. Statist. Comput. , vol. 9, pp. 862-874,1988(该算法可称为Gilbert/Peierls算法，简称G/P算法)。最后一个步骤是数值分配，即将第k 列的数值计算结果(是一个列向量)中行号小于等于k的部分分配给U矩阵，行号大于等于k的部分分配给L矩阵。依次对k=l，2，…，η循环执行上述步骤，即完成了对整个矩阵A的LU分解。虽然目前已经有很多基于LU分解的软件包，比如SuperLU，KLU，PARDISO，UMFPack 等共几十个软件，但是在这些软件中，并行的版本非常少，而针对通用多核CPU平台的并行软件更是寥寥无几。由于LU分解过程中的高度数据依赖性，以及并行时多线程操作所带来的额外代价，造成了并行的LU分解软件数据如此之少。
技术实现思路
本专利技术的目的是提供一种针对电路仿真中电路矩阵的并行LU分解方法。本专利技术提出使用消去图来表示LU分解过程中的数据依赖性并用来进行并行的任务调度。该算法从对电路矩阵的符号分析结果中提取出消去图，然后基于该消去图的不同结构，提出两种不同的并行方式群并行与流水线并行。在LU数值分解过程中这两种并行方式动态的调用，以适应不同的数据依赖性。，其特征在于，是在计算机中按照以下步骤实现的步骤(1)，输入要仿真解析的电路的网单；步骤(2)，建立nXn 的电路稀疏矩阵，包括 add20、circuit_l、circuit_2、 add32、rajat03、coupled、circuit_3、onetonel、onetone2、cktll752_dc_l、circuit_4、 ASIC_100ks、ASIC_100k、del、trans4、G2_circuit、transient、ASIC_320ks、ASIC_320k、 rajat30>ASIC_680ks>ASIC_680k>G3_circuit>FreescaleUrajat31 \)JsR circuit5M^|^ 的任何一个电路稀疏矩阵；步骤(3)，选择对角块模式和非对角块模式中的任何一种模式，对步骤O)中所建立的电路稀疏矩阵进行预处理，得到预处理之后的电路稀疏矩阵A ；步骤，根据John R. GiIbert和Timothy Peierls提出的非零元符号分析方法， 对步骤C3)得到的所述nXn的电路稀疏矩阵A进行符号分析，从第1列到第η列依次完成对所述电路稀疏矩阵A的下三角矩阵L和上三角矩阵U的η列内的非零结构的计算，A = LU ；步骤(5)，利用步骤(4)得到的U矩阵的非零结构，计算消去图消去图共有η个顶点，每个分立的顶点i唯一地对应于所述U矩阵中的第i列，顶点之间不存在任何连线， 但每个顶点都有一个“所处的层”属性值，用于表示每个顶点以及该顶点在所述U矩阵中对应的列在消去图中处于哪一层，该属性值记为level (i本文档来自技高网...
【技术保护点】
１．电路仿真时电路稀疏矩阵的基于消去图的并行分解方法，其特征在于，是在计算机中按照以下步骤实现的：步骤（１），输入要仿真解析的电路的网单；步骤（２），建立ｎ×ｎ的电路稀疏矩阵，包括ａｄｄ２０、ｃｉｒｃｕｉｔ＿１、ｃｉｒｃｕｉｔ＿２、ａｄｄ３２、ｒａｊａｔ０３、ｃｏｕｐｌｅｄ、ｃｉｒｃｕｉｔ＿３、ｏｎｅｔｏｎｅ１、ｏｎｅｔｏｎｅ２、ｃｋｔ１１７５２＿ｄｃ＿１、ｃｉｒｃｕｉｔ＿４、ＡＳＩＣ＿１００ｋｓ、ＡＳＩＣ＿１００ｋ、ｄｃ１、ｔｒａｎｓ４、Ｇ２＿ｃｉｒｃｕｉｔ、ｔｒａｎｓｉｅｎｔ、ＡＳＩＣ＿３２０ｋｓ、ＡＳＩＣ＿３２０ｋ、ｒａｊａｔ３０、ＡＳＩＣ＿６８０ｋｓ、ＡＳＩＣ＿６８０ｋ、Ｇ３＿ｃｉｒｃｕｉｔ、Ｆｒｅｅｓｃａｌｅ１、ｒａｊａｔ３１以及ｃｉｒｃｕｉｔ５Ｍ在内的任何一个电路稀疏矩阵；步骤（３），选择对角块模式和非对角块模式中的任何一种模式，对步骤（２）中所建立的电路稀疏矩阵进行预处理，得到预处理之后的电路稀疏矩阵Ａ；步骤（４），根据Ｊｏｈｎ Ｒ．Ｇｉｌｂｅｒｔ和Ｔｉｍｏｔｈｙ Ｐｅｉｅｒｌｓ提出的非零元符号分析方法，对步骤（３）得到的所述ｎ×ｎ的电路稀疏矩阵Ａ进行符号分析，从第１列到第ｎ列依次完成对所述电路稀疏矩阵Ａ的下三角矩阵Ｌ和上三角矩阵Ｕ的ｎ列内的非零结构的计算，Ａ＝ＬＵ；步骤（５），利用步骤（４）得到的Ｕ矩阵的非零结构，计算消去图：消去图共有ｎ个顶点，每个分立的顶点ｉ唯一地对应于所述Ｕ矩阵中的第ｉ列，顶点之间不存在任何连线，但每个顶点都有一个“所处的层”属性值，用于表示每个顶点以及该顶点在所述Ｕ矩阵中对应的列在消去图中处于哪一层，该属性值记为ｌｅｖｅｌ（ｉ），建立所述消去图的步骤如下：步骤（５．１），从第１列到第ｎ列循环所有列，对其中的每一列ｋ，取出所述Ｕ矩阵中第ｋ列中不包括对角线的上三角部分的所有非零元，共ｚ个，所述ｚ个非零元所在行的行号分别标记为ｊ１，ｊ２，…，ｊｚ；步骤（５．２），若：ｚ＝０，则ｌｅｖｅｌ（ｋ）＝１；步骤（５．３），若：ｚ＞０，则取出所述各行号对应的相同序号的列所处的层的属性值ｌｅｖｅｌ（ｊ１），ｌｅｖｅｌ（ｊ２），…，ｌｅｖｅｌ（ｊｚ），求出其中的最大值ｍ＝ｍａｘ（ｌｅｖｅｌ（ｊ１），ｌｅｖｅｌ（ｊ２），…，ｌｅｖｅｌ（ｊｚ）），于是，所述第ｋ列在所述消去图中所处的层的属性值为ｌｅｖｅｌ（ｋ）＝ｍ＋１；步骤（５．４），得到所述ｎ列中每一列所处的层的属性值之后，就得到了整个消去图，从该消去图中所有顶点所处的层的属性值中求出最大值，即为整个消去图的总层数，记为Ｌｔｏｔａｌ；步骤（６），根据步骤（４）得到的Ｌ矩阵和Ｕ矩阵的非零结构以及步骤（５）得到的消去图，按以下步骤进行并行的ＬＵ数值分解：步骤（６．１），输入线程数目Ｔｈｒｅａｄ，并创建每个线程，所述线程数目Ｔｈｒｅａｄ的值应大于或等于１，小于或等于所使用的计算机上的ＣＰＵ的核心数；步骤（６．２），输入一个用于考察所述消去图中各层顶点总数Ｖ（ｌｖ）的阈值Ｖｔｈ，该阈值Ｖｔｈ在所述线程数目Ｔｈｒｅａｄ的１倍－１０倍之间取值；步骤（６．３），从第１层到第Ｌｔｏｔａｌ层依次考察所述消去图中每一层ｌｖ的顶点总数Ｖ（ｌｖ）：若：Ｖ（ｌｖ）＜Ｖｔｈ，则对第ｌｖ层采用流水线并行方法，若：Ｖ（ｌｖ）≥Ｖｔｈ，则对第ｌｖ层采用群并行方法，把当前任何一种所述并行方法的开始层记为Ｌｓｔａｒｔ，结束层记为Ｌｓｔｏｐ，Ｌｓｔａｒｔ＝ｌｖ，按以下方法从开始层Ｌｓｔａｒｔ开始向后考察各层，第ｌｖ＋１层，…，直到第Ｌｔｏｔａｌ层为止，寻找所述的结束层Ｌｓｔｏｐ：在采用流水线并行方法时，若有一层ｃ满足Ｖ（ｃ）≥Ｖｔｈ，则结束层Ｌｓｔｏｐ＝ｃ－１，在采用群并行的方法时，若有一层ｃ’满足Ｖ（ｃ’）＜Ｖｔｈ，则结束层Ｌｓｔｏｐ＝ｃ’－１，步骤（６．４），把从步骤（６．３）中得到的从所述消去图中的第Ｌｓｔａｒｔ到第Ｌｓｔｏｐ层这些层中的所有顶点号，根据步骤（６．３）中得到的相应的一种并行方法，按以下步骤进行多线程方式的并行ＬＵ数值分解：步骤（６．４．１），当从开始层Ｌｓｔａｒｔ到结束层Ｌｓｔｏｐ采用群并行方法时，在循环进行到其中任何一层时，把该层上的所有顶点平均分配给所述各个线程，所述各个线程都采用Ｊｏｈｎ Ｒ．Ｇｉｌｂｅｒｔ和Ｔｉｍｏｔｈｙ Ｐｅｉｅｒｌｓ提出的ＬＵ数值分解方法对分配到的各列进行并行的ＬＵ数值分解，等待所有线程都计算完成后再进入下一层，循环往复，遍历从开始层Ｌｓｔａｒｔ到结束层Ｌｓｔｏｐ的所有层，步骤（６．４．２），当从开始层Ｌｓｔａｒｔ到结束层Ｌｓｔｏｐ采用流水线并行方法时，把从开始层Ｌｓｔａｒｔ到结束层Ｌｓｔｏｐ的所有层中的列号按层的顺序排成一个序列，并给这个序列中每一个列都置一个标志，初始时所有标志都是“未完成”，然后每个线程从该序列的头部开始依次取出一个列号，并从...
【技术特征摘要】
【专利技术属性】
技术研发人员：汪玉，陈晓明，武伟，杨华中，
申请(专利权)人：清华大学，
类型：发明
国别省市：11
还没有人留言评论。发表了对其他浏览者有用的留言会获得科技券。