本文将考虑二维分数阶发展型方程的正式的二阶BDF ADI紧致差分格式^{[5, 11, 19-20]}

$ u_t(x, y, t)-I^{(\alpha)}\Delta u(x, y, t)=f(x, y, t), \, \, (x, y, t)\in \Omega\times (0, T], $

(1.1)

其中 $I^{(\alpha)}g(t)$表示函数 $g(t)$的 $\alpha$阶Riemann-Liouville分数阶积分, 定义为^{[14, 22, 26]}

$ I^{(\alpha)}g(t)=\frac{1}{\Gamma(\alpha)}\int_0^{t}(t-s)^{\alpha-1}g(s){\rm d}s, \qquad t>0, \qquad 0< \alpha <1, $

区域 $\Omega=(0, x_R)\times(0, y_R)$, $\partial\Omega$为 $\Omega$的边界, $\Delta=\frac{\partial^2}{\partial x^2}+\frac{\partial^2}{\partial y^2}$为二维拉普拉斯算子.边界条件为

$ u(x, y, t)=\phi (x, y, t), \qquad (x, y, t)\in \partial\Omega\times (0, T], $

(1.2)

初值条件为

$ u(x, y, 0)=u^0(x, y), \qquad (x, y)\in \Omega, $

(1.3)

其中 $f(x, y, t), \phi (x, y, t)$和 $u^0(x, y)$为已知的光滑函数.此模型广泛应用于科学领域, 如黏弹性力的波传播、带有记忆材料的热传导、多孔黏弹性介质的压缩、非牛顿流体学等(参见文献[6-8, 14, 22]).

近年来, 国内外有很多人研究了一维方程(1.1)-(1.3)的数值方法.例如, Mclean, Thomee时间方向使用了向后Euler、Crank-Nicolson和二阶BDF, 空间方向用Galerkin有限元方法, 并给出了齐次方程的正则性估计和各数值方法的稳定性和收敛性^[18]. Kim等使用Galerkin谱和谱配置方法进行时间离散, 得到最优阶误差估计^[10]. Pani等考虑了空间方向的 $H^1$-Galerkin方法, 给出了时间方向的向后Euler、Crank-Nicolson和二阶BDF的误差分析^[21].也有许多研究者考虑了该方程的有限差分方法, Lopez-Marcos^[15], 汤涛^[27]和陈红斌等^[1-3]分别考虑了向后Euler差分格式, 改正的Crank-Nicolson型差分格式, 二阶BDF差分格式, 改正的Crank-Nicolson型和二阶BDF紧致差分格式, 给出了格式的稳定性和收敛性.

对于高维分数阶微分方程的有限差分方法有诸多困难, 其一分数阶导数是非局部的, 具有历史依赖性的特点, 这就意味过去层的解不得不保存以用来计算当前层, 它将需要更多的存储空间.其二对于传统的隐式差分格式, 这个计算复杂性和CPU的使用时间是相当可观的^[29].近年来, 有不少研究者考虑用交替方向隐式差分格式求解高维分数阶微分方程.如二维次扩散方程的ADI差分格式^[30]和二维分数阶扩散-双曲方程的ADI紧致差分格式^[31].二维时间分数阶发展型方程(1.1)-(1.3)的ADI样条配置方法^[19-20]和ADI有限元方法^{[9, 12]}、向后Euler ADI差分格式^[11]、分数阶梯形型ADI差分格式^[4]、二阶BDF ADI差分格式^[5]等. ADI格式是一种既无条件稳定又可用追赶法求解的格式, 它最大的优点是将高维问题分解成一系列一维问题来求解.

本文的主要目的是建立方程(1.1)-(1.3)高效的数值算法, 并推导相应的误差分析.在时间方向上, 利用二阶向后微分公式离散一阶时间导数, 积分项用二阶卷积求积公式近似, 用四阶精度的紧致差分离散二阶空间导数得到全离散紧致差分格式.基于论文^{[23, 31]}的思想, 我们定义了新的能量模, 记为 $\|\cdot\|_A$, 并证明了它与标准的 $L^2$模等价.基于与卷积求积相对应的实二次型的非负性, 我们证明了差分格式在 $\|\cdot\|_A$模下无条件稳定和收敛.其收敛阶为 $O(k^{\alpha+1}+h_1^4+h_2^4)$, 其中 $k$为时间步长, $h_1$和 $h_2$分别是空间 $x$和 $y$方向的步长.在时间方向上相比于向后Euler ADI差分格式^[11], 在空间上相比于二阶BDF ADI差分格式^[5], 我们讨论的数值方法在收敛阶上都有明显提高.

本文结构如下:第2节我们将介绍一些记号和有用的引理.第3节我们将构建正式的二阶BDF ADI紧致差分格式.第4节我们研究差分格式的稳定性和收敛性.和理论相吻合的数值算例放在第5节.最后一节给出本文的简要总结.

给定正整数 $N$, 时间步长 $k=T/N$, 那么 $t_n=nk(k=0, 1, \cdots, N)$.时间区间 $[0, T]$剖分为 $\Omega_n=\{t_n| 0\leq n\leq N\}$.对于给定 $\Omega_n$上的网格函数 $w=\{w^n|0\leq n\leq N\}$, 记

$ \delta_t w^n=\frac{1}{k}(w^n-w^{n-1}), \qquad D_t^{(2)}w^n=\frac{3}{2}\delta_tw^n-\frac{1}{2}\delta_tw^{n-1}. $

为了空间离散, 对正整数 $M_1, M_2$, 记 $h_1=x_R/M_1, h_2=y_R/M_2$, $h=\max\{h_1, h_2\}$, $x_i=ih_1 (i=0, 1, \cdots, M_1)$和 $y_j=jh_2 (j=0, 1, \cdots, M_2)$.记 $\overline{\Omega}_h=\{(x_i, y_j)|0\leq i \leq M_1, 0\leq j \leq M_2\}$, $\Omega_h=\overline{\Omega}_h\bigcap\Omega$, $\partial\Omega_h=\Omega_h\bigcap\partial\Omega $.对于 $\overline{\Omega}_h$上的任意网格函数 $v=\{v_{ij}|0\leq i\leq M_1, 0\leq j\leq M_2\}$, 引进如下记号

$ \delta_x v_{i+\frac{1}{2}, j}=\frac{1}{h_1}(v_{i+1, j}-v_{ij}), \qquad \delta_x^2 v_{ij}=\frac{1}{h_1}(\delta_x v_{i+\frac{1}{2}, j}-\delta_x v_{i-\frac{1}{2}, j}), \\ \delta_y\delta_x v_{i+\frac{1}{2}, j+\frac{1}{2}}=\frac{1}{h_2}(\delta_x v_{i+\frac{1}{2}, j+1}-\delta_x v_{i+\frac{1}{2}, j}), \qquad \delta_y\delta_x^2 v_{i, j+\frac{1}{2}}=\frac{1}{h_2}(\delta_x^2 v_{i, j+1}-\delta_x^2 v_{i, j}). $

按如上类似定义记号 $\delta_y v_{i, j+\frac{1}{2}}, $ $\delta_y^2 v_{ij}, \delta_x\delta_y v_{i+\frac{1}{2}, j+\frac{1}{2}}, $ $\delta_x\delta_y^2 v_{i+\frac{1}{2}, j}, \delta_x^2\delta_y^2 v_{ij}$, 离散的拉普拉斯算子定义为 $\Delta_h v_{ij}=\delta_x^2 v_{ij}+\delta_y^2 v_{ij}$.

对于任意的网格函数 $v=\{v_{ij}|0\leq i\leq M_1, 0\leq j\leq M_2\}$和 $w=\{w_{ij}|0\leq i\leq M_1, $ $0\leq j\leq M_2\}$, 记

$ \langle v, w \rangle=\sum\limits_{i=1}^{M_1-1}\sum\limits_{j=1}^{M_2-1}h_1h_2v_{ij}w_{ij}, \qquad \|v\|=\sqrt{\langle v, v \rangle}, \\ \|\delta_x v\|=\sqrt{\sum\limits_{i=0}^{M_1-1}\sum\limits_{j=1}^{M_2-1}h_1h_2(\delta_x v_{i+\frac{1}{2}, j})^2}, \quad \|\delta_y v\|=\sqrt{\sum\limits_{i=1}^{M_1-1}\sum\limits_{j=0}^{M_2-1}h_1h_2(\delta_y v_{i, j+\frac{1}{2}})^2}, \\ \|\delta_x\delta_y v\|=\sqrt{\sum\limits_{i=0}^{M_1-1}\sum\limits_{j=0}^{M_2-1}h_1h_2(\delta_x\delta_y v_{i+\frac{1}{2}, j+\frac{1}{2}})^2}, \quad \|\delta_x\delta_y^2 v\|=\sqrt{\sum\limits_{i=0}^{M_1-1}\sum\limits_{j=1}^{M_2-1}h_1h_2(\delta_x\delta_y^2 v_{i+\frac{1}{2}, j})^2}.\\ \|\Delta_hv\|=\sqrt{\sum\limits_{i=1}^{M_1-1}\sum\limits_{j=1}^{M_2-1}h_1h_2|\Delta_hv_{ij}|^2}. $

设 $V_h=\{v|v=\{v_{ij}|(x_i, y_j)\in\overline{\Omega}_h\}, \, \, $如 $\, \, (x_i, y_j)\in\partial{\Omega}_h\, \, $则 $ \, \, v_{ij}=0\}$.对任意 $\Omega_h$上的网格函数 $v=\{v_{ij}|0\leq i\leq M_1, 0\leq j\leq M_2\}$, 引进如下平均算子

$ {\cal H}_xv_{ij}= \left\{\begin{array}{ll} \frac{1}{12}(v_{i-1, j}+10v_{ij}+v_{i+1, j}), ~~& 1\leq i\leq M_1-1, \, \, 1\leq j\leq M_2, \\[2mm] v_{ij},&i=0 \, \, \mbox{或} \, \, M_1, \, \, 1\leq j\leq M_2. \end{array} \right.\\ {\cal H}_yv_{ij}= \left\{\begin{array}{ll} \frac{1}{12}(v_{i, j-1}+10v_{ij}+v_{i, j+1}), ~~ & 1\leq j\leq M_2-1, \, \, 1\leq i\leq M_1, \\[2mm] v_{ij}, &j=0\, \, \mbox{或} \, \, M_2, \, \, 1\leq i\leq M_1. \end{array} \right. $

显然

$ {\cal H}_xv_{ij}=(I+\frac{h_1^2}{12}\delta_x^2)v_{ij}, \qquad {\cal H}_yv_{ij}=(I+\frac{h_2^2}{12}\delta_y^2)v_{ij}, \qquad (x_i, y_j)\in \Omega_h, $

其中 $I$表示恒等算子.为了方便叙述, 我们介绍下面的记号

$ {\cal H}v_{ij}={\cal H}_x{\cal H}_yv_{ij}, \qquad \Lambda_hv_{ij}={\cal H}\Delta_hv_{ij}=({\cal H}_y\delta_x^2+{\cal H}_x\delta_y^2)v_{ij}, \qquad (x_i, y_j)\in \Omega_h. $

为了我们后面主要结果的分析, 我们给出如下有用的辅助引理.

引理2.1^[13]  对任意的网格函数 $w, v\in V_h$, 有

$ -\sum\limits_{i=1}^{M_1-1}\sum\limits_{j=1}^{M_2-1}h_1h_2(\delta_x^2v_{ij})w_{ij} =\sum\limits_{i=0}^{M_1-1}\sum\limits_{j=1}^{M_2-1}h_1h_2(\delta_xw_{i+\frac{1}{2}, j})(\delta_xv_{i+\frac{1}{2}, j}), \\ -\sum\limits_{i=1}^{M_1-1}\sum\limits_{j=1}^{M_2-1}h_1h_2(\delta_y^2v_{ij})w_{ij} =\sum\limits_{i=1}^{M_1-1}\sum\limits_{j=0}^{M_2-1}h_1h_2(\delta_yw_{i, j+\frac{1}{2}})(\delta_yv_{i, j+\frac{1}{2}}). $

引理2.2^[25-26]  对任意的网格函数 $w, v\in V_h$, 有

$ \sum\limits_{i=1}^{M_1-1}\sum\limits_{j=1}^{M_2-1}h_1h_2(\delta_x^2\delta_y^2w_{ij})v_{ij} =\sum\limits_{i=0}^{M_1-1}\sum\limits_{j=0}^{M_2-1}h_1h_2(\delta_x\delta_yw_{i+\frac{1}{2}, j+\frac{1}{2}})(\delta_x\delta_yv_{i+\frac{1}{2}, j+\frac{1}{2}}). $

引理2.3^[2-3]  对任意的网格函数 $w, v\in V_h$, 有

$ |(\delta_x^2v^n, v^m)|\leq \frac{4}{h_1^2}\|v^n\|\|v^m\|, \qquad |(\delta_y^2v^n, v^m)|\leq \frac{4}{h_2^2}\|v^n\|\|v^m\|. $

引理2.4^{[29, 32]}  设 $g(x)\in C^6[x_{j-1}, x_{j+1}]$, $\zeta(s)=5(1-s)^3-3(1-s)^5$, 则

$ \begin{eqnarray*} && \frac{1}{12}[g''(x_{i+1})+10g''(x_i)+g''(x_{i-1})]\\ &=&\frac{1}{h_1^2}[g(x_{i+1})-2g(x_i)+g(x_{i-1})] +\frac{h_1^4}{360}\int_0^1[g^{(6)}(x_i-sh)+g^{(6)}(x_i+sh)]\zeta(s){\rm d}s. \end{eqnarray*} $

本章中 $C$表示一个正常数, 在不同的地方可能不相同, 它们仅仅依赖于方程的解和给定的边界函数, 与时间步长 $k$和空间步长 $h_1, h_2$无关.

本节将构建方程(1.1)-(1.3)的正式的二阶BDF ADI紧致差分格式.

定义网格函数

$ U_{ij}^n=u(x_i, y_j, t_n), \qquad f_{ij}^n=f(x_i, y_j, t_n), \qquad (x_i, y_j)\in \Omega_h, \qquad 1\leq n\leq N. $

首先, 我们介绍下列二阶卷积求积公式(参见文献[16])

$ Q_n^{(\alpha)}(\varphi)=k^{\alpha}\sum\limits_{p=1}^nq_{n-p}^{(\alpha)}\varphi^p+k^{\alpha}\tilde{q}_n^{(\alpha)}\varphi^0, $

(3.1)

其中求积权系数 $q_p^{(\alpha)}$由它的发生函数所确定

$ K(\delta(z))=[\delta(z)]^{-\alpha}=\sum\limits_{p=0}^{\infty}q_p^{(\alpha)}z^p. $

(3.2)

这里 $K(s)$表示卷积积分核的Laplace变换, $\delta(z)$是线性多步法发生多项式的商所确定的一个有理函数.具体地, 对于二阶向后微分公式, 那么

$ \delta(z)=(1-z)+\frac{1}{2}(1-z)^2=\frac{1}{2}(3-4z+z^2). $

(3.3)

为了使积分真正达到二阶, 取校正的权系数 $\tilde{q}_n^{(\alpha)}$使求积公式对多项式 $\varphi=1$准确成立, 也就是说

$ k^{\alpha}\sum\limits_{p=1}^nq_{n-p}^{(\alpha)}+k^{\alpha}\tilde{q}_n^{(\alpha)}=\frac{1}{\Gamma(\alpha)}\int_0^{t_n}(t_n-s)^{\alpha-1}{\rm d}s =\frac{1}{\Gamma(\alpha+1)} t_n^{\alpha}, $

即

$ \sum\limits_{p=1}^nq_{n-p}^{(\alpha)}+\tilde{q}_n^{(\alpha)}=\frac{n^{\alpha}}{\Gamma(\alpha+1)}. $

(3.4)

下面给出求积误差 $E(\varphi)(t_n)=I^{(\alpha)}\varphi (t_n)-Q_n^{(\alpha)}(\varphi)$, 其中 $Q_n^{(\alpha)}(\varphi)$由(3.1)式定义.

引理3.1^{[3, 6, 16]}  设 $\varphi$是在 $0< t \leq T$上实的, 连续可微的函数, 且 $\varphi_{tt}$在 $ 0<t<T $上连续可积.那么基于二阶向后微分公式(3.3)的卷积求积近似(3.1)和(3.2)式的误差界为

$ |E(\varphi)(t_n)|=|I^{(\alpha)}\varphi (t_n)-Q_n^{(\alpha)}(\varphi)|\leq Ck^{\alpha+1}, $

对于固定的 $T< \infty$, 常数 $C$不依赖与 $k$和 $t_n\in (0, T]$.

利用带积分余项的Taylor级数展开, 我们有下列引理.

引理3.2^{[13, 30]}  设 $u(x, y, t)\in C_{x, y}^{4, 4}([0, x_R]\times [0, y_R])$, 那么成立

$ \begin{eqnarray*} && \frac{\partial^2u}{\partial x^2}(x_i, y_j, t_n)\\ &=&\delta_x^2U_{ij}^n -\frac{h_1^2}{6}\int_0^1\bigg[\frac{\partial^4u}{\partial x^4}(x_i+sh_1, y_j, t_n)+\frac{\partial^4u}{\partial x^4}(x_i-sh_1, y_j, t_n)\bigg](1-s)^3{\rm d}s, \\[3mm] && \frac{\partial^2u}{\partial y^2}(x_i, y_j, t_n)\\ &=&\delta_y^2U_{ij}^n -\frac{h_2^2}{6}\int_0^1\bigg[\frac{\partial^4u}{\partial y^4}(x_i, y_j+sh_2, t_n)+\frac{\partial^4u}{\partial y^4}(x_i, y_j-sh_2, t_n)\bigg](1-s)^3{\rm d}s. \end{eqnarray*} $

下面将给出 $(R_1)_{ij}^n=I^{(\alpha)}\Delta u(x_i, y_j, t_n)-Q_n^{(\alpha)}(\Delta_h u_{ij})$的界.

引理3.3^[5]  设 $u(x, y, t)\in C_{x, y, t}^{4, 4, 2}$ $([0, x_R]\times [0, y_R]\times (0, T])$, 那么成立

$ |(R_1)_{ij}^n|=|I^{(\alpha)}\Delta u(x_i, y_j, t_n)-Q_n^{(\alpha)}(\Delta_hU_{ij}^n)| \leq C(k^{\alpha+1}+h_1^2+h_2^2), \quad 1\leq n \leq N. $

引理3.4^[3]  设 $ g(t) \in C^3[t_{n-2}, t_n]$.对 $n\geq 2$, 成立

$ \begin{eqnarray*} && g'(t_n)-\frac{1}{2k}[3g(t_n)-4g(t_{n-1})+g(t_{n-2})]\\ &=&-k^2\int_0^1[g^{(3)}(t_n-sk)-2g^{(3)}(t_n-2sk)](1-s)^2{\rm d}s. \end{eqnarray*} $

下面将构建方程(1.1)-(1.3)的正式的二阶BDF ADI紧致差分格式.

在结点 $(x_i, y_j, t_n)$处考虑方程(1.1), 有

$ u_t(x_i, y_j, t_n)-I^{(\alpha)}\Delta u(x_i, y_j, t_n)=f(x_i, y_j, t_n), \, \, (x_i, y_j)\in\Omega_h, \, \, 1\leq n\leq N. $

(3.5)

由引理3.3和引理3.4, 有

$ I^{(\alpha)}\Delta u(x_i, y_j, t_n)=Q_n^{(\alpha)}(\Delta_hU_{ij})+(R_1)_{ij}^n, \, \, (x_i, y_j)\in \Omega_h, \, \, 1\leq n\leq N, $

(3.6)

$ u_t(x_i, y_j, t_1)=\delta_tU_{ij}^1+(R_2)_{ij}^1, \quad (x_i, y_j)\in \Omega_h, $

(3.7)

$ u_t(x_i, y_j, t_n)=D_t^{(2)}U_{ij}^n+(R_2)_{ij}^n, \quad (x_i, y_j)\in \Omega_h, \quad 2\leq n\leq N, $

(3.8)

其中

$ (R_2)_{ij}^1=-k\int_0^1\frac{\partial^2U}{\partial t^2}(x_i, y_j, t_0+sk)(1-s){\rm d}s, \\ (R_2)_{ij}^n=-k^2\int_0^1 \bigg[\frac{\partial^3U}{\partial t^3}(x_i, y_j, t_n-sk)-2\frac{\partial^3U} {\partial t^3}(x_i, y_j, t_n-2sk)\bigg](1-s)^2{\rm d}s, \, \, 2\leq n\leq N. $

将(3.6)-(3.8)式代入(3.5)式, 得

$ \delta_tU_{ij}^1-k^{\alpha}(q_0^{(\alpha)}\Delta_hU_{ij}^1+\tilde{q}_1^{(\alpha)}\Delta_hU_{ij}^0)=f_{ij}^1+(R_1)_{ij}^1-(R_2)_{ij}^1, \, \, (x_i, y_j)\in \Omega_h. $

(3.9)

$ D_t^{(2)}U_{ij}^n-k^{\alpha} \bigg(\sum\limits_{p=1}^nq_{n-p}^{(\alpha)}\Delta_hU_{ij}^p+\tilde{q}_n^{(\alpha)}\Delta_hU_{ij}^0\bigg) \\ = f_{ij}^n+(R_1)_{ij}^n-(R_2)_{ij}^n, \quad (x_i, y_j)\in \Omega_h, \quad 2\leq n\leq N. $

(3.10)

利用带积分余项的Taylor展开公式, 我们有(参见文献[13])

$ \begin{equation} 0=\delta_x^2\delta_y^2\delta_tU_{ij}^n-(R_3)_{ij}^n, \quad (x_i, y_j)\in \Omega_h, \quad 1\leq n\leq N, \end{equation} $

(3.11)

其中

$ (R_3)_{ij}^n=\frac{1}{2}\int_0^1\delta_x^2\delta_y^2\bigg[\frac{\partial U}{\partial t}(x_i, y_j, t_n-\frac{sk}{2}) +\frac{\partial U}{\partial t}(x_i, y_j, t_n+\frac{sk}{2})\bigg]{\rm d}s. $

在(3.9)和(3.10)式的左端分别增加小量项 $(k^{\alpha+1}q_0^{(\alpha)})^2\delta_x^2\delta_y^2\delta_tU_{ij}^1$和 $\frac{2}{3}(k^{\alpha+1}q_0^{(\alpha)})^2\delta_x^2\delta_y^2\delta_tU_{ij}^n$, 得

$ \begin{eqnarray} &&\delta_tU_{ij}^1-k^{\alpha}(q_0^{(\alpha)}\Delta_hU_{ij}^1+\tilde{q}_1^{(\alpha)}\Delta_hU_{ij}^0)+(k^{\alpha+1}q_0^{(\alpha)})^2\delta_x^2\delta_y^2\delta_tU_{ij}^1 \\ & =&f_{ij}^1+R_{ij}^1, \quad (x_i, y_j)\in \Omega_h, \end{eqnarray} $

(3.12)

$ \begin{eqnarray} && D_t^{(2)}U_{ij}^n-k^{\alpha} \bigg(\sum\limits_{p=1}^nq_{n-p}^{(\alpha)}\Delta_hU_{ij}^p+\tilde{q}_n^{(\alpha)}\Delta_hU_{ij}^0\bigg)+ \frac{2}{3}(k^{\alpha+1}q_0^{(\alpha)})^2\delta_x^2\delta_y^2\delta_tU_{ij}^n \\ &=&f_{ij}^n+R_{ij}^n, \quad (x_i, y_j)\in \Omega_h, \quad 2\leq n\leq N, \end{eqnarray} $

(3.13)

其中

$ R_{ij}^n=(R_1)_{ij}^n-(R_2)_{ij}^n+(R_3)_{ij}^n, \quad (x_i, y_j)\in \Omega_h, \qquad 1\leq n\leq N. $

用平均算子 ${\cal H}={\cal H}_x{\cal H}_y$作用在(3.12)-(3.13)式上, 可得

$ \begin{eqnarray} && {\cal H}\delta_tU_{ij}^1-k^{\alpha}(q_0^{(\alpha)}\Lambda_hU_{ij}^1+\tilde{q}_1^{(\alpha)}\Lambda_hU_{ij}^0)+(k^{\alpha+1}q_0^{(\alpha)})^2\delta_x^2\delta_y^2\delta_tU_{ij}^1 \\ & =&{\cal H}f_{ij}^1+{\cal H}R_{ij}^1, \quad (x_i, y_j)\in \Omega_h, \end{eqnarray} $

(3.14)

$ \begin{eqnarray} && {\cal H}D_t^{(2)}U_{ij}^n-k^{\alpha} \bigg(\sum\limits_{p=1}^nq_{n-p}^{(\alpha)}\Lambda_hU_{ij}^p+\tilde{q}_n^{(\alpha)}\Lambda_hU_{ij}^0\bigg)+\frac{2}{3}(k^{\alpha+1}q_0^{(\alpha)})^2\delta_x^2\delta_y^2\delta_tU_{ij}^n \\ &=&{\cal H}f_{ij}^n+{\cal H}R_{ij}^n, \quad (x_i, y_j)\in \Omega_h, \quad 2\leq n\leq N, \end{eqnarray} $

(3.15)

其中

$ {\cal H}R_{ij}^n={\cal H}(R_1)_{ij}^n-{\cal H}(R_2)_{ij}^n+{\cal H}(R_3)_{ij}^n, \qquad (x_i, y_j)\in \Omega_h, \qquad 1\leq n\leq N. $

因此, 利用引理2.4, 存在一个不依赖与 $h_1, h_2$和 $k$的常数 $C$满足

$ \begin{equation} |{\cal H}R_{ij}^n|\leq C(k^{\alpha+1}+h_1^4+h_2^4), \qquad (x_i, y_j)\in \Omega_h, \qquad 1\leq n\leq N. \end{equation} $

(3.16)

另外, 成立如下的初值和边值条件

$ \begin{equation} U_{ij}^n=\phi(x_i, y_j, t_n), \qquad (x_i, y_j)\in\partial\Omega_h, \qquad 1\leq n\leq N. \end{equation} $

(3.17)

$ \begin{equation} U_{ij}^0=u^0(x_i, y_j), \qquad (x_i, y_j)\in\Omega_h. \end{equation} $

(3.18)

略去(3.14)和(3.15)式右端的微小项 ${\cal H}R_{ij}^n$, 并用 $u_{ij}^n$替代它的近似解 $U_{ij}^n$, 我们得到如下正式的二阶BDF ADI紧致差分格式

$ \begin{equation} {\cal H}\delta_tu_{ij}^1-k^{\alpha}(q_0^{(\alpha)}\Lambda_hu_{ij}^1+\tilde{q}_1^{(\alpha)}\Lambda_hu_{ij}^0)+(k^{\alpha+1}q_0^{(\alpha)})^2\delta_x^2\delta_y^2\delta_tu_{ij}^1\\ ={\cal H}f_{ij}^1, \quad (x_i, y_j)\in \Omega_h. \end{equation} $

(3.19)

$ \begin{equation} \begin{array}[b]{rl} & {\cal H}D_t^{(2)}u_{ij}^n-k^{\alpha} \bigg(\sum\limits_{p=1}^nq_{n-p}^{(\alpha)}\Lambda_hu_{ij}^p+\tilde{q}_n^{(\alpha)}\Lambda_hu_{ij}^0 \bigg) +\frac{2}{3}(k^{\alpha+1}q_0^{(\alpha)})^2\delta_x^2\delta_y^2\delta_tu_{ij}^n \\[2mm] =& {\cal H}f_{ij}^n, \quad (x_i, y_j)\in \Omega_h, \quad 2\leq n\leq N. \end{array} \end{equation} $

(3.20)

$ u_{ij}^n=\phi(x_i, y_j, t_n), \qquad (x_i, y_j)\in \partial\Omega_h, \qquad 1\leq n\leq N. $

(3.21)

$ u_{ij}^0=u^0(x_i, y_j), \qquad (x_i, y_j)\in\Omega_h. $

(3.22)

将(3.19)-(3.20)式写成熟悉的ADI形式. (3.19)-(3.20)式分别等价于下列表达式

$ {\cal H}u_{ij}^1-k^{\alpha+1}q_0^{(\alpha)}\Lambda_hu_{ij}^1+(k^{\alpha+1}q_0^{(\alpha)})^2\delta_x^2\delta_y^2u_{ij}^1 ={\cal H}F_{ij}^1, \quad (x_i, y_j)\in\Omega_h, $

(3.23)

$ {\cal H}u_{ij}^n-\frac{2}{3}k^{\alpha+1}q_0^{(\alpha)}\Lambda_hu_{ij}^n+(\frac{2}{3}k^{\alpha+1}q_0^{(\alpha)})^2\delta_x^2\delta_y^2u_{ij}^n ={\cal H}F_{ij}^n,\\ \quad (x_i, y_j)\in\Omega_h, \ 2\leq n\leq N, $

(3.24)

其中

$ {\cal H}F_{ij}^1={\cal H}u_{ij}^0+k^{\alpha+1}\tilde{q}_1^{(\alpha)}\Lambda_hu_{ij}^0+(k^{\alpha+1}q_0^{(\alpha)})^2\delta_x^2\delta_y^2u_{ij}^0+k{\cal H}f_{ij}^1, \\ {\cal H}F_{ij}^n={\cal H}u_{ij}^{n-1}+\frac{1}{3}k{\cal H}\delta_tu_{ij}^{n-1}+\frac{2}{3}k^{\alpha+1} \bigg[\sum\limits_{p=1}^{n-1}q_{n-p}^{(\alpha)}\Lambda_hu_{ij}^p+\tilde{q}_n^{(\alpha)}\Lambda_hu_{ij}^0\bigg]\\ +(\frac{2}{3}k^{\alpha+1}q_0^{(\alpha)})^2\delta_x^2\delta_y^2u_{ij}^{n-1}+\frac{2}{3}k{\cal H}f_{ij}^n, \qquad 2\leq n\leq N. $

方程(3.23)-(3.24)可写成

$ ({\cal H}_x-k^{\alpha+1}q_0^{(\alpha)}\delta_x^2)({\cal H}_y-k^{\alpha+1}q_0^{(\alpha)}\delta_y^2)u_{ij}^1={\cal H}F_{ij}^1, \quad (x_i, y_j)\in \Omega_h.\\ ({\cal H}_x-\frac{2}{3}k^{\alpha+1}q_0^{(\alpha)}\delta_x^2)({\cal H}_y-\frac{2}{3}k^{\alpha+1}q_0^{(\alpha)}\delta_y^2)u_{ij}^n={\cal H}F_{ij}^n, \quad (x_i, y_j)\in \Omega_h, \quad 2\leq n\leq N. $

通过求解两个独立的一维问题得到 $\{u_{ij}^n\}$.引入中间变量

$ \bar{u}_{ij}^1=({\cal H}_y-k^{\alpha+1}q_0^{(\alpha)}\delta_y^2)u_{ij}^1, \quad (x_i, y_j)\in \Omega_h, \\ \bar{u}_{ij}^n=({\cal H}_y-\frac{2}{3}k^{\alpha+1}q_0^{(\alpha)}\delta_y^2)u_{ij}^n, \quad (x_i, y_j)\in \Omega_h, \qquad 2\leq n\leq N. $

对于固定的 $j\in \{1, 2, \cdots, M_2-1\}$, 求解下列方程组得到 $\{\bar{u}_{ij}^n\}$

$ ({\cal H}_x-k^{\alpha+1}q_0^{(\alpha)}\delta_x^2)\bar{u}_{ij}^1= {\cal H}F_{ij}^1, \quad (x_i, y_j)\in \Omega_h, \\ ({\cal H}_x-\frac{2}{3}k^{\alpha+1}q_0^{(\alpha)}\delta_x^2)\bar{u}_{ij}^n={\cal H}F_{ij}^n, \quad (x_i, y_j)\in \Omega_h, \quad n\geq 2, \\ \bar{u}_{0j}^n=({\cal H}_y-k^{\alpha+1}q_0^{(\alpha)}\delta_y^2)u_{0j}^n, \quad \bar{u}_{Mj}^n=({\cal H}_y-k^{\alpha+1}q_0^{(\alpha)}\delta_y^2)u_{Mj}^n, \quad 1\leq n\leq N. $

一旦 $\{\bar{u}_{ij}^n\}$被确定, 固定 $i\in \{1, 2, \cdots, M_1-1\}$, 求解下列方程组

$ ({\cal H}_y-k^{\alpha+1}q_0^{(\alpha)}\delta_y^2)u_{ij}^1=\bar{u}_{ij}^1, \qquad (x_i, y_j)\in \Omega_h.\\ ({\cal H}_y-\frac{2}{3}k^{\alpha+1}q_0^{(\alpha)}\delta_y^2)u_{ij}^n=\bar{u}_{ij}^n, \qquad (x_i, y_j)\in \Omega_h, \qquad n\geq 2.\\ u_{i0}^n=\phi(x_i, y_0, t_n), \qquad u_{iM}^n=\phi(x_i, y_M, t_n), \qquad 1\leq n\leq N, $

得到我们想要的解 $\{u_{ij}^n\}$.

为了以后研究稳定性和收敛性, 我们将介绍新的模.

定义 $V_h=\{v|v=\{v_{ij}|(x_i, y_j)\in\overline{\Omega}_h\}, $如 $ (x_i, y_j)\in\partial{\Omega}_h$则 $ v_{ij}=0\}$是网格函数 $\overline{\Omega}_h$.对任意的网格函数 $v, w \in V_h$, 记

$ \langle v, w \rangle_A=\langle{\cal H}v, w \rangle, \\ \|v\|_A^2=\langle v, v \rangle_A =\|v\|^2-\frac{h_1^2}{12}\|\delta_xv\|^2-\frac{h_2^2}{12}\|\delta_yv\|^2+\frac{h_1^2h_2^2}{144}\|\delta_x\delta_yv\|^2. $

下面我们将证明 $\|\cdot\|_A$与标准的 $L^2$模等价.

引理4.1  对 $\forall \, \, v\in V_h$, 有

$ \frac{1}{3}\|v\|^2 \leq \|v\|_A^2 \leq \frac{10}{9}\|v\|^2. $

证  由 $\|\delta_xv\|$, $\|\delta_yv\|$和 $\|\delta_x\delta_yv\|$的定义 $\|\delta_xv\|$, 有

$ \begin{eqnarray*} h_1^2\|\delta_xv\|^2&=&h_1^2\sum\limits_{i=0}^{M_1-1}\sum\limits_{j=1}^{M_2-1}h_1h_2(\delta_xv_{i+\frac{1}{2}, j})^2 \\ & \leq& 2\sum\limits_{i=0}^{M_1-1}\sum\limits_{j=1}^{M_2-1}h_1h_2[(v_{i+1, j})^2+(v_{ij})^2]\\ &\leq& 4\|v\|^2, \\[3mm] h_2^2\|\delta_yv\|^2&=&h_2^2\sum\limits_{i=1}^{M_1-1}\sum\limits_{j=0}^{M_2-1}h_1h_2(\delta_yv_{i, j+\frac{1}{2}})^2 \\ &\leq& 2\sum\limits_{i=1}^{M_1-1}\sum\limits_{j=0}^{M_2-1}h_1h_2[(v_{i, j+1})^2+(v_{ij})^2]\\ &\leq& 4\|v\|^2 \end{eqnarray*} $

和

$ h_1^2h_2^2\|\delta_x\delta_yv\|^2\leq 4h_1^2\|\delta_xv\|^2\leq 16\|v\|^2. $

结合 $\langle v, v \rangle_A$内积的定义, 有

$ \|v\|_A^2=\|v\|^2-\frac{h_1^2}{12}\|\delta_xv\|^2-\frac{h_2^2}{12}\|\delta_yv\|^2+\frac{h_1^2h_2^2}{144}\|\delta_x\delta_yv\|^2. $

证毕.

下面给出卷积结构下实二次型的非负性质的一般结果.

引理4.2^{[3, 16]}  若 $\{a_0, a_1, \cdots, a_n, \cdots\}$为在 $D=\{z\in {\cal C}:|z|\leq 1\}$上 $\hat{a}(z)=\sum\limits_{n=0}^{\infty}a_nz^n$解释的实值序列, 那么对任意的正整数 $N$和任意的 $(V^1, V^2, \cdots, V^N)\in {\Bbb R}^N$, 有

$ \sum\limits_{n=1}^N\bigg(\sum\limits_{p=1}^na_{n-p}V^p\bigg)V^n\geq 0, $

当且仅当

$ {\cal \Re}\hat{a}(z)\geq 0, \qquad \mbox{对} \qquad z\in D, $

其中 ${\cal \Re}(z)$表示复数 $z$的实部.

因为积分核 $\beta(t)=\frac{1}{\Gamma{(\alpha)}}t^{\alpha-1}$是正类型, 也就是说, 当 ${\cal \Re}(s)\geq 0$, ${\cal \Re}\tilde{\beta}(s)\geq 0$, 且当 $|z|<1$, ${\cal \Re}[\frac{1}{2}(3-4z+z^2)]>0$ (参见文献[17, 引理4.1]).因此对 $|z|<1$有 ${\cal \Re}\tilde{\beta}(\frac{(1-z)(3-z)}{2})=$ ${\cal \Re}[\frac{1}{2}(3-4z+z^2)]>0$.故发生函数(3.2)满足引理4.2的条件.

通过能量方法, 我们将推导正式的二阶BDF ADI紧致差分格式(3.19)-(3.22)的稳定性.

定理4.3  设 $\{u_{ij}^n|(x_i, y_j)\in \Omega_h, 1\leq n\leq N\}$是正式的二阶BDF ADI紧致差分格式(3.19)-(3.22)的解, 在 $\partial\Omega_h$上 $u_{ij}^n=0$, 那么成立

$ \|u^n\|_A \leq C(T)\|u^0\|_A+2(k^{\alpha}q_0^{(\alpha)})^2\|\delta_x^2\delta_y^2u^0\|+2k\sum\limits_{l=1}^n\|{\cal H}f^l\|. $

证  当 $n=1$时, 在(3.19)式两端乘以 $h_1h_2u_{ij}^1$并对 $i, j$在 $(x_i, y_j)\in\Omega_h$上求和, 可得

$ \begin{eqnarray} && \langle{\cal H}\delta_tu^1, u^1\rangle-k^{\alpha}(q_0^{(\alpha)}\langle\Lambda_hu^1, u^1\rangle+\tilde{q}_1^{(\alpha)}\langle\Lambda_hu^0, u^1\rangle) \\ &&+(k^{\alpha}q_0^{(\alpha)})^2\langle\delta_x^2\delta_y^2\delta_tu^1, u^1\rangle =\langle{\cal H}f^1, u^1\rangle. \end{eqnarray} $

(4.1)

当 $n\geq 2$时, 在(3.20)式两端乘以 $h_1h_2u_{ij}^n$并对 $i, j$在 $(x_i, y_j)\in\Omega_h$求和, 可得

$ \begin{eqnarray} && \langle{\cal H}D_t^{(2)}u^n, u^n\rangle-k^{\alpha} \bigg(\sum\limits_{p=1}^nq_{n-p}^{(\alpha)}\langle\Lambda_hu^p, u^n\rangle+\tilde{q}_n^{(\alpha)}\langle\Lambda_hu^0, u^n\rangle\bigg) \\ &&+\frac{2}{3}(k^{\alpha}q_0^{(\alpha)})^2\langle\delta_x^2\delta_y^2\delta_tu^n, u^n\rangle =\langle{\cal H}f^n, u^n\rangle. \end{eqnarray} $

(4.2)

当 $n\geq 2$时, 对(4.1)和(4.2)式求和, 有

$ \begin{eqnarray} && k\langle{\cal H}\delta_tu^1, u^1\rangle+k\sum\limits_{n=2}^N\langle{\cal H}D_t^{(2)}u^n, u^n\rangle +k\langle(k^{\alpha}q_0^{(\alpha)})^2\delta_x^2\delta_y^2\delta_tu^1, u^1\rangle \\ &&+\sum\limits_{n=2}^Nk\langle\frac{2}{3}(k^{\alpha}q_0^{(\alpha)})^2\delta_x^2 \delta_y^2\delta_tu^n, u^n\rangle\\ &=&k^{\alpha+1}\sum\limits_{n=1}^N\sum\limits_{p=1}^nq_{n-p}^{(\alpha)}\langle\Lambda_hu^p, u^n\rangle +k^{\alpha+1}\sum\limits_{n=1}^N\tilde{q}_n^{(\alpha)}\langle\Lambda_hu^0, u^n\rangle+k\sum\limits_{n=1}^N\langle{\cal H}f^n, u^n\rangle. \end{eqnarray} $

(4.3)

现逐项估计(4.3)式.首先, 令 $\triangle_mu_{ij}^n=u_{ij}^n-u_{ij}^{n-m}, $ $m=1, 2$, 有

$ kD_t^{(2)}u_{ij}^n=2\triangle_1u_{ij}^n-\frac{1}{2}\triangle_2u_{ij}^n. $

利用引理2.1, 有

$ \begin{eqnarray*} && 2\langle{\cal H}\triangle_mu^n, u^n\rangle \\ &=&2\sum\limits_{i=1}^{M_1-1}\sum\limits_{j=1}^{M_2-1}h_1h_2{\cal H}\triangle_mu_{ij}^nu_{ij}^n \\ &=&2\sum\limits_{i=1}^{M_1-1}\sum\limits_{j=1}^{M_2-1}h_1h_2(1+\frac{h_1^2}{12}\delta_x^2)(1+\frac{h_2^2}{12}\delta_y^2)\triangle_mu_{ij}^nu_{ij}^n\\ &=&2\sum\limits_{i=1}^{M_1-1}\sum\limits_{j=1}^{M_2-1}h_1h_2(u_{ij}^n -u_{ij}^{n-m})u_{ij}^n +2\frac{h_1^2}{12} \sum\limits_{i=1}^{M_1-1}\sum\limits_{j=1}^{M_2-1}h_1h_2\delta_x^2(u_{ij}^n-u_{ij}^{n-m})u_{ij}^n\\ &&+2\frac{h_2^2}{12}\sum\limits_{i=1}^{M_1-1}\sum\limits_{j=1}^{M_2-1} h_1h_2\delta_y^2(u_{ij}^n-u_{ij}^{n-m})u_{ij}^n +2\frac{h_1^2h_2^2}{144}\sum\limits_{i=1}^{M_1-1}\\ &&\sum\limits_{j=1}^{M_2-1}h_1h_2\delta_x^2\delta_y^2(u_{ij}^n-u_{ij}^{n-m})u_{ij}^n\\ &=&2\sum\limits_{i=1}^{M_1-1}\sum\limits_{j=1}^{M_2-1}h_1h_2(u_{ij}^n -u_{ij}^{n-m})u_{ij}^n -2\frac{h_1^2}{12}\sum\limits_{i=0}^{M_1-1}\sum\limits_{j=1}^{M_2-1}h_1h_2\delta_x(u_{ij}^n-u_{ij}^{n-m})\delta_xu_{ij}^n\\ &&-2\frac{h_2^2}{12}\sum\limits_{i=1}^{M_1-1}\sum\limits_{j=0}^{M_2-1}h_1 h_2\delta_y (u_{ij}^n-u_{ij}^{n-m})\delta_yu_{ij}^n\\ &&+2\frac{h_1^2h_2^2}{144}\sum\limits_{i=0}^{M_1-1}\sum\limits_{j=0}^{M_2-1}h_1h_2\delta_x\delta_y(u_{ij}^n-u_{ij}^{n-m})\delta_x\delta_yu_{ij}^n\\ &=&\sum\limits_{i=1}^{M_1-1}\sum\limits_{j=1}^{M_2-1}h_1h_2[(u_{ij}^n)^2-(u_{ij}^{n-m})^2+(u_{ij}^n-u_{ij}^{n-m})^2]\\ &&-\frac{h_1^2}{12}\sum\limits_{i=0}^{M_1-1}\sum\limits_{j=1}^{M_2-1}h_1h_2[(\delta_xu_{ij}^n)^2-(\delta_xu_{ij}^{n-m})^2+(\delta_xu_{ij}^n-\delta_xu_{ij}^{n-m})^2]\\ &&-\frac{h_2^2}{12}\sum\limits_{i=1}^{M_1-1}\sum\limits_{j=0}^{M_2-1}h_1h_2[(\delta_yu_{ij}^n)^2-(\delta_yu_{ij}^{n-m})^2+(\delta_yu_{ij}^n-\delta_yu_{ij}^{n-m})^2]\\ &&+\frac{h_1^2h_2^2}{144}\sum\limits_{i=0}^{M_1-1}\sum\limits_{j=0}^{M_2-1}h_1h_2[(\delta_x\delta_yu_{ij}^n)^2-(\delta_x\delta_yu_{ij}^{n-m})^2+(\delta_x\delta_yu_{ij}^n-\delta_x\delta_yu_{ij}^{n-m})^2]\\ &=&\|u^n\|_A^2-\|u^{n-m}\|_A^2+\|u^n-u^{n-m}\|_A^2=\triangle_m\|u^n\|_A^2+\|\triangle_mu^n\|_A^2, \end{eqnarray*} $

因此可得

$ \begin{eqnarray} k\sum\limits_{n=2}^N\langle{\cal H}D_t^{(2)}u^n, u^n\rangle &=&\sum\limits_{n=2}^N\langle 2{\cal H}\triangle_1u^n-\frac{1}{2}{\cal H}\triangle_2u^n, u^n\rangle \\ &=&\sum\limits_{n=2}^N(\triangle_1\|u^n\|_A^2-\frac{1}{4}\triangle_2\|u^n\|_A^2) \\ && +\sum\limits_{n=2}^N(\|\triangle_1u^n\|_A^2-\frac{1}{4}\|\triangle_2u^n\|_A^2). \end{eqnarray} $

(4.4)

方程(4.4)右边的第一项变成

$ \sum\limits_{n=2}^N(\triangle_1\|u^n\|_A^2-\frac{1}{4}\triangle_2\|u^n\|_A^2)\\ =\frac{3}{4}\|u^N\|_A^2-\frac{1}{4}\|u^{N-1}\|_A^2-\frac{3}{4}\|u^1\|_A^2+\frac{1}{4}\|u^0\|_A^2. $

(4.5)

因为 $\triangle_2u^n=\triangle_1u^n+\triangle_1u^{n-1}$和 $(a+b)^2\leq 2(a^2+b^2), \, \, a, b\in {\Bbb R}$, 由下列的方法, 可得(4.4)式右边的第二项的界

$ \sum\limits_{n=2}^N(\|\triangle_1u^n\|_A^2-\frac{1}{4}\|\triangle_2u^n\|_A^2)\\ \geq\sum\limits_{n=2}^N[\|\triangle_1u^n\|_A^2-\frac{1}{4}(\|\triangle_1u^n\|_A+\|\triangle_1u^{n-1}\|_A)^2] \\ \geq\sum\limits_{n=2}^N[\|\triangle_1u^n\|_A^2-\frac{1}{2}(\|\triangle_1u^n\|_A^2+\|\triangle_1u^{n-1}\|_A^2)] \\ =\frac{1}{2}\sum\limits_{n=2}^N(\|\triangle_1u^n\|_A^2-\|\triangle_1u^{n-1}\|_A^2) \\ =\frac{1}{2}(\|\triangle_1u^N\|_A^2-\|\triangle_1u^1\|_A^2). $

(4.6)

因此

$ \begin{eqnarray} k\langle{\cal H}\delta_tu^1, u^1\rangle=\langle{\cal H}\triangle_1u^1, u^1\rangle=\frac{1}{2}(\|u^1\|_A^2-\|u^0\|_A^2+\|\triangle_1u^1\|_A^2). \end{eqnarray} $

(4.7)

结合(4.5)-(4.7)式, 得到

$ \begin{eqnarray} && k\langle{\cal H}\delta_tu^1, u^1\rangle+k\sum\limits_{n=2}^N\langle{\cal H}D_t^{(2)}u^n, u^n\rangle \\ & =&\frac{1}{2}(\|u^1\|_A^2-\|u^0\|_A^2+\|\triangle_1u^1\|_A^2)+\frac{1}{2}(\|\triangle_1u^N\|_A^2-\|\triangle_1u^1\|_A^2) \\ &&+(\frac{3}{4}\|u^N\|_A^2-\frac{1}{4}\|u^{N-1}\|_A^2-\frac{3}{4}\|u^1\|_A^2+\frac{1}{4}\|u^0\|_A^2) \\ & \geq&\frac{3}{4}\|u^N\|_A^2-\frac{1}{4}(\|u^{N-1}\|_A^2+\|u^1\|_A^2+\|u^0\|_A^2). \end{eqnarray} $

(4.8)

其次, 易得

$ \begin{eqnarray} \langle\delta_tu^n, u^n\rangle&=&\langle\frac{1}{k}(u^n-u^{n-1}), \frac{1}{2}[(u^n-u^{n-1})+(u^n+u^{n-1})]\rangle \\ &\geq &\frac{1}{2k}[\|u^n\|^2-\|u^{n-1}\|^2]=\frac{1}{2}\delta_t\|u^n\|^2. \end{eqnarray} $

(4.9)

利用引理2.2, 有

$ \begin{eqnarray} \langle\delta_x^2\delta_y^2\delta_tu^n, u^n\rangle=\delta_t\langle\delta_x^2\delta_y^2u^n, u^n\rangle=\delta_t\langle\delta_x\delta_yu^n, \delta_x\delta_yu^n\rangle=\delta_t\|\delta_x\delta_yu^n\|^2. \end{eqnarray} $

(4.10)

结合(4.9)-(4.10)式, 可得

$ \begin{eqnarray} && k\langle(k^{\alpha}q_0^{(\alpha)})^2\delta_x^2\delta_y^2\delta_tu^1, u^1\rangle+\sum\limits_{n=2}^Nk\langle\frac{2}{3}(k^{\alpha}q_0^{(\alpha)})^2\delta_x^2\delta_y^2\delta_tu^n, u^n\rangle \\ &=& (k^{\alpha}q_0^{(\alpha)})^2(\|\delta_x\delta_yu^1\|^2-\|\delta_x\delta_yu^0\|^2)+\frac{2}{3}(k^{\alpha}q_0^{(\alpha)})^2(\|\delta_x\delta_yu^N\|^2-\|\delta_x\delta_yu^1\|^2) \\ & =& \frac{2}{3}(k^{\alpha}q_0^{(\alpha)})^2\|\delta_x\delta_yu^N\|^2+\frac{1}{3}(k^{\alpha}q_0^{(\alpha)})^2\|\delta_x\delta_yu^1\|^2-(k^{\alpha}q_0^{(\alpha)})^2\|\delta_x\delta_yu^0\|^2 \\ &\geq& -(k^{\alpha}q_0^{(\alpha)})^2\|\delta_x\delta_yu^0\|^2=-(k^{\alpha}q_0^{(\alpha)})^2\langle\delta_x^2\delta_y^2u^0, u^0\rangle. \end{eqnarray} $

(4.11)

再者, 首先利用引理2.1, 然后对于固定的 $i, j$, 交换求和指标, 再利用引理4.2可得(4.3)式右边的第一项非负

$ \sum\limits_{n=1}^N\sum\limits_{p=1}^nq_{n-p}^{(\alpha)} \langle\Lambda_hu^p, u^n\rangle\\ =\sum\limits_{n=1}^N\sum\limits_{p=1}^nq_{n-p}^{(\alpha)}\langle({\cal H}_y\delta_x^2+{\cal H}_x\delta_y^2)u^p, u^n\rangle \\ =-h_1h_2\sum\limits_{n=1}^N\sum\limits_{p=1}^nq_{n-p}^{(\alpha)}\bigg(\sum\limits_{i=0}^{M_1-1}\sum\limits_{j=1}^{M_2-1}\delta_xu_{ij}^p\delta_xu_{ij}^n+\frac{h_2^2}{12}\\ \sum\limits_{i=0}^{M_1-1}\sum\limits_{j=0}^{M_2-1}\delta_y\delta_xu_{ij}^p\delta_y\delta_xu_{ij}^n \\ +\sum\limits_{i=1}^{M_1-1}\sum\limits_{j=0}^{M_2-1}\delta_yu_{ij}^p\delta_yu_{ij}^n+\frac{h_1^2}{12}\sum\limits_{i=0}^{M_1-1}\sum\limits_{j=0}^{M_2-1}\delta_x\delta_yu_{ij}^p\delta_x\delta_yu_{ij}^n\bigg) \\ =-h_1h_2\bigg(\sum\limits_{i=0}^{M_1-1}\sum\limits_{j=1}^{M_2-1}\sum\limits_{n=1}^N\sum\limits_{p=1}^nq_{n-p}^{(\alpha)}{\cal H}_y\delta_xu_{ij}^p{\cal H}_y\delta_xu_{ij}^n \\ +\sum\limits_{i=1}^{M_1-1}\sum\limits_{j=0}^{M_2-1}\sum\limits_{n=1}^N\sum\limits_{p=1}^nq_{n-p}^{(\alpha)}{\cal H}_x\delta_yu_{ij}^p{\cal H}_x\delta_yu_{ij}^n\bigg) \\ \leq 0. $

(4.12)

因此, 利用(4.8), (4.11)-(4.12)式, 有

$ \begin{eqnarray} && \frac{3}{4}\|u^N\|_A^2-\frac{1}{4}(\|u^{N-1}\|_A^2+\|u^1\|_A^2+\|u^0\|_A^2) -(k^{\alpha}q_0^{(\alpha)})^2\langle\delta_x^2\delta_y^2u^0, u^0\rangle \\ &\leq &k^{\alpha+1}\sum\limits_{n=1}^N \bigg[\sum\limits_{p=1}^nq_{n-p}^{(\alpha)}\langle\Lambda_hu^p, u^n\rangle+\tilde{q}_n^{(\alpha)}\langle\Lambda_hu^0, u^n\rangle\bigg]+k\sum\limits_{n=1}^N\langle{\cal H}f^n, u^n\rangle \\ &\leq& k^{\alpha+1}\sum\limits_{n=1}^N\tilde{q}_n^{(\alpha)}\langle\Lambda_hu^0, u^n\rangle+k\sum\limits_{n=1}^N\langle{\cal H}f^n, u^n\rangle. \end{eqnarray} $

(4.13)

选取满足 $\|u^M\|_A=\max\limits_{0\leq n\leq N}{\|u^n\|_A}$的 $M$, 利用引理2.3, 引理4.1和Cauchy-Schwarz不等式, 那么

$ \begin{eqnarray*} \|u^M\|_A^2&\leq&\frac{1}{3}(\|u^{M-1}\|_A^2+\|u^1\|_A^2+\|u^0\|_A^2)+\frac{4}{3}(k^{\alpha}q_0^{(\alpha)})^2\|\delta_x^2\delta_y^2u^0\|\|u^M\|_A\\ &&+\frac{2k^{\alpha+1}}{3}\sum\limits_{n=1}^N|\tilde{q}_n^{(\alpha)}|(\frac{1}{h_1^2}+\frac{1}{h_2^2})\|u^0\|_A\|u^M\|_A+\frac{4}{3}k\sum\limits_{n=1}^N\|{\cal H}f^n\|\|u^M\|_A. \end{eqnarray*} $

可得

$ \begin{eqnarray*} \|u^M\|_A&\leq&\frac{1}{2}(\|u^1\|_A+\|u^0\|_A)+2(k^{\alpha}q_0^{(\alpha)})^2\|\delta_x^2\delta_y^2u^0\| \\ && +8k^{\alpha+1}\sum\limits_{n=1}^N|\tilde{q}_n^{(\alpha)}|(\frac{1}{h_1^2}+\frac{1}{h_2^2})\|u^0\|_A+2k\sum\limits_{n=1}^N\|{\cal H}f^n\|, \end{eqnarray*} $

所以, 对 $n\geq 2$, 有

$ \begin{eqnarray} \|u^n\|_A&\leq&\frac{1}{2}(\|u^1\|_A+\|u^0\|_A)+2(k^{\alpha}q_0^{(\alpha)})^2\|\delta_x^2\delta_y^2u^0\| \\ && +8k^{\alpha+1}\sum\limits_{n=1}^N|\tilde{q}_n^{(\alpha)}|(\frac{1}{h_1^2}+\frac{1}{h_2^2})\|u^0\|_A+2k\sum\limits_{n=1}^N\|{\cal H}f^n\|. \end{eqnarray} $

(4.14)

由(4.1)式, 利用引理4.2, 有

$ \begin{eqnarray*} \frac{1}{2}(\|u^1\|_A^2-\|u^0\|_A^2)&\leq& k^{\alpha+1}[q_0^{(\alpha)}\langle\Lambda_hu^1, u^1\rangle+\tilde{q}_1^{(\alpha)}\langle\Lambda_hu^0, u^1\rangle]+k\langle{\cal H}f^1, u^1\rangle\\ &\leq& k^{\alpha+1}\tilde{q}_1^{(\alpha)}\langle\Lambda_hu^0, u^1\rangle+k\langle{\cal H}f^1, u^1\rangle. \end{eqnarray*} $

因为 $2ab\leq a^2+b^2, \, \, a, b\in {\Bbb R}$, 利用引理2.3和引理4.1, 有

$ \begin{equation} \|u^1\|_A\leq\|u^0\|_A+8\frac{k^{\alpha+1}}{h^2}|\tilde{q}_1^{(\alpha)}|(\frac{1}{h_1^2}+\frac{1}{h_2^2})\|u^0\|_A+2k\|{\cal H}f^1\|. \end{equation} $

(4.15)

对 $1\leq n\leq N$, $\tilde{q}_n^{(\alpha)}=O(k^{\alpha}n^{\alpha-1})$ (参见文献[16])和利用文献[15, (4.7式]), 得

$ \begin{equation} \sum\limits_{n=1}^N\mid \tilde{q}_n^{(\alpha)}\mid\leq C\sum\limits_{n=1}^N(k^{\alpha}n^{\alpha-1})\leq Ck^{\alpha}NN^{\alpha-1}\leq C(Nk)^{\alpha}\leq C(T). \end{equation} $

(4.16)

利用(4.14)-(4.16)式, 我们完成了证明.

通过能量方法, 我们推导出正式的二阶BDF ADI紧致差分格式(3.19)-(3.22)的收敛性.

令

$ e_{ij}^n=u(x_i, y_j, t_n)-u_{ij}^n, \qquad (x_i, y_j)\in\Omega_h, \quad 0\leq n\leq N. $

将(3.5)式减去(3.14)-(3.15)式, 可得误差方程

$ \begin{eqnarray} && {\cal H}\delta_te_{ij}^1-k^{\alpha}(q_0^{(\alpha)}\Lambda_he_{ij}^1+\tilde{q}_1^{(\alpha)}\Lambda_he_{ij}^0)+(k^{\alpha+1}q_0^{(\alpha)})^2\delta_x^2\delta_y^2\delta_te_{ij}^1 \\ & =&{\cal H}R_{ij}^1, \quad (x_i, y_j)\in \Omega_h. \end{eqnarray} $

(4.17)

$ \begin{eqnarray} && {\cal H}D_t^{(2)}e_{ij}^n-k^{\alpha} \bigg(\sum\limits_{p=1}^nq_{n-p}^{(\alpha)}\Lambda_he_{ij}^p+\tilde{q}_n^{(\alpha)}\Lambda_he_{ij}^0\bigg)+\frac{2}{3}(k^{\alpha+1}q_0^{(\alpha)})^2\delta_x^2\delta_y^2\delta_te_{ij}^n \\ &=&{\cal H}R_{ij}^n, \quad (x_i, y_j)\in \Omega_h, \quad 2\leq n\leq N. \end{eqnarray} $

(4.18)

$ e_{ij}^n=0, \qquad (x_i, y_j)\in \partial\Omega_h, \qquad 1\leq n\leq N. $

(4.19)

$ e_{ij}^0=0, \qquad (x_i, y_j)\in\Omega_h. $

(4.20)

由定理4.3, 有

$ \|e^n\|_A\leq C(T)\|e^0\|_A+(k^{\alpha}q_0^{(\alpha)})^2\|\delta_x^2\delta_y^2e^0\|+2k\sum\limits_{l=1}^n\|{\cal H}R^l\|. $

将(3.16)式代入上面的不等式, 可得

$ \|e^n\|_A\leq C(k^{\alpha+1}+h_1^4+h_2^4). $

(4.21)

利用引理4.1, 定理4.3和(4.21)式, 我们有如下收敛性结果.

定理4.4  设定解问题(1.1)-(1.3)在区域 $\Omega\times (0, T]$上有光滑解, $\{u_{ij}^n|(x_i, y_j)\in\Omega_h, 1\leq n\leq N\}$是正式的二阶BDF ADI紧致差分格式(3.19)-(3.22)的解, 那么成立

$ \max\limits_{1\leq n\leq N}\|u(x_i, y_j, t_n)-u_{ij}^n\|\leq C(k^{\alpha+1}+h_1^4+h_2^4). $

本节中取 $x_R=y_R=\pi$和 $T=1$, 取相同的空间步长, 即: $h_1=h_2=h$.分别用向后Euler ADI差分格式^[11]、二阶BDF ADI差分格式^[5]和正式的二阶BDF ADI紧致差分格式(3.19)-(3.22)离散定解问题(1.1)-(1.3).用 $u_{BE}$, $u_{BDF}$和 $u_{CBDF}$分别表示它们的数值解.

在时空域 $\Omega\times (0, T]$, 记

$ E_{BE}(k, h)=\max\limits_{1\leq n\leq N}\|u(x_i, y_j, t_n)-{u_{BE}}_{ij}^n\|, $

(5.1)

$ rate_{BE}^t=\log_2\bigg(\frac{E_{BE}(2k, h)}{E_{BE}(k, h)}\bigg), \qquad rate_{BE}^x=\log_2\bigg(\frac{E_{BE}(k, 2h)}{E_{BE}(k, h)}\bigg). $

(5.2)

记号 $E_{BDF}(k, h)$, $rate_{BDF}^t$, $rate_{BDF}^x$, $E_{CBDF}(k, h)$, $rate_{CBDF}^t$, $rate_{CBDF}^x$如上类似定义.

算例  在我们的算例中, 精确解为

$ u(x, y, t)=\sin x\sin y-\frac{t^{\alpha+1}}{\Gamma (\alpha+2)}\sin 2 x\sin 2 y, $

(5.3)

初始条件 $u^0(x, y)=\sin x\sin y$和非齐次项为

$ f(x, y, t)=\frac{t^{\alpha}}{\Gamma (\alpha+1)}(2\sin x\sin y-\sin 2x\sin 2y)-\frac{8t^{2\alpha+1}}{\Gamma (2\alpha+2)}\sin 2x\sin 2y. $

(5.4)

固定空间步数 $M_1=M_2=1024$, 表 1给出了 $\alpha$=0.25, 0.5, 0.75, 0.95时向后Euler ADI差分格式和正式的二阶BDF ADI紧致差分格式的最大模误差和时间误差阶.数值结果表明向后Euler ADI差分格式的时间收敛阶约为一阶, 正式的二阶BDF ADI紧致差分格式的时间收敛阶 $\approx 1+\alpha$.由于 $u_{tt}(x, t)$在 $t=0$处奇异, 在一致时间步长 $k$下, 正式的二阶BDF ADI紧致差分格式的二阶精度不能达到.然而相比于向后Euler ADI差分格式^[11], 时间方向的精度明显提高.

固定时间步数 $N=2048$, 表 2给出了 $\alpha$=0.25, 0.5, 0.75, 0.95时二阶BDF ADI差分格式^[5]和紧致差分格式(3.19)-(3.22)的最大模误差和相应的空间误差阶.数值结果表明差分格式的空间收敛阶为二阶, 而紧致差分格式的空间收敛阶为四阶.这与定理4.4的结论相吻合.

本文建立和分析了二维分数阶发展型方程的正式的二阶BDF ADI紧致差分格式.基于与卷积求积相对应的实二次型的非负性, 利用能量方法证明了差分格式在一类等价于标准的 $L^2$的新模下是稳定性和收敛性.理论结果表明紧致差分格式的收敛阶为 $O(k^{\alpha+1}+h_1^4+h_2^4)$, 其中 $k$为时间步长, $h_1$和 $h_2$分别是空间 $x$和 $y$方向的步长.最后, 数值算例验证了理论分析的正确性.