九对角线性方程组的追赶法

郑茂波

(成都工业学院 信息与计算科学系，成都 610031)

九对角线性方程组的追赶法

郑茂波

(成都工业学院 信息与计算科学系，成都 610031)

通过矩阵的分解，利用三对角线性方程组追赶法思想，推导出九对角线性方程组的追赶法。理论推导表明：对于九对角线性方程组的求解，该算法运算量和存储量均为线性量级。数值实验表明: 该算法比高斯消去法和其他一些迭代法有明显的速度和内存优势，极大地提高了线性方程的求解速度。

九对角矩阵;线性方程组;追赶法;对角占优

在现代科学与工程计算领域，线性方程组的求解极为重要[1-2]。几乎有效的数值方法，最终都归结于求解线性方程组。通过矩阵分解引出的追赶法[3-7]是求解线性方程组的有效方法。九对角线性方程组当其阶数较大时，虽然用其他解法(比如常见的高斯消元法，LU分解等)也可实现，但是精度和速度都不够理想。本文在三对角线性方程组追赶法的研究基础上，推导出针对九对角线性方程组的快速追赶算法，速度和精度都有较大提高，并用数值算例进行了验证。

1 算法的推导

设有九对角线性方程组

(1)

简记为Ax=h。其中A为非奇异的九对角带状矩阵，系数满足丨i-j丨gt;4时，aij=0，特别地，当

时，A为对角占优矩阵，利用三角分解，将矩阵A分解为2个三角矩阵的乘积。

A=LU[3]

(2)

其中:L为下三角阵，U为上三角阵，具体如下：

(3)

根据式(1)～(3)，利用矩阵乘法，比较等式两边的对应元素有：

(4)

容易证明矩阵A非奇异当且仅当αi≠0，i=1,2,…,n。因此，当A非奇异时可以实现矩阵A的LU分解。从而原方程(1)等价于下面2个方程关于x求解：

Ly=h

(5)

Ux=y

(6)

首先，用以下的递推式求出待定系数：

α1=c1,β1=d1/α1,q1=e1/α1,n1=f1/α1

v1=ki/αi(i=1,2,…,n-4)

(7)

γ2=b2,α2=c2-γ2β1,β2=(d2-γ2q1)/α2,q2=(e2-γ2n1)/α2

ni=(fi-γivi-1)/αi(i=2,3,…,n-3)

(8)

z3=α3,γ3=b3-z3β1,α3=c3-z3q1-γ3β2,β3=(d3-z3n1-γ3q3)/α3

qi=(ei-zivi-2-γini-1)/αi(i=3,4,…,n-2)

(9)

m4=g4,z4=a4-m4β1,γ4=b4-m4q1-z4β2,α4=c4-m4n1-z4q2-γ4β3

βi=(di-mivi-3-zini-2-γiqi-1)/αi(i=4,5…,n-1)

(10)

ui=ji,mi=gi-uiβi-4,zi=αi-uiqi-4-miβi-3,

γi=bi-uini-4-miqi-3-ziβi-2

αi=ci-uivi-4-mini-3-ziqi-2-γiβi-1(i=5,6,…,n)

(11)

1)设Ly=h，求y。算法如下：

y1=h1/α1,

y2=(h2-γ2y1)/α2,

y3=(h3-z3y1-γ3y2)/α3,

“沽源优质蔬菜专家工作站”的试验基地设在裕农公司沽源生产基地，是北京市农林科学院在京津冀地区合作共建的首批专家工作站之一。该工作站旨在结合张家口坝上地区农业生产现状，以生菜等叶类蔬菜为重点，以规模化绿色生产为主题，通过产学研合作新模式，推行水肥高效管理、生态调控、生物防治、土壤修复等方面亟需的相关科技支撑，提升沽源蔬菜产业的竞争力，保障环京地区集约化种植业的健康发展。授牌仪式后，鲁西集团特种肥北方公司毕景龙总监做了关于液体肥使用及发展的相关学术报告，与参会人员就相关技术、资源在沽源地区规模化叶菜绿色生产的应用开展了深入地交流和讨论。

y4=(h4-m4y1-z4y2-γ4y3)/α4,

yi=(hi-uiyi-4-miyi-3-ziyi-2-γiyi-1)/αi(i=5,6,…,n)

(12)

2)设Ux=y，求x。算法如下：

xn=yn

xn-1=yn-1-βn-1xn,

xn-2=yn-2-βn-2xn-1-qn-2xn,

xn-3=yn-3-βn-3xn-2-qn-3xn-1-nn-3xn,

(13)

2 运算量和内存的估计

2.1 运算量的估计

表1 用追赶法求解的运算量估计

表2 求解n阶九对角(对角占优但不对称)线性方程组示例

表3 求解n阶九对角(对角不占优)线性方程组示例

如表1所示，用追赶法求解n阶九对角线性方程组大约共需要23n个乘法和19n个加法，其运算量级完全是线性的，为O(23n)。而大家熟知的高斯消元法的运算量是O(n3/3)。可见当n较大时，追赶法的运算量比高斯消元法少。

2.2 内存的估计

由式(1)，存储n阶九对角线性方程组的系数向量和右端向量需要10n个存储单元，由式(7)～(11)，计算系数mi,zi,γi,αi,βi,qi,ni,vi需要8n个存储单元，由式(12)～(13)计算中间变量y和最终结果x需要n个存储单元。所以，求解n阶九对角线性方程组总共需要19n个存储单元。从这个结果来看，用追赶法解对角方程所需的存储单元也为线性的。

3 数值试验及与其他算法比较

本文用Matlab编程实现九对角线性方程组的求解，用2个例子调用该函数进行了试验，并且将本文的追赶法(febs)与LU分解法(LU)，qr分解法(qr)，最小二乘法(lsqr)，稳定双共轭梯度法(bicgstab)，拟最小残余法(qmr)，双共轭梯度法(bicg)，共轭梯度法(cgs)等算法进行了比较。

例1：设n阶九对角线性方程组(1)中的系数为ai=4，bi=6，ci=30，di=3，ei=5，fi=1，gi=2，ji=1，ki=7，hi=1 000。这个对角方程组是对角占优的但是不对称，用追赶法及以上方法对比试验，迭代条件为：容忍误差1e-10或者迭代次数30次，计算结果见表2。

例2：设n阶九对角线性方程组(1)中的系数为ai=0，bi=-1，ci=5，di=2，fi=0，gi=1，ji=2，ki=1，hi=1 000。此方程组对角不占优，用追赶法及以上各种方法对比试验，迭代条件都为：容忍误差1e-10或者迭代次数30次，计算结果见表2。

4 结语

从理论上看，用追赶法求解n阶九对角线性方程组的运算量是线性量级，比其他算法的运算量大为减少，精度也比较好。追赶法对内存的需求也很少。数值试验表明运算的速度与理论估计基本一致。同时，追赶法有较好的适应性，当迭代法效果较差时追赶法仍然有比较良好的计算性能。此外，算法也比较容易实现。

[1] 蔡大用,白峰杉.高等数值分析[M].北京:清华大学出版社,1997.

[2] 李庆扬,王能超,易大义.数值分析[M].北京:施普格林出版社,2001.

[3] 徐树方,高立,张平文.数值线性代数[M].北京:北京大学出版社，2009.

[4] 王礼广,蔡放,熊岳山.五对角方程组追赶法[J].南华大学学报:自然科学版,2008,22(1):1-4.

[5] 叶新涛,张志.七对角线性方程组的追赶法[J].绍兴文理学院学报,2010,30(10):1-5.

[6] 夏爱生,李长国,胡宝安,等.求解五对角方程组追赶法[J].军事交通学院学报,2008,10(4):87-89.

[7] 李欣.一类顺序主子式全为正的三对角矩阵[J].攀枝花学院学报,2004,21(4):101-104.

Catch-upMethodofNineDiagonalLinearEquations

ZHENGMaobo

(Department of Information and Computing Science, Chengdu Technological University, Chengdu 610031, China)

A catch-up algorithm is derived for solutions of linear equation system with nine diagonal matrix using ones with triune diagonal matrix.It is deduced theoretically that the operational amount and the memory amount of the algorithm are linear. It is shown in the numerical experiments that this method has some advantages in computational cost and memory need evidently. It improves the calculational rates.

nine diagonal matrix;linear equation;catch-up method; diagonal domination

2013-04-02

成都工业学院科研项目“模具寿命预估研究”(KY1211007B)

郑茂波(1984- )，男(汉族)，四川营山人，讲师，硕士，研究方向：微分方程数值解。

O241.6

A

2095-5383(2013)02-0026-03