修正RLW方程的一种混合有限元方法的数值分析与模拟

The Numerical Analysis and Simulation Based on A Mixed Finite Element Method for the MRLW Equation

摘要: 本文针对修正的RLW方程提出并讨论了一种二阶向后差分的混合有限元方法。在空间方向上使用混合Galerkin有限元方法来近似，在时间上采用向后差q = u_x分二阶离散格式来近似。并且得到了方程的未知解u在L²模和H¹模下的最优误差估计以及在L²模下的最优误差估计。为了对数值理论进行有效性验证，我们通过一些数值算例给出一些数值模拟结果。

Abstract: A second-order backward-difference mixed finite element (MFE) method for modified regularized long wave (MRLW) equation is proposed and discussed in this paper. The spatial direction is approximated by the mixed Galerkin method using mixed linear space finite elements, and the time direction is considered by backward difference scheme with second-order convergence rate. The optimal error estimates for u in L² and H1-norms and its flux q = u_x and in L²-norm are derived. Some numerical results are given to test our theoretical analysis and illustrate the efficiency of the studied method.

1. 引言

本文考虑带有初边值条件的修正RLW方程问题：

$\{\begin{array}{l}{u}_t+u_x+6u^2{u}_x-\mu u_{xxt}=0,\left(x,t\right)\in I\times J\\ u\left(a,t\right)=u\left(b,t\right)=0,t\in \bar{J},\\ u\left(x,0\right)=u_0\left(x\right),x\in \bar{I},\end{array}$ (1.1)

其中 $I=\left(a,b\right)$， $J=\left[0,T\right]\left(0<T<\infty \right)$。初始函数 $u_0\left(x\right)$ 为已知函数，系数 $\mu$ 为正常数。

近几年来，物理、生物、流体力学等研究领域均出现了很多的非线性现象，而这些现象都可以由一些非线性发展方程的数学模型来描述。特别地是，越来越多的学者将关注点放在了非线性色散波研究领域中较为重要的非线性RLW方程中 [1]。孤立波在非线性分散介质中是波包或脉冲，由于非线性和色散影响之间的动态平衡，这些波的波形会保持不变，而孤波是一种非常特殊的类型，在与其他波进行碰撞后波形依然保持不变。RLW方程对于更为一般的KdV方程是一种非线性分散波的替代描述。我们考虑文献 [1] - [23] 中的数学理论以及数值方法。文献 [5] 中讨论了RLW方程解的存在性和唯一性，因为只有在有限的初边值条件下方程才有解析解，因此他们将焦点放在了数值的角度。在一些文献中已经有很多针对RLW方程的数值方法，比如差分方法 [5] [6] [15]，基于Galerkin和排列原则的有限元方法 [7] [8] [9] [10] [11]，混合有限元方法 [12] [13] [14] [23]，无网格法 [16] 等。

Mei和Chen [18] 针对RLW方程研究了基于Galerkin方法的显式两步方法。在本文中，我们的目的是对修正RLW方程研究一种显式两步混合元法，利用混合Galerkin方法逼近空间方向，向后二阶差分格式逼近时间方向。我们在全离散显式两步混合格式下，得到了未知函数u在 $L^2$ 模和 $H^1$ 模下最优误差估计和其辅助变量 $q=u_x$ 基于 $L^2$ 模最优误差估计，并与其他的数值方法的准确性进行了比较。通过与文献 [18] [24] [25] 中的数值方法进行比较，我们得到了u与 $q=u_x$ 的近似值。

文章的整体结构如下所示，第二部分中，给出了一个显式多步混合方法以及数值过程。第三部分中，证明了u与 $q=u_x$ 分别在 $L^2$ 模和 $H^1$ 模下基于二阶向后差分混合方法的最优误差估计。第四部分中，给出一些数值算例，对理论的结果进行数值验证，最后对整篇论文做了总结。在整篇文章中，C表示一般的正常数，与剖分步长h与 $\Delta t$ 无关。

2. 混合有限元数值方法

引入辅助变量 $q=u_x$，我们可以将方程(1.1)整理为如下形式：

$u_x=q$ (2.1)

$u_t+q+6u^{2}q-\mu q_{xt}=0.$ (2.2)

要得到方程(2.1)~(2.2)的混合弱形式，也就是找到 $\left\{u,q\right\}:\left[0,T\right]\mapsto H_0^1\ast L^2\left(\Omega \right)$ 使得如下系统成立：

$\left(u_x,w\right)=\left(q,w\right),\forall w\in L^2\left(\Omega \right)$ (2.3)

$\left(u_t,v\right)+\left(q,v\right)+6\left(u^{2}q,v\right)+\mu \left(q_t,v_x\right)=0,\forall v\in H_0^1$ .(2.4)

使 $V_h$ 与 $W_h$ 分别为 $H_0^1$ 和 $L^2\left(\Omega \right)$ 的有限维子空间，其定义如下：

$V_h=\left\{v_h|v_h\in C^0\left(\bar{I}\right),{v_h|}_{I_j}\in P_k\left(I_j\right),\forall I_j\in T_h,v_h\left(a\right)=v_h\left(b\right)=0\right\}\subset H_0^1,$

$W_h=\left\{w_h|w_h\in L^2\left(\Omega \right),{w_h|}_{I_j}\in P_r\left(I_j\right),\forall I_j\in T_h\right\}\subset L^2\left(\Omega \right),$

其中 $T_h$ 为小区间的集合，即将 $\bar{I}=\left[a,b\right]$ 剖分为N个子区间，分别为 $I_j=\left[x_j,x_{j+1}\right]$， $j=0,1,2,\cdots ,\left(N-1\right)$， $h_j=x_{j+1}-x_j$， $h=\underset{0\le j\le N-1}{\max }{h}_j$， $P_m\left(I_j\right)$ 表示在 $I_j$ 中小于等于m的多项式。方程(2.3)和(2.4)的半离散混合有限元方法就是找到 $\left\{u_h,q_h\right\}:\left[0,T\right]\mapsto V_h\times W_h$，使得：

$\left(u_{hx},w_h\right)=\left(q_h,w_h\right),\forall w_h\in W_h,$ (2.5)

$\left(u_{ht},v_h\right)+\left(q_h,v_h\right)+6\left({\left(u_h\right)}^2{q}_h,v_h\right)+\mu \left(q_{ht},v_{hx}\right)=0,\forall v_h\in V_h.$ (2.6)

3. 二阶向后差分混合方法和最优误差估计

3.1. 二阶向后差分方法和相关定理

在这部分，我们将对修正的RLW方程在二阶向后差分混合有限元下讨论先验误差估计。对于整个离散过程，我们对时间区间 $\left[0,T\right]$ 进行剖分， $0<t_0<t_1<\cdots <t_N=T$，其中步长为 $\Delta t=T/N$，N为正整数。 $\varphi$ 为 $\left[0,T\right]$ 上的光滑函数，定义 ${\varphi }^n=\varphi \left(t_n\right)$，以及 ${\partial }_t^2{u}^{n+1}=\frac{3u^{n+1}-4u^n+u^{n-1}}{2\Delta t}$。在时间 $t=t_{n+1}$ 处，方程(2.3)-(2.4)有如下形式；

$\left(u_x^{n+1},w\right)=\left(q^{n+1},w\right),\forall w\in L_2\left(\Omega \right),$ (3.1)

$\left(u_t^{n+1},v\right)+\left(q^{n+1},v\right)+6\left({\left(u^{n+1}\right)}^2{q}^{n+1},v\right)+\mu \left(q_t^{n+1},v_x\right)=0,\forall v\in H_0^1.$ (3.2)

将方程(3.1)~(3.2)在时间方向上利用向后二阶差分近似，我们得到线性化的时间离散形式：

$\left(u_x^{n+1},w\right)=\left(q^{n+1},w\right),\forall w\in L_2\left(I\right),$ (3.3)

$\begin{array}{l}\left({\partial }_t^2{u}^{n+1},v\right)+\mu \left({\partial }_t^2{q}^{n+1},v_x\right)+\left(q^{n+1},v\right)+6\left(2{\left(u^n\right)}^2{q}^n-{\left(u^{n-1}\right)}^2{q}^{n-1},v\right)\\ =\left({\tau }_1^{n+1},v\right)+\left({\tau }_2^{n+1},v_x\right)+6\left({\tau }_3^{n+1},v\right),\forall v\in H_0^1,\end{array}$ (3.4)

其中

$\begin{array}{l}{\tau }_1^{n+1}={\partial }_t^2{u}^{n+1}-u_t^{n+1},\\ {\tau }_2^{n+1}=\mu \left({\partial }_t^2{q}_x^{n+1}-q_{xt}^{n+1}\right),\\ {\tau }_3^{n+1}=\left(2{\left(u^n\right)}^2{q}^n-{\left(u^{n-1}\right)}^2{q}^{n-1}\right)-{\left(u^{n+1}\right)}^2{q}^{n+1}.\end{array}$ (3.5)

引理3.1：(文献 [26] )对于 ${\tau }_i^{n+1},i=1,2,3$，下列估计成立：

$\Vert {\tau }_1^{n+1}\Vert +\Vert {\tau }_2^{n+1}\Vert +\Vert {\tau }_3^{n+1}\Vert \le C\Delta t^2.$ (3.6)

引理3.2：(文献 [26] )对于数列 $\left\{w^n\right\}$，下列不等式成立：

$\left({\partial }_t^2{w}^{n+1},w^{n+1}\right)\ge \frac{1}{4\Delta t}\left[{\Vert w^{n+1}\Vert }^2+{\Vert 2w^{n+1}-w^n\Vert }^2-{\Vert w^n\Vert }^2-{\Vert 2w^n-w^{n-1}\Vert }^2\right].$ (3.7)

现在我们需要在 $V_h\times W_h$ 中找到一组 $\left\{u_h^{n+1},q_h^{n+1}\right\}$ 满足：

$\left(u_{hx}^{n+1},w_h\right)=\left(q_h^{n+1},w_h\right),\forall w_h\in W_h,$ (3.8)

和

$\left({\partial }_t^2{u}_h^{n+1},v_h\right)+\mu \left({\partial }_t^2{q}_h^{n+1},v_{hx}\right)+\left(q_h^{n+1},v_h\right)+6\left(2{\left(u_h^n\right)}^2{q}_h^n-{\left(u_h^{n-1}\right)}^2{q}_h^{n-1},v_h\right)=0,\forall v_h\in V_h.$ (3.9)

为了先验误差估计的理论分析，针对我们的方程引入了一些投影。

引理3.3：(文献 [12] [17] [27] )若存在一个非线性算子 ${\Pi }_h:L^2\left(\Omega \right)\to W_h$ 使得下式成立

$\left(q-{\Pi }_{h}q,v_{hx}\right)=0,v_h\in V_h,$ (3.10)

则下面的估计式成立：

$\Vert q-{\Pi }_{h}q\Vert \le C{h}^{k+1}{\Vert q\Vert }_{k+1},$ (3.11)

引理3.4：(文献 [12] [17] [27] )对于任意的 $u\in H_0^1$ 成立

$\left(u_x-P_h{u}_x,w_h\right)=0,w_h\in W_h,$ (3.12)

则有

$\Vert u_x-P_h{u}_x\Vert +h{\Vert u-P_{h}u\Vert }_1\le C{h}^{k+1}{\Vert u\Vert }_{k+1}$ (3.13)

为了全离散误差估计，我们将误差可以写为：

$u^n-u_h^n=\left(u^n-P_h{u}^n\right)+\left(P_h{u}^n-u_h^n\right)={\eta }^n+{\varsigma }^n,$

$q^n-q_h^n=\left(q^n-{\Pi }_h{q}^n\right)+\left({\Pi }_h{q}^n-q_h^n\right)={\rho }^n+{\xi }^n.$

联立方程(3.10)，(3.13)，(3.3)，(3.4)，(3.8)和(3.9)，在 $t=t_{n+1}$ 处可以得到如下的误差方程：

$\left({\xi }^{n+1}+{\rho }^{n+1},w_h\right)=\left({\varsigma }_x^{n+1},w_h\right),\forall w_h\in W_h,$ (3.14)

$\begin{array}{l}\left({\partial }_t^2{\varsigma }^{n+1},v_h\right)+\mu \left({\partial }_t^2{\xi }^{n+1},v_{hx}\right)+6\left(2\left({\left(u^n\right)}^2{q}^n-{\left(u_h^n\right)}^2{q}_h^n\right)-\left({\left(u^{n-1}\right)}^2{q}^{n-1}-{\left(u_h^{n-1}\right)}^2{q}_h^{n-1}\right),v_h\right)\\ +\left({\rho }^{n+1}+{\xi }^{n+1},v_h\right)=-\left({\partial }_t^2{\eta }^{n+1},v_h\right)+\left({\tau }_1^{n+1},v_h\right)+\left({\tau }_2^{n+1},v_{hx}\right)+6\left({\tau }_3^{n+1},v_h\right),v_h\in V_h.\end{array}$ (3.15)

3.2. 最优误差估计

这部分中将得到全离散最优误差估计，并且获得了如下的定理3.5

定理3.5：已知 $u_h^0\in V_h,q_h^0\in W_h$，则

$\begin{array}{l}\Vert u^n-u_h^n\Vert +\Vert q^n-q_h^n\Vert \le C\left(h^{k+1}+\Delta t^2\right),\\ {\Vert u^n-u_h^n\Vert }_1\le C\left(h^k+\Delta t^2\right).\end{array}$

证明：由(3.14)可得

$\left({\partial }_t^2{\xi }^{n+1},w_h\right)=\left({\partial }_t^2{\varsigma }_x^{n+1},w_h\right)-\left({\partial }_t^2{\rho }^{n+1},w_h\right),\forall w_h\in W_h,$ (3.16)

在(3.15)和(3.16)中分别令 $v_h={\varsigma }^{n+1}$， $w_h=\mu {\varsigma }_x^{n+1}$，可以得到下面方程

$\begin{array}{l}\left({\partial }_t^2{\varsigma }^{n+1},{\varsigma }^{n+1}\right)+\mu \left({\partial }_t^2{\varsigma }_x^{n+1},{\varsigma }_x^{n+1}\right)\\ =\mu \left({\partial }_t^2{\rho }^{n+1},{\varsigma }_x^{n+1}\right)-6\left(2\left({\left(u^n\right)}^2{q}^n-{\left(u_h^n\right)}^2{q}_h^n\right)-\left({\left(u^{n-1}\right)}^2{q}^{n-1}-{\left(u_h^{n-1}\right)}^2{q}_h^{n-1}\right),{\varsigma }^{n+1}\right)\\ +\left({\rho }^{n+1}+{\xi }^{n+1},{\varsigma }^{n+1}\right)+\left({\tau }_1^{n+1},{\varsigma }^{n+1}\right)+\left({\tau }_2^{n+1},{\varsigma }_x^{n+1}\right)+6\left({\tau }_3^{n+1},{\varsigma }^{n+1}\right)\\ \le C({\Vert 2\left({\left(u^n\right)}^2{q}^n-{\left(u_h^n\right)}^2{q}_h^n\right)-\left({\left(u^{n-1}\right)}^2{q}^{n-1}-{\left(u_h^{n-1}\right)}^2{q}_h^{n-1}\right)\Vert }^2\\ \begin{array}{c}\\ \end{array}+{\Vert {\rho }^{n+1}\Vert }^2+{\Vert {\xi }^{n+1}\Vert }^2+{\Vert {\tau }_1^{n+1}\Vert }^2+{\Vert {\tau }_2^{n+1}\Vert }^2+{\Vert {\tau }_3^{n+1}\Vert }^2+{\Vert {\varsigma }^{n+1}\Vert }^2+{\Vert {\varsigma }_x^{n+1}\Vert }^2).\end{array}$ (3.17)

现在我们估计不等式可以得到

$\begin{array}{l}{\Vert 2\left({\left(u^n\right)}^2{q}^n-{\left(u_h^n\right)}^2{q}_h^n\right)-\left({\left(u^{n-1}\right)}^2{q}^{n-1}-{\left(u_h^{n-1}\right)}^2{q}_h^{n-1}\right)\Vert }^2\\ \le {\Vert 2{\left(u^n\right)}^2\left(q^n-q_h^n\right)+2\left({\left(u^n\right)}^2-{\left(u_h^n\right)}^2\right)q_h^n-\left({\left(u^{n-1}\right)}^2\left(q^{n-1}-q_h^{n-1}\right)+\left[{\left(u^{n-1}\right)}^2-{\left(u_h^{n-1}\right)}^2\right]q_h^{n-1}\right)\Vert }^2\\ \le C{\displaystyle \underset{j=n-1}{\overset{n}{\sum }}\left[{\Vert {\left(u^j\right)}^2\Vert }_{\infty }^2{\left(\Vert {\xi }^j\Vert +\Vert {\rho }^j\Vert \right)}^2+{\left({\Vert u^j\Vert }_{\infty }+{\Vert u_h^j\Vert }_{\infty }\right)}^2{\Vert q_h^j\Vert }_{\infty }^2{\left(\Vert {\eta }^j\Vert +\Vert {\varsigma }^j\Vert \right)}^2\right]}\\ \le C\left({\Vert {\rho }^n\Vert }^2+{\Vert {\rho }^{n-1}\Vert }^2+{\Vert {\eta }^n\Vert }^2+{\Vert {\eta }^{n-1}\Vert }^2+{\Vert {\xi }^n\Vert }^2+{\Vert {\xi }^{n-1}\Vert }^2+{\Vert {\varsigma }^n\Vert }^2+{\Vert {\varsigma }^{n-1}\Vert }^2\right).\end{array}$ (3.18)

在方程(3.14)中，我们令 $w_h={\xi }^{n+1}$，使用Cauchy-Schwarz不等式可以得到

${\Vert {\xi }^{n+1}\Vert }^2\le C\left({\Vert {\rho }^{n+1}\Vert }^2+{\Vert {\varsigma }_x^{n+1}\Vert }^2\right).$ (3.19)

现在联合方程(3.17)，(3.18)，(3.19)可以得到

$\begin{array}{l}\left({\partial }_t^2{\varsigma }^{n+1},{\varsigma }^{n+1}\right)+\mu \left({\partial }_t^2{\varsigma }_x^{n+1},{\varsigma }_x^{n+1}\right)\\ =\mu \left({\partial }_t^2{\rho }^{n+1},{\varsigma }_x^{n+1}\right)-6\left(2\left({\left(u^n\right)}^2{q}^n-{\left(u_h^n\right)}^2{q}_h^n\right)-\left({\left(u^{n-1}\right)}^2{q}^{n-1}-{\left(u_h^{n-1}\right)}^2{q}_h^{n-1}\right),{\varsigma }^{n+1}\right)\\ +\left({\rho }^{n+1}+{\xi }^{n+1},{\varsigma }^{n+1}\right)+\left({\tau }_1^{n+1},{\varsigma }^{n+1}\right)+\left({\tau }_2^{n+1},{\varsigma }_x^{n+1}\right)+6\left({\tau }_3^{n+1},{\varsigma }^{n+1}\right)\\ \le C({\Vert {\rho }^{n+1}\Vert }^2+{\Vert {\tau }_1^{n+1}\Vert }^2+{\Vert {\tau }_2^{n+1}\Vert }^2+{\Vert {\tau }_3^{n+1}\Vert }^2+{\displaystyle \underset{i=n-1}{\overset{n+1}{\sum }}\left({\Vert {\varsigma }^i\Vert }^2+{\Vert {\varsigma }_x^i\Vert }^2\right)}\\ \begin{array}{c}\\ \end{array}+{\Vert {\partial }_t^2{\rho }^{n+1}\Vert }^2+{\Vert {\rho }^n\Vert }^2+{\Vert {\rho }^{n-1}\Vert }^2+{\Vert {\eta }^n\Vert }^2+{\Vert {\eta }^{n-1}\Vert }^2).\end{array}$ (3.20)

利用引理3.2，并且用 $4\Delta t$ 乘方程(3.20)，我们得到

$\begin{array}{l}{\Vert {\varsigma }^{n+1}\Vert }_1^2+{\Vert 2{\varsigma }^{n+1}-{\varsigma }^n\Vert }_1^2-{\Vert {\varsigma }^n\Vert }_1^2-{\Vert 2{\varsigma }^n-{\varsigma }^{n-1}\Vert }_1^2\\ \le C({\Vert {\rho }^{n+1}\Vert }^2+{\Vert {\tau }_1^{n+1}\Vert }^2+{\Vert {\tau }_2^{n+1}\Vert }^2+{\Vert {\tau }_3^{n+1}\Vert }^2+{\displaystyle \underset{i=n-1}{\overset{n+1}{\sum }}\left({\Vert {\varsigma }^i\Vert }^2+{\Vert {\varsigma }_x^i\Vert }^2\right)}\\ \begin{array}{c}\\ \end{array}+{\Vert {\partial }_t^2{\rho }^{n+1}\Vert }^2+{\Vert {\rho }^n\Vert }^2+{\Vert {\rho }^{n-1}\Vert }^2+{\Vert {\eta }^n\Vert }^2+{\Vert {\eta }^{n-1}\Vert }^2).\end{array}$ (3.21)

对n求和并使用离散Gronwall引理，可以得到：

$\begin{array}{l}{\Vert {\varsigma }^{n+1}\Vert }_1^2+{\Vert 2{\varsigma }^{n+1}-{\varsigma }^n\Vert }_1^2-{\Vert {\varsigma }^n\Vert }_1^2-{\Vert 2{\varsigma }^n-{\varsigma }^{n-1}\Vert }_1^2\\ \le C({\Vert {\rho }^{n+1}\Vert }^2+{\Vert {\tau }_1^{n+1}\Vert }^2+{\Vert {\tau }_2^{n+1}\Vert }^2+{\Vert {\tau }_3^{n+1}\Vert }^2\\ +{\Vert {\partial }_t^2{\rho }^{n+1}\Vert }^2+{\Vert {\rho }^n\Vert }^2+{\Vert {\rho }^{n-1}\Vert }^2+{\Vert {\eta }^n\Vert }^2+{\Vert {\eta }^{n-1}\Vert }^2).\end{array}$ (3.22)

联立方程(3.20)和(3.19)可以得到

$\begin{array}{c}{\Vert {\xi }^{n+1}\Vert }^2\le C({\Vert {\rho }^{n+1}\Vert }^2+{\Vert {\tau }_1^{n+1}\Vert }^2+{\Vert {\tau }_2^{n+1}\Vert }^2+{\Vert {\tau }_3^{n+1}\Vert }^2\\ +{\Vert {\partial }_t^2{\rho }^{n+1}\Vert }^2+{\Vert {\rho }^n\Vert }^2+{\Vert {\rho }^{n-1}\Vert }^2+{\Vert {\eta }^n\Vert }^2+{\Vert {\eta }^{n-1}\Vert }^2).\end{array}$ (3.23)

联立方程(3.22)~(3.23)，引理3.3~3.4以及三角不等式可以得到定理中的结论。

4. 数值实验

在这部分中，为了检验显式两步方法的有效性和准确性，我们给出了不同时间层上单波传播以及双波碰撞数值算例。为了检验数值格式的守恒，我们可以考虑如下三个不变量：

$\begin{array}{l}{I}_1={\displaystyle {\int }_a^{b}u\text{d}x}\approx h{\displaystyle \underset{j=1}{\overset{N}{\sum }}{u}_j^n},\\ I_2={\displaystyle {\int }_a^b\left(u^2+\mu {\left(u_x\right)}^2\right)\text{d}x}\approx h{\displaystyle \underset{j=1}{\overset{N}{\sum }}\left[{\left(u_j^n\right)}^2+\mu {\left[{\left(u_x\right)}_j^n\right]}^2\right]},\\ I_3={\displaystyle {\int }_a^b\left(u^4-\mu {\left(u_x\right)}^2\right)\text{d}x}\approx h{\displaystyle \underset{j=1}{\overset{N}{\sum }}\left[{\left(u_j^n\right)}^4-\mu {\left[{\left(u_x\right)}_j^n\right]}^2\right]}.\end{array}$

其中 $I_1,I_2,I_3$ 通常被称为质量，动量和能量。

我们考虑MRLW方程(1.1)，其单波解为

$u\left(x,t\right)=\sqrt{c}\text{sech}\left(p\left(x-\left(c+1\right)t-x_0\right)\right)$ (4.1)

其中 $p=\sqrt{\frac{c}{\mu \left(1+c\right)}}$， $x_0$ 和c是任意的常数。变量的解析值为

$I_1=\frac{\pi \sqrt{c}}{p},I_2=\frac{2c}{p}+\frac{2\mu pc}{3},I_3=\frac{4c^2}{3p}-\frac{2\mu pc}{3}.$

我们讨论单孤立波的运动，并且令初始条件为

$u\left(x,0\right)=\sqrt{c}\text{sech}\left(p\left(x-x_0\right)\right)$ (4.2)

其中 $c=1,x_0=20$。在整个过程中，我们令 $0\le x\le 100,0\le t\le 10$。

表1中，令空间剖分的步长为 $h=0.125$，时间剖分的步长为 $\Delta t=0.0125$，表中列出了三个不变量 $I_1,I_2,I_3$ 在不同时间 $t=2,4,6,8,10$ 的具体数值。从表中的结果可以看出我们的方法对不变量保持的是比较好的。

表2中，令空间剖分步长固定为 $h=0.125$，得到了u在不同时间层 $t=2,4,6,8,10$ 上，时间剖分步长为 $\Delta t=0.025,0.05,0.1$ 时，在 $L^2$ 模下的最优误差结果。从表2中我们容易看到时间的收敛阶为2。同理，q的相关结果在表3中体现。

表4中，我们固定时间剖分步长为 $\Delta t=0.01$，得到了u在不同时间层 $t=2,4,6,8,10$ 上，空间剖分步长为 $h=0.2,0.4,0.8$ 下的 $L^2$ 模最优误差结果。从表中易知，在空间方向上的收敛阶为2。同理，q的相关结果在表5中体现。

Table 1. Invariants of single solitary wave. h = 0 , 125 , Δ t = 0.0125 , c 1 = 1 , c 2 = 0.5

表1. 孤立波的不变量。 $h=0,125,\Delta t=0.0125,c_1=1,c_2=0.5$

Table 2. Convergence order and error in L2-norm for u of time. h = 0.125 , c = 1

表2. u的L²误差和时间收敛阶。 $h=0.125,c=1$

Table 3. Convergence order and error in L2-norm for q of time. h = 0. 125 , c = 1

表3. q的L²误差和时间收敛阶。 $h=0.\text{125},c=\text{1}$

Table 4. Convergence order and error in L2-norm for u of space. Δ t = 0.01 , c = 1

表4. u的L²误差和空间收敛阶。 $\Delta t=0.01,c=1$

Table 5. Convergence order and error in L2-norm for q of space. Δ t = 0.01 , c = 1

表5. q的L²误差和空间收敛阶。 $\Delta t=0.01,c=1$

Figure 1. Surface for exact solution u

图1. 精确解u

Figure 2. Surface for numerical solution u_h

图2. 数值解u_h

图1为精确解u在 $\left(x,t\right)\in \left[0,100\right]\times \left[0,10\right]$ 下的三维图像，图2为 $h=0.2,\Delta t=0.02$ 时，数值解u_h的三维图像。从图1和图2可以看出，数值解u_h很好的近似了精确解u。图3和图4分别为q的精确值与 $h=0.2,\Delta t=0.02$ 下的数值解q_h的三维对比图，我们可以看出数值解q_h很好的近似了精确解q。

Figure 3. Surface for exact solution q

图3. 精确解q

Figure 4. Surface for numerical solution q_h

图4. 数值解q_h

Figure 5. Comparison between u and u_h

图5. u和u_h对比图

Figure 6. Comparison between q and q_h

图6. q和q_h对比图

Figure 7. Numerical solution u_h

图7. 数值解u_h

Figure 8. Numerical solution q_h

图8. 数值解q_h

图5表示在 $h=0.2,\Delta t=0.02$ 时，不同时间层 $t=2,4,6,8,10$ 上数值解u_h和精确值u的对比图。同理，图6为q的相关对比图。图5和图6表示的单孤立波朝右运动时u和q保持波形和速度不变。可以看出数值解u_h和q_h可以很好的近似精确值u和q。

图7和图8表示的是单孤立波在不同时间层 $t=2,4,6,8,10$ 的三维线图。

从上述表1~5和图1~8的分析，我们知道数值结果与定理3.6中的理论结果是比较吻合的，并且也验证了的我们的方法对于MRLW方程的有效性。

5. 结论

所本文中，我们对于MRLW方程提出并分析了一种显式二阶向后差分的混合有限元方法，即在空间方向上利用混合有限元方法，在时间方向上使用二阶向后差分方法。我们讨论了方法的数值过程，对该方法的全离散格式证明了相关的理论结果，并且与其他数值相比较，我们获得了很好的误差收敛阶。从数值算例的结果可以看到，我们的方法对于MRLW方程是比较有效的。

基金项目

本文由内蒙古自然科学基金项目(2016MS0102)资助。