續(xù)小磊，馬 丁 (寧夏大學(xué)數(shù)學(xué)計(jì)算機(jī)學(xué)院，寧夏 銀川 750021)

求解五對(duì)角和九對(duì)角線性方程組的追趕法

續(xù)小磊，馬 丁 (寧夏大學(xué)數(shù)學(xué)計(jì)算機(jī)學(xué)院，寧夏 銀川 750021)

利用追趕法求解三對(duì)角線性方程組的思想，推導(dǎo)出求解五對(duì)角和九對(duì)角線性方程組的追趕法。此方法不必選主元、計(jì)算量小、存儲(chǔ)量小、避免了中間結(jié)果數(shù)量級(jí)的巨大增長(zhǎng)和舍入誤差的嚴(yán)重積累、運(yùn)算速度快而且Matlab程序編寫也較為簡(jiǎn)單。

追趕法；稀疏矩陣；五對(duì)角矩陣；九對(duì)角矩陣

求解偏微分方程經(jīng)常用到差分法，比如五點(diǎn)差分和九點(diǎn)差分，這樣的差分格式最終是以線性方程組的形式出現(xiàn)，而這些方程組的系數(shù)矩陣就會(huì)是五對(duì)角或者九對(duì)角矩陣[1-4]。對(duì)于線性方程組的求解主要有直接法和迭代法，最基本的直接法是高斯消去法，但是高斯消去法計(jì)算量比較大，而且由于對(duì)角矩陣含有大量的0元素，如果使用高斯消去法將會(huì)對(duì)計(jì)算機(jī)的存儲(chǔ)造成極大的浪費(fèi)，而且需要較長(zhǎng)的計(jì)算時(shí)間。而迭代法一般也需要足夠的迭代次數(shù)才能達(dá)到所需要的精度，甚至有時(shí)候發(fā)散，而且程序編寫復(fù)雜，不易實(shí)現(xiàn)。在這種情況下，使用追趕法求解線性方程組就比較理想[5-6]。下面，筆者將利用追趕法的思想來(lái)求解五對(duì)角和九對(duì)角線性方程組。

1 五對(duì)角線性方程組追趕法的推導(dǎo)

設(shè)有如下五對(duì)角線性方程組：

(1)

將其簡(jiǎn)記為：

Ax=f

(2)

若五對(duì)角矩陣A滿足如下條件：

(1)ai,ei≠0,i=1,…,n-2;bi,di≠0,i=1,…,n-1；

(2)|c1|>|d1|+|e1|; |c2|>|b1|+|d2|+|e2|; |ci| > |ai-2| + |bi-1| + |di| + |ei|,i=3,…,n-2;|cn-1|>|an-3|+|bn-2|+|dn-1|; |cn|>|an-2|+|bn-1|；

那么矩陣A存在唯一的Crout分解，其分解形式為：

A=LU

(3)

式中，L為下三角矩陣；U為上三角矩陣，具體形式如下：

(4)

利用矩陣的乘法，對(duì)照式(3)容易得出：

si=ai,i=1,2,…,n-2;m1=b1;mi=bi-si-1pi-1,i=2,3,…,n-1;

l1=c1;l2=c2-m1p1;li=ci-si-2qi-2-mi-1pi-1,i=3,4,…,n;

將方程組Ax=f改寫為L(zhǎng)Ux=f，新的方程可分解成下面的形式：

(5)

綜合分解與求解的過程，可得到解五對(duì)角方程組的追趕法。類似于三對(duì)角方程組，它也分為追和趕的過程，其中計(jì)算yi的過程稱之為追的過程，計(jì)算xi的過程稱之為趕的過程。

首先通過解方程組Ly=f求得yi的值，即追的過程。將矩陣形式轉(zhuǎn)化為方程組可得：

將si、mi、li的值代入即可求得yi的值。

現(xiàn)在可以通過解方程組Ux=y而求出xi的值，即趕的過程，由方程組Ux=y，可得：

最后將pi、qi的值代入即可求得xi的值。

2 九對(duì)角線性方程組追趕法的推導(dǎo)

設(shè)有如下九對(duì)角線性方程組：

(6)

簡(jiǎn)記為：

Ax=f

(7)

設(shè)此九對(duì)角矩陣A滿足下列條件：

(1)ai≠0,i=1,…,n-4;zi≠0,i=1,…,n-3;hi≠0,i=1,…,n-2;bi≠0,i=1,…,n-1;di≠0,i=1,…,n-1;wi≠0,i=1,…,n-2;ui≠0,i=1,…,n-3;ei≠0,i=1,…,n-4；

(2)|c1|>|d1|+|wi|+|ui|+|e1|; |c2|>|b1|+|d2|+|w2|+|u2|+|e2|,|c3|>|h1|+|b2|+|d3|+|w3|+|u3|+|e3|; |c4|>|z1|+|h2|+|b3|+|d4|+|w4|+|u4|+|e4|,|ci|>|ai-4|+|zi-3|+|hi-2|+|bi-1|+|di|+|wi|+|ui|+|ei|,i=5,…,n-4,|cn-3|>|an-7|+|zn-6|+|hn-5|+|bn-4|+|dn-3|+|wn-3|+|un-3||cn-2|>|an-6|+|zn-5|+|hn-4|+|bn-3|+|dn-2|+|wn-2|,|cn-1|>|an-5|+|zn-4|+|hn-3|+|bn-2|+|dn-1|; |cn|>|an-4|+|zn-3|+|hn-2|+|bn-1|；

則A存在唯一的Crout分解，其分解形式為：

A=LU

(8)

式中，L為下三角矩陣；U為上三角矩陣，具體形式如下：

(9)

(10)

利用矩陣的乘法，通過式(9)可得：

si=ai,i=1,…,n-4;g1=z1;gi=zi-si-1pi-1,i=2,…,n-3;

t1=h1;

t2=h2-g1p1;ti=hi-si-2oi-2-gi-1pi-1,i=3,…,n-2;

m1=b1;m2=b2-t1p1;

m3=b3-g1o1-t2p2;mi=bi-si-3vi-3-gi-2oi-2-ti-1pi-1,i= 4,…,n-1;

l1=c1;

l2=c2-m1p1;l3=c3-t1o1-m2p2;l4=c4-g1v1-t2o2-m3p3;

li=ci-si-4qi-4-gi-3vi-3-ti-2oi-2-mi-1pi-1,i=5,…,n;

將方程組Ax=f改寫為L(zhǎng)Ux=f，同前面求解五對(duì)角方程組一樣，用追趕法求解九對(duì)角方程組也是分為追的過程和趕的過程，其中計(jì)算yi的過程稱為追的過程，計(jì)算xi的過程稱為趕的過程。

先解方程組Ly=f，即追的過程,可得：

(11)

則：

(12)

將si、gi、ti、mi和li的值代入式(12)即可求得yi的值。

再解方程組Ux=y，即趕的過程，可得：

(13)

則：

(14)

將pi、oi、vi和qi的值代入式(14)即可求得xi的值。

3 數(shù)值試驗(yàn)

例1[5]以如下一個(gè)n維九對(duì)角線性方程組為例，分別用高斯消去法和筆者推導(dǎo)的追趕法進(jìn)行求解，比較2種方法的運(yùn)算時(shí)間。在Matlab編程中，對(duì)于高斯消去法中的矩陣采用稀疏矩陣進(jìn)行存儲(chǔ)。

如表1所示，在取n=30,40,50,60時(shí)，分別用高斯消去法和追趕法求解了這個(gè)線性方程組，很明顯能夠看出，追趕法所需要的運(yùn)算時(shí)間遠(yuǎn)小于消去法，而且這是在對(duì)高斯消去法中的矩陣?yán)孟∈璐鎯?chǔ)之后的結(jié)果。

表1 消去法和追趕法運(yùn)算時(shí)間的比較

4 結(jié) 語(yǔ)

將求解三對(duì)角方程組的追趕法推廣到求解五對(duì)角和九對(duì)角線性方程組，這一推廣提供了現(xiàn)成的求解五對(duì)角和九對(duì)角方程組的方法和程序，方便直接利用。通過一個(gè)九對(duì)角線性方程組的數(shù)值算例將追趕法和高斯消去法進(jìn)行對(duì)比，數(shù)值試驗(yàn)結(jié)果表明，追趕法在運(yùn)算量、存儲(chǔ)量以及計(jì)算效率上均具有較大的優(yōu)勢(shì)。

[1]王禮廣,蔡放,熊岳山.五對(duì)角線性方程組追趕法[J]. 南華大學(xué)學(xué)報(bào)(自然科學(xué)版),2008, 22(1): 1-4.

[2]王晶昕,薛靜. 用五參數(shù)法求解擬五對(duì)角方程組[J]. 吉林師范大學(xué)學(xué)報(bào)(自然科學(xué)版), 2008, 29(2): 5-9.

[3]李文強(qiáng), 劉曉. 追趕法并行求解循環(huán)三對(duì)角方程組[J]. 科技導(dǎo)報(bào), 2009, 27(18): 90-93.

[4]陳志,高旅端. 求解大規(guī)模稀疏線性方程組的算法[J]. 北京工業(yè)大學(xué)學(xué)報(bào), 2001(9): 262-265.

[5]李慶揚(yáng), 王能超, 易大義. 數(shù)值分析[M]. 北京: 清華大學(xué)出版社, 2008.

[6]蔡大用, 白峰杉. 高等數(shù)值分析[M]. 北京: 清華大學(xué)出版社, 1997.

[編輯] 洪云飛

O241.6

A

1673-1409(2013)25-0005-05

2013-06-17

續(xù)小磊(1986-)，男，碩士生，現(xiàn)主要從事偏微分方程數(shù)值解方面的研究工作。