求解Poisson方程改進(jìn)的交替方向隱式迭代格式

葛志昊, 劉富豪

(河南大學(xué)數(shù)學(xué)與統(tǒng)計(jì)學(xué)院，開(kāi)封 475004)

1 引言

Poisson 方程是一類非?；A(chǔ)且重要的物理方程，在圖像修復(fù)[1]、計(jì)算流體力學(xué)[2]、熱源識(shí)別[3]等領(lǐng)域均有重要應(yīng)用．因此，構(gòu)造一種簡(jiǎn)便而且快速的求解Poisson 方程的算法是很有必要的．交替方向隱式迭代法具有比一般迭代法更高的收斂速度[4,5]，歷史上第一個(gè)交替方向隱式迭代法是由Peaceman 和Rachford 于1955 年在文獻(xiàn)[6]中提出的．1962 年，Varga 將交替方向隱式迭代法用于求解橢圓方程[4]，并將其命名為Peaceman-Rachford iterative method，簡(jiǎn)稱PR 迭代法．該方法廣泛應(yīng)用于各種凸優(yōu)化問(wèn)題[7,8]．

經(jīng)典的交替方向隱式迭代格式[5]在計(jì)算大規(guī)模問(wèn)題時(shí)，會(huì)產(chǎn)生諸多計(jì)算上的不便與缺點(diǎn)，如計(jì)算程序編寫(xiě)困難、計(jì)算耗時(shí)較長(zhǎng)、計(jì)算過(guò)程中占用內(nèi)存過(guò)大等．為了克服上述缺點(diǎn)，本文給出了改進(jìn)的交替方向隱式迭代格式，相比于經(jīng)典的交替方向隱式迭代格式，改進(jìn)的交替方向隱式迭代格式通過(guò)求解較低維數(shù)的矩陣方程，降低了程序?qū)崿F(xiàn)難度，大大減少了計(jì)算量，并進(jìn)一步證明了改進(jìn)的交替方向隱式迭代格式與經(jīng)典的交替方向隱式迭代格式的等價(jià)性．因此，改進(jìn)的交替方向隱式迭代格式更加便于存儲(chǔ)，同時(shí)克服了經(jīng)典的交替方向隱式迭代格式在計(jì)算上的缺點(diǎn)．

本文其余部分結(jié)構(gòu)如下：在第2 部分的2.1 節(jié)，我們介紹經(jīng)典的交替方向隱式迭代格式；在2.2 節(jié)，我們給出改進(jìn)的交替方向隱式迭代格式；在2.3 節(jié)，我們證明經(jīng)典的交替方向隱式迭代格式和改進(jìn)的交替方向隱式迭代格式之間的等價(jià)性．在第3 部分，我們給出數(shù)值算例驗(yàn)證改進(jìn)的交替方向隱式迭代格式在計(jì)算量和程序?qū)崿F(xiàn)等方面的優(yōu)勢(shì)．最后，給出結(jié)論總結(jié)改進(jìn)的交替方向隱式迭代格式的優(yōu)勢(shì)和不足．

2 改進(jìn)的交替方向隱式迭代格式

考慮Poisson 方程Dirichlet 邊值問(wèn)題

其中

下面，對(duì)? 進(jìn)行矩形網(wǎng)格剖分．取定沿x 軸，y 軸的步長(zhǎng)h1和h2，作兩族與坐標(biāo)軸平行的直線

兩族直線的交點(diǎn)(ih1,jh2)稱為網(wǎng)點(diǎn)或節(jié)點(diǎn)，記為(xi,yj)．

問(wèn)題(1)和(2)的五點(diǎn)差分格式為

其中uij表示解函數(shù)u(x,y)在節(jié)點(diǎn)(xi,yj)處的數(shù)值解，fij=f(xi,yj)．

2.1 經(jīng)典的交替方向隱式迭代格式

首先，我們給出矩陣的自然排列的定義．

定義1給定U =(uij)m×n∈ Rm×n，令

稱向量u 為矩陣U 的自然排列，并記作u = {uij}，其中U(:,j)T表示矩陣U 的第j 列的轉(zhuǎn)置．

下面，我們介紹經(jīng)典的交替方向隱式迭代格式[4,5]：考慮問(wèn)題(1)和(2)的五點(diǎn)差分格式(3)和(4)，定義矩陣L1和L2分別為

由(5)和(3)，得

于是，經(jīng)典的交替方向隱式迭代格式(以下簡(jiǎn)稱：經(jīng)典的迭代格式)為

其中τk是迭代參數(shù)．按層合并，可得

2.2 改進(jìn)的交替方向隱式迭代格式

接下來(lái)，我們給出改進(jìn)的交替方向隱式迭代格式：考慮問(wèn)題(1)和(2)的五點(diǎn)差分格式(3)和(4)，定義矩陣分別為

由(10)，可將(3)寫(xiě)成

于是，改進(jìn)的交替方向隱式迭代格式(以下簡(jiǎn)稱：改進(jìn)的迭代格式)為

其中τk是迭代參數(shù)．按層合并，得

2.3 改進(jìn)的迭代格式與經(jīng)典的迭代格式的等價(jià)性

定理1改進(jìn)的迭代格式(13)和(14)與經(jīng)典的迭代格式(8)和(9)等價(jià)．

證明 只需證明(8)與(13)的等價(jià)性，可以類似證明(9)與(14)的等價(jià)性，這里不再贅述．

考慮內(nèi)點(diǎn)(xi,yj)處的五點(diǎn)差分格式．在(8)中u 的排列方式為自然排列．因此，點(diǎn)(xi,yj)處的數(shù)值解uij對(duì)應(yīng)于向量u 中的第m(j?1)+i 個(gè)分量，而對(duì)應(yīng)于矩陣U 第i 行第j 列的元素．

不失一般性，選取(8)中第m(j ?1)+i 個(gè)方程和(13)第i 行，第j 列的方程．

由 (8)，可知矩陣 I + τkL1的第 m(j ? 1)+i 行為

矩陣 I ? τkL2的第 m(j ? 1)+i 行為

由(15)和(16)，知(8)的第m(j ?1)+i 個(gè)方程為

同理，可得(13)第i 行，第j 列的方程為

經(jīng)化簡(jiǎn)，可以看出(17)與(18)相同．因此，在正則內(nèi)點(diǎn)處，兩種方法建立的方程是等價(jià)的．

對(duì)于非正則內(nèi)點(diǎn)處，我們可以在(17)和(18)中取i = 1,m (j = 1,2,··· ,n)和j =1,n (i=1,2,··· ,m)，并將邊界條件 u0,j=um+1,j=ui,0=ui,n+1=0 (i=1,2,··· ,m,j =1,2,··· ,n)代入其中得到．于是，兩種迭代格式建立的方程是等價(jià)的．

2.4 改進(jìn)的迭代格式的優(yōu)勢(shì)

雖然，經(jīng)典的迭代格式在每一次迭代過(guò)程中，只需交替解兩個(gè)具有三對(duì)角系數(shù)矩陣的方程(8)和(9)，但L1, L2均為mn×mn 的三對(duì)角矩陣，而且L1L2．這不僅會(huì)增大問(wèn)題的計(jì)算難度，而且會(huì)增加許多計(jì)算量．

改進(jìn)的迭代格式在每次迭代過(guò)程中，需要交替求解兩個(gè)較低維數(shù)的具有三對(duì)角系數(shù)矩陣的矩陣方程(13)和(14)，這等價(jià)于求解n 個(gè)m×m 和m 個(gè)n×n 具有三對(duì)角系數(shù)矩陣的線性方程組，而且這兩個(gè)矩陣方程的三對(duì)角系數(shù)矩陣在U 為方陣的情況下是相等的，區(qū)別僅在于：在(13)中，系數(shù)矩陣左乘未知矩陣，在(14)中，系數(shù)矩陣右乘未知矩陣，相比于經(jīng)典的交替方向隱式迭代格式，這是比較容易處理的．

為了比較兩種迭代格式在實(shí)際計(jì)算中的不同，現(xiàn)基于Gauss 消去法對(duì)求解(8)和(9)與(13)和(14)的過(guò)程，進(jìn)行計(jì)算量的預(yù)估．若使用追趕法[9]對(duì)上面四個(gè)方程進(jìn)行求解，則求解(8)和(9)的算法總計(jì)算量為

求解(13)和(14)的算法總計(jì)算量為

對(duì)比(19)和(20)可以發(fā)現(xiàn)，在盡可能避開(kāi)零參與運(yùn)算的情況下，相比于經(jīng)典的迭代格式，改進(jìn)的迭代格式在計(jì)算量上減少了很多，表1 和表2 也證實(shí)了這一點(diǎn)．事實(shí)上，經(jīng)典的迭代格式在實(shí)際計(jì)算中，若是盡可能避開(kāi)零參與運(yùn)算和存儲(chǔ)的情況，則會(huì)大大增加程序編寫(xiě)的難度和長(zhǎng)度．而改進(jìn)的迭代格式在這方面，是遠(yuǎn)遠(yuǎn)優(yōu)于經(jīng)典的迭代格式的．另外，由于改進(jìn)的迭代格式和經(jīng)典的迭代格式是等價(jià)的，這個(gè)結(jié)論已在第2.3 節(jié)給出證明，所以改進(jìn)的迭代格式在計(jì)算精度方面并沒(méi)有提升．但改進(jìn)的迭代格式在數(shù)據(jù)存儲(chǔ)方面，相比于經(jīng)典的的迭代格式要少得多，而且要比后者方便的多．

3 數(shù)值實(shí)驗(yàn)

在問(wèn)題(1)和(2)中，取

此時(shí)，問(wèn)題(1)和(2)的精確解為

記

取步長(zhǎng)

表1: 經(jīng)典的迭代格式的誤差及運(yùn)行時(shí)間

考慮迭代格式(13)和(14)，取步長(zhǎng)

表2: 改進(jìn)的迭代格式的誤差及運(yùn)行時(shí)間

從表1 和表2 可以直觀地看出，利用改進(jìn)的迭代格式的計(jì)算誤差與利用經(jīng)典的迭代格式計(jì)算誤差并沒(méi)有明顯的差異，但在相同的步長(zhǎng)下，經(jīng)典的迭代格式用時(shí)遠(yuǎn)遠(yuǎn)超過(guò)改進(jìn)的迭代格式，而且兩者之間的時(shí)間比也隨著步長(zhǎng)h 的縮小而增大，如：當(dāng)時(shí)，時(shí)間比為0.2188/0.1164 ≈ 1.9；當(dāng)時(shí)，時(shí)間比為2.1853/0.4411 ≈5.0．對(duì)于經(jīng)典的交替方向隱式迭代格式，當(dāng)步長(zhǎng)時(shí)，程序運(yùn)行時(shí)間過(guò)長(zhǎng)，而改進(jìn)的迭代格式仍然很好．從圖1 和圖2 可以看出，步長(zhǎng)越小，改進(jìn)的迭代格式對(duì)真解的擬合效果越好．

圖1: 經(jīng)典的迭代格式，(a)和(b)中的左圖為真解，右圖為數(shù)值解

圖2: 改進(jìn)的迭代格式，(a)和(b)中的左圖為真解，右圖為數(shù)值解

由表1 和表2 中的數(shù)據(jù)可知，兩種迭代格式的誤差的差值在10?10范圍內(nèi)，并且在以相同的步長(zhǎng)和迭代停止條件下，兩種迭代格式的迭代步數(shù)相等，這說(shuō)明在不考慮由求解器不同而產(chǎn)生的誤差時(shí)，每次迭代過(guò)程中的有效數(shù)據(jù)也是等價(jià)的，這與定理1 的結(jié)論一致．

4 結(jié)論

交替方向隱式迭代法是一種具有比一般迭代格式更高的收斂速度的求解Poisson 方程的方法，但經(jīng)典的交替方向隱式迭代格式在當(dāng)求解大規(guī)模問(wèn)題時(shí)，會(huì)產(chǎn)生計(jì)算程序編寫(xiě)困難，計(jì)算耗時(shí)較長(zhǎng)，計(jì)算過(guò)程中占用內(nèi)存過(guò)大等問(wèn)題．本文給出了改進(jìn)的交替方向隱式迭代格式(13)和(14)，通過(guò)求解較低維數(shù)的矩陣方程，降低了程序?qū)崿F(xiàn)難度，大大減少了計(jì)算量．在計(jì)算量和程序?qū)崿F(xiàn)等方面，改進(jìn)的迭代格式(13)和(14)遠(yuǎn)遠(yuǎn)優(yōu)于經(jīng)典的迭代格式(8)和(9)．同時(shí)，本文證明了兩種迭代格式的等價(jià)性．因此，改進(jìn)的迭代格式不僅具有經(jīng)典的迭代格式在理論分析層面的收斂性，而且在計(jì)算上更加簡(jiǎn)便快速．

然而，與經(jīng)典的交替方向隱式迭代格式一樣，改進(jìn)的交替方向隱式迭代格式僅適用求解區(qū)域?yàn)榫匦螀^(qū)域或可以轉(zhuǎn)化為矩形求解區(qū)域的線性問(wèn)題和半線性問(wèn)題．本文以五點(diǎn)差分格式為基礎(chǔ)，導(dǎo)出改進(jìn)的交替方向隱式迭代格式，事實(shí)上，仿照文中的過(guò)程，可以以九點(diǎn)差分格式導(dǎo)出更高階的交替方向隱式迭代格式，九點(diǎn)差分格式下的改進(jìn)的交替方向隱式迭代格式與五點(diǎn)差分格式下的迭代格式在計(jì)算過(guò)程方面的區(qū)別僅在于：五點(diǎn)差分格式得到的迭代格式需要求解兩個(gè)系數(shù)矩陣為三對(duì)角矩陣的線性矩陣方程，而九點(diǎn)差分格式下的迭代格式需要求解兩個(gè)系數(shù)矩陣為五對(duì)角矩陣的線性矩陣方程，在兩個(gè)方向取的節(jié)點(diǎn)數(shù)相同的情況下，同一差分格式對(duì)應(yīng)的矩陣方程系數(shù)矩陣是相等的．