2001年, Dullin, Gottwald和Holm从Euler方程出发, 考虑浅水层自由表面的无旋不可压缩无粘水波运动时, 得到了一类带线性和非线性色散项的新型浅水波方程, 即DGH方程[1]
其中$\omega, \gamma, \alpha$是常数. DGH方程被提出以后引起了广泛的关注.方程(1.1)具有一类奇异的孤立波解即尖峰孤立波$u(x, t)=c{\rm e}^{-|x-ct|}$, 尖峰孤立波在许多模型中相继被发现. DGH方程包含了两类相对独立的可积孤立波方程.当$\alpha\rightarrow 0$时, DGH方程变为KdV方程
KdV方程(1.2)具有光滑孤立波解
当$\gamma=0$时, DGH方程转化为另一类重要的可积方程, 即C-H方程
近几十年, 许多学者对方程(1.1)-(1.3)已经做了许多研究, 其中, Camassa和Holm在研究浅水波运动规律时, 推导出一类非线性水波方程(1.3).对任意的$ \omega $, C-H方程具有一个Lax对, 具有双哈密顿结构和无穷个守恒量[2-3].对$\omega=0 $, C-H方程有尖峰孤立波解, 更进一步研究表明, C-H方程具有多重尖峰孤波解.在文献[4-5]中, 作者对方程(1.3)进行了数值模拟, 且得到方程(1.3)的守恒性质.在文献[6]中, 作者对方程(1.3)进行了对称性和可积性研究.在文献[7]中, 作者对方程(1.3)的孤立子解和双孤立子解进行了研究.在文献[8]中, 作者对方程(1.1)的哈密尔顿结构和解的整体存在性及Blow-up现象进行了研究. Tian, Gui和Liu研究了DGH方程的Cauchy问题的局部适定性理论, 整体适定性理论, 解的极限行为, 孤立波的轨道稳定性, 散射理论, 新型尖峰孤立波解, 同时给出了DGH方程的散射数据.郭柏林和刘正荣通过使用平面的自治系统的定性分析研究方法, 研究了DGH方程的尖峰孤立波解.
由于问题(1.1)的非线性, 在实际应用中要求得方程(1.1)的初边值问题的精确解几乎不可能, 大部分情况下只能用数值方法来模拟方程(1.1), 但是对方程(1.1)的数值模拟的研究却很少.在文献[9]中, 给出了方程(1.1)数值方法的一些理论研究, 但是没有给出具体的数值方法.鉴于此, 本文的目的就是用数值方法来模拟方程(1.1)的初边值问题.如果采用一般的数值方法(如差分法、有限元法、谱方法), 这些方法都不是保结构的算法, 都会人为代入一定的耗散性, 都会破坏问题本来的一些物理性质, 在进行长期数值模拟时都会失真.由冯康先生开创的Hamilton系统的辛几何算法[10]在理论和应用上已日趋成熟并得到了蓬勃发展.现在, 一定程度上保持原问题的能力是评价一个数值算法优劣的标准.但用此方法对偏微分方程进行辛离散时具有局限性, 具体表现在此方法是个全局性的概念, 依赖于边界条件.为克服此局限性, Bridges和Reich引入了一个保持守恒型偏微分方程局部多辛结构的多辛积分的概念[11-13].本文的目的就是利用多辛算法来模拟方程(1.1), 这种算法可以进行长时间的数值模拟.由于这种算法不会破坏本身的物理性质, 所以本文的算法可以模拟孤立波解.
本文利用多辛Fourier拟谱方法对满足周期边界条件的DGH方程(1.1)进行了数值模拟.在第二节里通过引入正则动量, 验证了DGH方程(1.1)具有多辛哈密尔顿结构, 并证实此格式具有多辛守恒律、局部能量守恒和动量守恒.第三节给出了DGH方程(1.1)的离散多辛Fourier拟谱方法, 并证实此格式在离散格式下仍保持多辛守恒律.在第四节给出了数值模拟, 验证格式有长时间的稳定性.
大量的偏微分方程可以表示为下列多辛哈密尔顿偏微分方程的形式[14-18]
其中$M, K\in {\Bbb R}^{n\times n}(n\geq 3)$是反对称矩阵.$S:{\Bbb R}^n\rightarrow {\Bbb R}$是光滑函数, 称为哈密尔顿函数, $\nabla_zS(Z)$为函数$S(z)$的梯度. Bridges和Reich证明多辛哈密尔顿系统具有三个守恒律, 即多辛守恒律、局部能量守恒律、局部动量守恒律.
定理2.1 多辛哈密尔顿系统(2.1)满足多辛守恒律
其中
分别对应于$t$和$x$方向上的辛结构.
定理2.2 多辛哈密尔顿系统(2.1)具有局部能量守恒律
局部动量守恒律
$E$为能量密度; $F$为能量流; $I$为动量密度; $G$为动量流.
若引入矩阵分解
可以得到多辛守恒律、局部能量守恒律、局部动量守恒律的另一种表示, 其中
与辛结构相比, 多辛结构具有的守恒律, 在任意时空都成立, 不依赖于边界条件.大量的守恒方程都具有多辛结构, 例如KdV方程、KP方程、薛定谔方程等.下面我们建立DGH方程(1.1)的多辛结构.
引入正则动量
系统(1.1)可以表示为下面一阶方程组
定义状态变量
系统(2.5)可以表示为多辛哈密尔顿偏微分方程的形式(2.1), 其中
Hamilton函数
系统(2.5)满足多辛守恒律(2.2), 其中
系统(2.5)具有能量守恒律(2.3), 其中
系统(2.5)具有动量守恒律(2.4), 其中
多辛是哈密尔顿偏微分方程的一个几何特征, 用数值方法模拟多辛哈密尔顿偏微分方程时, 自然希望能保持连续系统的这一特征.基于这个想法, 称一种能保持多辛守恒律的离散数值方法为多辛算法.
利用某种数值方法对多辛哈密尔顿系统(2.1)进行离散, 可以得到下面离散格式
其中$z_j^n$为离散点$x=x_j, t=t_n$处$z$的近似值, $ \partial_t^{j, n}, \partial_x^{j, n}$分别表示微分算子$\partial_t, \partial_x$在这种数值方法下离散点$x=x_j, t=t_n$处的离散值.
定义3.1 若一个数值方法满足离散多辛守恒律
称这样的算法为多辛算法.
Bridges和Reich结合传统谱方法在Fourier空间上提出了多辛Fourier谱离散格式, 基于这种思想, 陈景波和秦孟兆构造了发展方程的多辛Fourier拟谱格式.下面我们建立方程(1.1)的多辛Fourier拟谱方法.设$u(x, t)$是空间$[0, L]$上的周期函数, 即$u(0, t)=u(L, t)$.用$I_Nu(x, t)$表示为$u(x, t)$ 的插值函数, 配置点为$x_j=h\cdot j, j=0, 1, 2, \cdots N-1$, 其中$h =\frac{L}{N}$, N为偶数, 具体表达式为
其中$u_j=u(x_j, t)$, $g_j(x_k)=\delta_j^k$.$g_j(x)$是$\frac{N}{2}$次正交三角多项式
其中$C_l=1(|l|\neq \frac{N}{2}), C_{\frac{N}{2}}=C_{-\frac{N}{2}}=2, \mu=\frac{2\pi}{L}$.通过直接计算, 得到
为了建立DGH方程(1.1)的Fourier拟谱多辛算法, 引入双线性定义
应用插值公式(3.3)得
Fourier拟谱方法的关键是计算$\frac{{\rm d}^k}{{\rm d}x^k}I_N(x, t)$ 在$x=x_i$处的导数值, 对(3.5)式直接微分得到
其中$D_k$为谱微分矩阵
定义3.2 为了进行数值分析, 引入向量的Hadamard内积
其中$u=(u_0, u_1, u_2, \cdots u_{N-1})^T$, $v=(v_0, v_1, v_2, \cdots v_{N-1})^T$.
定理3.1 为了进行数值分析, 引入向量的Hadamard内积微分形式为
其中$u=(u_0, u_1, u_2, \cdots u_{N-1})^T$和$v=(v_0, v_1, v_2, \cdots v_{N-1})^T$, $u_i=u(x_i, t), v_i=v(x_i, t), $$i=0, 1, 2 \cdots N-1.$
利用谱微分矩阵, 对$x$方向进行离散, 得到系统(2.5)的半离散系统
定理3.2 系统(2.5)的多辛Fourier拟谱半离散格式(3.7)具有$N$个半离散多辛守恒律
和整体辛守恒律
用隐式中点辛格式对半离散格式(3.7)时间方向进行离散, 得到系统(2.5)的全离散系统
定理3.3 系统(2.5)的多辛Fourier拟谱全离散格式(3.9)具有$N$个全离散多辛守恒律
和整体辛守恒
多辛Fourier拟谱方法不能满足离散局部能量守恒律, 为此我们考虑守恒律误差, 记
$(R_E)_i^j$称为局部能量守恒律误差.
在周期边界条件下, 我们可以定义离散能量为
定理3.4 利用文献[19]类似的方法, 我们可以得到局部能量守恒律误差(3.12)和能量误差(3.13)满足
利用文献[19]类似的方法, 我们可得到如果$\Delta t$足够小, 则存在两个不依赖于$\Delta t$, $\Delta x$的两个常数$C_1, C_2, $且有
为了说明多辛算法的诸多优点, 本文用多辛Fourier拟谱方法离散DGH方程(1.1)对应的多辛Hamilton方程(2.5), 消去辅助变量$\phi, w, \Phi, \psi$, 得到DGH方程(1.1)的多辛Fourier拟谱三层格式为
在利用上面三层格式进行数值模拟时, 可用下面二层格式计算上面三层格式第二层上的初值
其中, $U=[u_0, u_1, u_2, \cdots u_{N-1}]^T, U^{\frac{1}{2}}=\frac{U^j+U^{j+1}}{2}.$
取参数$\alpha=1, \omega=1, \gamma=1$, 考虑下面DGH方程(1.1)的孤立波的初边值问题(周期边界条件)
由文献[20], 我们可以得到问题(4.1)有孤立波解
考虑时间步长$\Delta t=0.001, $空间步长$\Delta x=0.01$.利用多辛Fourier拟谱方法对DGH方程(1.1)的初边值问题(周期边界条件) (4.1)在区间$x\in[0, 2\pi], t\in [0, 100]$内模拟时, 得到DGH方程(1.1)的初边值问题(4.1)的孤立波的数值解.图 1给出了DGH方程的初边值问题(4.1)随时间的演化图, 表 1给出了精确解、多辛Fourier拟谱格式、有限差分法的比较.
取参数$\alpha=1, \omega=1, \gamma=1$, 考虑下面DGH方程(1.1)的尖峰孤立波的初边值问题(周期边界条件)
考虑时间步长$\Delta t=0.001, $空间步长$\Delta x=0.01$.利用多辛Fourier拟谱方法对DGH方程(1.1)的初边值问题(周期边界条件) (4.2)在区间$x\in[-10, 10], t\in[0, 100]$内模拟时, 得到DGH方程(4.2)的孤立波的数值解.图 2给出了DGH方程的初边值问题(4.2)随时间的演化图, 图 3给出了离散能量守恒律误差.