苏保金 姜子文

( 山东师范大学数学与统计学院,250358,济南 )

1 引 言

本文考虑如下的二维粘性波动方程

Dutt-BΔut-Δu+Aut=f(x,y,t,u),

其中,D,B,A是常系数,f为外力.

粘性波动方程描述了很多物理问题,如粘性介质中声波的传播,微尺度热量传播,流体中颗粒的随机移动等问题.目前,粘性波动方程的数值求解方法已经有不少,最具代表性的是差分方法和有限体积元方法,其目的都是建立适用于该问题的高效的数值求解格式.Dai[1]对于一维热传导方程给出了紧致差分格式.谢建强[2]则对一维粘性波动方程给出了一种三层紧致差分格式.Zhang J和Harfash A J[3,4]则针对类似的高维问题给出了高精度有限差分方法.王同科和高理平[5,6]则对二维问题分别给出了有限体积元和有限元方法.

本文通过Taylor展式得到二维拟线性粘性波动方程关于时间二阶的格式，再通过作用紧算子得到该方程关于空间四阶的紧致差分格式.数值实验验证了方法的精度和有效性.

2 格式构造

本文考虑如下的二维粘性波动方程的初边值问题

(1)

其中,Δ表示Laplace算子,Ω表示二维有界轴平行区域,∂Ω表示区域Ω的边界,T,D,B,A为正常数,u为待求函数,f,u0,u1为已知函数.

引入下列记号

在点(xi,yj,tn)处运用中心差分法离散时间导数,由Taylor展开可知

(2)

其中

在剖分网格上定义两个紧差分算子

则由Taylor展开可知

其中

在等式(2)两端同时作用算子Ax和Ay,且两个算子是可交换的,最终得到

(3)

其中

(4)

3 数值算例

本节将给出具体算例说明格式(4)的有效性.

例1系数A,B,D都取1，空间剖分步长hx=hy=h.精确解u=t10sin(πx)sin(πy),源……