变带宽存储形式下对称正定稀疏方程组的解法

董延亮

【摘要】稀疏方程组一般都是大型系数矩阵方程组，如果按照直接法求解，往往需要耗费大量的内存资源和计算时间.本文则给出了一种使用变带宽压缩存储的方式求解对称正定稀疏方程组的高效解法.

【关键词】变宽带存储；稀疏方程组；解法

稀疏方程组一般都是大型系数矩阵方程组，如果按照直接解法的一般原理，往往需要耗费大量的内存资源和计算时间.但是对于一些特殊形式的系数矩阵，我们可以根据其自身特点寻找相应的矩阵存储形式，进而得到一种更为节省内存、提高运算速度的新解法.这里我们以对称正定稀疏方程组为例，应用变带宽存储的方法压缩储存系数矩阵，给出相应方程组的高效算法.

下用变带宽压缩存储法求解对称正定稀疏线性方程组

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

设矩阵A各行的第一个非零元素是aimi，i=1，2，…，n.记mij=max（mi，mj），则由式（1-2）和（1-3），根据矩阵乘法规则，推出矩阵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）的求解就归结为方程组（1-4）的求解.

以下给出对称正定稀疏方程组（1-1）的算法，其中A∈Rn×n是对称正定稀疏阵，X，B∈Rn×m（m≥1）.

1.采用变带宽压缩存储方式

用一维数组a（1：p）和一维整型数组d（0：n）存放对称矩阵A.

这里a（1：p）=（a11，a2m2，a22，a3m3，…，a33，…，aimi，…，aii，…，anmn，…，ann），其中p=∑ni=1（i-mi+1）是数组a中元素的个数.另外用d（0）=0，d（i），i=1，2，…，n来标记对角元素aii在数组a中的序号.矩阵A的元素aij处于数组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存放在二维数组b（1：n，1：m）的位置.最终解X也存放在b（1：n，1：m）的位置.

2.计算矩阵L，D和B～

L的元素存放在数组a（1：p）中A元素的相应位置，L的对角线元素不存.D的元素di存放在a（1：p）中aii相应的位置.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，则转（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.