趙 鑫, 呂毅斌

(昆明理工大學(xué) 理學(xué)院, 昆明 650500)

0 引 言

共形映射可將復(fù)雜區(qū)域映射到正則狹縫域,同時(shí)保持曲線之間的夾角在大小和方向上不變,這一性質(zhì)常用于復(fù)雜區(qū)域中研究的問(wèn)題,例如應(yīng)用于流體力學(xué)、水波和彈性力學(xué)等領(lǐng)域內(nèi)的許多實(shí)際問(wèn)題[1-7]．但是,將復(fù)雜區(qū)域共形映射到正則狹縫域的解析式難以求得,針對(duì)這一問(wèn)題,許多研究者以計(jì)算機(jī)為工具,提出了有效的數(shù)值方法和思想[8-19]．例如:1966年,Symm[8-10]首次使用計(jì)算機(jī)計(jì)算了共形映射的數(shù)值解,該方法利用位勢(shì)理論提出了使用第一類Fredholm積分方程,計(jì)算了單連通區(qū)域和雙連通區(qū)域的數(shù)值共形映射映射函數(shù)．1988年,Amano[11]簡(jiǎn)化了Symm的方法,提出了基于模擬電荷法的數(shù)值共形映射計(jì)算法(Amano法),用于將任意單連通區(qū)域共形映射到單位圓盤(pán)．后經(jīng)不斷發(fā)展,其也可用于計(jì)算多連通區(qū)域到多種正則域的共形映射函數(shù)[12-15]．2009年,Nasser等[16-17]將求解共形映射函數(shù)的問(wèn)題轉(zhuǎn)化為求解一個(gè)Riemann-Hilbert問(wèn)題,提出了帶廣義Neumann核的第二類Fredholm積分方程法,并將其應(yīng)用于碳納米管增強(qiáng)復(fù)合材料內(nèi)穩(wěn)態(tài)熱傳導(dǎo)的仿真模擬．

1916年,Koebe[20]提出了39種具有廣泛應(yīng)用前景的正則狹縫域,并將其歸為5類．目前,雖然存在許多計(jì)算多連通區(qū)域共形映射的數(shù)值方法,但是其研究?jī)?nèi)容大都針對(duì)于將低連通度區(qū)域共形映射到較為簡(jiǎn)單的第一類正則狹縫域[14-16,21-25]．區(qū)域的連通度以及正則狹縫域復(fù)雜度的提高會(huì)導(dǎo)致共形映射的計(jì)算量和誤差大大增加,而模擬電荷法在計(jì)算共形映射的過(guò)程中不計(jì)算數(shù)值積分,并且使用最大模原理評(píng)價(jià)誤差,所以具有計(jì)算時(shí)間短、計(jì)算精度高、避免計(jì)算奇異積分等優(yōu)勢(shì)．因此,本文基于模擬電荷法,研究了將有界高連通度區(qū)域映射到更為復(fù)雜的,對(duì)數(shù)螺旋狹縫單位圓環(huán)[26]的共形映射計(jì)算法．首先,我們給出了正規(guī)化條件和共形映射的函數(shù)形式,并使用Laplace方程基本解的線性組合近似函數(shù)形式中的待定函數(shù)．其次,根據(jù)Dirichlet邊界條件建立未知量應(yīng)滿足的約束方程組,并提出利用BiCR算法[27]求解病態(tài)約束方程組,得到模擬電荷,進(jìn)而構(gòu)造出高精度的近似共形映射函數(shù)．最后,在數(shù)值實(shí)驗(yàn)中,通過(guò)對(duì)有界10連通區(qū)域和有界16連通區(qū)域進(jìn)行數(shù)值實(shí)驗(yàn),驗(yàn)證了本文算法的有效性,并模擬了有界高連通度區(qū)域內(nèi)螺旋點(diǎn)渦的繞流．

本文主要部分結(jié)構(gòu)如下:第1節(jié)給出了將有界高連通度區(qū)域映射到對(duì)數(shù)螺旋狹縫單位圓環(huán)(S)的共形映射函數(shù)與調(diào)和函數(shù)(a(z))的近似形式,并根據(jù)邊界條件構(gòu)造出約束方程組．第2節(jié)提出了利用BiCR算法計(jì)算約束方程的矩陣形式,得到了精度更高的模擬電荷．第3節(jié)簡(jiǎn)述了將本文共形映射計(jì)算法應(yīng)用于繞流模擬的原理．第4節(jié)通過(guò)數(shù)值實(shí)驗(yàn)驗(yàn)證了該算法的有效性,并給出了有界高連通度區(qū)域內(nèi)螺旋點(diǎn)渦的繞流仿真結(jié)果．第5節(jié)對(duì)全文進(jìn)行了總結(jié)和展望．

1 有界高連通度區(qū)域到對(duì)數(shù)螺旋狹縫單位圓環(huán)域的共形映射

1.1 共形映射函數(shù)及邊界條件

在z平面上,區(qū)域D是一個(gè)連通度為n(≥10)的有界高連通度區(qū)域,其邊界?D=C由n條封閉的Jordan曲線C1,C2,…,Cn組成,即C=C1∪C2∪…∪Cn．共形映射函數(shù)w=f(z)將曲線C1映射到單位圓|w|=R1=1,將曲線C2映射到圓盤(pán)|w|=R2,將曲線Cm映射到單位圓環(huán)內(nèi)斜角為θm的對(duì)數(shù)螺旋狹縫Sm(m=3,4,…,n),見(jiàn)圖1,圖中“·”代表約束點(diǎn),“+”代表模擬電荷點(diǎn)(后同)．斜角是指對(duì)數(shù)螺旋狹縫與從坐標(biāo)系原點(diǎn)發(fā)出的射線相交所形成的夾角,始終是一個(gè)固定值θ．值得注意的是, 當(dāng)斜角θ=0和θ=π/2時(shí), 對(duì)數(shù)螺旋狹縫將分別退化為從原點(diǎn)發(fā)出的射線和以原點(diǎn)為圓心的圓?。C上所述, 共形映射函數(shù)w=f(z)的邊界值應(yīng)滿足

(1)

圖1 基于模擬電荷法的有界高連通度區(qū)域共形映射Fig. 1 Conformal mappings of bounded high connectivity regions based on the charge simulation method

其中,z∈Cm;R2,…,Rn為待定實(shí)常數(shù);θ1=θ2=π/2,θ3,…,θn為給定實(shí)常數(shù)．

對(duì)于有界多連通區(qū)域D,滿足正規(guī)化條件f(0)>0的共形映射函數(shù)w=f(z)是唯一確定的,函數(shù)形式可定義為

(2)

其中,c=f(0)為正實(shí)常數(shù),v是邊界C2內(nèi)側(cè)一定點(diǎn),將式(2)代入式(1),得

z∈Cm,m=1,2,…,n．

(3)

1.2 模擬電荷法的原理

根據(jù)模擬電荷法的原理,a(z)可由Laplace方程基本解的線性組合近似,即

(4)

其中,Q0為未知復(fù)常數(shù),電荷Qlj為待定實(shí)數(shù),ζlj是配置在曲線Cl外的模擬電荷點(diǎn),見(jiàn)圖1．為了將有界多連通區(qū)域映射到區(qū)域S,式(2)需滿足以下條件[15]:

即

(5)

② 保持近似函數(shù)在問(wèn)題區(qū)域D上的坐標(biāo)系放縮不變性[15,28],有

(6)

③ 正規(guī)化條件,由f(0)=c可得

(7)

結(jié)合式(4),消去Q0,可得到

(8)

(9)

即可保證A(z)在區(qū)域D中連續(xù),進(jìn)而能夠在問(wèn)題區(qū)域上構(gòu)造一個(gè)連續(xù)的近似映射函數(shù)．

根據(jù)上述條件①、②和③,通過(guò)限制近似函數(shù)(9)中的z取Cm上的約束點(diǎn)zmk,令其滿足邊界條件(1),則可得到關(guān)于電荷Qlj和Hm滿足的N1+N2+…+Nn+n維約束方程組,即

zmk∈Cm;m=1,2,…,n;k=1,2,…,Nm．

(10)

將求得的電荷Qlj代入式(9),可解得a(z)的近似函數(shù)A(z)．隨后將A(z)代入式(2),便可得到共形映射函數(shù)f(z)的近似函數(shù),記為F(z)．

2 求解模擬電荷的BiCR法

現(xiàn)將約束方程組(10)寫(xiě)成標(biāo)準(zhǔn)線性方程組Ax=b的形式:

A∈(N1+N2+…+Nn+n)×(N1+N2+…+Nn+n),x∈(N1+N2+…+Nn+n),b∈(N1+N2+…+Nn+n),

其中,A是非對(duì)稱的,Nm(1,2,…,n)表示每條邊界上的模擬電荷點(diǎn)數(shù)．

在基于模擬電荷法的共形映射計(jì)算法中,線性方程組的求解精度對(duì)最終的映射結(jié)果起著關(guān)鍵性作用．?dāng)?shù)值共形映射的精度要求越高,所取的模擬電荷點(diǎn)數(shù)就需要越多,但線性方程組的系數(shù)矩陣隨著模擬電荷點(diǎn)數(shù)的增加,條件數(shù)cond(A)也隨之增大．所以,針對(duì)病態(tài)矩陣A,本文使用BiCR算法求解線性方程組．BiCR算法是由Sogabe在conjugate residual(CR)算法的基礎(chǔ)上提出的,可以高效求解非對(duì)稱線性方程組,從而提高所求模擬電荷的精度,進(jìn)而能夠得到高精度的近似共形映射函數(shù)．求解模擬電荷的BiCR法標(biāo)準(zhǔn)步驟如下．

算法1 求解模擬電荷的BiCR法

1 輸入A,b,x0,ε,nmax

4 forn=0,1,2,…

5 if |rn|>ε,n

6pn=rn+βn-1pn-1;

9xn+1=xn+αnpn;

10rn+1=rn-αnApn;

13 else break; end if

14 end for

15 輸出xn+1．

上述算法中,x0是給定的初始近似值,取值為零向量,r0=b-Ax0是初始向量,ε是給定誤差,nmax為最大迭代次數(shù)．

3 螺旋點(diǎn)渦的繞流模擬

在流體力學(xué)領(lǐng)域中,繞流是指繞過(guò)置于無(wú)限流體中的物體的流動(dòng),或物體無(wú)限流體中運(yùn)動(dòng),是自然界和工程中常見(jiàn)的黏性流體運(yùn)動(dòng)形式．其在能源動(dòng)力工程、海洋工程和航空航天工程等領(lǐng)域中具有重要意義[29-31]．通過(guò)使用共形映射將復(fù)雜的多連通區(qū)域映射為簡(jiǎn)單的正則狹縫域后,根據(jù)問(wèn)題域內(nèi)的復(fù)勢(shì)可以由正則域上的復(fù)勢(shì)經(jīng)過(guò)共形映射得到,可計(jì)算出平面無(wú)旋流在有界多連通區(qū)域內(nèi)各點(diǎn)的復(fù)勢(shì),進(jìn)而模擬出該平面無(wú)旋流的流線和等勢(shì)線．使用共形映射計(jì)算多連通區(qū)域內(nèi)的復(fù)勢(shì)的原理如下:

設(shè)D(z)=φ(x,y)+iψ(x,y)和S(w)分別是區(qū)域D和區(qū)域S內(nèi)的復(fù)勢(shì),其中φ表示速度勢(shì),ψ表示流函數(shù)．共形映射函數(shù)w=f(z)將區(qū)域D映射到區(qū)域S,同時(shí)也將區(qū)域D內(nèi)的流線族映射到區(qū)域S內(nèi)的流線上．原因是定義在區(qū)域D內(nèi)的調(diào)和函數(shù)經(jīng)過(guò)共形映射后仍是調(diào)和函數(shù),因此,當(dāng)S(w)已知時(shí),區(qū)域D內(nèi)的復(fù)勢(shì)可表示為

S(w)=S(f(z))=D(z),z∈D．

(11)

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

在MATLAB 2018b的環(huán)境下,對(duì)連通度為n(≥10)的有界高連通度區(qū)域共形映射到對(duì)數(shù)螺旋狹縫單位圓環(huán)進(jìn)行數(shù)值實(shí)驗(yàn),檢驗(yàn)本文算法的有效性和精度．共形映射誤差由

(12)

步驟1 設(shè)定正整數(shù)Nl,實(shí)常數(shù)θ1=π/2,θ2=π/2,θ3,…,θn,合適的模擬電荷點(diǎn)ζlj(l=1,2,…,n;j=1,2,…,Nl)與對(duì)應(yīng)數(shù)量的約束點(diǎn)zmk(m=1,2,…,n;k=1,2,…,Nm),邊界Cl內(nèi)側(cè)的相應(yīng)定點(diǎn)ζl0(l=1,2,…,n)和v．

步驟2 將ζlj,zmk,ζl0代入式(10),得到將有界高階連通區(qū)域映射到區(qū)域S的約束方程組,求解約束方程組,得到實(shí)常數(shù)Hl、正實(shí)常數(shù)c和電荷Qlj．

在步驟2求解約束方程組時(shí),本文使用BiCR算法進(jìn)行求解,并與同是迭代法的GMRES(m)(the generalized minimal residual)算法[12,34]進(jìn)行比較,將前者記為method 1,后者記為method 2．廣義極小殘差法GMRES(m)是在Krylov子空間討論求解線性方程組的方法,可以直接應(yīng)用于求解系數(shù)矩陣密集且不對(duì)稱的大型線性方程組．

例1 將z平面上的有界10連通區(qū)域(圖2)共形映射到w平面上具有8條對(duì)數(shù)螺旋狹縫的單位圓環(huán)(圖3)．給定圖3中對(duì)數(shù)螺旋狹縫的斜角為

圖2 有界10連通區(qū)域網(wǎng)格圖 圖3 帶有對(duì)數(shù)螺旋狹縫的單位圓環(huán) Fig. 2 The grid diagram of bounded 10 connected domains Fig. 3 The unit circular ring with logarithmic spiral slits

圖4 例1中共形映射誤差曲線 圖5 例1中螺旋點(diǎn)渦的繞流模擬 Fig. 4 The conformal mapping error curves in example 1 Fig. 5 Simulation of flow over the spiral point vortex in example 1

圖6 有界16連通區(qū)域網(wǎng)格圖 圖7 帶有對(duì)數(shù)螺旋狹縫的單位圓環(huán) Fig. 6 The grid diagram of bounded 16 connected domains Fig. 7 The unit circular ring with logarithmic spiral slits

圖8 例2中共形映射誤差曲線 圖9 例2中螺旋點(diǎn)渦的繞流模擬 Fig. 8 The conformal mapping errors curves in example 2 Fig. 9 Simulation of flow over the spiral point vortex in example 2

5 結(jié) 論

本文基于模擬電荷法,提出了將有界高連通度區(qū)域映射到對(duì)數(shù)螺旋狹縫單位圓環(huán)的共形映射計(jì)算法,模擬了有界高連通度區(qū)域內(nèi)螺旋點(diǎn)渦的繞流．在約束方程組的求解中,本文提出利用BiCR算法構(gòu)造高精度的共形映射函數(shù),得到了比基于GMRES(m)算法改進(jìn)的共形映射計(jì)算法更好的新算法,即基于BiCR算法改進(jìn)的共形映射計(jì)算法,并通過(guò)數(shù)值實(shí)驗(yàn)驗(yàn)證了本文算法的有效性．在今后的研究中,我們將繼續(xù)研究具有廣泛應(yīng)用前景的正則狹縫域,其次我們也會(huì)將共形映射用于求解熱傳導(dǎo)和圖像處理等領(lǐng)域的實(shí)際問(wèn)題．