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

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

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

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

追趕法；稀疏矩陣；五對角矩陣；九對角矩陣

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

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

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

(1)

將其簡記為：

Ax=f

(2)

若五對角矩陣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)

利用矩陣的乘法，對照式(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改寫為LUx=f，新的方程可分解成下面的形式：

(5)

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

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

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

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

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

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

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

(6)

簡記為：

Ax=f

(7)

設(shè)此九對角矩陣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改寫為LUx=f，同前面求解五對角方程組一樣，用追趕法求解九對角方程組也是分為追的過程和趕的過程，其中計算yi的過程稱為追的過程，計算xi的過程稱為趕的過程。

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

(11)

則：

(12)

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

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

(13)

則：

(14)

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

3 數(shù)值試驗

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

如表1所示，在取n=30,40,50,60時，分別用高斯消去法和追趕法求解了這個線性方程組，很明顯能夠看出，追趕法所需要的運算時間遠小于消去法，而且這是在對高斯消去法中的矩陣?yán)孟∈璐鎯χ蟮慕Y(jié)果。

表1 消去法和追趕法運算時間的比較

4 結(jié) 語

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

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

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

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

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

[5]李慶揚, 王能超, 易大義. 數(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ù)值解方面的研究工作。