病態(tài)線性方程組的一類迭代改進(jìn)法

關(guān)晉瑞，溫瑞萍

(太原師范學(xué)院 數(shù)學(xué)系，山西 晉中 030619)

考慮線性方程組的求解:

Ax=b，

(1)

其中A∈n×n是非奇異的矩陣，b∈n是給定的向量.線性方程組廣泛出現(xiàn)在科學(xué)計(jì)算和工程應(yīng)用的許多領(lǐng)域中，如大地測量、結(jié)構(gòu)分析、網(wǎng)絡(luò)分析、數(shù)據(jù)分析、最優(yōu)化、非線性方程組以及微分方程數(shù)值解中的許多問題最終都?xì)w結(jié)為線性方程組的求解[1-3].如何快速有效地求解線性方程組是數(shù)值代數(shù)的核心問題[4-8]，也是目前仍在繼續(xù)研究的重要課題之一[9-12].

線性方程組的數(shù)值方法大體上可分為直接法和迭代法兩大類.國內(nèi)外許多學(xué)者對此進(jìn)行了深入的研究，提出了眾多的有效方法，如經(jīng)典的Gauss消元法、平方根法、Jacobi迭代法、Gauss-Seidel迭代法、SOR迭代法、共軛梯度法和Krylov子空間方法[1,5,8]，以及近年來著名的HSS迭代法[9]等.這些都是求解線性方程組的有效的數(shù)值方法.

本文研究病態(tài)線性方程組的求解問題.當(dāng)線性方程組系數(shù)矩陣的條件數(shù)較大時(shí)或者系數(shù)矩陣接近奇異時(shí)，方程組比較病態(tài)，這時(shí)直接求解得到解的精度較低甚至是完全錯(cuò)誤的.對此，本文提出了一種迭代改進(jìn)方法，對原方程組的解進(jìn)行迭代改進(jìn)，在一定情況下可以有效提高所求得解的精度.

1 預(yù)備知識

本節(jié)介紹文中的符號記法以及后面將要用到的一些預(yù)備知識.

求解線性方程組(1)最簡單的一類迭代方法是分裂迭代法，它的基本思路是，給定系數(shù)矩陣的一個(gè)分裂:

A=M-N，

其中M是可逆的，代入式(1)得到原方程組的等價(jià)形式Mx=Nx+b，進(jìn)而可以得到迭代法:

xk+1=M-1Nxk+M-1b，k=0,1,2,….

(2)

迭代法(2)稱為單步線性定常迭代法，其中M-1N為迭代矩陣.

如果對于任意給定的初始向量x0∈n，由式(2)生成的序列{xk}都是適定的且收斂于方程組(1)的解x*，則稱(2)是收斂的.下述定理給出了單步線性定常迭代法收斂的充要條件.

定理1[13]單步線性定常迭代法(2)收斂當(dāng)且僅當(dāng)其迭代矩陣M-1N的譜半徑ρ(M-1N)<1.

2 一類迭代改進(jìn)法

本節(jié)提出一類迭代法以求解線性方程組(1)，并給出新方法的收斂性分析.

設(shè)α>0是一個(gè)給定的參數(shù)，利用分裂A=(αI+A)-αI，可以得到迭代格式:

(αI+A)xk+1=αxk+b，k=0,1,2,….

(3)

或者:

xk+1=(αI+A)-1αxk+(αI+A)-1b，k=0,1,2,…

其中x0∈n是給定的初始向量.

迭代法(3)作為一種迭代方法似乎沒有太大的應(yīng)用價(jià)值，這是因?yàn)榕c原方程組(1)相比，迭代法(3)需要求解一系列的方程組才能得到原方程組的解.但是，同原方程組相比，迭代法(3)中系數(shù)矩陣更加對角占優(yōu)，當(dāng)原方程組很病態(tài)或者接近奇異時(shí)，適當(dāng)選取參數(shù)可以使得式(3)的系數(shù)矩陣病態(tài)降低從而易于求解.此外，如果利用某種方法已得到原方程組的一個(gè)近似解，那么可以利用迭代(3)進(jìn)行改進(jìn)，得到精度更高的一個(gè)解，這是迭代法(3)的理論意義及價(jià)值.

以下分析迭代法(3)的收斂性.可以證明，當(dāng)系數(shù)矩陣為正穩(wěn)定矩陣時(shí)，迭代法總是無條件收斂的.

定義1[14]設(shè)A∈n×n，若A的每個(gè)特征值都具有正實(shí)部，則稱A為正穩(wěn)定矩陣.

定理2 設(shè)線性方程組(1)的系數(shù)矩陣A∈n×n為正穩(wěn)定矩陣，則對任意的參數(shù)α>0，迭代法(3)是收斂的.

證明式(3)的迭代矩陣為Lα=(αI+A)-1α，由定理1可知迭代法(3)收斂當(dāng)且僅當(dāng)ρ(Lα)<1.因此，只需要驗(yàn)證Lα的每個(gè)特征值的模都嚴(yán)格小于1.

正穩(wěn)定矩陣包含了很多常見的矩陣，如實(shí)對稱正定矩陣、非對稱正定矩陣[9]、非奇異M-矩陣[15]等，因此迭代法(3)的應(yīng)用范圍是較為廣泛的.

容易發(fā)現(xiàn)，迭代法(3)的數(shù)值效果與參數(shù)α的選取密切相關(guān).一方面，參數(shù)α越大，迭代(3)的系數(shù)矩陣越對角占優(yōu)，從而越容易求解；另一方面，參數(shù)α又不能太大，否則迭代(3)會(huì)收斂很慢.因此，選取合適的參數(shù)很重要.在應(yīng)用中，一般應(yīng)該把參數(shù)α盡量選擇小些.

另外，在迭代法(3)的每步計(jì)算中，如何有效求解方程組也是比較重要的，這時(shí)可以利用直接法去求解，如LU分解，也可以用迭代法如SOR或者共軛梯度法等去求解.

當(dāng)線性方程組(1)的系數(shù)矩陣為一般的方陣時(shí)，可以考慮其對應(yīng)的正則化方程組:

ATAx=ATb.

應(yīng)用迭代法(3)于上述正則化方程組可得迭代:

(αI+ATA)xk+1=αxk+ATb，k=0,1,2,…

(4)

由于ATA是實(shí)對稱正定的，根據(jù)定理2上述迭代是收斂的.迭代法(4)的缺陷是增大了條件數(shù)，但是擴(kuò)大了迭代法(3)的應(yīng)用范圍.

3 數(shù)值算例

本節(jié)中，通過兩個(gè)數(shù)值例子來驗(yàn)證所提出的新方法的可行性和有效性,并分別給出利用Matlab中基本命令“x=A”直接求解以及利用新方法計(jì)算的結(jié)果.實(shí)驗(yàn)中，所有的程序都用MATLAB(R2012a)編寫并在筆記本(2 G CPU和8 G內(nèi)存)上運(yùn)行.

例1 考慮線性方程組Ax=b，其中A為n階Hilbert矩陣，b為n維列向量，且選取b使得方程組的精確解為x*=(1,1,…,1)T.隨著n的增大，系數(shù)矩陣的條件數(shù)急劇增長，方程組變得非常病態(tài)[13].

表1 兩種方法求得解的精度Tab.1 The accuracy of solution by the two methods

此外，當(dāng)n=20時(shí)，對于參數(shù)α不同的值，表2中給出了迭代法(3)要達(dá)到同樣精度所需要的迭代次數(shù).從表2中可知，參數(shù)α選取太大，迭代法(3)收斂很慢，只有參數(shù)選取較小些，方法才是比較有效的.

表2 參數(shù)取不同值時(shí)迭代法(3)所需的迭代次數(shù)Tab.2 The number of iterations by iteration method (3) with different parameters

例2 考慮線性方程組Ax=b，其中：

選取b使得方程組的精確解為x*=(1,1，…,1)T.隨著n的增大，系數(shù)矩陣的條件數(shù)增長很快，方程組變得非常病態(tài)[13].

表3 兩種方法求得解的精度Ta.3 The accuracy of solution by the two methods

表4 迭代法(4)所得解的精度與迭代次數(shù)的關(guān)系Tab.4 The precision of the solution by iteration method (4) with the number of iterations

由上述兩個(gè)數(shù)值算例可見，新方法是可行的，而且也是有效的.

4 結(jié)語

本文提出了一類求解病態(tài)線性方程組的迭代法.理論分析和數(shù)值例子表明，新方法是可行的，而且在一定情況下可以求得較高精度的解.但是新方法中的參數(shù)選擇比較困難，如何選擇合適的參數(shù)是一個(gè)較為復(fù)雜的問題，這有待進(jìn)一步的研究.