劉 穎，謝偉松**，何振鵬

（1. 天津大學(xué)數(shù)學(xué)學(xué)院，天津 300354；2. 中國(guó)民航大學(xué)航空工程學(xué)院，天津 300300）

0 引 言

現(xiàn)代科學(xué)、技術(shù)和工程中的大量問(wèn)題的數(shù)學(xué)模型都可以用偏微分方程來(lái)表示，大多數(shù)自然科學(xué)問(wèn)題的基本方程本身就是偏微分方程.但是，絕大多數(shù)偏微分方程定解問(wèn)題的解不能以實(shí)用的解析形式來(lái)表示.如果方程的形式復(fù)雜，或者求解區(qū)域的形狀復(fù)雜，其解析解就更難以得到.這就需要尋找偏微分方程定解問(wèn)題的近似解，也就是數(shù)值解.

隨著計(jì)算機(jī)技術(shù)的飛速發(fā)展，微分方程的數(shù)值解法得到了前所未有的發(fā)展和應(yīng)用.有限差分法是求解偏微分方程定解問(wèn)題的主要數(shù)值方法之一［1-3］，此方法已經(jīng)有了比較成熟的算法格式和軟件，并能對(duì)其進(jìn)行比較完善的理論分析.目前，利用有限差分法求解雙曲型方程的研究成果比較豐富. Bou‐jaoui 等［4］研究了一維波動(dòng)方程邊界反饋問(wèn)題的有限差分離散方法，并證明了解的能量是呈指數(shù)衰減的；Boujaoui［5］研究了具有內(nèi)阻尼的一維波動(dòng)方程的有限差分全離散形式，并利用可觀測(cè)不等式證明了阻尼系統(tǒng)能量的衰減率一致依賴于網(wǎng)格長(zhǎng)度；Wang 等［6］采用有限差分法對(duì)時(shí)變偏微分方程尤其是二階波動(dòng)方程進(jìn)行數(shù)值求解，采用具有守恒性質(zhì)的有限差分算子對(duì)內(nèi)點(diǎn)進(jìn)行數(shù)值離散，同時(shí)采用相關(guān)的離散方法對(duì)邊界和網(wǎng)格界面條件進(jìn)行處理，得到了精度較高的數(shù)值實(shí)驗(yàn)結(jié)果.現(xiàn)有文獻(xiàn)多數(shù)都是在非移動(dòng)邊界域［4-20］，也就是固定區(qū)域上進(jìn)行的研究.針對(duì)移動(dòng)邊界域上的研究相對(duì)較少，本文研究?jī)?nèi)容相對(duì)于在非移動(dòng)邊界域上的研究難度加大，是實(shí)際工程問(wèn)題中的一個(gè)難題.本文將采用有限差分方法在移動(dòng)邊界域上對(duì)一類退化波動(dòng)方程的數(shù)值解法進(jìn)行研究，并給出2 種數(shù)值求解該類方程的差分格式.這2 種差分格式對(duì)時(shí)間的離散均采用二階中心差分格式，對(duì)空間的離散分別采用二階差分和四階差分格式，其截?cái)嗾`差分別為O(τ2+h2)和O(τ2+h4).同時(shí)進(jìn)行相應(yīng)的數(shù)值模擬，給出具有一定精度的數(shù)值實(shí)驗(yàn)結(jié)果.

1 問(wèn)題陳述

本文主要利用有限差分法研究移動(dòng)邊界域上一類退化波動(dòng)方程的數(shù)值解法，考慮下述非齊次初邊值問(wèn)題

式中：xp為非線性系數(shù)；0

2 差分格式

為了使用有限差分法求解問(wèn)題（1）～（3），將求解區(qū)域

進(jìn)行剖分.取正整數(shù)m和n，并記

本文后面涉及到第j層節(jié)點(diǎn)都以此為標(biāo)準(zhǔn).

定義Ωτh上的網(wǎng)格函數(shù)

圖1 網(wǎng)格剖分

式中

2.1 二階精度差分格式

在內(nèi)部矩形網(wǎng)格節(jié)點(diǎn)處考慮方程（1），有

代入（4）并忽略高階無(wú)窮小量O(τ2+h2)，得到

由初值條件（3），有

由方程（1），可知

根據(jù)帶積分余項(xiàng)的Taylor 展開式以及式（3），可得

在邊界x=lk(t)上的非矩形網(wǎng)格處采用非等距差分格式（圖2）.如果邊界點(diǎn)左側(cè)的節(jié)點(diǎn)記作(i，j)，右側(cè)節(jié)點(diǎn)記作(i+ 1，j)，點(diǎn)(i+ 1，j)實(shí)際不在求解域上.因?yàn)檫吔缬覀?cè)節(jié)點(diǎn)不在求解域上，對(duì)于求解域外矩形網(wǎng)格點(diǎn)不編號(hào)，將邊界點(diǎn)編號(hào)記作(i+1，j).

圖2 非等距差分格式

忽略高階無(wú)窮小量O(τ2+h2)，得到方程（1）的差分格式：

式中h′如下：

2.2 四階精度差分格式

首先將方程（1）改寫為

為保證四階精度，利用待定系數(shù)法得到uxx和ux的計(jì)算公式

式（14）即為求解退化波動(dòng)方程（1）的時(shí)間二階精度、空間四階精度的差分格式.此格式為一個(gè)7 點(diǎn)格式，對(duì)于空間方向第2 層以及距離斜邊1 個(gè)點(diǎn)和2 個(gè)點(diǎn)處無(wú)法利用此格式，還需要單獨(dú)建立其差分格式.

對(duì)于空間第2 層處差分格式，同樣采用待定系數(shù)法求得uxx和ux的計(jì)算公式

將上式和式（13）代入式（10）并忽略高階無(wú)窮小量O(τ2+h4)，整理得

對(duì)于距離斜邊2 個(gè)點(diǎn)處差分格式，同理可得

對(duì)于距離斜邊1 個(gè)點(diǎn)處的值，則不易采用待定系數(shù)法求得差分格式，本文采用相應(yīng)精度的插值法進(jìn)行計(jì)算.

3 格式檢驗(yàn)

通過(guò)數(shù)值模擬方法分別給出2 種差分格式下離散問(wèn)題的解，并進(jìn)行數(shù)值計(jì)算.文中所涉及到的誤差定義如下：最大絕對(duì)誤差

初值條件u(x，0) =φ(x) = ex，ut(x，0) =ψ(x) = ex，邊值條件u(0，t)=α(t)=et，u(lk(t)，t)=β(t)=elk(t)+t，即

該問(wèn)題的精確解為u(x，t) = ex+t.

令p=0.5，k=0.2，T=3，通過(guò)數(shù)值模擬可以得到問(wèn)題的數(shù)值解，同時(shí)可以得到精確解的圖像（圖3）.

圖3 精確解的圖像

當(dāng)h=0.05，τ=0.025、0.012 5、0.006 25 時(shí)，使用空間二階精度差分格式和空間四階精度差分格式所得誤差結(jié)果如圖4 所示.結(jié)果表明：用這2 種差分格式得到的數(shù)值解與精確解非常接近，絕對(duì)誤差也比較小.因此，該差分格式符合實(shí)際需要，驗(yàn)證了2 種差分格式的精確性.

圖4 不同時(shí)間步長(zhǎng)時(shí)2 種差分格式得到的數(shù)值解與精確解的誤差圖像

當(dāng)t=T= 3，空間步長(zhǎng)為0.05，時(shí)間步長(zhǎng)分別為0.025、0.012 5 和0.006 25 時(shí)，采用空間二階和四階精度差分格式所得數(shù)值解的絕對(duì)誤差對(duì)比曲線如圖5 所示.空間四階精度差分格式的絕對(duì)誤差比空間二階精度差分格式的絕對(duì)誤差更小，說(shuō)明空間四階精度差分格式優(yōu)于空間二階精度差分格式.

圖5 2 種差分格式誤差比較

當(dāng)h= 0.05，τ= 0.025 時(shí)，空間二階和四階精度差分格式所得數(shù)值解的最大絕對(duì)誤差分別為0.041 5和6.930 3×10?4；當(dāng)h=0.05，τ=0.012 5時(shí)，2種差分格式對(duì)應(yīng)的最大絕對(duì)誤差分別為0.025 8 和4.931 3×10?4；當(dāng)h= 0.05，τ= 0.006 25 時(shí)，2 種差分格式對(duì)應(yīng)的最大絕對(duì)誤差分別為0.026 3 和5.168 1×10?4.可知，步長(zhǎng)相同時(shí)，空間四階精度差分格式最大絕對(duì)誤差比二階的更小，并且空間二階精度差分格式的最大絕對(duì)誤差達(dá)到了10?2，而空間四階精度差分格式的最大絕對(duì)誤差能達(dá)到10?4，更加充分地說(shuō)明空間四階精度差分格式比空間二階精度差分格式具有更高的精確性.

4 結(jié)束語(yǔ)

本文利用有限差分法研究了移動(dòng)邊界域上一類退化波動(dòng)方程的數(shù)值解法，針對(duì)此類退化波動(dòng)方程，對(duì)時(shí)間方向的離散采用二階中心差分格式，對(duì)空間方向的離散分別采用二階和四階差分格式.最后通過(guò)數(shù)值模擬驗(yàn)證了有限差分法在解決此類問(wèn)題時(shí)是有效可行的，且空間四階精度差分格式比空間二階精度差分格式具有更高的精確性.