賴智超　羅曉群　張其林

摘要：采用列壓縮稀疏（Compressed Sparse Column， CSC）矩陣存儲策略對矩陣LDL分解前進(jìn)行填充元優(yōu)化排序；基于消去樹進(jìn)行LDL符號分解，使之獨立于數(shù)值分解，避免多余的內(nèi)存消耗，減少不必要的數(shù)值運算.利用矩陣非零元的分布特性分析并實現(xiàn)超節(jié)點LDL分解算法，將稀疏矩陣的分解運算變?yōu)橐幌盗谐砻芫仃囘\算，并使用優(yōu)化的BLAS函數(shù)庫加速分解.測試表明：算法在成倍地提高計算速度的同時進(jìn)一步降低內(nèi)存消耗，適用于大規(guī)模的結(jié)構(gòu)計算.

關(guān)鍵詞：稀疏矩陣； LDL分解； 消去樹； 符號分解； 超節(jié)點

中圖分類號： TP301

文獻(xiàn)標(biāo)志碼：B

現(xiàn)代結(jié)構(gòu)工程的計算經(jīng)有限元法處理后往往歸結(jié)為一個規(guī)模龐大的線性代數(shù)方程組的求解.大規(guī)模的矩陣存儲和快速的求解是制約大型結(jié)構(gòu)計算的關(guān)鍵.

稀疏矩陣是指矩陣中大多數(shù)數(shù)值為零的矩陣.對于較大規(guī)模的稀疏矩陣，非零元占總元素個數(shù)比例非常小.結(jié)構(gòu)分析中得到的矩陣就是一個具有很好特性的對稱正定（Symmetric Positive Definite，SPD）稀疏矩陣.本文的關(guān)注點是結(jié)構(gòu)工程中的SPD稀疏矩陣的直接法求解計算.目前，已經(jīng)有不少穩(wěn)定高效的求解庫，如csparse，taucs，cholmod和superLU等，它們主要是使用cholesky分解和LU分解將稀疏矩陣A分解為LLT或

LU的形式.由于分解后的L中非零元相對于A會增多，所以稱之為“填充元”.[1]因此，一般將稀疏矩陣的分解分為2部分：計算填充元的位置，稱之為“符號分解”；進(jìn)行實際的數(shù)值分解計算，稱之為“數(shù)值分解”.

第一個通用快速SPD稀疏矩陣求解器實現(xiàn)庫是SPARSPAK庫[2]，具有里程碑的意義.隨后，稀疏矩陣快速直接解法的算法研究與圖論算法緊密結(jié)合，如消去樹[3]性質(zhì)研究、最小度算法[4]和圖剖分算法[5]等，同時算法的高效性還充分考慮到計算機(jī)硬件的內(nèi)存布局和cache命中率等，因此稀疏矩陣算法得到迅速發(fā)展，已經(jīng)是一個集工程結(jié)構(gòu)、計算數(shù)學(xué)和計算機(jī)應(yīng)用科學(xué)于一體的綜合性學(xué)術(shù)研究領(lǐng)域.

本文采用列壓縮稀疏（Compressed Sparse Column， CSC）矩陣存儲[6]策略，對矩陣LDL分解前進(jìn)行優(yōu)化排序減少填充元，基于消去樹進(jìn)行LDL符號分解，使之獨立于LDL數(shù)值分解.由于非零元相同的稀疏矩陣符號分解結(jié)果相同，因此可以重復(fù)使用，避免多余的內(nèi)存消耗，減少不必要的數(shù)值運算.進(jìn)一步充分利用矩陣非零元的分布特性，分析并實現(xiàn)超節(jié)點LDL分解算法，使稀疏矩陣的分解運算轉(zhuǎn)變成一系列稠密矩陣運算，并使用優(yōu)化的BLAS函數(shù)庫加速分解，在成倍地提高計算速度的同時進(jìn)一步降低內(nèi)存消耗，適用于大規(guī)模結(jié)構(gòu)計算問題.

1 稀疏三角系統(tǒng)Lx=b的求解

本節(jié)討論稀疏三角系統(tǒng)Lx=b的求解，其中：

L為稀疏下三角矩陣，b為稀疏向量.該算法是矩陣分解算法內(nèi)部的基本運算子單元，擁有高效的稀疏右端項求解算法，直接影響矩陣分解算法的效率，因而具有重要意義.

1.1 求解原理

因為CSC格式存儲矩陣的特點是能夠快速地讀取一列中的元素，因此算法以列的方式對矩陣操作，提高求解速度.式（1）～（3）可以滿足要求.

由于右端項b是稀疏的，因此解x大多也是稀疏的.所以，首先要快速地判斷求解結(jié)果x中非零元的位置模式，之后采用式（1）～（3）求解x.

1.2 非零元位置模式求解算法

該問題可以抽象為一個圖論數(shù)學(xué)模型——有向圖的可達(dá)性遍歷[7]，其中，

有向圖G構(gòu)建方式為

摘要：采用列壓縮稀疏（Compressed Sparse Column， CSC）矩陣存儲策略對矩陣LDL分解前進(jìn)行填充元優(yōu)化排序；基于消去樹進(jìn)行LDL符號分解，使之獨立于數(shù)值分解，避免多余的內(nèi)存消耗，減少不必要的數(shù)值運算.利用矩陣非零元的分布特性分析并實現(xiàn)超節(jié)點LDL分解算法，將稀疏矩陣的分解運算變?yōu)橐幌盗谐砻芫仃囘\算，并使用優(yōu)化的BLAS函數(shù)庫加速分解.測試表明：算法在成倍地提高計算速度的同時進(jìn)一步降低內(nèi)存消耗，適用于大規(guī)模的結(jié)構(gòu)計算.

關(guān)鍵詞：稀疏矩陣； LDL分解； 消去樹； 符號分解； 超節(jié)點

中圖分類號： TP301

文獻(xiàn)標(biāo)志碼：B

現(xiàn)代結(jié)構(gòu)工程的計算經(jīng)有限元法處理后往往歸結(jié)為一個規(guī)模龐大的線性代數(shù)方程組的求解.大規(guī)模的矩陣存儲和快速的求解是制約大型結(jié)構(gòu)計算的關(guān)鍵.

稀疏矩陣是指矩陣中大多數(shù)數(shù)值為零的矩陣.對于較大規(guī)模的稀疏矩陣，非零元占總元素個數(shù)比例非常小.結(jié)構(gòu)分析中得到的矩陣就是一個具有很好特性的對稱正定（Symmetric Positive Definite，SPD）稀疏矩陣.本文的關(guān)注點是結(jié)構(gòu)工程中的SPD稀疏矩陣的直接法求解計算.目前，已經(jīng)有不少穩(wěn)定高效的求解庫，如csparse，taucs，cholmod和superLU等，它們主要是使用cholesky分解和LU分解將稀疏矩陣A分解為LLT或

LU的形式.由于分解后的L中非零元相對于A會增多，所以稱之為“填充元”.[1]因此，一般將稀疏矩陣的分解分為2部分：計算填充元的位置，稱之為“符號分解”；進(jìn)行實際的數(shù)值分解計算，稱之為“數(shù)值分解”.

第一個通用快速SPD稀疏矩陣求解器實現(xiàn)庫是SPARSPAK庫[2]，具有里程碑的意義.隨后，稀疏矩陣快速直接解法的算法研究與圖論算法緊密結(jié)合，如消去樹[3]性質(zhì)研究、最小度算法[4]和圖剖分算法[5]等，同時算法的高效性還充分考慮到計算機(jī)硬件的內(nèi)存布局和cache命中率等，因此稀疏矩陣算法得到迅速發(fā)展，已經(jīng)是一個集工程結(jié)構(gòu)、計算數(shù)學(xué)和計算機(jī)應(yīng)用科學(xué)于一體的綜合性學(xué)術(shù)研究領(lǐng)域.

本文采用列壓縮稀疏（Compressed Sparse Column， CSC）矩陣存儲[6]策略，對矩陣LDL分解前進(jìn)行優(yōu)化排序減少填充元，基于消去樹進(jìn)行LDL符號分解，使之獨立于LDL數(shù)值分解.由于非零元相同的稀疏矩陣符號分解結(jié)果相同，因此可以重復(fù)使用，避免多余的內(nèi)存消耗，減少不必要的數(shù)值運算.進(jìn)一步充分利用矩陣非零元的分布特性，分析并實現(xiàn)超節(jié)點LDL分解算法，使稀疏矩陣的分解運算轉(zhuǎn)變成一系列稠密矩陣運算，并使用優(yōu)化的BLAS函數(shù)庫加速分解，在成倍地提高計算速度的同時進(jìn)一步降低內(nèi)存消耗，適用于大規(guī)模結(jié)構(gòu)計算問題.

1 稀疏三角系統(tǒng)Lx=b的求解

本節(jié)討論稀疏三角系統(tǒng)Lx=b的求解，其中：

L為稀疏下三角矩陣，b為稀疏向量.該算法是矩陣分解算法內(nèi)部的基本運算子單元，擁有高效的稀疏右端項求解算法，直接影響矩陣分解算法的效率，因而具有重要意義.

1.1 求解原理

因為CSC格式存儲矩陣的特點是能夠快速地讀取一列中的元素，因此算法以列的方式對矩陣操作，提高求解速度.式（1）～（3）可以滿足要求.

由于右端項b是稀疏的，因此解x大多也是稀疏的.所以，首先要快速地判斷求解結(jié)果x中非零元的位置模式，之后采用式（1）～（3）求解x.

1.2 非零元位置模式求解算法

該問題可以抽象為一個圖論數(shù)學(xué)模型——有向圖的可達(dá)性遍歷[7]，其中，

有向圖G構(gòu)建方式為

摘要：采用列壓縮稀疏（Compressed Sparse Column， CSC）矩陣存儲策略對矩陣LDL分解前進(jìn)行填充元優(yōu)化排序；基于消去樹進(jìn)行LDL符號分解，使之獨立于數(shù)值分解，避免多余的內(nèi)存消耗，減少不必要的數(shù)值運算.利用矩陣非零元的分布特性分析并實現(xiàn)超節(jié)點LDL分解算法，將稀疏矩陣的分解運算變?yōu)橐幌盗谐砻芫仃囘\算，并使用優(yōu)化的BLAS函數(shù)庫加速分解.測試表明：算法在成倍地提高計算速度的同時進(jìn)一步降低內(nèi)存消耗，適用于大規(guī)模的結(jié)構(gòu)計算.

關(guān)鍵詞：稀疏矩陣； LDL分解； 消去樹； 符號分解； 超節(jié)點

中圖分類號： TP301

文獻(xiàn)標(biāo)志碼：B

現(xiàn)代結(jié)構(gòu)工程的計算經(jīng)有限元法處理后往往歸結(jié)為一個規(guī)模龐大的線性代數(shù)方程組的求解.大規(guī)模的矩陣存儲和快速的求解是制約大型結(jié)構(gòu)計算的關(guān)鍵.

稀疏矩陣是指矩陣中大多數(shù)數(shù)值為零的矩陣.對于較大規(guī)模的稀疏矩陣，非零元占總元素個數(shù)比例非常小.結(jié)構(gòu)分析中得到的矩陣就是一個具有很好特性的對稱正定（Symmetric Positive Definite，SPD）稀疏矩陣.本文的關(guān)注點是結(jié)構(gòu)工程中的SPD稀疏矩陣的直接法求解計算.目前，已經(jīng)有不少穩(wěn)定高效的求解庫，如csparse，taucs，cholmod和superLU等，它們主要是使用cholesky分解和LU分解將稀疏矩陣A分解為LLT或

LU的形式.由于分解后的L中非零元相對于A會增多，所以稱之為“填充元”.[1]因此，一般將稀疏矩陣的分解分為2部分：計算填充元的位置，稱之為“符號分解”；進(jìn)行實際的數(shù)值分解計算，稱之為“數(shù)值分解”.

第一個通用快速SPD稀疏矩陣求解器實現(xiàn)庫是SPARSPAK庫[2]，具有里程碑的意義.隨后，稀疏矩陣快速直接解法的算法研究與圖論算法緊密結(jié)合，如消去樹[3]性質(zhì)研究、最小度算法[4]和圖剖分算法[5]等，同時算法的高效性還充分考慮到計算機(jī)硬件的內(nèi)存布局和cache命中率等，因此稀疏矩陣算法得到迅速發(fā)展，已經(jīng)是一個集工程結(jié)構(gòu)、計算數(shù)學(xué)和計算機(jī)應(yīng)用科學(xué)于一體的綜合性學(xué)術(shù)研究領(lǐng)域.

本文采用列壓縮稀疏（Compressed Sparse Column， CSC）矩陣存儲[6]策略，對矩陣LDL分解前進(jìn)行優(yōu)化排序減少填充元，基于消去樹進(jìn)行LDL符號分解，使之獨立于LDL數(shù)值分解.由于非零元相同的稀疏矩陣符號分解結(jié)果相同，因此可以重復(fù)使用，避免多余的內(nèi)存消耗，減少不必要的數(shù)值運算.進(jìn)一步充分利用矩陣非零元的分布特性，分析并實現(xiàn)超節(jié)點LDL分解算法，使稀疏矩陣的分解運算轉(zhuǎn)變成一系列稠密矩陣運算，并使用優(yōu)化的BLAS函數(shù)庫加速分解，在成倍地提高計算速度的同時進(jìn)一步降低內(nèi)存消耗，適用于大規(guī)模結(jié)構(gòu)計算問題.

1 稀疏三角系統(tǒng)Lx=b的求解

本節(jié)討論稀疏三角系統(tǒng)Lx=b的求解，其中：

L為稀疏下三角矩陣，b為稀疏向量.該算法是矩陣分解算法內(nèi)部的基本運算子單元，擁有高效的稀疏右端項求解算法，直接影響矩陣分解算法的效率，因而具有重要意義.

1.1 求解原理

因為CSC格式存儲矩陣的特點是能夠快速地讀取一列中的元素，因此算法以列的方式對矩陣操作，提高求解速度.式（1）～（3）可以滿足要求.

由于右端項b是稀疏的，因此解x大多也是稀疏的.所以，首先要快速地判斷求解結(jié)果x中非零元的位置模式，之后采用式（1）～（3）求解x.

1.2 非零元位置模式求解算法

該問題可以抽象為一個圖論數(shù)學(xué)模型——有向圖的可達(dá)性遍歷[7]，其中，

有向圖G構(gòu)建方式為