由冯康先生开创的Hamilton系统的辛几何算法发展至今, 在理论上已较为完善.由于这一算法具有长时间的数值稳定性, 能够很好地保持Hamilton系统的辛几何结构的性质, 在现代物理学和力学研究中发挥着重要作用.在实际研究与实践中, 有许多问题需要进行长时间的数值模拟计算, 因此, 对具有长时间数值行为的辛算法的研究具有重要的理论与实际意义.近年来, 孤立波方程定解问题的研究一直是科学研究的一个热门课题, 取得了丰硕的研究成果.在许多物理领域中具有重要意义.本文研究一类具有重要意义的方程(ZK-BBM方程)

$ \begin{equation} u_t+u_x-a(u^2)_x-(bu_{xt}+ku_{yt})_x=0. \end{equation} $

(1.1)

在文献[1]对系统(1.1) 的解的定态问题做了一些研究, 但是这些研究都是给出系统(1.1) 的部分精确解.本文在第二节里验证了ZK-BBM方程具有Hamilton多辛格式, 并证实此格式具有多辛守恒律、局部能量守恒律.第三节给出了ZK-BBM方程的离散多辛Preissmann格式, 并证实此格式在离散格式下忍保持多辛守恒律.在第四节里我们给出了ZK-BBM方程组的离散多辛Preissmann格式的误差分析, 此格式具有误差 $\circ(\Delta t+\Delta x^2+\Delta y^2)$.在第五节给出了一个数值模拟, 验证了本文的算法不仅简单, 而且有长时间的稳定性.

关于多辛的定义, 一切耗散效应可以忽略的方程都可以写成下列哈密顿偏微分方程的形式^[2-8]

$ \begin{equation} Mz_t+Kz_x+Lz_y=\nabla_zS(z), \end{equation} $

(2.1)

其中 $M, K, L\in R^{n\times n}(n\geq 3)$是反对称矩阵, $S:R^n->R$是光滑函数, 称为Hamilton函数. $\nabla_zS(Z)$为函数 $S(z)$的梯度.

方程(2.1) 满足多辛守恒律

$ \begin{eqnarray}&& \frac{\partial}{\partial t}w+\frac{\partial}{\partial x}k+\frac{\partial}{\partial y}l=0, \end{eqnarray} $

(2.2)

其中 $ w=\frac{1}{2}dz\wedge Mdz, k=\frac{1}{2}dz\wedge Kdz, l=\frac{1}{2}dz\wedge Ldz.$

方程(2.1) 具有能量守恒律

$ \begin{eqnarray}&&\frac{\partial }{\partial t} E+\frac{\partial }{\partial x} F+\frac{\partial }{\partial y} G=0, \end{eqnarray} $

(2.3)

其中 $ E=S(z)-\frac{1}{2}z^TKz_x-\frac{1}{2}z^TLz_y, F=\frac{1}{2}z^TKz_t, G=\frac{1}{2}z^TLz_t. $

对系统(1.1), 引入正则动量

$ u_t=-2p_x, u_t=w, bu_x+ku_y=v, u=\varphi_x, $

方程(1.1) 可以变为等价的哈密顿偏微分方程的形式

$ \begin{equation} \displaystyle \left\{ \begin{array}{lcl} -\frac{1}{2}u_t-p_x=0, \\ \frac{1}{2}\varphi_t-\frac{1}{2}v_t-\frac{1}{2}bw_x-\frac{1}{2}kw_y=p-u+au^2, \\ \frac{1}{2}u_t=\frac{1}{2}w, \\ \frac{1}{2}bu_x+\frac{1}{2}ku_y=\frac{1}{2}v, \\ \varphi_x=u. \end{array} \right.\end{equation} $

(2.4)

定义状态变量

$ z=(\varphi, u, v, w, p), $

可以把方程(2.4) 写成哈密顿偏微分方程的形式(2.1), 其中

$ \begin{eqnarray*}&& M=\left[\begin{array}{ccccc} 0&-\frac{1}{2}&0&0&0\\ \frac{1}{2}&0&-\frac{1}{2}&0&0\\ 0&\frac{1}{2}&0&0&0\\ 0&0&0&0&0\\ 0&0&0&0&0 \end{array}\right], K=\left[\begin{array}{ccccc} 0&0&0&0&-1\\ 0&0&0&-\frac{1}{2}b&0\\ 0&0&0&0&0\\ 0&\frac{1}{2}b&0&0&0\\ 1&0&0&0&0 \end{array}\right], \\ &&L=\left[\begin{array}{ccccc} 0&0&0&0&0\\ 0&0&0&-\frac{1}{2}k&0\\ 0&0&0&0&0\\ 0&\frac{1}{2}k&0&0&0\\ 0&0&0&0&0 \end{array}\right], \\ && S(z)=pu-\frac{1}{2}u^2+\frac{a}{3}u^3+\frac{1}{2}vw. \end{eqnarray*} $

方程(2.4) 满足多辛守恒律(2.2), 其中

$ \begin{eqnarray*}&&w=\frac{1}{2}(du\wedge d\varphi+dv\wedge du), \\ &&k=\frac{1}{2}(2dp\wedge d\varphi+bdw\wedge bu), \\ &&l=\frac{1}{2}kdw\wedge du.\end{eqnarray*} $

方程(2.4) 具有能量守恒律为(2.3), 其中

$ \displaystyle \left\{ \begin{array}{lcl} \displaystyle E=pu-\frac{1}{2}u^2+\frac{a}{3}u^3+\frac{1}{2}vw-\frac{1}{2}(pu+\frac{1}{2}bwu_x-\frac{1}{2}buw_x-\varphi p_x)-\frac{1}{4}(kw u_y-uw_y), \\ \displaystyle F=\frac{1}{2}(p\varphi_t +\frac{1}{2}bwu_t-\frac{1}{2}buw_t-\varphi p_t), \\ \displaystyle G=\frac{1}{4}(kw u_t-uw_t). \end{array} \right. $

多辛形式的一个重要性质是:它的局部守恒的概念.多辛是Hamilton偏微分方程的一个几何性质, 我们用数值方法模拟多辛偏微分方程时, 我们自然希望能反映这个性质.基于这个想法, Bridges和Reich引入了多辛积分的概念, 即一种能保持多辛守恒律的离散数值方法.

求解多辛Hamilton系统(2.1) 的Preissmann格式可表示为

$ \begin{equation} MD_tz_i^{j+1/2, l+1/2}+KD_xz_{i+1/2}^{j, l+1/2}+LD_yz_{i+1/2}^{j+1/2, l}=\nabla_zS(z_{i+1/2}^{j+1/2, l+1/2})_{i+1/2}^{j+1/2, l+1/2}, \end{equation} $

(3.1)

其中

$ D_tu_i=\frac{u_{i+1}-u_i}{\Delta t}, ~D_xu^j=\frac{u^{j+1}-u^j}{\Delta x}, ~D_yu^l=\frac{u^{l+1}-u^l}{\Delta y}, $

定义3.1  称数值方法(3.1) 为多辛积分, 如果其满足如下离散多辛守恒律

$ \begin{equation} D_tw_{i}^{j+1/2, l+1/2}+D_xk_{i+1/2}^{j, l+1/2}+D_yl_{i+1/2}^{j+1/2, l}=0, \end{equation} $

(3.2)

其中

$ \displaystyle \left\{\begin{array}{lcl} w_{i}^{j+1/2, l+1/2}=\frac{1}{2}dz_{i}^{j+1/2, l+1/2}\Lambda Mdz_{i}^{j+1/2, l+1/2}, \\ \displaystyle k_{i+1/2}^{j, l+1/2}=\frac{1}{2}dz_{i+1/2}^{j, l+1/2}\Lambda Kdz_{i+1/2}^{j, l+1/2}, \\ \displaystyle l_{i+1/2}^{j+1/2, l}=\frac{1}{2}dz_{i+1/2}^{j+1/2, l}\Lambda Ldz_{i+1/2}^{j+1/2, l}. \end{array} \right. $

在进行数值求解偏微分方程组(2.4) 我们希望构造的数值方法严格满足上述守恒律, 即具有多辛性质.本文采用Preissmann格式对偏微分方程组(2.4) 进行离散

$ \begin{equation} \displaystyle \left\{ \begin{array}{lcl} & \displaystyle-\frac{1}{2}\frac{u_{i+1}^{j+1/2, l+1/2}-u_{i}^{j+1/2, l+1/2}}{\Delta t}-\frac{p_{i+1/2}^{j+1, l+1/2}-p_{i+1/2}^{j, l+1/2}}{\Delta x}=0, \\ & \displaystyle\frac{1}{2}\frac{\varphi_{i+1}^{j+1/2, l+1/2}-\varphi_{i}^{j+1/2, l+1/2}}{\Delta t}-\frac{1}{2}\frac{v_{i+1}^{j+1/2, l+1/2}-v_{i}^{j+1/2, l+1/2}}{\Delta t}\\& - \displaystyle \displaystyle\frac{1}{2}b\frac{w_{i+1/2}^{j+1, l+1/2}-u_{i+1/2}^{j, l+1/2}}{\Delta x}-\frac{1}{2}k\frac{w_{i+1/2}^{j+1/2, l+1}-w_{i+1/2}^{j+1/2, l}}{\Delta y}\\ &= p_{i+1/2}^{j+1/2, l+1/2}-u_{i+1/2}^{j+1/2, l+1/2}+a(u_{i+1/2}^{j+1/2, l+1/2})^2, \\ &\displaystyle \frac{1}{2}\frac{u_{i+1}^{j+1/2, l+1/2}-u_{i}^{j+1/2, l+1/2}}{\Delta t}=\frac{1}{2}w_{i+1/2}^{j+1/2, l+1/2}, \\ &\displaystyle\frac{1}{2}bu_x+\frac{1}{2}ku_y=\frac{1}{2}v_{i+1/2}^{j+1/2, l+1/2}, \\ & \displaystyle\frac{\varphi_{i+1/2}^{j+1, l+1/2}-\varphi_{i+1/2}^{j, l+1/2}}{\Delta x}=u_{i+1/2}^{j+1/2, l+1/2}. \end{array} \right.\end{equation} %\[ %\frac{1}{\Delta t}\left[\begin{array}{ccccc} %\frac{1}{2}(\varphi_{i+1}^{j+1/2, l+1/2}-\varphi_{i}^{j+1/2, l+1/2})\\ %0\\ %0\\ %-\frac{1}{2}(u_{i+1}^{j+1/2, l+1/2}-u_{i}^{j+1/2, l+1/2})\\ %0 % \end{array}\right]+\frac{1}{\Delta x}\left[\begin{array}{ccccc} %\delta(v_{i+1/2}^{j+1, l+1/2}-v_{i+1/2}^{j, l+1/2})\\ %-\delta(u_{i+1/2}^{j+1, l+1/2}-u_{i+1/2}^{j, l+1/2})\\ %\varphi_{i+1/2}^{j+1, l+1/2}-\varphi_{i+1/2}^{j, l+1/2}\\ %0\\ %0 %\end{array}\right]\]\[=\frac{1}{\Delta y}\left[\begin{array}{ccccc} %\gamma(u_{i+1/2}^{j+1/2, l+1}-u_{i+1/2}^{j+1/2, l})+\rho(w_{i+1/2}^{j+1/2, l+1}-w_{i+1/2}^{j+1/2, l})\\ %0\\ %0\\ %0\\ %-\rho(u_{i+1/2}^{j+1/2, l+1}-u_{i+1/2}^{j+1/2, l}) % \end{array}\right]=\left[\begin{array}{ccccc} %\delta(v_{i+1/2}^{j+1, l+1/2}-v_{i+1/2}^{j, l+1/2})\\ %-\delta(u_{i+1/2}^{j+1, l+1/2}-u_{i+1/2}^{j, l+1/2})\\ %\varphi_{i+1/2}^{j+1, l+1/2}-\varphi_{i+1/2}^{j, l+1/2}\\ %0\\ %0 % \end{array}\righ $

(3.3)

格式(3.3) 是多辛格式, 因为具有下面多辛守恒律.

定理3.1  离散格式(3.3) 是多辛格式, 且保持下面的离散多辛守恒律

$ \begin{equation} D_tw_{i}^{j+1/2, l+1/2}+D_xk_{i+1/2}^{j, l+1/2}+D_yl_{i+1/2}^{j+1/2, l}=0，\end{equation} $

(3.4)

其中

$ \begin{equation} \displaystyle \left\{ \begin{array}{lcl} \displaystyle w_{i}^{j+1/2, l+1/2}=du_{i}^{j+1/2, l+1/2}\Lambda d \varphi_{i}^{j+1/2, l+1/2}, \\ \displaystyle k_{i+1/2}^{j, l+1/2}=\delta du_{i+1/2}^{j, l+1/2}\Lambda dv_{i+1/2}^{j, l+1/2}+dp_{i+1/2}^{j, l+1/2} \Lambda d\varphi_{i+1/2}^{j, l+1/2}, \\ \displaystyle l_{i+1/2}^{j+1/2, l}=\rho du_{i+1/2}^{j+1/2, l} \Lambda dw_{i+1/2}^{j+1/2, l}. \end{array} \right. \end{equation} $

(3.5)

证  对方程(3.1) 变分可以得到(3.1) 的变分形式为

$ \begin{eqnarray} && MD_tdz_i^{j+1/2, l+1/2}+KD_xdz_{i+1/2}^{j, l+1/2}+D_ydz_{i+1/2}^{j+1/2, l}\nonumber\\ &=& \nabla_{zz}S(z_{i+1/2}^{j+1/2, l+1/2})_{i+1/2}^{j+1/2, l+1/2}dz_{i+1/2}^{j+1/2, l+1/2}. \end{eqnarray} $

(3.6)

对方程(3.6) 与 $dz_{i+1/2}^{j+1/2, l+1/2}$做外积, 注意到

$ dz_{i+1/2}^{j+1/2, l+1/2}\Lambda\nabla_{zz}S(z_{i+1/2}^{j+1/2, l+1/2})_{i+1/2}^{j+1/2, l+1/2}dz_{i+1/2}^{j+1/2, l+1/2}=0. $

我们可以得到离散守恒律(3.4).

对非线性Hamilton系统, 离散局部能量和动量守恒律不能精确满足, 但是可以定义下面局部能量和动量误差.

定义3.2  记

$ \begin{equation} R_E=D_tE_i^{j+1/2, l+1/2}+D_xF_{1+1/2}^{j, l+1/2}+D_yG_{1+1/2}^{j+1/2, l}, \end{equation} $

(3.7)

称 $R_E$是局部能量守恒律的误差, 其中

$ \begin{eqnarray*}&& E_{i+1/2}^{j+1/2, l+1/2}=p_{i+1/2}^{j+1/2, l+1/2}u_{i+1/2}^{j+1/2, l+1/2}-\frac{1}{2}(u_{i+1/2}^{j+1/2, l+1/2})^2+ \frac{a}{3}(u_{i+1/2}^{j+1/2, l+1/2})^3\\ && +\frac{1}{2}v_{i+1/2}^{j+1/2, l+1/2}w_{i+1/2}^{j+1/2, l+1/2}- \frac{1}{2}(p_{i+1/2}^{j+1/2, l+1/2}u_{i+1/2}^{j+1/2, l+1/2}+\frac{1}{2}bw_{i+1/2}^{j+1/2, l+1/2}D_xu_{i+1/2}^{j, l+1/2}\\ &&-\frac{1}{2}bu_{i+1/2}^{j+1/2, l+1/2} D_xw_{i+1/2}^{j, l+1/2}-\varphi_{i+1/2}^{j+1/2, l+1/2} D_xp_{i+1/2}^{j, l+1/2})\\ &&-\frac{1}{4}(kw_{i+1/2}^{j+1/2, l+1/2} D_yu_{i+1/2}^{j+1/2, l}-u_{i+1/2}^{j+1/2, l+1/2}D_yw_{i+1/2}^{j+1/2, l}), \\ && F_{i+1/2}^{j+1/2, l+1/2}=\frac{1}{2}(p_{i+1/2}^{j+1/2, l+1/2}D_t\varphi_{i}^{j+1/2, l+1/2} +\frac{1}{2}bw_{i+1/2}^{j+1/2, l+1/2}D_tu_{i}^{j+1/2, l+1/2}\\ && -\frac{1}{2}bu_{i+1/2}^{j+1/2, l+1/2}D_tw_{i}^{j+1/2, l+1/2}-\varphi_{i+1/2}^{j+1/2, l+1/2} D_tp_{i}^{j+1/2, l+1/2}), \\ && G_{i+1/2}^{j+1/2, l+1/2}=\frac{1}{4}(kw_{i+1/2}^{j+1/2, l+1/2} D_tu_{i}^{j+1/2, l+1/2}-u_{i+1/2}^{j+1/2, l+1/2}D_tw_{i}^{j+1/2, l+1/2}). \end{eqnarray*} $

假设 $z$是充分光滑函数, 将函数 $z$在离散点 $(t_i, x_j, y_l)$处分别关于 $t, x, y$进行泰勒展开的,

$ \begin{equation} \displaystyle \left\{ \begin{array}{lcl}\displaystyle z_{i+1}^{j+1/2, \ l+1/2}=z_{i}^{j+1/2, \ l+1/2}+\Delta t (z_t)_{i+1/2}^{j+1/2, \ l+1/2}+\frac{1}{2}(\Delta t)^2 (z_{tt})_{i+1/2}^{j+1/2, \ l+1/2}+\cdots, \\ \displaystyle z_{i+1/2}^{j+1, \ l+1/2}=z_{i+1/2}^{j, \ l+1/2}+\Delta x (z_x)_{i+1/2}^{j+1/2, \ l+1/2}+\frac{1}{2}(\Delta x)^2 (z_{xx})_{i+1/2}^{j+1/2, \ l+1/2}+\cdots, \\ \displaystyle z_{i+1/2}^{j+1/2, \ l+1}=z_{i+1/2}^{j+1/2, \ l}+\Delta y (z_y)_{i+1/2}^{j+1/2, \ l+1/2}+\frac{1}{2}(\Delta y)^2 (z_{yy})_{i+1/2}^{j+1/2, \ l+1/2}+\cdots. \end{array} \right. \end{equation} $

(4.1)

把上面的形式可以写成下面的等价形式

$ \begin{equation} \displaystyle \left\{ \begin{array}{lcl} \displaystyle\frac{z_{i+1}^{j+1/2, \ l+1/2}-z_{i}^{j+1/2, \ l+1/2}}{\Delta t}= (z_t)_{i+1/2}^{j+1/2, \ l+1/2}+\frac{1}{2}(\Delta t) (z_{tt})_{i+1/2}^{j+1/2, \ l+1/2}+\circ(\Delta t), \\ \displaystyle \frac{z_{i+1/2}^{j+1, \ l+1/2}-z_{i+1/2}^{j, \ l+1/2}}{\Delta x}= (z_x)_{i+1/2}^{j+1/2, \ l+1/2}+\frac{1}{2}(\Delta x) (z_{xx})_{i+1/2}^{j+1/2, \ l+1/2}+\circ(\Delta x), \\ \displaystyle \frac{z_{i+1/2}^{j+1/2, \ l+1}-z_{i+1/2}^{j+1/2, \ l}}{\Delta y }=(z_y)_{i+1/2}^{j+1/2, \ l+1/2}+\frac{1}{2}(\Delta y) (z_{yy})_{i+1/2}^{j+1/2, \ l+1/2}+\circ(\Delta y). \end{array} \right. \end{equation} $

(4.2)

为了误差, 记 $z_t=(z_t)_{i+1/2}^{j+1/2, \ l+1/2}$, 则离散格式(3.3) 可以改写成

$ \begin{equation} \displaystyle \left\{\begin{array}{lcl} -\frac{1}{2}(u_t+\frac{1}{2}u_{tt}\Delta t)-(p_x+\frac{1}{2}p_{xx}\Delta x)=0, \\ \frac{1}{2}(\varphi_t+\frac{1}{2}\varphi_{tt}\Delta t)-\frac{1}{2}(v_t+\frac{1}{2}v_{tt}\Delta t)\\ -\frac{1}{2}b(w_x+\frac{1}{2} w_{xx}\Delta x)-\frac{1}{2}k(w_y+\frac{1}{2}w_{yy}\Delta y)=p-u-au^2, \\ \frac{1}{2}(u_t+\frac{1}{2}u_{tt}\Delta t)=\frac{1}{2}w, \\ \frac{1}{2}b(u_x+\frac{1}{2}u_{xx}\Delta x)+\frac{1}{2}k(u_y+\frac{1}{2}u_{yy}\Delta y)=\frac{1}{2}v, \\ (\varphi_x+\frac{1}{2}\varphi_{xx}\Delta x)=u, \end{array} \right.\end{equation} $

(4.3)

把方程组(4.3) 的第三个公式和第四个公式代入第二个公式得

$ \begin{eqnarray} &&\frac{1}{2}(\varphi_t+\frac{1}{2}\varphi_{tt}\Delta t)-\frac{1}{2}((b(u_{xt}+\frac{1}{2}u_{xxt}\Delta x)+k(u_{yt}+\frac{1}{2}u_{yyt}\Delta y))\nonumber\\&&+ \frac{1}{2}(b(u_{xtt}+\frac{1}{2}u_{xxtt}\Delta x)+k(u_{ytt}+\frac{1}{2}u_{yytt}\Delta y))\Delta t)\nonumber\\&&- \frac{1}{2}b((u_{xt}+\frac{1}{2}u_{xtt}\Delta t)+\frac{1}{2}(u_{xxt}+\frac{1}{2}u_{xxtt}\Delta t)\Delta x)\nonumber\\&&- \frac{1}{2}k((u_{yt}+\frac{1}{2}u_{ytt}\Delta t)+\frac{1}{2}(u_{yyt}+\frac{1}{2}u_{yytt}\Delta t)\Delta y)=p-u-au^2, \end{eqnarray} $

(4.4)

对(4.4) 关于 $x$求导数可以得到

$ \begin{eqnarray} p_x&=&\frac{1}{2}(\varphi_{xt}+\frac{1}{2}\varphi_{xtt}\Delta t)-\frac{1}{2}((b(u_{xxt}+\frac{1}{2}u_{xxxt}\Delta x)+k(u_{xyt}+\frac{1}{2}u_{xyyt}\Delta y))\nonumber\\&&+ \frac{1}{2}(b(u_{xxtt}+\frac{1}{2}u_{xxxtt}\Delta x)+k(u_{xytt}+\frac{1}{2}u_{xyytt}\Delta y))\Delta t)\nonumber\\&&- \frac{1}{2}b((u_{xxt}+\frac{1}{2}u_{xxtt}\Delta t)+\frac{1}{2}(u_{xxt}+\frac{1}{2}u_{xxxtt}\Delta t)\Delta x)\nonumber\\&&- \frac{1}{2}k((u_{xyt}+\frac{1}{2}u_{xytt}\Delta t)+\frac{1}{2}(u_{xyyt}+\frac{1}{2}u_{xyytt}\Delta t)\Delta y)+u_x+2auu_x, \end{eqnarray} $

(4.5)

把(4.5) 代入(4.3) 的第一个公式得

$ \begin{eqnarray} &&-\frac{1}{2}(u_t+\frac{1}{2}u_{tt}\Delta t)-(\frac{1}{2}(\varphi_{xt}+\frac{1}{2}\varphi_{xtt}\Delta t)\\&&-\frac{1}{2}((b(u_{xxt}+\frac{1}{2}u_{xxxt}\Delta x)+k(u_{xyt}+\frac{1}{2}u_{xyyt}\Delta y))\nonumber\\&&+ \frac{1}{2}(b(u_{xxtt}+\frac{1}{2}u_{xxxtt}\Delta x)+k(u_{xytt}+\frac{1}{2}u_{xyytt}\Delta y))\Delta t)\nonumber\\&&- \frac{1}{2}b((u_{xxt}+\frac{1}{2}u_{xxtt}\Delta t)+\frac{1}{2}(u_{xxt}+\frac{1}{2}u_{xxxtt}\Delta t)\Delta x)\nonumber\\&&- \frac{1}{2}k((u_{xyt}+\frac{1}{2}u_{xytt}\Delta t)+\frac{1}{2}(u_{xyyt}+\frac{1}{2}u_{xyytt}\Delta t)\Delta y)\\&&+u_x+2auu_x+\frac{1}{2}p_{xx}\Delta x)=0, \end{eqnarray} $

(4.6)

化简可得

$ \begin{eqnarray} &&u_t+u_x-a(u^2)_x-(bu_{xt}+ku_{yt})_x\nonumber\\&=& -\frac{1}{4}u_{tt}\Delta t-(\frac{1}{4}\varphi_{xtt}\Delta t-\frac{1}{2}(\frac{1}{2}u_{xxxt}\Delta x+k\frac{1}{2}u_{xyyt}\Delta y)\nonumber\\&&+ \frac{1}{2}(b(u_{xxtt}+\frac{1}{2}u_{xxxtt}\Delta x)+k(u_{xytt}+\frac{1}{2}u_{xyytt}\Delta y))\Delta t)\nonumber\\&&- \frac{1}{2}b((u_{xxt}+\frac{1}{2}u_{xxtt}\Delta t)+\frac{1}{2}(u_{xxt}+\frac{1}{2}u_{xxxtt}\Delta t)\Delta x)\nonumber\\&&- \frac{1}{2}k((u_{xyt}+\frac{1}{2}u_{xytt}\Delta t)+\frac{1}{2}(u_{xyyt}+\frac{1}{2}u_{xyytt}\Delta t)\Delta y)+\frac{1}{2}p_{xx}\Delta x). \end{eqnarray} $

(4.7)

所以本文的中心Preissmann格式具有精度 $\circ(\Delta t+\Delta x^2+\Delta y^2)$.

虽然中心Preissmann格式是一个多辛格式, 但在实际数值计算中需要计算辅助变量 $\varphi, v, w, p$, 计算量大大增加, 因此我们对式(3.3) 消去中间变量并整理, 可得一个等价于中心Preissmann格式的新的多辛格式

$ \begin{align} & \frac{u_{i+1}^{j+\frac{1}{2}, l+\frac{1}{2}}-u_{i}^{j+\frac{1}{2}, l+\frac{1}{2}}}{\Delta t} \\ & -b\frac{u_{i+1}^{j+\frac{3}{2}, l+\frac{1}{2}}-2u_{i+1}^{j+\frac{1}{2}, l+\frac{1}{2}}+u_{i+1}^{j-\frac{1}{2}, l+\frac{1}{2}}-u_{i}^{j+\frac{3}{2}, l+\frac{1}{2}}+2u_{i}^{j+\frac{1}{2}, l+\frac{1}{2}}-u_{i}^{j-\frac{1}{2}, l+\frac{1}{2}}}{\Delta {{x}^{2}}\Delta t} \\ & -k\frac{u_{i+1}^{j+1, l+1}-u_{i+1}^{j, l+1}-u_{i+1}^{j+1, l}+u_{i+1}^{j, l}-u_{i}^{j+1, l+1}+u_{i}^{j, l+1}+u_{i}^{j+1, l}-u_{i}^{j, l}}{\Delta x\Delta y\Delta t} \\ & +\frac{u_{i+\frac{1}{2}}^{j+1, l+\frac{1}{2}}-u_{i+\frac{1}{2}}^{j, l+\frac{1}{2}}}{\Delta x} \\ & -a\frac{{{(u_{i+\frac{1}{2}}^{j+1, l+\frac{1}{2}})}^{2}}-{{(u_{i+\frac{1}{2}}^{j, l+\frac{1}{2}})}^{2}}}{\Delta x}=0, \\ \end{align} $

其中

$ \begin{eqnarray*}&& u_{i+1}^{j+\frac{1}{2}, l+\frac{1}{2}}=\frac{u_{i+1}^{j+1, l+1}+u_{i+1}^{j, l+1}+u_{i+1}^{j+1, l}+u_{i+1}^{j, l}}{4}, \\ && u_{i}^{j+\frac{1}{2}, l+\frac{1}{2}}=\frac{u_{i}^{j+1, l+1}+u_{i}^{j, l+1}+u_{i}^{j+1, l}+u_{i}^{j, l}}{4}, \\ && u_{i+1}^{j+\frac{3}{2}, l+\frac{1}{2}}=\frac{u_{i+1}^{j+2, l+1}+u_{i+1}^{j+1, l+1}+u_{i+1}^{j+2, l}+u_{i+1}^{j+1, l}}{4}, \\ && u_{i+1}^{j-\frac{1}{2}, l+\frac{1}{2}}=\frac{u_{i+1}^{j, l+1}+u_{i+1}^{j-1, l+1}+u_{i+1}^{j, l}+u_{i+1}^{j-1, l}}{4}, \\ && u_{i}^{j+\frac{3}{2}, l+\frac{1}{2}}=\frac{u_{i}^{j+2, l+1}+u_{i}^{j+1, l+1}+u_{i}^{j+2, l}+u_{i}^{j+1, l}}{4}, \\ && u_{i}^{j-\frac{1}{2}, l+\frac{1}{2}}=\frac{u_{i}^{j, l+1}+u_{i}^{j-1, l+1}+u_{i}^{j, l}+u_{i}^{j-1, l}}{4}, \\ && u_{i+\frac{1}{2}}^{j+1, l+\frac{1}{2}}=\frac{u_{i+1}^{j+1, l+1}+u_{i}^{j+1, l+1}+u_{i+1}^{j+1, l}+u_{i}^{j+1, l}}{4}, \\ && u_{i+\frac{1}{2}}^{j, l+\frac{1}{2}}=\frac{u_{i+1}^{j, l+1}+u_{i}^{j, l+1}+u_{i+1}^{j, l}+u_{i}^{j, l}}{4}. \end{eqnarray*} $

为了说明多辛方法的诸多优点, 本文用Preissmann多辛格式, 并以此离散ZK-BBM方程(1.1).

考虑下面ZK-BBM方程初值问题：

$ \begin{equation} \displaystyle \left\{ \begin{array}{lcl} u_t+u_x-a(u^2)_x-(bu_{xt}+ku_{yt})_x=0, \\ u_0=f(x, y). \end{array} \right.\end{equation} $

(5.1)

情形1  当取定初值函数

$ f(x, y)=\frac{3}{4a}(\sec^2(\frac{1}{2}\sqrt{\frac{1}{b+k}}(x+y)+2)), $

则初值问题(5.1) 具有如下孤子解

$ \begin{equation} u(x, y, t)=\frac{3}{4a}(\sec^2(\frac{1}{2}\sqrt{\frac{1}{b+k}}(x+y-\frac{1}{2}t)+2)), b+k>0, \end{equation} $

(5.2)

取时间步长 $\Delta t=0.01, $空间步长 $\Delta x=0.01, \Delta y=0.01$, 利用Preissmann多辛格式, 在区间 $x\in[-20,20], y\in[-20,20]$内模拟孤子解(5.2), 得到ZK-BBM方程(1.1) 的数值解.数值解(5.2) 的演化过程如图 1, 2, 同时图 3给出了 $t\in[-20,20]$时段内的局部能量误差和局部动量误差.

从以上数值结果我们发现, 利用本文构造的多辛中心Preissmann多辛格式模拟的孤子解(5.2), 得到的波形和波速都不随时间变化而变化, 这说明多辛格式能够很好的保持孤子解的基本几何性质; 并具有良好的长时间数值行为.