三對角線性方程組演算法特點論文
三對角線性方程組演算法特點論文
一、概述
三對角線性方程組的求解是許多科學和工程計算中最重要也是最基本的問題之一。在核物理、流體力學、油藏工程、石油地震資料處理及數值天氣預報等許多領域的大規模科學工程和數值處理中都會遇到三對角系統的求解問題。很多三對角線性方程組的演算法可以直接推廣到求解塊三對角及帶狀線性方程組。由於在理論和實際應用上的重要性,近20年來三對角方程組的並行演算法研究十分活躍。
大規模科學計算需要高效能的平行計算機。隨著軟硬體技術的發展,高效能的平行計算機日新月異。現今,SMP可構成每秒幾十億次運算的系統,PVP和COW可構成每秒幾百億次運算的系統,而MPP和DSM可構成每秒萬億次運算或更高的系統。
高效能平行計算機只是給大型科學計算提供了計算工具。如何發揮平行計算機的潛在效能和對三對角系統進行有效求解,其關鍵在於抓住平行計算的特點進行並行演算法的研究和程式的設計與實現。另外,對處理機個數較多的平行計算系統,在設計並行演算法時必須解決演算法的可擴充套件性,並對可擴充套件性進行研究和分析。
二、問題的提出
設三對角線性方程組為
AX=Y(1)
式中:A∈Rn×n非奇異,αij=0,。X=(x1,x2,…xn)TY=(y1,y2,…yn)T。
此係統在許多演算法中被提出,因此研究其高效能並行演算法是很有理論和實際意義的。
三、並行求解三對角系統的直接解法
關於三對角線性方程組的直接求解已經有大量並行演算法,其中Wang的分裂法是最早針對實際硬體環境,基於分治策略提出的並行演算法。它不僅通訊結構簡單,容易推廣到一般帶狀線性方程組的並行求解,而且為相繼出現的許多其它並行演算法提供了可行的區域性分解策略。
近20年來求解三對角方程組的並行演算法都是基於分治策略,即透過將三對角方程組分解成P個小規模問題,求解這P個小規模問題,再將這些解結合起來得到原三對角方程組的解。一般求解三對角方程組的分治方法的計算過程可分為3個階段:一是消去,每臺處理機對子系統消元;二是求解縮減系統(需要通訊);三是回代,將縮減系統的解回代到每個子系統,求出最終結果。具體可分為以下幾類:
(一)遞推耦合演算法(RecursiveDoubling)
由Stone於1975年提出,演算法巧妙地把LU分解方法的時序性很強的遞推計算轉化為遞推倍增平行計算。D.J.Evans對此方法做了大量研究。P.Dubois和G.Rodrigue的研究表明Stone演算法是不穩定的。
(二)迴圈約化方法(CyclicReduction)
迴圈約化方法由Hockey和G.Golub在1965年提出,其基本思想是每次迭代將偶數編號方程中的奇變數消去,只剩下偶變數,問題轉變成求解僅由偶變數組成的規模減半的新三對角方程組。求解該新方程組,得到所有的偶變數後,再回代求解所有的奇變數。即約化和回代過程。由於其基本的算術操作可以向量化,適合於向量機。此方法有大量學者進行研究,提出了許多改進的方法。例如,Heller針對最後幾步的`短向量操作提出了不完全迴圈約化方法;R.Reulter結合IBM3090VF向量機的特點提出了局部迴圈約化法;P.Amodio針對分散式系統的特點改進了迴圈約化方法;最近針對此方法又提出對三對角方程組進行更大約化步的交替迭代策略。
(三)基於矩陣乘分解演算法
將係數矩陣A分解成A=FT,方程Ax=b化為Fy=b和Tx=y兩個方程組的並行求解。這種演算法又可以分為兩類:
1.分解。
如Wang的分裂法及其改進演算法就屬於這一類。P.Amodio在1993年對這類演算法進行了很好的總結,用本地LU、本地LUD和本地迴圈約化法求解,並在1995年提出基於矩陣乘分解的並行QR演算法。H.Michielse和A.VanderVorst改變Wang演算法的消元次序,提出了通訊量減少的演算法。李曉梅等將H.Michielse和A.VanderVorst演算法中的通訊模式從單向序列改為雙向並行,提出DPP演算法,是目前最好的三對角方程組分散式演算法之一。2000年駱志剛等中依據DPP演算法,利用計算與通訊重疊技術,減少處理機空閒時間取得了更好的並行效果。此類演算法要求解P-1階縮減系統。
2.不重疊分解。
例如Lawrie&Sameh演算法、Johsoon演算法、Baron演算法、Chawla在1991年提出的WZ分解演算法以及Mattor在1995年提出的演算法都屬於這一類。此類演算法要求解2P-2階縮減系統。
(四)基於矩陣和分解演算法
將係數矩陣分解成A=Ao+△A,這類演算法的共同特點是利用Sherman&Morrison公式將和的逆化為子矩陣逆的和。按矩陣分解方法,這種演算法又可分為兩類:
1.重疊分解。
這類演算法首先由Mehrmann在1990年提出,透過選擇好的分解在計算過程中保持原方程組係數矩陣的結構特性,具有好的數值穩定性,需要求解P-1階縮減系統。
2.不重疊分解。
Sun等在1992年提出的並行劃分LU演算法PPT演算法和並行對角佔優演算法PDD演算法均屬於這一類。需要求解2P-2階縮減系統。其中PDD演算法的通訊時間不隨處理機的變化而變化,具有很好的可擴充套件性。X.H.Sun和W.Zhang在2002年提出了兩層混合並行方法PTH,其基本思想是在PDD中嵌入一個內層三對角解法以形成一個兩層的並行,基本演算法是PDD,三對角系統首先基於PDD分解。PTH演算法也具有很好的可擴充套件性。
四、並行求解三對角系統的迭代解法
當稀疏線性方程組的係數矩陣不規則時,直接法在求解過程中會帶來大量非零元素,增加了計算量、通訊量和儲存量,並且直接法不易並行,不能滿足求解大規模問題的需要。因此通常使用迭代法來求解一般係數線性方程組和含零元素較多三對角線性方程組。迭代法包括古典迭代法和Krylov子空間迭代法。
古典迭代法包括Jacobi、Gauss-Seidel、SOR、SSOR等方法。通常採用紅黑排序、多色排序和多分裂等技術進行平行計算。
由於古典迭代法有收斂速度慢、並行效果不好等缺點,目前已較少用於直接求解大型稀疏線性方程組,而是作為預條件子和其它方法(如Krylov子空間方法)相結合使用。
Krylov子空間方法具有儲存量小,計算量小且易於並行等優點,非常適合於並行求解大型稀疏線性方程組。結合預條件子的Krylov子空間迭代法是目前並行求解大型稀疏線性方程組的最主要方法。
給定初值X0,求解稀疏線性方程組AX=Y。設Km為維子空間,一般投影方法是從m維仿射子空間X0+Km中尋找近似解Xm使之滿足Petrov-Galerkin條件:
Y-AXm┻Lm
其中Lm為另一個維子空間。如果Km是Krylov子空間,則上述投影方法稱為Krylov子空間方法。Krylov子空間Km(A,r0)定義為:
Km(A,r0)=span{r0,Ar0,A2r0,…,Am-1r0}
選取不同的Km和Lm就得到不同的Krylov子空間方法。主要演算法包括四類:基於正交投影方法、基於正交化方法、基於雙正交化方法、基於正規方程方法。
Krylov子空間迭代法的收斂速度依賴於係數矩陣特徵值的分佈,對於很多問題,直接使用迭代法的收斂速度特別慢,或者根本不收斂。因此使用預條件改變其收斂性,使中斷問題可解,並加速收斂速度是需要的。目前人們研究的預條件技術可分為四類:採用基於矩陣分裂的古典迭代法作為預條件子、採用不完全LU分解作預條件子、基於係數矩陣近似逆的預條件子、結合實際問題用多重網格或區域分解作預條件子。對Krylov子空間和預條件Krylov子空間方法有詳細的討論。
預條件Krylov子空間方法的平行計算問題一直是研究熱點,已提出了一系列好的並行演算法。目前預條件Krylov子空間方法的計算量主要集中在矩陣向量乘上。雖然學者們做了大量的研究工作,但是還沒找到效果好,又易於並行的預條件子。
需要特別指出的是,對於一般線性代數方程組的並行求解,其可擴充套件平行計算的研究已相對成熟,並已形成相應的並行軟體庫,如美國田納西亞州立大學和橡樹嶺國家實驗室研製的基於訊息傳遞計算平臺的可擴充套件線性代數程式庫ScaLAPACK和得克薩斯大學開發的介面更加友好的並行線性代數庫PLAPACK。我們借鑑其研究成果和研究方法,對三對角線性方程組並行演算法的研究是有幫助的。
五、結語
三對角線性方程組的直接解法,演算法豐富,程式較容易實現。但計算過程要增加計算量,並且大部分演算法都對係數矩陣的要求比較高。迭代解法適合於非零元素較多的情況,特別是結合預條件子的Krylov子空間迭代法已成為當前研究的熱點。
儘管三對角系統並行演算法的研究取得了很多成果。但是還存在一些問題:直接法中,分治策略帶來計算量和通訊量的增加,如何減少計算量和通訊量有待於進一步的研究;目前直接演算法均基於分治策略,如何把其它並行演算法設計技術,如平衡樹和流水線等技術應用到三對角系統的並行求解中也是需要引起重視的方向;對於非對稱系統還沒找到一種通用的Krylov子空間方法;Krylov子空間方法的並行實現時僅考慮係數矩陣與向量乘,對其它問題考慮不夠;以往設計的並行演算法缺乏對演算法可擴充套件性的考慮和分析。
【參考文獻】
[1]駱志剛,李曉梅,王正華.三對角線性方程組的一種有效分散式並行演算法[J].計算機研究與發展,2000,(7)