變帶寬存儲形式下對稱正定稀疏方程組的解法

董延亮

【摘要】稀疏方程組一般都是大型系數(shù)矩陣方程組，如果按照直接法求解，往往需要耗費大量的內(nèi)存資源和計算時間.本文則給出了一種使用變帶寬壓縮存儲的方式求解對稱正定稀疏方程組的高效解法.

【關(guān)鍵詞】變寬帶存儲；稀疏方程組；解法

稀疏方程組一般都是大型系數(shù)矩陣方程組，如果按照直接解法的一般原理，往往需要耗費大量的內(nèi)存資源和計算時間.但是對于一些特殊形式的系數(shù)矩陣，我們可以根據(jù)其自身特點尋找相應(yīng)的矩陣存儲形式，進而得到一種更為節(jié)省內(nèi)存、提高運算速度的新解法.這里我們以對稱正定稀疏方程組為例，應(yīng)用變帶寬存儲的方法壓縮儲存系數(shù)矩陣，給出相應(yīng)方程組的高效算法.

下用變帶寬壓縮存儲法求解對稱正定稀疏線性方程組

AX=B，（1-1）

其中A∈Rn×n是對稱正定稀疏陣，X，B∈Rn×m（m≥1）.

先對A，再對B分別作如下的分解：

A=LDLT，（1-2）

A=LDB～.（1-3）

其中L為單位下三角矩陣（lii=1，i=1，2，…，n，lij=0，i

設(shè)矩陣A各行的第一個非零元素是aimi，i=1，2，…，n.記mij=max（mi，mj），則由式（1-2）和（1-3），根據(jù)矩陣乘法規(guī)則，推出矩陣L，D，B～的元素的下列計算公式.

lij=（aij-∑j-1k=mijlikdkljk）/dj，j=mi，mi+1，…，i-1，i=2，3，…，n，

di=aii-∑i-1k=mil2ikdk，i=1，2，…，n，

b～ik=（bik-∑i-1j=milijdjb～jk）/di，i=1，2，…，n，k=1，2，…，m.

將式（1-2）和（1-3）代入式（1-1），得到

LTX=B～.（1-4）

這樣，方程組（1-1）的求解就歸結(jié)為方程組（1-4）的求解.

以下給出對稱正定稀疏方程組（1-1）的算法，其中A∈Rn×n是對稱正定稀疏陣，X，B∈Rn×m（m≥1）.

1.采用變帶寬壓縮存儲方式

用一維數(shù)組a（1：p）和一維整型數(shù)組d（0：n）存放對稱矩陣A.

這里a（1：p）=（a11，a2m2，a22，a3m3，…，a33，…，aimi，…，aii，…，anmn，…，ann），其中p=∑ni=1（i-mi+1）是數(shù)組a中元素的個數(shù).另外用d（0）=0，d（i），i=1，2，…，n來標記對角元素aii在數(shù)組a中的序號.矩陣A的元素aij處于數(shù)組a中的第d（i）-i+j個位置上，即

aij=a（d（i）-i+j），

aimi中的列下標mi的計算公式為

mi=i-（d（i）-d（i-1））+1，i=1，2，…，n.

矩陣B存放在二維數(shù)組b（1：n，1：m）的位置.最終解X也存放在b（1：n，1：m）的位置.

2.計算矩陣L，D和B～

L的元素存放在數(shù)組a（1：p）中A元素的相應(yīng)位置，L的對角線元素不存.D的元素di存放在a（1：p）中aii相應(yīng)的位置.B～存放在b（1：n，1：m）的位置.

對于i=1，2，…，n

（1）計算I=d（i）-i；MI←mi=i-（d（i）-d（i-1））+1.

（2）對于j=MI，MI+1，…，i

J=d（j）-j；MJ←mj=j-（d（j）-d（j-1））+1；

M←mij=max（MI，MJ）=max（mi，mj）

a（I+j）：=a（I+j）-∑j-1k=Ma（I+k）·a（d（k））·a（J+k）

如果j=i，則轉(zhuǎn)（3），否則a（I+j）：=a（I+j）/a（d（j））；

b（i，k）：=b（i，k）-a（I+j）·a（d（j））·b（j，k），k=1，2，…，m.

（3）如果a（I+i）=0，則計算停止（這時A非正定）

否則b（i，k）：=b（i，k）/a（d（i）），k=1，2，…，m.

3.求解方程組LTX=B～

對于i=n，n-1，…，1

（1）計算I=d（i）-i；

MI←mi=i-（d（i）-d（i-1））+1.

（2）對于j=MI，MI+1，…，i-1

b（j，k）：=b（j，k）-a（I+j）·b（i，k），k=1，2，…，m.