王小花,潘文峰
熱傳導(dǎo)方程的完美匹配層公式及其穩(wěn)定性分析
王小花,潘文峰
(武漢理工大學(xué)理學(xué)院,湖北武漢430070)
為了求解無(wú)界空間中的熱傳導(dǎo)方程,基于Laplace變換,引入若干個(gè)輔助變量,提出了一個(gè)熱傳導(dǎo)方程的完美匹配層(PML)公式。通過分析偏微分算子特征值實(shí)部的符號(hào)和特征向量的完備性,得到了PML方程的穩(wěn)定性。在二維空間中,常系數(shù)PML方程的柯西問題是弱穩(wěn)定;在三維空間中,常系數(shù)PML方程的柯西問題是強(qiáng)穩(wěn)定。數(shù)值實(shí)驗(yàn)結(jié)果表明:熱傳導(dǎo)方程PML公式的絕對(duì)誤差最大值大約是1.5×10-3,經(jīng)典Dirichlet邊界條件和Neumann邊界條件的絕對(duì)誤差最大值大約是2.5×10-2和3.0×10-2。因此,熱傳導(dǎo)方程PML公式可以顯著提高數(shù)值解的準(zhǔn)確性。
熱傳導(dǎo)方程;完美匹配層;穩(wěn)定性;輔助變量;Laplace變換
熱量在無(wú)界介質(zhì)中的傳播是一個(gè)重要而又基本的問題。對(duì)于波的傳播問題,完美匹配層(perfectly matched layer,PML)算法[1]已經(jīng)被證明是一種有效且精確的方法?;诙S傅里葉變換,文獻(xiàn)[2]推導(dǎo)出了一個(gè)關(guān)于對(duì)流擴(kuò)散方程的PML公式,并且證明了在吸收層的外邊界,反射是非常小的。同時(shí),數(shù)值實(shí)驗(yàn)也驗(yàn)證了完美匹配層算法的有效性和準(zhǔn)確性。對(duì)于波方程,目前有多種非常有效的PML算法,例如輔助微分方程完美匹配層(auxiliary differential equation perfectly matched layer,ADE-PML)算法[3]和復(fù)頻移完美匹配層(complex frequency shifted perfectly matched layer,CFS-PML)算法[4]。對(duì)于依賴時(shí)間的熱傳導(dǎo)方程,方程的穩(wěn)定性是非常重要的,不穩(wěn)定的數(shù)值方法可能會(huì)導(dǎo)致計(jì)算區(qū)域內(nèi)的數(shù)值解被污染。已經(jīng)有很多學(xué)者研究了偏微分方程的穩(wěn)定性[5-8],然而,關(guān)于熱傳導(dǎo)問題的完美匹配層算法的研究還比較少,特別是類似于波方程的ADE-PML算法。本文基于Laplace變換,引入若干個(gè)輔助變量,嘗試推導(dǎo)出一個(gè)適用于熱傳導(dǎo)方程的ADE-PML算法,以期提高數(shù)值解的準(zhǔn)確性,并針對(duì)熱傳導(dǎo)方程的PML公式進(jìn)行穩(wěn)定性分析。
考慮一個(gè)溫度場(chǎng)u在一個(gè)無(wú)界的三維空間中進(jìn)行傳導(dǎo),并假設(shè)所有的熱源初始分布集中在一個(gè)給定的區(qū)域Ω=[-a1,a1]×[-a2,a2]×[-a3,a3],a1,a2,a3>0,因此,在區(qū)域Ω內(nèi)部,溫度場(chǎng)u(x1,x2,x3,t)滿足如下方程:
本文希望截?cái)酂o(wú)界區(qū)域并將計(jì)算限制在有限的區(qū)域Ω內(nèi)部。因此,假設(shè)Ω被包含于一個(gè)厚度為L(zhǎng)i(i=1,2,3)的PML,其中Li(i=1,2,3)為在xi(i=1,2,3)方向的厚度。在吸收層的內(nèi)部,溫度場(chǎng)u滿足一個(gè)修正的無(wú)熱源的以指數(shù)速度衰減的熱傳導(dǎo)方程。
類似于文獻(xiàn)[5,9-10],設(shè)u^表示u的Laplace變換,定義如下:
引進(jìn)坐標(biāo)變換
在區(qū)域Ω之外,u^滿足如下方程:
其中:ζi在吸收層xi>ai是正的,在Ω內(nèi)部是0。若要求u^滿足如下的方程
顯然u在Ω的內(nèi)部保持不變,但是在吸收層內(nèi)呈指數(shù)衰減。
文獻(xiàn)[5]在Laplace變換域上給出了雙曲拋物系統(tǒng)標(biāo)準(zhǔn)的PML算法,但是涉及到了太多的輔助變量。本文嘗試使用較少的輔助變量,推導(dǎo)出一個(gè)新的關(guān)于熱傳導(dǎo)方程的PML公式。
由坐標(biāo)變換式(4)可知:
設(shè)γi=γi(ζi;s),i=1,2,3,則:
根據(jù)式(6)可得:
根據(jù)式(7),可以得到如下的恒等式:
結(jié)合式(8)和式(9),可得:
接下來,引進(jìn)一些輔助變量:
使用Laplace逆變換將方程(10)變換回時(shí)域,可得:
其中:
在Ω內(nèi)部,吸收函數(shù)ζi(i=1,2,3)和輔助變量φ取值為0,因此,方程(12)又簡(jiǎn)化成了方程(1)。可以看到:該P(yáng)ML公式只需要5個(gè)輔助變量φ1,φ2,φ3,ψ1,ψ2,并且沒有涉及到高階導(dǎo)數(shù),算法的實(shí)現(xiàn)更加簡(jiǎn)單。
類似的推導(dǎo)可知,在二維空間中不存在變量ψ2,φ3,ζ3,此時(shí)PML方程簡(jiǎn)化成如下形式:
其中:
本節(jié)中,將依次討論二維和三維情況下PML公式的穩(wěn)定性。對(duì)于拋物系統(tǒng)的穩(wěn)定性,文獻(xiàn)[11]有比較詳細(xì)的介紹,本文將引用文獻(xiàn)[11]中相關(guān)的理論證明方程(12)和方程(13)的穩(wěn)定性。
考慮一個(gè)一般的偏微分方程的柯西問題:
根據(jù)文獻(xiàn)[11]可知:柯西問題的適定性和穩(wěn)定性,可以通過算子P(iω)的特征值和特征向量來判斷。下面給出一個(gè)關(guān)于判別適定性和穩(wěn)定性的充要條件。
引理1[11]若存在不依賴于ω的常數(shù)C>0(resp.C=0)使得算子P(iω)的所有特征值λ滿足:
則柯西問題(14)是弱適定的(resp.弱穩(wěn)定的)。若對(duì)應(yīng)的特征向量是完備的,則柯西問題(14)是適定的(resp.穩(wěn)定的)。
利用以上的穩(wěn)定性理論,可以得到下面兩個(gè)穩(wěn)定性定理。
定理1在二維空間中,對(duì)于任意的ζ1,ζ2>0,若ζ1=ζ2,則PML公式的初值問題是弱穩(wěn)定的。
證明設(shè)U=(u,φ1,φ2,ψ)T,可以將PML方程(13)寫成如下形式:
其中:
設(shè)T=det(λI-P1(iω)),則:
由式(17)可知:只有ζ1=ζ2才能判斷PML方程的穩(wěn)定性。此外,當(dāng)ζ1=ζ2時(shí),矩陣E不可對(duì)角化,因此PML方程是弱穩(wěn)定的。
定理2在三維空間中,對(duì)于任意的ζ1,ζ2,ζ3>0,PML公式的初值問題是強(qiáng)穩(wěn)定的。
證明令U=(u,φ1,φ2,φ3,ψ1,ψ2),則可將式(12)寫成如下形式:
其中:
因此,P(iω)的所有特征值的實(shí)部小于等于0,同時(shí),矩陣E可對(duì)角化,故在三維空間中,對(duì)于任意的ζ1,ζ2,ζ3>0,PML公式的初值問題是強(qiáng)穩(wěn)定的。
本節(jié)通過具體的數(shù)值算例說明具有輔助微分方程的完美匹配層方法的有效性,并與經(jīng)典的狄利克雷(Dirichlet)邊界條件和紐曼(Neumann)邊界條件的結(jié)果相比較??紤]如下的初值問題:
其中:u0(x1,x2)=為初始的熱量分布,并且參數(shù)α=1。此時(shí),初值問題(20)具有解析解
對(duì)應(yīng)的PML設(shè)置如下:
其中:n=6;ζ=100;計(jì)算區(qū)域Ω=[-2,2]2;吸收層厚度L1=L2=0.3。對(duì)于對(duì)應(yīng)的PML方程,本文使用標(biāo)準(zhǔn)的有限元方法求解,采用線性基函數(shù),網(wǎng)格的最大尺度是0.1,在PML吸收層的外邊界,采用齊次的Dirichlet邊界條件。
當(dāng)t=0.7時(shí),不同的邊界條件所對(duì)應(yīng)的絕對(duì)誤差,如圖1所示。從圖1可以看出:PML邊界條件的絕對(duì)誤差最大值大約是1.5×10-3(見圖1a),經(jīng)典的Dirichlet邊界條件和Neumann邊界條件的絕對(duì)誤差最大值分別約為2.5×10-2和3.0×10-2(見圖1b和圖1c),因此,PML方法可以顯著提高數(shù)值解的準(zhǔn)確性。為了進(jìn)一步說明PML邊界條件對(duì)原方程的影響,取t=0.7,x1=0,此時(shí),數(shù)值解和解析解的差別如圖2所示。
圖1 不同邊界條件的絕對(duì)誤差
從圖2可以看出:在內(nèi)部區(qū)域,數(shù)值解和解析解幾乎完全重合;在吸收層內(nèi)部,u迅速衰減,到達(dá)吸收層的外部時(shí)幾乎為0,與波方程的結(jié)果類似。
本文推導(dǎo)出了一個(gè)熱傳導(dǎo)方程的PML公式,與文獻(xiàn)[5]中關(guān)于雙曲拋物系統(tǒng)的PML公式相比,本文的PML公式使用了較少的輔助變量。此外,證明了PML方程的柯西問題的穩(wěn)定性。在二維空間中,常系數(shù)的完美匹配層方程的柯西問題是弱穩(wěn)定的;在三維空間中,常系數(shù)的完美匹配層方程的柯西問題是強(qiáng)穩(wěn)定的。通過數(shù)值實(shí)驗(yàn)驗(yàn)證了熱傳導(dǎo)方程PML公式可以顯著提高數(shù)值解的準(zhǔn)確性。
[1]DRUSKIN V,GTTEL S,KNIZHNERMAN L.Near-optimal perfectly matched layers for indefinite Helmholtz problems[J].Mathematics,2015,58(1):686-710.
[2]LANTOS N,NATAF F.Perfectly matched layers for the heat and advection-diffusion equations[J].Journal of computational physics,2010,229(24):9042-9052.
[3]ATLURI S N,MARTIN R,KOMATITSCH D,et al.A high-order time and space formulation of the unsplit perfectly matched layer for the seismic wave equation using auxiliary dfferential equations(ADE-PML)[J].Computer modeling in egineering&sciences,2010,56(1):17-40.
[4]MA Y,CHEN S W,LU F,et al.The CFS-PML for 2-D WLP-FDTD method of dispersive materials[J].International journal of applied electromagnetics and mechanics,2016,51(4):349-361.
[5]APPELD,HAGSTROM T,KREISS G.Perfectly matched layers for hyperbolic systems:general formulation,wellposedness,and stability[J].SIAM journal on applied mathematics,2006,67(1):1-23.
[6]楊彩虹,胡志興.一類具有飽和發(fā)生率的SEIR模型的穩(wěn)定性[J].河南科技大學(xué)學(xué)報(bào)(自然科學(xué)版),2017,38(1): 78-83.
[7]COHEN G C,LIU Q H.Higher-order numerical methods for transient wave equations[J].Springer,2002:349-354.
[8]劉俊利,劉璐菊.具有媒體報(bào)道的傳染病模型穩(wěn)定性[J].河南科技大學(xué)學(xué)報(bào)(自然科學(xué)版),2016,37(2):88-91.[9]SJGREEN B,PETERSSON N A.Perfectly matched layers for Maxwell’s equations in second order formulation[J].Journal of computational physics,2004,209(1):19-46.
[10]ABARBANEL S,GOTTLIEB D,HESTHAVEN J S.Long time behavior of the perfectly matched layer equations in computational electromagnetic[J].Journal of scientific computing,2002,17(1):405-422.
[11]KREISS H O,LORENZ J.Initial-boundary value problems and the Navier-Stokes equations[M].London:Academic Press
Inc,1989:44-62.
O193
A
1672-6871(2017)05-0070-06
10.15926/j.cnki.issn1672-6871.2017.05.015
國(guó)家自然科學(xué)基金項(xiàng)目(11601402);湖北省自然科學(xué)基金項(xiàng)目(2014CFB865)
王小花(1990-),女,河南信陽(yáng)人,碩士生;潘文峰(1964-),男,湖北天門人,教授,博士,碩士生導(dǎo)師,主要研究方向?yàn)槠⒎址匠虜?shù)值解.
2016-12-21