基于LSQR法的外部數值保角逆變換計算法

呂毅斌,王 堅,王櫻子,吳 爽

(1.昆明理工大學理學院,云南昆明 650500)

(2.昆明理工大學計算中心,云南昆明 650500)

1 引言

保角變換是復變函數中非常重要的理論之一,廣泛應用于物理學和工學領域.特別是在電磁理論、膜和板的振動、彈性理論、熱傳輸、流體力學等方面有很多應用[1?5].通常將保角變換的求解方法分為解析法和數值法.解析法只能在一些特殊區域給出變換函數表達式,對于復雜區域問題沒有解決方法.因此,對于很多實際中的復雜問題必須采用數值法求解保角變換函數.數值保角變換計算法的主要有:積分方程式法[6]、正交多項式法[7?8]和有限差分法[9?10]等.德國人Steinbigler[11]首次提出用若干個虛設電荷來模擬電極表面上電荷分布的電場計算法,形成了模擬電荷法的基本思想;日本的天野要等數學學者從20世紀80年代開始對模擬電荷法和數值保角變換作了大量研究工作,并提出了基于模擬電荷法的數值保角變換計算法(天野法)[12?16].天野法適用于單連通區域及多連通區域的數值保角變換問題[15,17,18,19,20].

單連通區域的數值保角變換,分為內部數值保角變換和外部數值保角變換.本文通過對基于模擬電荷法的雙方向的內部數值保角變換計算法[13,15,16]的研究,提出了外部數值保角逆變換計算法.該方法的原理是基于模擬電荷法來求Laplace方程的Dirichlet問題的解,并通過預先建立的邊界對應關系構造從標準區域到問題區域的近似保角逆變換函數,誤差用正則函數的最大值原理進行評價.計算數值保角逆變換最主要的就是確定模擬電荷的位置和數量,以及約束方程組的求解.

LSQR方法[21?24]是Paige和Saunders提出的一種適用于求解系數矩陣為大型、稀疏矩陣線性方程組的方法.LSQR方法求解的思路是把任意系數矩陣方程化為系數矩陣為方陣的方程,然后利用Lanczos方法,求解最小二乘解.由于在求解過程中應用到QR分解,因此稱為LSQR(Least Square QR-factorization)方法.本文利用LSQR方法求解出了外部數值保角變換模擬電荷法中的約束方程組,得到電荷量和逆變換半徑,從而構造出近似逆保角變換函數,最后利用數值實驗驗證了所提計算法的有效性.

2 基于模擬電荷法的外部數值保角變換計算法

本節主要闡述了基于模擬電荷法的外部區域數值保角變換計算法[12,15,17,18].如圖1所示,C對于z平面上的任意Jordan曲線,圍繞C的外部區域為D,通過數值保角變換將區域D映射成w平面上的單位圓的外部|w|>1,=D+C.

圖1:基于模擬電荷法的外部區域數值保角變換（+代表模擬電荷點,?代表約束點）

保角變換函數w=f(z),f(z)滿足正規化條件f(∞)=∞,f(∞)>0時,可以表示如下

g(z)是Dirichlet型場勢問題

的解,其中h(z)是g(z)的共軛調和函數,且h(∞)=0.在下面的敘述中,F,G,H,Γ表示f,g,h,γ的近似值.由模擬電荷法,g(z)可以用C圍繞的區域內部里配置的N 個電荷點ξj作為極的對數勢場的1次結合

高度近似g(z).此時h(z)的高度近似函數為

未知電荷Qj通過滿足下面邊界條件進行求解,即

另外,由條件g(∞)=0,h(∞)=0,可得

因此,通過(2.5)式和(2.6)式能推導出以Qj(j=1,2,···,N)和logΓ作為未知數的N+1維線性方程組的構成如下

式中,

通過式(2.3),(2.4)和(2.7)得到近似保角變換函數

最后利用zi,Qj,Γ,ζj計算雙連通保角變換.

3 外部區域數值保角逆變換計算法

根據上節內容,本節提出基于模擬電荷法的外部區域數值逆保角變換計算法.如圖2在w平面上,單位圓圍成的外部區域|w|>1,通過數值保角逆變換將數值保角正變換映射成的單位圓的邊界及外部區域變換成z平面上封閉的Jordan曲線C及所圍成的外部區域D[16].

圖2:基于模擬電荷法的外部區域數值保角逆變換（+代表模擬電荷點,?代表約束點）

在不失一般性的情況下,假定映射函數z=f?(w)滿足正規化條件f?(∞)= ∞,f?(∞)>0 時是正則的,即

式中,γ?是變換半徑,g?(w)是Dirichlet型勢場問題

的解,其中 h?(w)是 g?(w)的共軛調和函數,且 h?(∞)=0.分別用 F?,G?,H?,Γ?表示 f?,g?,h?,γ?的近似值.

根據模擬電荷法,g?(w)可以用單位圓的內部區域配置的N?個電荷點ζ?j作為極的對數勢場的1次結合

高度近似,這里h?(w)可以用

高度近似.

同時,由條件 g?(∞)=0,h?(∞)=0 可得

其中zi是外部的數值正保角變換的約束點,wi是經過zi數值正保角變換得到的映射結果,通過wi來確定.因此,由(3.5)和(3.6)式可得以(1≤j≤N?)和logΓ?作為未知數的(N?+1)維線性方程組如下

其中

通過(3.3),(3.4)和(3.7)式可得近似保角逆變換函數

最后利用zi,wi,,Γ?,計算外部區域數值保角逆變換.

4 基于LSQR方法的數值保角逆變換模擬電荷求解

將約束方程組(3.7)式寫成標準線性方程組的

(3)在今后的研究中可以繼續聯合實地監測數據，除植被因素外，將景觀要素和土壤要素以及周邊居民滿意度等要素，在生態重建效果評價中的重要性考慮進去。另外下一步工作中可以進一步結合多種評價方法，例如和層次分析法、灰色關聯度法、聚類分析法、模糊綜合評價法等做對比，對研究區的生態重建效果進行全面評價比較和分析。

形式,其中

約束方程的系數矩陣A是非對稱的且病態的,LSQR方法[20?23]是求解系數矩陣為病態的大型稀疏矩陣線性方程組的有效算法之一.利用Lanczos雙對角化方法來求解方程的最小二乘解 minkAx?bk2. 假定 Uk=[u1,···,uk]和 Vk=[v1,···,vk]是正交陣且 Lk為如下的(k+1)×k的下雙對角陣

用下列迭代方法可以實現A矩陣的雙對角分解

其中αi≥0,βi≥0.使上式(3.21)可寫成

可以確定

在滿足給定精度時停止迭代.我們希望krkk2盡量小,且Uk+1理論上是正交陣,取yk使ktk+1k2最小,解最小二乘問題minkβ1e1?Lkykk2.得到LSQR算法如下:

Algorithm 1 LSQR Algorithm Input:A,b,x0,ε.Initialize β1u1=b,α1v1=ATu1,h1=v1,e?1=β1,eρ1=α1.for i=1,2,3,···while stopping criterion is not satisfied do βi+1ui+1=Avi?αiui;αi+1ui+1=ATi?βi+1Vi;ρi=eρ2i+β2i+1,ci=eρi/ρi,si=βi+1/ρi;θi+1=siαi+1,eρi+1=?ciαi+1,?i=cie?i,e?i+1=sie?i;xi=xi?1+(e?i/eρi)hi;hi+1=vi+1? (θi+1/ρi)hi;if minkAxi?bk<ε;end if end while end for Output xi.q

上述算法中,ε是給定精度.

這里給出基于LSQR法的外部區域數值保角逆變換計算法的具體步驟如下.

步驟1通過外部區域數值保角正變換(2.9)得到映射點F(zi),將F(zi)的位置作為數值保角逆變換的約束點wi的位置.

步驟2根據約束點wi的位置配置外部區域數值保角逆變換模擬電荷點的位置.

步驟3 通過LSQR 方法求解約束方程組(3.7)得到模擬電荷,,···,和逆變換半徑 logΓ?.

步驟4對單位圓的邊界及外部區域的每一個點通過(3.3)和(3.4)式計算得到G?(w)和H?(w)后,構造近似保角逆變換函數(3.9),然后計算對應的變換點.

5 數值實驗

針對橙形為邊界的外部區域,在MATLAB 13b環境下,檢驗雙方向的外部區域數值保角變換計算方法的有效性.基于模擬電荷法的單連通區域的外部區域數值正保角變換的誤差由Ez=max(||f(z)|?1|)確定,外部數值保角逆變換的誤差由Ew=max(|f?(w)?z|)確定.

例1橙形邊界及外部區域的雙方向數值保角變換.邊界

約束點的位置由

確定,其中i=0,1,···,N?1.約束點分布在邊界上,邊界由粗實線表示.

保角正變換模擬電荷點的位置由下面公式給出

保角逆變換模擬電荷點的位置由下面公式給出

其中a=21/16,rz=rw=3,r>0是確定模擬電荷點位置的參數,模擬電荷點的分布在區域D的外部,如圖3所示為數值保角正變換的模擬電荷點分布.如圖4所示為數值保角逆變換的模擬電荷點分布.約束點和模擬電荷點一一對應,數量都為N.

圖3:橙形保角變換邊界及電荷點位置

圖4:橙形保角逆變換邊界及電荷點位置

圖5:橙形保角逆變換誤差曲線

圖6:橙形保角逆變換誤差曲線

分別用Amano和lsqr分別表示天野法和基于LSQR方法的單連通外部區域數值保角逆變換計算法.圖5給出的是當a=21/16,rz=rw=3兩種數值保角逆變換方法的誤差曲線,圖5說明誤差隨著電荷量的增大而減小,同時可看出lsqr的誤差值一直小于Amano的誤差值,在N=61時,Amano誤差為2.4243×10?5而lsqr誤差為1.3025×10?6,說明了本文采用的算法可以得到更高的誤差精度,驗證了算法的有效性.

圖6是當a=21/14,rz=rw=2時兩種方法的誤差曲線,由圖6可看出電荷點越多誤差值越小,且各個數量的電荷量上,lsqr的誤差值均比Amano的誤差值小.電荷點數為180時,Amano誤差為1.0366×10?7,lsqr誤差為1.6639×10?8,因此數值實驗再次驗證了外部數值保角逆變換計算法的有效性.

圖7:橙形邊界和外部區域及等高線

圖8:圖7的保角變換

圖9:圖10的保角逆變換

圖10:圖8的邊界和其外部區域及其等高線

圖7–10中的粗實線表示邊界,細實線表示等高線.圖7是橙形的邊界及外部區域等高線,圖8是圖7通過近似保角變換函數F(z)映射后得到的結果.由圖8可知,近似保角變換函數F(z)將橙形的邊界映射成了單位圓.圖10表示的是單位圓的邊界和外部區域|w|>1及其等高線,以及模擬電荷點的配置位置.圖9是圖10通過近似保角逆變換函數F?(w)映射后得到的結果.由圖9和圖10可知,近似逆保角變換函數F?(w)將橙單位圓邊界映射成了橙形邊界,實現了數值保角逆變換.

6 結束語

本文利用LSQR方法求解出了基于模擬電荷法的外部數值保角逆變換中的約束方程組,提出了LSQR方法的外部數值保角逆變換計算法.并通過數值實驗驗證了所提計算法的有效性并用等高線模擬了外部數值保角逆變換的計算結果.在今后的研究中,本方法同樣可以應用于多連通區域的數值保角逆變換問題.