求解帶有不連續(xù)波數(shù)的二維變系數(shù)Helmholtz 方程的一種高精度緊致差分方法

王 芳, 馮秀芳

(寧夏大學(xué)數(shù)學(xué)統(tǒng)計(jì)學(xué)院，銀川 750021)

0 引言

Helmholtz 方程是一類描述時(shí)諧波、航空聲學(xué)、水下聲學(xué)以及電磁波散射等現(xiàn)象的物理模型。求解該方程常用的數(shù)值方法有有限元方法[1]、邊界元方法[2]、有限差分方法[3–7]等。由于Helmholtz 方程本身的復(fù)雜性，研究其數(shù)值計(jì)算比較困難，特別是對(duì)于變系數(shù)、變波數(shù)和大波數(shù)問題仍需深入研究。因此，研究Helmholtz 方程的數(shù)值計(jì)算方法具有重要意義。

近年來(lái)，國(guó)內(nèi)外許多學(xué)者采用有限差分方法在求解Helmholtz 方程方面做了大量的研究。1995 年，Harari 和Turkel[3]建立了均勻網(wǎng)格下的四階格式和非均勻網(wǎng)格下的三階格式；1998 年，Singer 和Turkel[4]構(gòu)造了兩種具有四階精度的格式，一種格式是基于pad′e 近似的推廣，另一種格式是基于Helmholtz 方程本身，計(jì)算高階導(dǎo)數(shù)對(duì)其進(jìn)行修正；2007 年，Nabavi 等人[5]基于方程本身對(duì)高階導(dǎo)數(shù)近似，構(gòu)造了一種新的六階九點(diǎn)緊致差分格式。以上工作研究的都是波數(shù)為常數(shù)的情況。2011 年，Wong 和Li[6]提出了任意波數(shù)的Helmholtz 方程的有限差分格式；2013 年，Turkel 等人[7]構(gòu)造了二維和三維變波數(shù)Helmholtz 方程的六階緊致差分格式；Feng 等人[8–10]利用浸入界面方法，構(gòu)造了常系數(shù)下帶有不連續(xù)波數(shù)的Helmholtz 方程的三階、四階緊致差分格式。

本文針對(duì)帶有不連續(xù)波數(shù)的二維變系數(shù)Helmholtz 方程構(gòu)造了四階緊致差分格式。首先，介紹了二維變系數(shù)Helmholtz 方程的四階緊致差分格式[11–12]。其次，研究波數(shù)k是分段常數(shù)的界面問題，界面問題生活中隨處可見，例如不同的兩種物質(zhì)：水和油；或者相同的物質(zhì)，但處于不同的狀態(tài)：水和冰。浸入界面方法[13–14]是一種處理界面問題的方法，該方法適用于求解區(qū)域是笛卡爾坐標(biāo)系、極坐標(biāo)系或球面坐標(biāo)系下的均勻網(wǎng)格、自適應(yīng)網(wǎng)格或者三角剖分，在遠(yuǎn)離界面的地方使用標(biāo)準(zhǔn)的有限差分方法或者有限元方法；經(jīng)過(guò)界面時(shí)，根據(jù)界面關(guān)系對(duì)有限元或者有限差分格式進(jìn)行局部修改，使其在整個(gè)求解區(qū)域內(nèi)達(dá)到高精度。本文選取在網(wǎng)格線上的垂直界面，在文獻(xiàn)[12]的基礎(chǔ)上構(gòu)造了一種新的四階九點(diǎn)差分格式，該格式在界面處可以達(dá)到四階精度。最后，通過(guò)數(shù)值實(shí)驗(yàn)驗(yàn)證了該格式的有效性和可行性。

1 四階緊致差分格式的建立

考慮下列二維變系數(shù)Helmholtz 方程

其中各系數(shù)為

2 波數(shù)k(x,y)為分段常數(shù)的四階緊致差分格式

考慮(1)式中的波數(shù)k(x,y)為分段常數(shù)的二維變系數(shù)Helmholtz 方程，即

矩形區(qū)域?被劃分為??, ?+兩個(gè)區(qū)域，即??={(x,y)|e

圖1 k 為分段常數(shù)時(shí)，求解區(qū)域示意圖

若格式中網(wǎng)格節(jié)點(diǎn)全部在界面的一側(cè)時(shí)，如圖1 中A、C 部分所示，則使用格式(2)進(jìn)行求解；若格式中網(wǎng)格節(jié)點(diǎn)涉及界面區(qū)域兩側(cè)時(shí)，如圖1 中B 部分所示，根據(jù)浸入界面方法和有限差分方法可推導(dǎo)界面處四階緊致差分格式。

方程在界面x=α上滿足跳躍條件

引理1 假設(shè)(1)式的k(x,y)為分段常數(shù)，則f(x,y)也是分段的，可以得到如下跳躍條件

將(13)式代入(6)式，并整理得到界面處的四階緊致差分格式

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

為了驗(yàn)證本文構(gòu)造的四階緊致差分格式(14)式的精確性和有效性，考慮以下兩個(gè)具有精確解的問題，均采用Dirichlet 邊界條件，其中精度

表1 當(dāng)[k2]=8 與[k2]=24 時(shí)，差分格式的誤差與精度

表2 當(dāng)[k2]=99 與[k2]=399 時(shí)，差分格式的誤差與精度

表3 當(dāng)[k2]=75 與[k2]=375 時(shí)，差分格式的誤差與精度

表4 當(dāng)[k2]=875 與[k2]=300 時(shí)，差分格式的誤差與精度

表5 當(dāng)k2?=1, k2+ =9 與k2+ =25 時(shí)，差分格式的誤差與精度

表6 當(dāng)k2?=1, k2+ =100 與k2+ =400 時(shí)，差分格式的誤差與精度

表7 當(dāng)k2?=25, k2+ =100 與k2+ =400 時(shí)，差分格式的誤差與精度

4 結(jié)論

本文利用有限差分方法和浸入界面方法構(gòu)造了帶有不連續(xù)波數(shù)的二維變系數(shù)Helmholtz 方程的四階緊致差分格式。若格式的網(wǎng)格節(jié)點(diǎn)全部在界面的一側(cè)時(shí)，使用基礎(chǔ)文獻(xiàn)的四階緊致差分格式(2)，若格式的網(wǎng)格節(jié)點(diǎn)涉及界面兩側(cè)使用本文構(gòu)造的四階緊致差分格式，為了使整體能夠達(dá)到四階精度，在正側(cè)點(diǎn)處添加修正項(xiàng)，進(jìn)行修正，將其轉(zhuǎn)化為負(fù)側(cè)點(diǎn)，從而構(gòu)造了一種新的四階九點(diǎn)緊致差分格式(14)。最后，通過(guò)理論分析和數(shù)值驗(yàn)證可以得到該格式是可行的、可靠的。

表8 當(dāng)k2?=25, k2+ =900 與k2?=100, k2+ =400 時(shí)，差分格式的誤差與精度