数学物理学报, 2026, 46(5): 1857-1883

带 Delta 势的非线性 Schrödinger 方程的局部保结构算法

汪佳玲,*, 邹淑慧

南京信息工程大学数学与统计学院 南京 210044; 南京信息工程大学江苏省应用数学中心 南京 210044; 南京信息工程大学江苏省系统建模与数据分析国际合作联合实验室 南京 210044

Local Structure-Preserving Algorithms for the Nonlinear Schrödinger Equation with Delta Potential

Wang Jialing,*, Zou Shuhui

School of Mathematics and Statistics, Nanjing University of Information Science and Technology, Nanjing 210044; Center for Applied Mathematics of Jiangsu Province, Nanjing University of Information Science and Technology, Nanjing 210044; Jiangsu International Joint Laboratory of System Modeling and Data Analysis, Nanjing University of Information Science and Technology, Nanjing 210044

通讯作者: * 汪佳玲, E-mail:wjl19900724@126.com

收稿日期: 2025-04-1   修回日期: 2025-10-21  

基金资助: 江苏省本科高校 "高质量数理类课程教材改革研究" 专项课题(2025JYSLKJ027)
国家自然科学基金(11801277)
国家自然科学基金(12426524)

Received: 2025-04-1   Revised: 2025-10-21  

Fund supported: Special Project on High-Quality Math & Physics Textbook Reform of Jiangsu Undergraduate Universities(2025JYSLKJ027)
NSFC(11801277)
NSFC(12426524)

摘要

该文聚焦于一类特殊的带 Delta 势的非线性 Schrödinger 方程. 作者基于弱多辛形式, 通过复合构造方法系统地给出了弱多辛意义下的局部保结构算法构造的统一框架. 作者成功开发了四类多辛算法和两类局部能量守恒算法, 并讨论了它们在弱多辛意义下的局部和全局守恒律. 这些算法不依赖于边界条件, 可适用于任何满足弱多辛形式的偏微分方程. 数值实验结果表明所提算法优异的计算性能. 相关理论框架与数值算例可直接用于偏微分方程数值解方向的研究生课程教学, 为计算数学专业人才培养提供典型教学案例.

关键词： 带 Delta 势的非线性 Schrödinger 方程; 弱多辛形式; 多辛算法; 局部能量守恒算法

Abstract

This paper focuses on a special class of nonlinear Schrödinger equation with Delta potential. Based on the weak multisymplectic form, we systematically establish a unified framework for constructing the locall structure-preserving algorithms in the weak multisymplectic sense using the concatenating method. We successfully develop four multi-symplectic algorithms and two local energy-preserving algorithms, whose local and global conservation laws in the weak sense are also discussed afterwards. These algorithms are independent of boundary conditions and can be applied to any partial differential equations satisfying the weak multi-symplectic form. Numerical experiments demonstrate the excellent computational performance of the proposed algorithm. The relevant theoretical framework and numerical examples can be directly applied to postgraduate teaching courses on numerical solutions of partial differential equations, providing typical teaching cases for talent cultivation in computational mathematics.

Keywords： nonlinear Schrödinger equation with Delta potential; weak multi-symplectic form; multi-symplectic algorithm; local energy-preserving algorithm

PDF (1289KB) 元数据 多维度评价 相关文章 导出 EndNote| Ris| Bibtex  收藏本文

本文引用格式

汪佳玲, 邹淑慧. 带 Delta 势的非线性 Schrödinger 方程的局部保结构算法[J]. 数学物理学报, 2026, 46(5): 1857-1883

Wang Jialing, Zou Shuhui. Local Structure-Preserving Algorithms for the Nonlinear Schrödinger Equation with Delta Potential[J]. Acta Mathematica Scientia, 2026, 46(5): 1857-1883

1 引言

偏微分方程 (PDEs) 是描述自然界中各种物理现象的基本工具, 广泛应用于流体力学 [1,2], 量子力学[3], 材料科学[4]等领域. 守恒型 PDEs 可以写成如下多辛形式

$\begin{equation}\label{multi-symplectic} \mathrm{\textbf{M}}z_t+\textbf{K}z_x=\nabla_zS(z),\quad z \in \mathbb{R}^d,\quad (x,t) \in \mathbb{R}^2, \end{equation} $

其中, $ \textbf{M} $ 和 $ \textbf{K} $ 为两个反对称矩阵, $ z \in \mathbb{R}^d $ 为状态变量, $ S:\mathbb{R}^d \rightarrow \mathbb{R} $ 为光滑的标量函数. 该系统具备以下三种局部守恒律

$ \bullet $ 多辛守恒律

$\begin{equation} \label{MSCL} \partial_t\omega + \partial_x\kappa = 0, \quad \text{其中} \quad \omega = \mathrm{d}z \wedge \textbf{{M}}_+ \mathrm{d}z, \quad \kappa = \mathrm{d}z \wedge \textbf{{K}}_+ \mathrm{d}z, \end{equation}$

$ \bullet $ 局部能量守恒律

$\begin{equation} \label{LECL} \partial_tE+\partial_xF=0,\quad {\text{其中}}\quad E=S(z)-z^{\rm T}{\textbf{{K}}_+}z_x, \quad F=z^{\rm T}{\textbf{{K}}_+}z_t, \end{equation} $

$ \bullet $ 局部动量守恒律

$\begin{equation}\label{LMCL} \partial_tI+\partial_xG=0,\quad {\text{其中}}\quad I=z^{\rm T}{\textbf{{M}}_+}z_x, \quad G=S(z)-z^{\rm T}{\textbf{{M}}_+}z_t, \end{equation}$

其中, $ \textbf{{M}}_+ $ 和 $ \textbf{{K}}_+ $ 分别由 $ \textbf{{M}} $ 和 $ \textbf{{K}} $ 按如下方式分裂

$\begin{equation} \notag \textbf{{M}}=\textbf{{M}}_+ + \textbf{{M}}_-,\quad \textbf{{M}}_-=-\textbf{{M}}_+^ \mathrm{T},\quad \textbf{{K}}=\textbf{{K}}_+ + \textbf{{K}}_-,\quad \textbf{{K}}_-=-\textbf{{K}}_+^ \mathrm{T}. \end{equation}$

此外, 对于某类特殊的 PDEs, 例如非线性 Schrödinger (NLS) 方程, 还具备局部质量守恒律

$\begin{equation} \label{LNCL} \partial_t N_d +\partial_x N_f=0,\quad {\text{其中}}\quad N_d = \frac14 (\textbf{{M}}z)^\mathrm{T}\textbf{{M}}z, \quad N_f = \frac12 (\textbf{{M}}z)^\mathrm{T}\textbf{{K}}z. \end{equation} $

在适当的边界条件, 例如, 周期边界条件或齐次 Dirichlet 边界条件下, 上述局部守恒律可沿空间方向积分, 分别导出辛守恒律, 全局能量守恒律, 全局动量守恒律以及全局质量守恒律. 对于 PDEs 的数值求解而言, 保持其物理属性与数学结构的完整性至关重要[5]. 多辛框架[6-8]在该方面展现出显著优势: 相较于传统保结构算法, 局部保结构算法不仅能够同等看待时间与空间变量, 还摆脱了对边界条件的依赖. 这使得数值算法能在每个点上保持结构, 且无需预设方程是否具有合适的边界条件. 这些特性能够有效地维持动力系统中的守恒律, 防止能量, 动量等物理量的虚假耗散[9-13], 从而更准确地刻画系统的真实动力学行为.

非线性 Schrödinger 方程

$\begin{equation}\label{NLSEs} i \hbar \partial_t u = -\frac{\hbar ^2}{2m}\Delta u+a|u|^2 u+V_{ext} u, \end{equation}$

是现代科学领域中最重要的非线性模型之一. 其中, $ i $ 是虚数单位, $ \hbar $ 为约化普朗克常数, $ u=u(x,t) $ 为时空波函数, $ m $ 为粒子质量, $ a $ 为源于 Gross-Pitaevskii 理论的耦合常数, $ V_{ext}\equiv V_{ext}(\bar{r},t) $ 为外势场, $ \Delta $ 为拉普拉斯算子. NLS 方程已被广泛研究, 其应用涵盖非线性量子场论, 凝聚态物理与等离子体物理, 非线性光学[14]等诸多学科领域, 对全球科技与经济的发展起到了重要推动作用.

在外势场 $ V_{ext} = 0 $ 的情形下, 对该方程理论解的研究在过去数十年间受到广泛关注[15–17]. 目前, 已采用有限差分法 [18,19], 有限元法 [20], 傅里叶拟谱法 [21]等多种数值方法对其进行数值求解. 随着计算科学的不断发展, 保结构算法 [22–25]的构造已成为研究热点. 随着多辛形式[8,26]的提出, 局部保结构算法应运而生 [22], 其将空间方向与时间方向置于同等地位的特点, 为 NLS 方程的数值求解提供了全新视角. 基于此, 众多学者致力于利用多辛形式 (1.1)为 NLS 方程构造局部保结构算法. 例如, Cai[23] 提出了一种基于多辛形式的算法, 为耦合 NLS 方程构造了局部能量与局部动量守恒格式; Wang[24] 以 NLS 方程为例, 基于多辛形式发展了一系列局部保结构算法构造的统一框架. 这些算法不仅针对特定方程有效, 更为 PDEs 的数值求解提供了统一而高效的途径.

当外势场 $ V_{ext} \neq 0 $ 时, 方程 (1.6) 呈现出独特的数学性质, 特别地, 当外势场取为 $ V_{ext} = \gamma \delta(x) $, 其中 \(\delta(x)\) 为支集位于 \(x=0\) 的 Delta 函数时, 该点势阱模型构成对理想量子阱的数学抽象, 在数学物理领域具有深远意义. 借助这一简化模型, 可更清晰地揭示量子束缚态的内在特性, 为极端条件下粒子行为的深入研究以及量子力学基本原理的探讨提供关键的理论价值 [27]. 对于带 Delta 函数的 PDEs 研究, 研究者已提出多种策略. 具体而言, 通过将 Delta 函数视为单点间断或有限支集函数, 文献[28] 构造了高精度的有限差分格式; 文献[29]发展了保守小波配置方法; Zhou[30] 则基于区域分解技术, 为含 Delta 势的 NLS 方程设计了高效的有限差分-Chebyshev 配置耦合算法. 此外, 亦可将该问题视为界面问题处理, 相关研究见文献[31,32]. 尽管 Delta 函数的引入使得解在 Delta 势所在位置上出现尖点型奇异性, 文献[33] 提出的函数变换方法可将原方程转化为含间断漂移项的新形式, 从而提升解的正则性. Bai[34] 针对一类含 Delta 势的一维定态 Schrödinger 方程, 进一步提出了显式跳跃浸入界面方法.

然而, 当外势场 $ V_{ext} \neq 0 $ 时, 方程 (1.6) 已不再满足经典的多辛形式 (1.1), 那是否还有某些守恒性质呢? 由此自然引发以下问题: 该方程是否仍满足某种广义结构? 若存在某类结构, 能否据此构造相应的数值格式以实现高精度求解? 针对上述问题, Bai 引入了弱多辛形式的概念, 并在此基础上发展了多辛 Runge-Kutta 格式与多辛 Runge-Kutta-Nyström 格式[35,36]. 此外, 基于哈密顿结构, Bai 进一步探讨了带 Delta 势的 NLS 方程的能量守恒算法 [37], 但在处理过程中仍摆脱不了对边界条件的依赖. 上述工作为带 Delta 势的 NLS 方程的数值求解与理论分析奠定了坚实基础, 并对构造其局部保结构算法具有重要的指导意义.

本文研究如下非光滑的 NLS 方程

$\begin{equation} \label{NLSE} iu_t=-\frac12 u_{xx} - \lvert u \rvert^2 u+\gamma \delta(x) u, \end{equation}$

即外势场 $ V_{ext} = \gamma \delta(x) $. 针对带 Delta 势的 NLS 方程, 尽管 Bai[35,36] 基于弱多辛形式提出了多辛 Runge-Kutta 格式与多辛 Runge-Kutta-Nyström 格式, 但其研究在系统性与普适性方面仍存在不足. 同时, 除多辛守恒以外, 其他局部物理量的守恒性的研究也值得进一步讨论. 为弥补该缺陷, 本文以带 Delta 势的 NLS 方程 (1.7) 为例, 成功给出了基于弱多辛形式的局部保结构算法构造的统一框架. 该框架涵盖四类多辛算法与两类局部能量守恒算法, 不仅适用于方程 (1.7), 更可推广至一切可化为弱多辛形式的 PDEs. 算法构造的核心思想源于复合构造方法与平均向量场 (AVF) 方法[24], 这促使我们在弱多辛形式下对空间与时间方向予以同等处理, 从而系统地推导出一系列局部保结构算法. 相较于标准多辛形式 (1.1), 弱多辛形式在算法构造中带来了新的困难: 一方面, 其复杂的变分结构与边界条件增加了离散过程中保持多辛精度的难度; 另一方面, 由于 Delta 势函数破坏了 NLS 方程在空间方向的平移不变性, 局部能量与动量守恒是否依然成立亦需深入探讨. 值得庆幸的是, 我们成功构造了四类多辛算法与两类局部能量守恒算法, 并严格证明了离散弱多辛守恒律与离散弱局部能量守恒律, 确保了算法在数值模拟中的有效性与稳定性.

本文结构安排如下. 第 2 节给出若干离散算子的定义及其性质. 第 3 节介绍带 Delta 势的 NLS 方程的弱多辛形式的理论基础, 随后利用复合构造方法, 从弱多辛形式出发构造四类多辛算法与两类局部能量守恒算法, 并严格证明其离散弱局部守恒律; 进一步地, 在适当的边界条件下, 探讨了相应的弱全局守恒律. 第 4 节数值实验展示了所提出算法的优良性能. 第 5 节给出结论与展望.

2 算子定义及其性质

为便于构造局部保结构算法, 我们引入如下有限差分算子

$ \delta _x^+=\frac {f^n_{j+1}-f^{n}_j} {\Delta x}, \quad \delta _x^-=\frac {f^n_j-f^{n}_{j-1}} {\Delta x}, \quad \delta_t^+=\frac{f^{n+1}_j-f^n_j}{\Delta t}, \quad \delta _t^-=\frac {f^n_j-f^{n-1}_j} {\Delta t}, $

与平均算子

$ A_x f_j^n = \frac{f_j^n+f^{n}_{j+1}}{2}, \quad A_t f_j^n = \frac{f_j^n+f^{n+1}_j}{2}. $

上述算子满足如下性质

$ \bullet $ 交换律

$ \delta _x\delta_t f_j^n=\delta_t\delta_xf_j^n, \quad A_xA_tf_j^n=A_tA_xf_j^n, \quad \delta_xA_tf_j^n=A_t\delta_xf_j^n, \quad A_x\delta_tf_j^n=\delta_tAxf_j^n. $

$ \bullet $ 离散 Leibnitz 准则

$ \delta_x(f\cdot g)_j^n=(af_{j+1}^n+(1-a)f_j^n)\cdot \delta_xg_j^n+\delta_xf_j^n\cdot ((1-a)g_{j+1}^n+ag_j^n), \quad \forall \:0 \leqslant a \leqslant 1. $

特别地,

$\begin{array}{ll} a=0, & \delta_{x}(f \cdot g)_{j}^{n}=f_{j}^{n} \cdot \delta_{x} g_{j}^{n}+\delta_{x} f_{j}^{n} \cdot g_{j+1}^{n} \\ a=\frac{1}{2}, & \delta_{x}(f \cdot g)_{j}^{n}=A_{x} f_{j}^{n} \cdot \delta_{x} g_{j}^{n}+\delta_{x} f_{j}^{n} \cdot A_{x} g_{j}^{n} \\ a=1, & \delta_{x}(f \cdot g)_{j}^{n}=f_{j+1}^{n} \cdot \delta_{x} g_{j}^{n}+\delta_{x} f_{j}^{n} \cdot g_{j}^{n} \end{array}$

取 $ f=g $ 时, 可得

$ \delta_x (\frac{1}{2}(f_j^n)^2)=\delta_x f_j^n\cdot A_x f_j^n, \quad \delta_x (\frac12 f_{j-1}^n\cdot f_j^n)=f_j^n\cdot A_x \delta_x f_{j-1}^n. $

注 2.1 上面的符号 $ \delta _x $ 表示 $ \delta _x^+ $ 和 $ \delta _x^- $. 类似地, 我们可以在时间方向上得到相应的离散 Leibnitz 准则. 离散 Leibnitz 准则在证明算法的弱局部守恒律方面起着重要作用.

3 带 Delta 势的 NLS 方程的弱多辛形式和局部保结构算法

在这一部分中, 我们将简要介绍带 Delta 势的 NLS 方程的弱多辛形式. 基于此, 再利用复合构造方法给出一系列局部保结构算法构造的统一框架, 包括四类多辛算法和两类局部能量守恒算法.

3.1 带 Delta 势的 NLS 方程的弱多辛形式

对于在 $ (\alpha,t) $ 的某个去心邻域内有定义的函数 $ u(x,t) $, 我们定义

$ u_{x=\alpha}^- = \lim\limits_{x\to \alpha^-}u(x,t), $$ u_{x=\alpha}^+ = \lim\limits_{x\to \alpha^+}u(x,t). $

鉴于本文专注于界面点 $ x = 0 $, 因此, 我们简化

$ u_{x=\alpha}^{-}=\lim _{x \rightarrow \alpha^{-}} u(x, t), \quad u_{x=\alpha}^{+}=\lim _{x \rightarrow \alpha^{+}} u(x, t). $

方程 (1.7) 中的波函数 $ u $ 是连续的, 但由于 Delta 势的存在, 其一阶导数在 $ x=0 $ 处间断, 即

$ \partial_x u^+ - \partial_x u^- = 2\gamma u(0,t). $

令 $ u=p+iq $, 共轭动量 $ v=p_x, w=q_x $, 引入状态变量 $ z={(p,q,v,w)}^\mathrm{T} $, 则带 Delta 势的 NLS 方程 (1.7) 等价于以下多辛形式, 其中势函数 $ \gamma \delta(x) $ 被一些跳跃条件所替代

$\begin{equation} \label{weak MSHS} \left\lbrace \begin{array}{ll} \textbf{M}z_t+ {\textbf{K}}z_x=\nabla_zS(z), \quad x \neq 0,\\ z^+ - z^- =2\gamma {(0,0,p^-,q^-)}^\mathrm{T}=\mathbf{{A_{jump}}} z^-, \end{array} \right. \end{equation}$

其中, $ z^\pm={(p^\pm,q^\pm,v^\pm,w^\pm)}^\mathrm{T} $, 矩阵

$ \mathbf{M}=\left(\begin{array}{cccc} 0 & -2 & 0 & 0 \\ 2 & 0 & 0 & 0 \\ 0 & 0 & 0 & 0 \\ 0 & 0 & 0 & 0 \end{array}\right), \quad \mathbf{K}=\left(\begin{array}{cccc} 0 & 0 & 1 & 0 \\ 0 & 0 & 0 & 1 \\ -1 & 0 & 0 & 0 \\ 0 & -1 & 0 & 0 \end{array}\right), \quad \mathbf{A}_{\text {jump }}=\left(\begin{array}{cccc} 0 & 0 & 0 & 0 \\ 0 & 0 & 0 & 0 \\ 2 \gamma & 0 & 0 & 0 \\ 0 & 2 \gamma & 0 & 0 \end{array}\right), $

哈密顿量 $ S(z)=\frac12 a(p^2+q^2)^2 - \frac12 (v^2+w^2) $.

系统 (3.1) 表明, $ v $ 和 $ w $ 关于 $ x $ 的一阶导数, 即 $ p $ 和 $ q $ 关于 $ x $ 的二阶导数, 在界面点 $ x=0 $ 处连续, 即

$\begin{equation} \partial _x v^+ = \partial _x v^-, \quad \partial _x w^+ = \partial _x w^-. \end{equation}$

从上面的观点来看, 当 $ x\neq0 $ 时, 系统 (3.1) 等价于一般的多辛形式 (1.1). 因此, 多辛守恒律, 局部能量守恒律, 局部动量守恒律和局部质量守恒律在 $ x\neq0 $ 时是成立的. 我们只需要关注这些守恒律在界面点 $ x=0 $ 处是否成立.

注 3.1 在微分几何的背景下, 特别从一阶微分形式的视角来看, 状态变量周围的局部线性化应保持原系统的内在性质. 例如, 方程 \(dp^+ = dp^-\) 反映了 $ dp $ 在界面点处的连续性, 即 $ p $ 在该点连续, 故 $ {dp}^+ $ 应等于 $ {dp}^- $. 类似地, 条件 $ dv^+ = dv^- + 2\gamma dp^- $ 表明 $ dv $ 在界面点处应呈现与 $ v $ 相同的跳跃间断特性, 这一性质同样体现在 $ dw^+ = dw^- + 2\gamma dq^- $ 中[36].

定义 3.1 文献[35] 对于多辛形式 (1.1) 所满足的某种局部守恒律, 如果积分形式

$\begin{equation} \label{weak sense form} \begin{array}{ll} & \displaystyle \int_{x_1}^{x_2} \displaystyle \int_{t_1}^{t_2} ( \partial_t C_d(x,t) + \partial_x C_f (x,t) ) \mathrm{d}x \mathrm{d}t\\ = & \displaystyle \int_{x_1}^{x_2} (C_d(x,t_2)-C_d(x,t_1))\mathrm{d}x + \displaystyle \int_{t_1}^{t_2} (C_f(x_2,t)-C_f(x_1,t))\mathrm{d}t\\ = & 0 \end{array} \end{equation} $

在任意矩形时空区间 $ \left[ t_1,t_2 \right] \times \left[ x_1,x_2\right] $ 内成立, 我们称该局部守恒律在某个界面点 $ (x_0,t_0) $ 处弱意义下成立, 其中 $ x_0 \in (x_1,x_2), t_0 \in (t_1,t_2) $. 这里, $ C_d $ 表示相应的密度项, $ C_f $ 表示相应的通量密度.

定理 3.1[弱多辛守恒律] 系统 (3.1) 在界面点 $ x=0 $ 处弱意义下满足多辛守恒律 (1.2).

证 对于时空域中包含界面点 $ x=0 $ 的任意矩形区域 $ [x_1,x_2]\times[t_1,t_2] $, 即 $ 0 \in (x_1,x_2) $, 我们将验证多辛守恒律 (1.2) 在弱意义下的有效性. 结合 (1.2) 和 (3.3) 式, 我们得到

$\begin{equation} \notag \begin{array}{ll} & \displaystyle \int_{x_1}^{x_2} \big[\omega(z(x,t_2)) - \omega(z(x,t_1)) \big] \mathrm{d}x + \displaystyle \int_{t_1}^{t_2} \big[\kappa(z(x_2,t)) - \kappa(z(x_1,t)) \big] \mathrm{d}t \\ = & \displaystyle \int_{x_1}^{-\epsilon_1} \big[\omega(z(x,t_2)) - \omega(z(x,t_1)) \big] \mathrm{d}x + \displaystyle \int_{t_1}^{t_2} \big[\kappa(z(-\epsilon_1,t)) - \kappa(z(x_1,t)) \big] \mathrm{d}t \\ &+ \displaystyle \int_{-\epsilon_1}^{\epsilon_2} \big[\omega(z(x,t_2)) - \omega(z(x,t_1)) \big] \mathrm{d}x + \displaystyle \int_{t_1}^{t_2} \big[\kappa(z(\epsilon_2,t)) - \kappa(z(-\epsilon_1,t)) \big] \mathrm{d}t \\ &+ \displaystyle \int_{\epsilon_2}^{x_2} \big[\omega(z(x,t_2)) - \omega(z(x,t_1)) \big] \mathrm{d}x + \displaystyle \int_{t_1}^{t_2} \big[\kappa(z(x_2,t)) - \kappa(z(\epsilon_2,t)) \big] \mathrm{d}t,\\ \end{array} \end{equation}$

其中正实数 $ \epsilon_1 $ 和 $ \epsilon_2 $ 足够小, 并且考虑到在任意点 $ x \neq 0 $ 处, 多辛守恒定律 (1.2) 在通常意义下成立, 上述表达式等号右侧的第一行和第三行自然为零. 因此, 上式等号右侧即为

$\begin{equation} \label{weak MSCL proof} \begin{array}{ll} & \displaystyle \int_{-\epsilon_1}^{\epsilon_2} \big[\omega(z(x,t_2)) - \omega(z(x,t_1))\big] \mathrm{d}x + \displaystyle \int_{t_1}^{t_2} \big[\kappa(z(\epsilon_2,t)) - \kappa(z(-\epsilon_1,t))\big] \mathrm{d}t \\ = & 2 \displaystyle \int_{-\epsilon_1}^{\epsilon_2} \big(dq \wedge dp \big)\Big|_{t_1}^{t_2} \mathrm{d}x + \displaystyle \int_{t_1}^{t_2} \big(dp\wedge dv + dq\wedge dw \big) \Big|_{-\epsilon_1}^{\epsilon_2} \mathrm{d}t. \end{array} \end{equation}$

在考虑 $ \epsilon_1 $ 和 $ \epsilon_2 $ 趋于零的极限时, 结合系统 (3.1) 的第二个方程中给出的跳跃条件, 即 $ dp $ 和 $ dq $ 在界面处的连续性, 我们可以得到式 (3.4) 等号右侧的第一项趋于零. 紧接着, 根据 $ dp^+=dp^-, dq^+=dq^- $ 和跳跃条件 $ \partial_t p^+ = \partial_t p^- $, $ \partial_t q^+ = \partial_t q^- $, (3.4) 式等号右侧第二项即为

$\begin{equation} \notag \begin{array}{ll} & \displaystyle \int_{t_1}^{t_2} (dp\wedge dv + dq\wedge dw) \Big|_{-\epsilon_1}^{\epsilon_2} \mathrm{d}t \\ = & \displaystyle \int_{t_1}^{t_2} (dp^+ \wedge dv^+ + dq^+ \wedge dw^+ - dp^- \wedge dv^- + dq^- \wedge dw^- ) \mathrm{d}t \\ = & \displaystyle \int_{t_1}^{t_2} (dp^- \wedge (dv^+ -dv^-) + dq^- \wedge (dw^+ - dw^-) ) \mathrm{d}t \\ = & \displaystyle \int_{t_1}^{t_2} (dp^- \wedge (2\gamma dp^-) + dq^- \wedge (2\gamma dq^-) ) \mathrm{d}t\\ = & 0. \end{array} \end{equation}$

因此, 可以得到

$ \displaystyle \int_{x_1}^{x_2} \big[\omega(z(x,t_2)) - \omega(z(x,t_1)) \big] \mathrm{d}x + \displaystyle \int_{t_1}^{t_2} \big[\kappa(z(x_2,t)) - \kappa(z(x_1,t)) \big] \mathrm{d}t=0. $

结合定义 3.1, 我们得到系统 (3.1) 在界面点 $ x=0 $ 处弱意义下满足多辛守恒律 (1.2).

综上, 定理得证.

虽然系统 (3.1) 由于界面点的出现而不同于标准的多辛形式 (1.1), 但系统 (3.1) 满足弱意义下的多辛守恒律, 因此可以认为它是弱多辛形式 [35].

定理 3.2[弱局部能量守恒律] 系统 (3.1) 在界面点 $ x=0 $ 处弱意义下满足局部能量守恒律 (1.3).

证 根据定理 3.1 及系统 (3.1), 我们可以得到

$\begin{equation}\notag \begin{array}{ll} & \displaystyle \int_{-\epsilon_1}^{\epsilon_2} \big[E(z(x,t_2)) - E(z(x,t_1)) \big] \mathrm{d}x + \displaystyle \int_{t_1}^{t_2} \big[F(z(\epsilon_2,t)) - F(z(-\epsilon_1,t)) \big] \mathrm{d}t \\ = & \dfrac12 \displaystyle \int_{-\epsilon_1}^{\epsilon_2} \Big( -\big(p^2+q^2\big)^2 - (p \partial_x v + q \partial_x w) \Big) \Big|_{t_1}^{t_2} \mathrm{d}x + \dfrac12 \displaystyle \int_{t_1}^{t_2}\Big(p\partial_t v + q \partial_t w -v \partial_t p - w \partial_t q\Big) \Big|_{-\epsilon_1}^{\epsilon_2} \mathrm{d}t \\ = & \dfrac12 \displaystyle \int_{t_1}^{t_2}\Big(p\partial_t v + q \partial_t w -v \partial_t p - w \partial_t q\Big) \Big|_{-\epsilon_1}^{\epsilon_2} \mathrm{d}t, \end{array} \end{equation}$

其中第二个等号成立是因为 $ p $ 和 $ q $ 在界面点 $ x = 0 $ 处连续. 根据连续性条件 $ p^+ = p^-, q^+ = q^- $ 和跳跃条件 $ \partial_t p^+ = \partial_t p^- $, $ \partial_t q^+ = \partial_t q^- $, 可以得到

$\begin{equation}\notag \begin{array}{ll} & \dfrac12 \displaystyle \int_{t_1}^{t_2} \Big(p\partial_t v + q \partial_t w -v \partial_t p - w \partial_t q\Big) \Big|_{-\epsilon_1}^{\epsilon_2} \mathrm{d}t\\ = & \dfrac12 \displaystyle \int_{t_1}^{t_2} \Big( p^-(\partial_t v^+ - \partial_t v^-) + q^-(\partial_t w^+ - \partial_t w^-) + \partial_t p^- (v^- - v^+) + \partial_t q^- ( w^- - w^+ ) \Big) \mathrm{d}t\\ = & \dfrac12 \displaystyle \int_{t_1}^{t_2} \Big(p^- (2\gamma\partial_t p^-) + q^- (2\gamma\partial_t q^-) - \partial_t p^- (2\gamma p^-) - \partial_t q^-(2\gamma q^-) \Big) \mathrm{d}t\\ = & 0. \end{array} \end{equation}$

于是, 我们有

$ \displaystyle \int_{x_1}^{x_2} \big[E(z(x,t_2)) - E(z(x,t_1)) \big] \mathrm{d}x + \displaystyle \int_{t_1}^{t_2} \big[F(z(x_2,t)) - F(z(x_1,t)) \big] \mathrm{d}t =0. $

综上, 定理得证.

定理 3.3(弱局部质量守恒律) 系统 (3.1) 在界面点 $ x=0 $ 处弱意义下满足局部质量守恒律 (1.5).

证 类似于定理 3.1 和 3.2 所示, 我们需要证明

$\begin{equation}\notag \displaystyle \int_{-\epsilon_1}^{\epsilon_2} \big[N_d(z(x,t_2)) - N_d(z(x,t_1)) \big] \mathrm{d}x + \displaystyle \int_{t_1}^{t_2} \big[N_f(z(\epsilon_2,t)) - N_f(z(-\epsilon_1,t)) \big] \mathrm{d}t =0. \end{equation}$

将 $ N_d $ 和 $ N_f $ 代入上述等式, 得到

$\begin{equation}\notag \displaystyle \int_{-\epsilon_1}^{\epsilon_2} \Big( p^2+q^2\Big) \Big|_{t_1}^{t_2} \mathrm{d}x + \displaystyle \int_{t_1}^{t_2}\Big(pw - q v\Big) \Big|_{-\epsilon_1}^{\epsilon_2} \mathrm{d}t=0. \end{equation}$

当 $ \epsilon_1 $ 与 $ \epsilon_2 $ 趋于零时, 可以得到

$\begin{equation}\notag \begin{array}{ll} & \displaystyle \int_{-\epsilon_1}^{\epsilon_2} \Big( p^2+q^2\Big) \Big|_{t_1}^{t_2} \mathrm{d}x + \displaystyle \int_{t_1}^{t_2}\Big(pw - q v\Big) \Big|_{-\epsilon_1}^{\epsilon_2} \mathrm{d}t\\ = &0 + \displaystyle \int_{t_1}^{t_2} \Big(p^- (w^+ - w^-) - q^- (v^+ -v^-)\Big) \mathrm{d}t\\ = & \displaystyle \int_{t_1}^{t_2} \Big(p^- (2\gamma q^-) - q^- (2\gamma p^-)\Big) \mathrm{d}t\\ =& 0. \end{array} \end{equation}$

综上, 定理得证.

注 3.2 由于 Delta 势函数 $ \delta(x) $ 的存在破坏了空间平移不变性, 因此局部动量守恒律在 $ x = 0 $ 处不成立.

基于弱多辛形式, 我们可以为带 Delta 势的 NLS 方程构造一系列局部保结构算法. 具体地, 我们采用不同于线方法与交替方向法的复合构造方法, 其基本思想源于 Runge-Kutta 方法, 即将PDEs按时间与空间变量分离处理的策略.

引入中间变量 $ m $ 或 $ n $, 弱多辛形式 (3.1) 可被写成

$\begin{equation} \left\lbrace \begin{array}{ll} \textbf{{K}}z_x=m, \\ \textbf{{M}} z_t=\nabla_zS(z) - m,\\ z^+ - z^- =2\gamma {(0,0,p^-,q^-)}^\mathrm{T}=\mathbf{{{A_{jump}}}} z^-, \end{array} \right. \end{equation}$

或者

$\begin{equation} \left\lbrace \begin{array}{ll} \textbf{{M}} z_t=n, \\ \textbf{{K}} z_x=\nabla_zS(z) - n,\\ z^+ - z^- =2\gamma {(0,0,p^-,q^-)}^\mathrm{T}=\mathbf{{A_{jump}}} z^-, \end{array} \right. \end{equation}$

即仅为传统的多辛形式增加了一个额外的跳跃条件. 通过离散常微分方程 (3.5a) 或者 (3.5b), 并消去中间变量 $ m $ 或 $ n $, 即可得到基于弱多辛形式的一系列局部保结构算法. 接下来, 我们以方程 (3.5a) 为例来进行详细介绍.

3.2 带 Delta 势的 NLS 方程的多辛算法

在这一部分中, 基于复合构造方法, 从弱多辛形式 (3.1) 出发, 我们将使用辛 Euler 方法, 隐式中点格式来构造带 Delta 势的 NLS 方程的多辛算法.

为实现数值离散, 需对时空区域进行网格划分. 在空间 $ [a,b]\, (0 \in (a, b)) $ 上采用等距网格, 定义节点坐标为 $ x_j=a+j\Delta x $, 其中$ j=0,1,\cdots,J $, 步长 $ \Delta x=\frac{b-a}{J} $, $ J $ 为正整数. 在时间维度上设置 $ t_n=n\Delta t $, 其中 $ n=1,\cdots,N $, 步长 $ \Delta t=T / N $, $ T $ 为总时间, $ N $ 为正整数. 因此, 界面点 $ x=0 $ 位于某个 $ 0<I<J $ 的子区间 $ [x_I,x_{I+1}) $ 中. 我们在该点处插入一个新网格点 $ x_{I'} $, 将原区间 $ [x_I,x_ {I+1}) $ 分割为两个子区间. 于是, 空间网格点序列变为 $ a = x_0 < x_1 < \cdots < x_I \leqslant x_{I'} =0 <x_{I+1} < \dots < x_J = b $, 两个新子区间的网格长度分别表示为 $ \Delta x_1 (\equiv x_{I'}-x_I) $ 和 $ \Delta x_2 (\equiv x_{I+1}-x_{I'}) $.

注 3.3 根据上述设置可知, $ \Delta x_1 + \Delta x_2 = \Delta x $ 且满足 $ 0\leq \Delta x_1 <\Delta x $. 当 $ \Delta x_1 = 0 $ 时, 意味着界面点 $ x = 0 $ 在初次等距剖分后已为网格点, 即网格点 $ x_{I} $ 与 $ x_{I'} $ 重合, 无需进一步划分.

3.2.1 多辛算法 I (MS I)

我们同时对方程组 (3.5a) 的前两个方程进行离散, 在时间和空间方向上均采用隐式中点格式. 由于界面点 $ x=0 $ 的存在, 对于节点$ (x_j,t_n) $ 且 $ j\notin \left\lbrace I,I' \right\rbrace $ 的情况, 可以得到

$\begin{equation}\label{MS I a} \left\lbrace \begin{array}{ll} \textbf{{K}}\delta_x^+ A_t z_j^n=A_tA_x m_j^n, \\ \textbf{{M}} \delta_t^+ A_x z_j^n=\nabla_zS(A_tA_x z_j^n) - A_tA_x m_j^n. \end{array} \right. \end{equation}$

而对于节点 $ (x_I, t_n) $ 和 $ (x_{I'}, t_n) $, 分别可以得到

$\begin{equation}\label{MS I b} \left\lbrace \begin{array}{ll} \textbf{{K}}\delta_{x_1}^+ A_t z_I^n=A_tA_x m_I^n, \\ \textbf{{M}} \delta_t^+ A_x z_I^n=\nabla_zS(A_tA_x z_I^n) - A_tA_x m_I^n \end{array} \right. \end{equation}$

和

$\begin{equation}\label{MS I c} \left\lbrace \begin{array}{ll} \textbf{{K}}\delta_{x_2}^+ A_t z_{I'}^n=A_tA_x m_{I'}^n, \\ \textbf{{M}} \delta_t^+ A_x z_{I'}^n=\nabla_zS(A_tA_x z_{I'}^n) - A_tA_x m_{I'}^n. \end{array} \right. \end{equation}$

消除辅助变量 $ m $, 得到

$\begin{equation}\label{MS I} \left\lbrace \begin{array}{ll} \textbf{{M}} \delta_t^+ A_x z_j^n + \textbf{{K}}\delta_x^+ A_t z_j^n=\nabla_zS(A_tA_x z_j^n), \quad j \notin \left\lbrace I,I' \right\rbrace, \\ \textbf{{M}} \delta_t^+ A_x z_I^n + \textbf{{K}}\delta_{x_1}^+ A_t z_I^n=\nabla_zS(A_tA_x z_I^n),\\ \textbf{{M}} \delta_t^+ A_x z_{I'}^n + \textbf{{K}}\delta_{x_2}^+ A_t z_{I'}^n=\nabla_zS(A_tA_x z_{I'}^n). \end{array} \right. \end{equation}$

在这里, 我们使用以下符号

$\begin{equation} \notag \delta_{x_1}^+ z_I^n=\frac{z_{I',-}^n - z_{I}^n}{\Delta x_1}, \quad \delta_{x_2}^+ z_{I'}^n=\frac{z_{I+1}^n - z_{I',+}^n}{\Delta x_2}, \end{equation}$

且在界面点 $ x_{I'} $ 处需满足跳跃条件

$\begin{equation}\label{jump conditions} z_{I',+}^n - z_{I',-}^n = \mathbf{{A_{jump}}} z_{I',-}^n. \end{equation}$

如注 3.3 所述, 需特别说明的是: 当 $ \Delta x_1 = 0 $ 时, 公式 (3.6b) 自然失效, 仅保留公式 (3.6c).

定理 3.4 算法 (3.6a) 是著名的多辛 Preissmann 格式, 而算法 (3.7) 因跳跃点 $ x=0 $ 的存在增加了两个特殊方程, 该算法满足如下离散弱多辛守恒律

(a) $ j \notin \left\lbrace I,I'\right\rbrace $,

$\begin{equation} \label{MSCL I a} \delta _t^+ (A_x dz_j^n \wedge \textbf{M}_+ A_xdz_j^n) +\delta _x^+(A_tdz_j^n \wedge \textbf{K}_+A_tdz_j^n) = 0; \end{equation}$

(b) $ j=I $,

$\begin{equation}\label{MSCL I b} \delta _t^+ (A_x dz_I^n \wedge\textbf{ M}_+ A_xdz_I^n) +\delta _{x_1}^+(A_tdz_I^n \wedge \textbf{K}_+A_tdz_I^n) = 0; \end{equation}$

(c) $ j=I' $,

$\begin{equation}\label{MSCL I c} \delta _t^+ (A_x dz_{I'}^n \wedge \textbf{M}_+ A_xdz_{I'}^n) +\delta _{x_2}^+(A_tdz_{I'}^n \wedge \textbf{K}_+A_tdz_{I'}^n) = 0. \end{equation}$

证 算法 (3.7) 的变分形式为

$\begin{equation}\label{MS I variational form} \left\lbrace \begin{array}{ll} \textbf{{M}} \delta_t^+ A_x dz_j^n + \textbf{{K}}\delta_x^+ A_t dz_j^n=S"(A_tA_x z_j^n)A_tA_x dz_j^n, \quad j \notin \left\lbrace I,I' \right\rbrace,\\ \textbf{{M}} \delta_t^+ A_x dz_I^n + \textbf{{K}}\delta_{x_1}^+ A_t dz_I^n=S"(A_tA_x z_I^n)A_tA_xdz_I^n,\\ \textbf{{M}} \delta_t^+ A_x dz_{I'}^n + \textbf{{K}}\delta_{x_2}^+ A_t dz_{I'}^n=S"(A_tA_x z_{I'}^n)A_tA_x dz_{I'}^n. \end{array} \right. \end{equation}$

对于 $ j \notin \left\lbrace I,I'\right\rbrace $ 的情况, 取式 (3.12) 中第一个方程与 $ A_tA_x dz_j^n $ 做外积, 可以得到

$\begin{equation}\notag \textbf{{M}} \delta_t^+ A_x dz_j^n \wedge A_tA_x dz_j^n+ \textbf{{K}}\delta_x^+ A_t dz_j^n \wedge A_tA_x dz_j^n=0. \end{equation}$

由于 Hessian 矩阵 $ S"(A_tA_x z_j^n) $ 的对称性, 利用矩阵分裂, 外积性质及离散 Leibnitz 准则, 可以推导出

$\begin{equation}\notag \begin{array}{rl} \textbf{{M}} \delta_t^+ A_x dz_j^n \wedge A_tA_x dz_j^n = & \textbf{{M}}_+ \delta_t^+ A_x dz_j^n \wedge A_tA_x dz_j^n - \textbf{{M}}_+ ^\mathrm{T} \delta_t^+ A_x dz_j^n \wedge A_tA_x dz_j^n\\ = & \textbf{{M}}_+ \delta_t^+ A_x dz_j^n \wedge A_tA_x dz_j^n - \delta_t^+ A_x dz_j^n \wedge \textbf{{M}}_+ A_tA_x dz_j^n\\ = & \textbf{{M}}_+ \delta_t^+ A_x dz_j^n \wedge A_tA_x dz_j^n + \textbf{{M}}_+ A_tA_x dz_j^n \wedge \delta_t^+ A_x dz_j^n\\ = & \delta_t^+(A_x dz_j^n \wedge \textbf{{M}}_+A_x dz_j^n). \end{array} \end{equation}$

同样地,

$\begin{equation}\notag \textbf{{K}}\delta_x^+ A_t dz_j^n \wedge A_tA_x dz_j^n=\delta_x^+(A_t dz_j^n \wedge \textbf{{K}}_+A_t dz_j^n). \end{equation}$

由此可以得到等式 (3.9) 成立. 类似地, 仅需将上述推导中的空间步长 $ \Delta x $ 分别替换为 $ \Delta x_1 $ 和 $ \Delta x_2 $, 即可分别得到结论 (b) 与 (c).

我们知道, 局部守恒律与边界条件无关, 但全局守恒律则不然. 我们假设解在整个区间内具有足够光滑性 (除界面点 $ x=0 $ 外). 原则上, $ a $ 和 $ b $ 可取无穷大, 但为方便起见限定为有限值. 为进一步推进研究, 我们通过规定以下合适条件来完善系统 (3.1), 即,

$ \bullet $ 初值条件

$ z|_{t=0} = z_0(x) \equiv {((p_0(x), q_0(x), v_0(x), w_0(x))} ^ \mathrm{T}, $

$ \bullet $ 边界条件

$\begin{equation} \label{IBVC} \left\lbrace \begin{array}{ll} z(a, t) = z(b, t), \quad \forall \; t \geq 0,\\ \partial _x z(a, t) = \partial _x z(b, t), \quad \forall \; t \geq 0, \end{array} \right. \end{equation}$

其中, 初值条件满足弱多辛形式 (3.1) 中最后一行的跳跃条件. (3.13) 式称为周期边界条件, 当 $ z(a, t) =0, z(b, t) =0 $ 时, 则称为齐次 Dirichlet 边界条件.

定理 3.5 在边界条件 (3.13) 下, 算法 (3.7) 在时间方向上保持如下离散弱辛守恒律

$\begin{equation} \begin{array}{rl} & \Delta x\sum \limits_{j \neq {I,I'}} \left( A_x dz_j^{n+1} \wedge \textbf{M}_+ A_x dz_j^{n+1} \right) \\ & + \Delta {x_1} \left( A_xdz_I^{n+1} \wedge \textbf{M}_+A_xdz_I^{n+1} \right) + \Delta {x_2} \left( A_xdz_{I'}^{n+1} \wedge \textbf{M}_+ A_xdz_{I'}^{n+1} \right) \\ = & \Delta x\sum \limits_{j \neq {I,I'}} \left( A_x dz_j^{n} \wedge \textbf{M}_+ A_x dz_j^{n} \right) \\ & + \Delta {x_1} \left( A_xdz_I^{n} \wedge \textbf{M}_+A_xdz_I^{n} \right) + \Delta {x_2} \left( A_xdz_{I'}^{n} \wedge \textbf{M}_+ A_xdz_{I'}^{n} \right). \end{array} \end{equation}$

证 首先, 对于 $ j \notin \left\lbrace I,I'\right\rbrace $ 的情形, 展开 (3.9) 式, 可以得到

$\begin{equation} \notag \begin{array}{rl} 0=& \dfrac{A_x dz_j^{n+1} \wedge \textbf{{M}}_+ A_xdz_j^{n+1} - A_x dz_j^n \wedge \textbf{{M}}_+A_xdz_j^n}{\Delta t} \\ &+ \dfrac{A_t dz_{j+1}^{n} \wedge \textbf{{K}}_+ A_t dz_{j+1}^{n} - A_t dz_j^n \wedge \textbf{{K}}_+A_t dz_j^n}{\Delta x}. \end{array} \end{equation}$

同样地, 对于 $ j = I $ 和 $ I' $, 我们得到

$\begin{equation} \notag \begin{array}{rl} 0= & \dfrac{A_x dz_I^{n+1} \wedge \textbf{{M}}_+ A_xdz_I^{n+1} - A_x dz_I^n \wedge \textbf{{M}}_+A_xdz_I^n}{\Delta t} \\ &+ \dfrac{A_t dz_{I',-}^{n} \wedge \textbf{{K}}_+ A_t dz_{I',-}^{n} - A_t dz_I^n \wedge\textbf{ {K}}_+A_t dz_I^n}{\Delta x_1} \end{array} \end{equation}$

和

$\begin{equation} \notag \begin{array}{rl} 0=& \dfrac{A_x dz_{I'}^{n+1} \wedge \textbf{{M}}_+ A_xdz_{I'}^{n+1} - A_x dz_{I'}^n \wedge \textbf{{M}}_+A_xdz_{I'}^n}{\Delta t} \\ &+ \dfrac{A_t dz_{I+1}^{n} \wedge \textbf{{K}}_+ A_t dz_{I+1}^{n} - A_t dz_{I',+}^n \wedge \textbf{{K}}_+A_t dz_{I',+}^n}{\Delta x_2}. \end{array} \end{equation}$

将上述三个方程的等号两边分别乘以 $ \Delta x \Delta t $, $ \Delta x_1 \Delta t $ 和 $ \Delta x_2 \Delta t $ 求和, 可以得到

$ \mathbf{I}+\mathbf{I} \mathbf{I}=0, $

其中,

$\begin{equation} \notag \begin{array}{rl} \mathbf{I} = & \Delta x \sum \limits_{j \neq {I,I'}} (A_x dz_j^{n+1} \wedge \textbf{{M}}_+A_xdz_j^{n+1} - A_xdz_j^n \wedge \textbf{{M}}_+A_xdz_j^n)\\ &+ \Delta x_1 (A_x dz_I^{n+1} \wedge \textbf{{M}}_+A_xdz_I^{n+1} - A_xdz_I^n \wedge \textbf{{M}}_+A_xdz_I^n) \\ &+ \Delta x_2 (A_x dz_{I'}^{n+1} \wedge \textbf{{M}}_+A_xdz_{I'}^{n+1} - A_xdz_{I'}^n \wedge \textbf{{M}}_+A_xdz_{I'}^n),\\ \mathbf{I} \mathbf{I} = & \Delta t \sum \limits_{j \neq {I,I'}} (A_tdz_{j+1}^n \wedge \textbf{{K}}_+A_t dz_{j+1}^n - A_tdz_j^n \wedge \textbf{{K}}_+A_tdz_j^n) \\ &+ \Delta t (A_tdz_{I',-}^n \wedge \textbf{{K}}_+A_t dz_{I',-}^n - A_tdz_{I}^n \wedge \textbf{{K}}_+A_tdz_I^n)\\ &+ \Delta t (A_tdz_{I+1}^n \wedge \textbf{{K}}_+A_t dz_{I+1}^n - A_tdz_{I',+}^n \wedge \textbf{{K}}_+A_tdz_{I',+}^n). \end{array} \end{equation}$

结合跳跃条件 (3.8) 和 (3.13) 式, $\textbf{II}$ 式可化简为

$\begin{equation}\notag \begin{array}{rl} & \Delta t (A_t dz_{I',-}^n \wedge \textbf{{K}}_+ A_t dz_{I',-}^n - A_t dz_{I',+}^n \wedge \textbf{{K}}_+A_t dz_{I',+}^n)\\ = & \Delta t (A_t dv_{I',-}^n \wedge A_t dp_{I',-}^n + A_t dw _{I',-}^n \wedge A_t dq_{I',-}^n - A_t dv_{I',+}^n \wedge A_t dp_{I',+}^n - A_t dw_{I',+} \wedge A_t dq_{I',+}^n)\\ = & \Delta t \left( (A_t dv _{I',-}^n - A_tdv_{I',+}^n) \wedge A_t dp _{I',-}^n + (A_t dw_{I',-}^n - A_t dw_{I',+}^n) \wedge A_tdq_{I',-}^n \right)\\ = &0. \end{array} \end{equation}$

至此, 定理得证.

3.2.2 多辛算法 II (MS II)

对于节点 $ (x_j,t_n) $, 其中 $ j\notin \left\lbrace I,I',I+1 \right\rbrace $, 使用辛 Euler 方法和隐式中点格式分别离散方程组 (3.5a) 的第一行和第二行, 可以得到

$\begin{equation}\label{MS II a} \left\lbrace \begin{array}{ll} \textbf{{K}}_+\delta_x^+ A_t z_j^n + \textbf{{K}}_-\delta_x^- A_t z_j^n=A_t m_j^n, \\ \textbf{{M}} \delta_t^+ z_j^n=\nabla_zS(A_t z_j^n) - A_t m_j^n. \end{array} \right. \end{equation}$

对于节点 $ (x_I,t_n) $, $ (x_{I'},t_n) $ 和 $ (x_{I+1},t_n) $, 分别可以得到

$\begin{equation}\label{MS II b} \left\lbrace \begin{array}{ll} \textbf{{K}}_+\delta_{x_1}^+ A_t z_I^n + \textbf{{K}}_-\delta_x^- A_t z_I^n=A_t m_I^n, \\ \textbf{{M}} \delta_t^+ z_I^n=\nabla_zS(A_t z_I^n) - A_t m_I^n, \end{array} \right. \end{equation}$
$\begin{equation}\label{MS II c} \left\lbrace \begin{array}{ll} \textbf{{K}}_+\delta_{x_2}^+ A_t z_{I'}^n + \textbf{{K}}_-\delta_{x_1}^- A_t z_{I'}^n=A_t m_{I'}^n, \\ \textbf{{M}} \delta_t^+ z_{I'}^n=\nabla_zS(A_t z_{I'}^n) - A_t m_{I'}^n \end{array} \right. \end{equation}$

和

$\begin{equation}\label{MS II d} \left\lbrace \begin{array}{ll} \textbf{{K}}_+\delta_x^+ A_t z_{I+1}^n + \textbf{{K}}_-\delta_{x_2}^- A_t z_{I+1}^n=A_t m_{I+1}^n, \\ \textbf{{M}} \delta_t^+ z_{I+1}^n=\nabla_z S(A_t z_{I+1}^n) - A_t m_{I+1}^n. \end{array} \right. \end{equation}$

消除辅助变量 $ m $, 得到

$\begin{equation}\label{MS II} \left\lbrace \begin{array}{ll} \textbf{{M}} \delta_t^+ z_j^n + \textbf{{K}}_+\delta_x^+ A_t z_j^n + \textbf{{K}}_-\delta_x^-A_tz_j^n=\nabla_zS(A_t z_j^n), \quad j \notin \left\lbrace I,I',I+1 \right\rbrace, \\ \textbf{{M}} \delta_t^+ z_I^n +\textbf{{K}}_+\delta_{x_1}^+ A_t z_I^n + \textbf{{K}}_-\delta_x^-A_tz_I^n=\nabla_zS(A_t z_I^n),\\ \textbf{{M}}\delta_t^+ z_{I'}^n + \textbf{{K}}_+\delta_{x_2}^+ A_t z_{I'}^n + \textbf{{K}}_-\delta_{x_1}^-A_tz_{I'}^n=\nabla_zS(A_t z_{I'}^n),\\ \textbf{{M}} \delta_t^+ z_{I+1}^n + \textbf{{K}}_+\delta_{x}^+ A_t z_{I+1}^n + \textbf{{K}}_-\delta_{x_2}^-A_tz_{I+1}^n=\nabla_zS(A_t z_{I+1}^n). \end{array} \right. \end{equation}$

在界面点 $ x_{I'} $ 处需满足跳跃条件

$\begin{equation}\notag z_{I',+}^n - z_{I',-}^n = \mathbf{{A_{jump}}}z_{I',-}^n. \end{equation}$

此外, 当 $ \Delta x_1 = 0 $ 时, 公式 (3.15b) 与 (3.15c) 自然失效, 仅公式 (3.15d) 成立.

定理 3.6 算法 (3.16) 满足如下离散弱多辛守恒律

(a) $ j \notin \left\lbrace I,I',I+1\right\rbrace $,

$\begin{equation}\label{MSCL II a} \delta _t^+ (\textbf{M}_+ dz_j^n \wedge dz_j^n) +\delta _x^+(\textbf{K}_+ A_tdz_j^n \wedge A_tdz_{j-1}^n) = 0; \end{equation}$

(b) $ j=I $,

$\begin{equation}\label{MSCL II b} \delta _t^+ (\textbf{M}_+ dz_I^n \wedge dz_I^n) +\delta _x^+(\textbf{K}_+ A_tdz_I^n \wedge A_tdz_{I-1}^n) = 0; \end{equation}$

(c) $ j=I' $,

$\begin{equation}\label{MSCL II c} \delta _t^+ (\textbf{M}_+ dz_{I'}^n \wedge dz_{I'}^n) +\delta _x^+(\textbf{K}_+ A_tdz_{I'}^n \wedge A_tdz_{I}^n) = 0; \end{equation}$

(d) $ j=I+1 $,

$\begin{equation}\label{MSCL II d} \delta _t^+ (\textbf{M}_+ dz_{I+1}^n \wedge dz_{I+1}^n) +\delta _x^+(\textbf{K}_+ A_tdz_{I+1}^n \wedge A_tdz_{I'}^n) = 0. \end{equation}$

证 算法 (3.16) 的变分形式为

$\begin{equation} \label{MS II variational form} \left\lbrace \begin{array}{ll} \textbf{{M}} \delta_t^+ dz_j^n + \textbf{{K}}_+\delta_x^+ A_t dz_j^n + \textbf{{K}}_-\delta_x^-A_t dz_j^n=S"(A_t z_j^n)A_t dz_j^n, \quad j \notin \left\lbrace I,I',I+1 \right\rbrace,\\ \textbf{{M}} \delta_t^+ dz_I^n + \textbf{{K}}_+\delta_{x_1}^+ A_t dz_I^n + \textbf{{K}}_-\delta_x^-A_t dz_I^n=S"(A_t z_I^n)A_t dz_I^n,\\ \textbf{{M}} \delta_t^+ dz_{I'}^n + \textbf{{K}}_+\delta_{x_2}^+ A_t dz_{I'}^n + \textbf{{K}}_-\delta_{x_1}^-A_t dz_{I'}^n=S"(A_t z_{I'}^n)A_t dz_{I'}^n,\\ \textbf{{M}} \delta_t^+ dz_{I+1}^n + \textbf{{K}}_+\delta_{x}^+ A_t dz_{I+1}^n + \textbf{{K}}_-\delta_{x_2}^-A_t dz_{I+1}^n=S"(A_t z_{I+1}^n)A_t dz_{I+1}^n. \end{array} \right. \end{equation}$

对于 $ j \notin \left\lbrace I,I',I+1 \right\rbrace $ 的情形, 取 (3.21) 式中第一个方程与 $ A_t dz_j^n $ 做外积, 可以得到

$\begin{equation}\notag \textbf{{M}} \delta_t^+ dz_j^n \wedge A_t dz_j^n + \textbf{{K}}_+\delta_x^+ A_t dz_j^n \wedge A_t dz_j^n + \mathbf{K_-}\delta_x^-A_t dz_j^n \wedge A_t dz_j^n=0. \end{equation}$

由于 Hessian 矩阵 $ S"(A_t z_j^n) $ 的对称性, 结合矩阵分裂, 外积性质及离散 Leibnitz 准则, 可以推导出

$\begin{equation}\notag \begin{array}{rl} \textbf{{M}} \delta_t^+ dz_j^n \wedge A_t dz_j^n = & \textbf{{M}}_+ \delta_t^+ dz_j^n \wedge A_tdz_j^n - \textbf{{M}}_+ ^\mathrm{T} \delta_t^+ dz_j^n \wedge A_t dz_j^n\\ = & \textbf{{M}}_+ \delta_t^+ dz_j^n \wedge A_tdz_j^n - \delta_t^+ dz_j^n \wedge \textbf{{M}}_+ A_t dz_j^n\\ = & \textbf{{M}}_+ \delta_t^+ dz_j^n \wedge A_tdz_j^n +\textbf{{M}}_+ A_t dz_j^n \wedge \delta_t^+ dz_j^n \\ = & \delta_t^+(\textbf{{M}}_+ dz_j^n \wedge dz_j^n). \end{array} \end{equation}$

此外,

$\begin{equation}\notag \begin{array}{rl} & \textbf{{K}}_+ \delta _x^+ A_t dz_j^n \wedge A_t dz_j^n + \textbf{{K}}_-\delta_x^-A_t dz_j^n \wedge A_t dz_j^n \\ = & \textbf{{K}}_+\delta _x^+A_t dz_j^n \wedge A_t dz_j^n - \textbf{{K}}_+^\mathrm{T} \delta _x^+A_t dz_{j-1}^n \wedge A_t dz_j^n\\ = & \textbf{{K}}_+\delta _x^+A_t dz_j^n \wedge A_t dz_j^n - \delta _x^+A_t dz_{j-1}^n \wedge \textbf{{K}}_+ A_tdz_j^n \\ = & \delta _x^+(\textbf{{K}}_+ A_tdz_j^n \wedge A_tdz_{j-1}^n). \end{array} \end{equation}$

由此可得等式 (3.17) 成立. 类似地, 仅需将上述推导中的空间步长 $ \Delta x $ 分别替换为 $ \Delta x_1 $ 和 $ \Delta x_2 $, 即可得到结论 (b), (c) 与 (d).

注 3.4 由于离散 Leibnitz 准则的特殊性, 在证明结论 (b), (c) 与 (d) 时要求空间步长保持一致, 即满足 $ \Delta x =\Delta x_1 = \Delta x_2 $ (对应 $ \Delta x = \Delta x_1, \Delta x_2 = 0 $ 或 $ \Delta x=\Delta x_2, \Delta x_1=0 $ 的情形). 此时界面点 $ x=0 $ 恰好与网格点重合.

定理 3.7 在边界条件 (3.13) 下, 算法 (3.16) 在时间方向上保持离散弱辛守恒律

$\begin{equation} \Delta x\sum \limits_{j = 0}^J \left(\textbf{ M}_+ dz_j^{n+1} \wedge dz_j^{n+1} \right) = \Delta x\sum \limits_{j =0}^J \left( \textbf{M}_+ dz_j^{n} \wedge dz_j^{n} \right). \end{equation} $

证 首先, 对于 $ j \notin \left\lbrace I,I',I+1 \right\rbrace $ 的情形, 展开式 (3.17), 可以得到

$\begin{equation} \notag \dfrac{\textbf{{M}}_+ dz_j^{n+1} \wedge dz_j^{n+1} - \textbf{{M}}_+ dz_j^n \wedge dz_j^n}{\Delta t} + \dfrac{\textbf{{K}}_+ A_t dz_{j+1}^{n} \wedge A_t dz_{j}^{n} -\textbf{{ K}}_+ A_t dz_j^n \wedge A_t dz_{j-1}^n}{\Delta x}=0. \end{equation}$

同样地, 对于 $ j = I,I' $ 和 $ I+1 $, 可以分别得到

$\begin{equation} \notag \dfrac{\textbf{{M}}_+ dz_I^{n+1} \wedge dz_I^{n+1} - \textbf{{M}}_+ dz_I^n \wedge dz_I^n}{\Delta t} + \dfrac{\textbf{{K}}_+ A_t dz_{I',-}^{n} \wedge A_t dz_{I}^{n} - \textbf{{K}}_+ A_t dz_I^n \wedge A_t dz_{I-1}^n}{\Delta x}=0, \end{equation}$
$\begin{equation} \notag \dfrac{\textbf{{M}}_+ dz_{I'}^{n+1} \wedge dz_{I'}^{n+1} - \textbf{{M}}_+ dz_{I'}^n \wedge dz_{I'}^n}{\Delta t} + \dfrac{\textbf{{K}}_+ A_t dz_{I+1}^{n} \wedge A_t dz_{I',+}^{n} - \textbf{{K}}_+ A_t dz_{I',-}^n \wedge A_t dz_{I}^n}{\Delta x}=0 \end{equation}$

和

$\begin{equation} \notag \dfrac{\textbf{{M}}_+ dz_{I+1}^{n+1} \wedge dz_{I+1}^{n+1} - \textbf{{M}}_+ dz_{I+1}^n \wedge dz_{I+1}^n}{\Delta t} + \dfrac{\textbf{{K}}_+ A_t dz_{I+2}^{n} \wedge A_t dz_{I+1}^{n} - \textbf{{K}}_+ A_t dz_{I+1}^n \wedge A_t dz_{I',+}^n}{\Delta x}=0. \end{equation}$

将上述四个方程的等号两侧分别乘以 $ \Delta x \Delta t $, 并对所有指标 $ j $ 求和, 可以得到

$ \mathbf{I I I}+\mathbf{I V}=0, $

其中,

$\begin{equation} \notag \begin{array}{rl} \textbf{III}= & \Delta x \sum \limits_{j \neq {I,I',I+1}} (\textbf{{M}}_+ dz_j^{n+1} \wedge dz_j^{n+1} - \textbf{{M}}_+ dz_j^n \wedge dz_j^n )\\ &+ \Delta x (\textbf{{M}}_+ dz_I^{n+1} \wedge dz_I^{n+1} - \textbf{{M}}_+ dz_I^n \wedge dz_I^n)\\ &+ \Delta x (\textbf{{M}}_+ dz_{I'}^{n+1} \wedge dz_{I'}^{n+1} - \textbf{{M}}_+ dz_{I'}^n \wedge dz_{I'}^n)\\ &+ \Delta x (\textbf{{M}}_+ dz_{I+1}^{n+1} \wedge dz_{I+1}^{n+1} - \textbf{{M}}_+ dz_{I+1}^n \wedge dz_{I+1}^n ),\\[0.3cm] \textbf{ IV }= & \Delta t \sum \limits_{j \neq {I,I',I+1}} (\textbf{{K}}_+ A_t dz_{j+1}^{n} \wedge A_t dz_{j}^{n} - \textbf{{K}}_+ A_t dz_j^n \wedge A_t dz_{j-1}^n) \\ &+ \Delta t (\textbf{{K}}_+ A_t dz_{I',-}^{n} \wedge A_t dz_{I}^{n} - \textbf{{K}}_+ A_t dz_I^n \wedge A_t dz_{I-1}^n)\\ &+ \Delta t (\textbf{{K}}_+ A_t dz_{I+1}^{n} \wedge A_t dz_{I',+}^{n} - \textbf{{K}}_+ A_t dz_{I',-}^n \wedge A_t dz_{I}^n)\\ &+ \Delta t (\textbf{{K}}_+ A_t dz_{I+2}^{n} \wedge A_t dz_{I+1}^{n} -\textbf{ {K}}_+ A_t dz_{I+1}^n \wedge A_t dz_{I',+}^n). \end{array} \end{equation}$

结合边界条件 (3.13), $\textbf{IV}$ 式满足

$\begin{equation}\notag \begin{array}{rl} & \Delta t (\textbf{{K}}_+ A_t dz_{I',-}^n \wedge A_t dz_{I}^n + \textbf{{K}}_+ A_t dz_{I+1}^n \wedge A_t dz_{I',+}^n)\\ &- \Delta t (\textbf{{K}}_+ A_t dz_{I',-}^n \wedge A_t dz_{I}^n + \textbf{{K}}_+ A_t dz_{I+1}^n \wedge A_t dz_{I',+}^n)\\ =&0. \end{array} \end{equation}$

综上所述, 定理得证.

注 3.5 定理 3.7 在证明过程中未利用跳跃条件 (3.8), 这可能是由辛 Euler 方法的特性所致, 即要求空间方向步长保持一致 $ \Delta x =\Delta x_1 = \Delta x_2 $. 因此该部分内容统一使用 $ \Delta x $ 表示空间步长.

3.2.3 多辛算法 III (MS III)

对于节点 $ (x_j,t_n) $, 其中 $ j\notin \left\lbrace I,I' \right\rbrace $, 采用隐式中点格式和辛 Euler 方法分别离散方程组 (3.5a) 的空间方向和时间方向, 可以得到

$\begin{equation}\label{MS III a} \left\lbrace \begin{array}{ll} \textbf{{K}}\delta_x^+ z_j^n = A_x m_j^n,\\ \textbf{{M}}_+\delta_t^+ A_xz_j^n + \textbf{{M}}_- \delta_t^-A_x z_j^n = \nabla_zS(A_x z_j^n) - A_x m_j^n. \end{array} \right. \end{equation}$

对于节点 $ (x_I,t_n) $ 和 $ (x_{I'},t_n) $, 分别可以得到

$\begin{equation}\label{MS III b} \left\lbrace \begin{array}{ll} \textbf{{K}}\delta_{x_1}^+ z_I^n = A_x m_I^n,\\ \textbf{{M}}_+\delta_t^+ A_xz_I^n + \textbf{{M}}_- \delta_t^-A_x z_I^n = \nabla_zS(A_x z_I^n) - A_x m_I^n \end{array} \right. \end{equation}$

和

$\begin{equation}\label{MS III c} \left\lbrace \begin{array}{ll} \textbf{{K}}\delta_{x_2}^+ z_{I'}^n = A_x m_{I'}^n,\\ \textbf{{M}}_+\delta_t^+ A_xz_{I'}^n + \textbf{{M}}_- \delta_t^-A_x z_{I'}^n = \nabla_zS(A_x z_{I'}^n) - A_x m_{I'}^n. \end{array} \right. \end{equation}$

消除辅助变量 $ m $, 得到

$\begin{equation}\label{MS III} \left\lbrace \begin{array}{ll} \textbf{{M}}_+\delta_t^+ A_xz_{j}^n + \textbf{{M}}_- \delta_t^-A_x z_{j}^n +\textbf{{K}}\delta_x^+ z_{j}^n = \nabla_zS(A_x z_{j}^n),\quad j \notin \left\lbrace I,I' \right\rbrace, \\ \textbf{{M}}_+\delta_t^+ A_xz_{I}^n + \textbf{{M}}_- \delta_t^-A_x z_{I}^n + \textbf{{K}}\delta_{x_1}^+ z_{I}^n = \nabla_zS(A_x z_{I}^n),\\ \textbf{{M}}_+\delta_t^+ A_xz_{I'}^n + \textbf{{M}}_- \delta_t^-A_x z_{I'}^n + \textbf{{K}}\delta_{x_2}^+ z_{I'}^n = \nabla_zS(A_x z_{I'}^n). \end{array} \right. \end{equation}$

同样地, 在界面点 $ x_{I'} $ 处需满足跳跃条件 (3.8). 如注 3.3 所述, 当 $ \Delta x_1 = 0 $ 时, 公式 (3.23b) 自然失效, 仅公式 (3.23c) 成立.

定理 3.8 算法 (3.24) 满足如下离散弱多辛守恒律

(a) $ j \notin \left\lbrace I,I'\right\rbrace $,

$\begin{equation}\label{MSCL III a} \delta _t^+ (\textbf{{M}}_+ A_x dz_j^n \wedge A_x dz_j^n) +\delta _x^+(\textbf{{K}}_+ dz_j^n \wedge dz_{j}^n) = 0; \end{equation}$

(b) $ j=I $,

$\begin{equation}\label{MSCL III b} \delta _t^+ (\textbf{{M}}_+ A_x dz_I^n \wedge A_x dz_I^n) +\delta _x^+(\textbf{{K}}_+ dz_I^n \wedge dz_{I}^n) = 0; \end{equation} $

(c) $ j=I' $,

$\begin{equation}\label{MSCL III c} \delta _t^+ (\textbf{{M}}_+ A_x dz_{I'}^n \wedge A_x dz_{I'}^n) +\delta _x^+(\textbf{{K}}_+ dz_{I'}^n \wedge dz_{I'}^n) = 0. \end{equation}$

证 算法 (3.24) 的变分形式为

$\begin{equation}\label{MS III variational form} \left\lbrace \begin{array}{ll} \textbf{{M}}_+\delta_t^+ A_x dz_{j}^n + \textbf{{M}}_- \delta_t^-A_x dz_{j}^n + \textbf{{K}}\delta_x^+ dz_{j}^n = S"(A_x z_{j}^n)A_x dz_j^n, \quad j \notin \left\lbrace I,I' \right\rbrace,\\ \textbf{{M}}_+\delta_t^+ A_x dz_{I}^n + \textbf{{M}}_- \delta_t^-A_x dz_{I}^n + \textbf{{K}}\delta_{x_1}^+ dz_{I}^n = S"(A_x z_{I}^n)A_x dz_{I}^n,\\ \textbf{{M}}_+\delta_t^+ A_x dz_{I'}^n + \textbf{{M}}_- \delta_t^-A_x dz_{I'}^n + \textbf{{K}}\delta_{x_2}^+ dz_{I'}^n = S"(A_x z_{I'}^n)A_x dz_{I'}^n. \end{array} \right. \end{equation}$

对于 $ j \notin \left\lbrace I,I' \right\rbrace $ 的情形, 取式 (3.28) 中第一个方程与 $ A_x dz_j^n $ 做外积, 可以得到

$\begin{equation}\notag \textbf{{M}}_+\delta_t^+ A_x dz_{j}^n \wedge A_xdz_j^n+ \textbf{{M}}_- \delta_t^-A_x dz_{j}^n \wedge A_xdz_j^n +\textbf{ {K}}\delta_x^+ dz_{j}^n \wedge A_xdz_j^n = 0. \end{equation}$

类似地, 利用离散 Leibnitz 准则, 可以得到

$\begin{equation}\notag \begin{array}{rl} & \textbf{{M}}_+\delta_t^+ A_x dz_{j}^n \wedge A_xdz_j^n+ \textbf{{M}}_- \delta_t^-A_x dz_{j}^n \wedge A_xdz_j^n\\ = & \textbf{{M}}_+\delta_t^+ A_x dz_{j}^n \wedge A_xdz_j^n + \textbf{{M}}_+ A_x dz_j^n \wedge \delta _t^+ A_xdz_j^{n-1}\\ = & \delta _t^+(\textbf{{M}}_+A_x dz_j^n \wedge A_xdz_j^{n-1}). \end{array} \end{equation}$

此外,

$\begin{equation}\notag \begin{array}{rl} \textbf{{K}}\delta_x^+ dz_j^n\wedge A_xdz_j^n = \delta_x^+(\textbf{{K}}_+ dz_j^n \wedge dz_j^n). \end{array} \end{equation}$

由此可得等式 (3.25) 成立. 类似地, 仅需将上述推导中的空间步长 $ \Delta x $ 分别替换为 $ \Delta x_1 $ 和 $ \Delta x_2 $, 即可得到结论 (b) 与 (c).

定理 3.9 在边界条件 (3.13) 下, 算法 (3.24) 在时间方向上保持离散弱辛守恒律

$\begin{equation} \begin{array}{rl} & \Delta x\sum \limits_{j \neq {I,I'}} \left( \textbf{M}_+A_x dz_j^{n+1} \wedge A_x dz_j^{n} \right) \\ & + \Delta {x_1} \left( \textbf{M}_+ A_xdz_I^{n+1} \wedge A_xdz_I^{n} \right) + \Delta {x_2} \left( \textbf{M}_+ A_xdz_{I'}^{n+1} \wedge A_xdz_{I'}^{n} \right) \\ = & \Delta x\sum \limits_{j \neq {I,I'}} \left( \textbf{M}_+ A_x dz_j^{n} \wedge A_x dz_j^{n-1} \right) \\ & + \Delta {x_1} \left( \textbf{M}_+ A_xdz_I^{n} \wedge A_xdz_I^{n-1} \right) + \Delta {x_2} \left( \textbf{M}_+A_xdz_{I'}^{n} \wedge A_xdz_{I'}^{n-1} \right). \end{array} \end{equation}$

3.2.4 多辛算法 IV (MS IV)

对于节点 $ (x_j,t_n) $, 其中 $ j\notin \left\lbrace I,I',I+1 \right\rbrace $, 采用辛 Euler 方法分别离散方程组 (3.5a) 的空间方向和时间方向, 可以得到

$\begin{equation}\label{MS IV a} \left\lbrace \begin{array}{ll} \textbf{{K}}_+ \delta_x^+ z_j^n + \textbf{{K}}_- \delta_x^- z_j^n= m_j^n,\\ \textbf{{M}}_+\delta_t^+ z_j^n + \textbf{{M}}_- \delta_t^- z_j^n = \nabla_zS(z_j^n) - m_j^n. \end{array} \right. \end{equation}$

对于节点 $ (x_I,t_n) $, $ (x_{I'},t_n) $ 和 $ (x_{I+1},t_n) $, 分别可以得到

$\begin{equation}\label{MS IV b} \left\lbrace \begin{array}{ll} \textbf{{K}}_+ \delta_{x_1}^+ z_I^n + \textbf{{K}}_- \delta_x^- z_I^n= m_I^n,\\ \textbf{{M}}_+\delta_t^+ z_I^n + \textbf{{M}}_- \delta_t^- z_I^n = \nabla_zS(z_I^n) - m_I^n, \end{array} \right. \end{equation}$
$\begin{equation}\label{MS IV c} \left\lbrace \begin{array}{ll} \textbf{{K}}_+ \delta_{x_2}^+ z_{I'}^n + \textbf{{K}}_- \delta_{x_1}^- z_{I'}^n= m_{I'}^n,\\ \textbf{{M}}_+\delta_t^+ z_{I'}^n + \textbf{{M}}_- \delta_t^- z_{I'}^n = \nabla_zS(z_{I'}^n) - m_{I'}^n \end{array} \right. \end{equation}$

和

$\begin{equation}\label{MS IV d} \left\lbrace \begin{array}{ll} \textbf{{K}}_+ \delta_x^+ z_{I+1}^n + \textbf{{K}}_- \delta_{x_2}^- z_{I+1}^n= m_{I+1}^n,\\ \textbf{{M}}_+\delta_t^+ z_{I+1}^n + \textbf{{M}}_- \delta_t^- z_{I+1}^n = \nabla_zS(z_{I+1}^n) - m_{I+1}^n. \end{array} \right. \end{equation}$

消除辅助变量 $ m $, 得到

$\begin{equation}\label{MS IV} \left\lbrace \begin{array}{ll} \textbf{{M}}_+\delta_t^+ z_j^n + \textbf{{M}}_-\delta _t^-z_j^n + \textbf{{K}}_+\delta_x^+ z_{j}^n + \textbf{{K}}_-\delta_x^- z_{j}^n = \nabla_zS(z_{j}^n),\quad j \notin \left\lbrace I,I',I+1 \right\rbrace, \\ \textbf{{M}}_+\delta_t^+ z_I^n + \textbf{{M}}_-\delta _t^-z_I^n + \textbf{{K}}_+\delta_{x_1}^+ z_{I}^n + \textbf{{K}}_-\delta_x^- z_{I}^n = \nabla_zS(z_{I}^n),\\ \textbf{{M}}_+\delta_t^+ z_{I'}^n + \textbf{{M}}_-\delta _t^-z_{I'}^n + \textbf{{K}}_+\delta_{x_2}^+ z_{I'}^n + \textbf{{K}}_-\delta_{x_1}^- z_{I'}^n = \nabla_zS(z_{I'}^n),\\ \textbf{{M}}_+\delta_t^+ z_{I+1}^n + \textbf{{M}}_-\delta _t^-z_{I+1}^n + \textbf{{K}}_+\delta_x^+ z_{I+1}^n + \textbf{{K}}_-\delta_{x_2}^- z_{I+1}^n = \nabla_zS(z_{I+1}^n). \end{array} \right. \end{equation}$

同样地, 在界面点 $ x_{I'} $ 处需满足跳跃条件 (3.8). 与算法 $\textbf{MS II}$ 的情况类似, 当 $ \Delta x_1 = 0 $ 时, 公式 (3.30b) 与 (3.30c) 自然失效, 仅公式 (3.30d) 成立.

定理 3.10 算法 (3.31) 满足如下离散弱多辛守恒律

(a) $ j \notin \left\lbrace I,I',I+1\right\rbrace $,

$\begin{equation}\label{MSCL IV a} \delta _t^+ (\textbf{M}_+ dz_j^n \wedge dz_j^{n-1}) +\delta _x^+(\textbf{K}_+ dz_j^n \wedge dz_{j-1}^n) = 0; \end{equation}$

(b) $ j=I $,

$\begin{equation}\label{MSCL IV b} \delta _t^+ (\textbf{M}_+ dz_I^n \wedge dz_{I}^{n-1}) +\delta _x^+(\textbf{K}_+ dz_I^n \wedge dz_{I-1}^n) = 0; \end{equation} $

(c) $ j=I' $,

$\begin{equation}\label{MSCL IV c} \delta _t^+ (\textbf{M}_+ dz_{I'}^n \wedge dz_{I'}^{n-1}) +\delta _x^+(\textbf{K}_+ dz_{I'}^n \wedge dz_{I}^n) = 0; \end{equation}$

(d) $ j=I+1 $,

$\begin{equation}\label{MSCL IV d} \delta _t^+ ( \textbf{M}_+ dz_{I+1}^n \wedge dz_{I+1}^{n-1}) +\delta _x^+(\textbf{K}_+ dz_{I+1}^n \wedge dz_{I'}^n) = 0. \end{equation}$

定理 3.11 在边界条件 (3.13) 下, 算法 (3.31) 在时间方向上保持离散弱辛守恒律

$\begin{equation} \Delta x\sum \limits_{j = 0}^J \left( \textbf{M}_+ dz_j^{n+1} \wedge dz_j^{n} \right) = \Delta x\sum \limits_{j =0}^J \left( \textbf{M}_+ dz_j^{n} \wedge dz_j^{n-1} \right). \end{equation} $

注 3.6 如注 3.5 所述, 步长统一记为 $ \Delta x $. 由于定理 3.9-定理 3.11 的证明过程分别与定理 3.5-定理 3.7 类似, 此处予以省略.

3.3 带 Delta 势的 NLS 方程的局部能量守恒算法

本节将采用 AVF 方法对时间方向进行离散去构造局部能量守恒算法. 当在空间方向应用辛 Euler 方法和隐式中点格式时, 鉴于界面点 $ x=0 $ 的存在, 需考虑远离界面点及界面点附近两种情形.

3.3.1 局部能量守恒算法 I (LEPS I)

首先, 对于节点 $ (x_j,t_n) $, 其中 $ j\notin \left\lbrace I,I' \right\rbrace $, 使用隐式中点格式和 AVF 方法分别离散方程组 (3.5a) 的空间方向和时间方向, 可以得到

$\begin{equation}\label{LEPS I a} \left\lbrace \begin{array}{ll} \textbf{{K}} \delta_x^+A_t z_j^n = A_tA_x m_j^n,\\ \textbf{{M}} \delta_t^+A_xz_j^n = \displaystyle \int_{0}^{1} \left( \nabla_zS((1-\xi)A_x z_j^n + \xi A_xz_j^{n+1}) - ((1-\xi)A_x m_j^n + \xi A_x m_j^{n+1})\right) \mathrm{d}\xi. \end{array} \right. \end{equation}$

对于节点 $ (x_I,t_n) $ 和 $ (x_{I'},t_n) $, 分别可以得到

$\begin{equation}\label{LEPS I b} \left\lbrace \begin{array}{ll} \textbf{{K}} \delta_{x_1}^+A_t z_I^n = A_tA_x m_I^n,\\ \textbf{{M}} \delta_t^+A_xz_I^n = \displaystyle \int_{0}^{1} \Big( \nabla_zS((1-\xi)A_x z_I^n + \xi A_xz_I^{n+1}) - ((1-\xi)A_x m_I^n + \xi A_x m_I^{n+1})\Big) \mathrm{d}\xi \end{array} \right. \end{equation} $

和

$\begin{equation}\label{LEPS I c} \left\lbrace \begin{array}{ll} \textbf{{K}} \delta_{x_2}^+A_t z_{I'}^n = A_tA_x m_{I'}^n,\\ \textbf{{M}} \delta_t^+A_xz_{I'}^n = \displaystyle \int_{0}^{1} \Big( \nabla_zS((1-\xi)A_x z_{I'}^n + \xi A_xz_{I'}^{n+1}) - ((1-\xi)A_x m_{I'}^n + \xi A_x m_{I'}^{n+1})\Big) \mathrm{d}\xi. \end{array} \right. \end{equation}$

消除辅助变量 $ m $, 得到

$\begin{equation}\label{LEPS I} \left\lbrace \begin{array}{ll} \textbf{{M}} \delta_t^+ A_x z_j^n +\textbf{{K}}\delta_x^+ A_t z_j^n=\displaystyle \int_{0}^{1} \left( \nabla_zS((1-\xi)A_x z_j^n + \xi A_xz_j^{n+1}) \right) \mathrm{d}\xi, \quad j \notin \left\lbrace I,I' \right\rbrace, \\ \textbf{{M}} \delta_t^+ A_x z_I^n + \textbf{{K}}\delta_{x_1}^+ A_t z_I^n=\displaystyle \int_{0}^{1} \Big( \nabla_zS((1-\xi)A_x z_I^n + \xi A_xz_I^{n+1}) \Big) \mathrm{d}\xi, \\ \textbf{{M}} \delta_t^+ A_x z_{I'}^n + \textbf{{K}}\delta_{x_2}^+ A_t z_{I'}^n=\displaystyle \int_{0}^{1} \Big( \nabla_zS((1-\xi)A_x z_{I'}^n + \xi A_xz_{I'}^{n+1}) \Big) \mathrm{d}\xi, \end{array} \right. \end{equation}$

其中在界面点 $ x_{I'} $ 处需满足跳跃条件

$\begin{equation}\notag z_{I',+}^n - z_{I',-}^n = \mathbf{{A_{jump}}} z_{I',-}^n. \end{equation}$

当 $ \Delta x_1 = 0 $ 时, 意味着公式 (3.37b) 自然失效, 仅保留公式 (3.37c).

定理 3.12 算法 (3.38) 满足如下离散弱局部能量守恒律

(a) $ j \notin \left\lbrace I,I'\right\rbrace $,

$\begin{equation} \label{LECL I a} \delta _t^+ \left( S(A_xz_j^n) + {(\delta _x^+ z_j^n)} ^\mathrm{T} \textbf{K}_+ A_xz_j^n \right) + \delta _x^+\left( -{(\delta _t^+z_j^n)}^\mathrm{T} \textbf{K}_+A_tz_j^n \right) = 0; \end{equation}$

(b) $ j=I $,

$\begin{equation}\label{LECL I b} \delta _t^+ \left( S(A_xz_I^n) + {(\delta _{x_1}^+ z_I^n)} ^\mathrm{T} \textbf{K}_+ A_xz_I^n \right) + \delta _{x_1}^+\left( -{(\delta _t^+z_I^n)}^\mathrm{T} \textbf{K}_+A_tz_I^n \right) = 0; \end{equation}$

(c) $ j=I' $,

$\begin{equation}\label{LECL I c} \delta _t^+ \left( S(A_xz_{T'}^n) + {(\delta _{x_2}^+ z_{I'}^n)} ^\mathrm{T} \textbf{K}_+ A_xz_{I'}^n \right) + \delta _{x_2}^+\left( -{(\delta _t^+z_{I'}^n)}^\mathrm{T} \textbf{K}_+A_tz_{I'}^n \right) = 0. \end{equation}$

证 对于 $ j \notin \left\lbrace I,I'\right\rbrace $ 的情况, 将算法 (3.38) 中的第一个表达式与 $ \delta_t^+ A_xz_j^n $ 做内积, 注意到

$\begin{equation}\notag {(\delta _t^+ A_xz_j^n)}^\mathrm{T} \textbf{{M}} \delta_tA_xz_j^n = 0, \end{equation}$

于是, 我们有

$\begin{equation}\label{1} {(\delta_t^+ A_xz_j^n)} ^\mathrm{T} \textbf{{K}} \delta_x^+ A_t z_j^n = {(\delta_t^+ A_xz_j^n)} ^\mathrm{T} \int_{0}^{1} \left( \nabla_zS((1-\xi)A_x z_j^n + \xi A_xz_j^{n+1}) \right) \mathrm{d}\xi. \end{equation}$

根据积分性质, 进一步推导可得

$\begin{equation}\label{2} \begin{array}{rl} & {\Big(\delta_t^+ A_xz_j^n\Big)} ^\mathrm{T} \displaystyle \int_{0}^{1} \left( \nabla_zS((1-\xi)A_x z_j^n + \xi A_xz_j^{n+1}) \right) \mathrm{d}\xi \\ = & \Big(\dfrac{A_xz_j^{n+1} - A_xz_j^n}{\Delta t} \Big)^\mathrm{T} \displaystyle \int_{0}^{1} \left( \nabla_zS((1-\xi)A_x z_j^n + \xi A_xz_j^{n+1}) \right) \mathrm{d}\xi \\ = & \dfrac 1 {\Delta t} \displaystyle \int_{0}^{1} \dfrac{\rm d}{{\rm d}\xi} S\left( (1-\xi)A_x z_j^n + \xi A_xz_j^{n+1}\right) \mathrm{d}\xi \\ = & \delta _t^+ S\Big(A_x z_j^n\Big). \end{array} \end{equation}$

结合 (3.42), (3.43) 式及矩阵分裂, 我们可以得到

$\begin{equation}\notag \delta _t^+S(A_xz_j^n) = {(\delta_t^+A_xz_j^n)}^ \mathrm{T} \textbf{{K}}_+\delta _x^+A_tz_j^n - {(\delta_x^+A_tz_j^n)}^\mathrm{T}\textbf{{K}}_+\delta_t^+A_xz_j^n. \end{equation}$

根据交换律与离散 Leibnitz 准则, 可以得到

$\begin{equation}\notag \begin{array}{rl} \delta _t^+ \left( {(\delta _x^+ z_j^n)} ^\mathrm{T} \textbf{{K}}_+ A_xz_j^n\right) = {(\delta _t^+ \delta_x^+z_j^n)}^\mathrm{T} \textbf{{K}}_+A_tA_xz_j^n + {(\delta _x^+A_tz_j^n)} ^\mathrm{T} \textbf{{K}}_+\delta_t^+A_xz_j^n,\\ \delta _x^+ \left( {(\delta _t^+z_j^n)}^\mathrm{T} \textbf{{K}}_+A_tz_j^n \right) = {(\delta _t^+ \delta_x^+z_j^n)}^\mathrm{T} \textbf{{K}}_+A_tA_xz_j^n + {(\delta _x^+A_tz_j^n)} ^\mathrm{T} \textbf{{K}}_+\delta_t^+A_xz_j^n. \end{array} \end{equation}$

通过上述两式, 可得等式 (3.39) 成立. 类似地, 只需将上述推导中的空间步长 $ \Delta x $ 分别替换为 $ \Delta x_1 $ 和 $ \Delta x_2 $, 即可得到结论 (b) 和 (c).

定理 3.13 在边界条件 (3.13) 下, 算法 (3.38) 在时间方向上保持如下离散弱全局能量守恒律

$\begin{equation} \begin{array}{rl} & \Delta x \sum \limits_{j \neq {I,I'}} \left( S(A_xz_j^{n+1}) + {(\delta _x^+ z_j^{n+1})} ^\mathrm{T} \textbf{{K}}_+ A_xz_j^{n+1} \right) \\ &+ \Delta x_1 \left( S(A_xz_I^{n+1}) + {(\delta _x^+ z_I^{n+1})} ^\mathrm{T} \textbf{{K}}_+ A_xz_I^{n+1} \right) + \Delta x_2 \left( S(A_xz_{I'}^{n+1}) + {(\delta _x^+ z_{I'}^{n+1})} ^\mathrm{T} \textbf{{K}}_+ A_xz_{I'}^{n+1} \right)\\ = & \Delta x \sum \limits_{j \neq {I,I'}} \left( S(A_xz_j^{n}) + {(\delta _x^+ z_j^{n})} ^\mathrm{T} \textbf{{K}}_+ A_xz_j^{n} \right) \\ &+ \Delta x_1 \left( S(A_xz_I^{n}) + {(\delta _x^+ z_I^{n})} ^\mathrm{T} \textbf{{K}}_+ A_xz_I^{n} \right) + \Delta x_2 \left( S(A_xz_{I'}^{n}) + {(\delta _x^+ z_{I'}^{n})} ^\mathrm{T} \textbf{{K}}_+ A_xz_{I'}^{n} \right). \end{array} \end{equation}$

证 首先, 对于 $ j \notin \left\lbrace I,I'\right\rbrace $ 的情况, 将 (3.39) 式展开, 可以得到

$\begin{equation} \notag \begin{array}{rl} 0= & \dfrac{S(A_xz_j^{n+1})-S(A_xz_j^n) + {(\delta _x^+ z_j^{n+1})} ^\mathrm{T} \textbf{{K}}_+ A_xz_j^{n+1} - {(\delta _x^+ z_j^n)} ^\mathrm{T} \textbf{{K}}_+ A_xz_j^n}{\Delta t} \\ &+ \dfrac{-{(\delta _t^+z_{j+1}^n)}^\mathrm{T} \textbf{{K}}_+A_tz_{j+1}^n + {(\delta _t^+z_j^n)}^\mathrm{T} \textbf{{K}}_+A_tz_j^n}{\Delta x}. \end{array} \end{equation}$

类似地, 对于 $ j = I $ 和 $ I' $, 分别可以得到

$\begin{equation} \notag \begin{array}{rl} 0= & \dfrac{S(A_xz_I^{n+1})-S(A_xz_I^n) + {(\delta _x^+ z_I^{n+1})} ^\mathrm{T} \textbf{{K}}_+ A_xz_I^{n+1} - {(\delta _x^+ z_I^n)} ^\mathrm{T} \textbf{{K}}_+ A_xz_I^n}{\Delta t} \\ &+ \dfrac{-{(\delta _t^+z_{I',-}^n)}^\mathrm{T} \textbf{{K}}_+A_tz_{I',-}^n + {(\delta _t^+z_I^n)}^\mathrm{T} \textbf{{K}}_+A_tz_I^n}{\Delta x} \end{array} \end{equation}$

和

$\begin{equation} \notag \begin{array}{rl} 0= & \dfrac{S(A_xz_{I'}^{n+1})-S(A_xz_{I'}^n) + {(\delta _x^+ z_{I'}^{n+1})} ^\mathrm{T} \textbf{{K}}_+ A_xz_{I'}^{n+1} - {(\delta _x^+ z_{I'}^n)} ^\mathrm{T} \textbf{{K}}_+ A_xz_{I'}^n}{\Delta t} \\ &+ \dfrac{-{(\delta _t^+z_{I+1}^n)}^\mathrm{T} \textbf{{K}}_+A_tz_{I+1}^n + {(\delta _t^+z_{I',+}^n)}^\mathrm{T} \textbf{{K}}_+A_tz_{I',+}^n}{\Delta x}. \end{array} \end{equation}$

将上述三个方程的等号两侧分别乘以 $ \Delta x \Delta t $, $ \Delta x_1 \Delta t $ 和 $ \Delta x_2 \Delta t $ 后, 对指标 $ j $ 进行求和, 可以得到

$ \mathbf{V}+\mathbf{V I}=0, $

其中,

$\begin{equation} \notag \begin{array}{rl} \textbf{V} = & \Delta x \sum \limits_{j \neq {I,I'}} \left( S(A_xz_j^{n+1}) - S(A_xz_j^{n}) + {(\delta _x^+ z_j^{n+1})} ^\mathrm{T} \textbf{{K}}_+ A_xz_j^{n+1} - {(\delta _x^+ z_j^{n})} ^\mathrm{T} \textbf{{K}}_+ A_xz_j^{n} \right) \\ &+ \Delta x_1 \left( S(A_xz_I^{n+1}) - S(A_xz_I^{n}) + {(\delta _x^+ z_I^{n+1})} ^\mathrm{T} \textbf{{K}}_+ A_xz_I^{n+1} - {(\delta _x^+ z_I^{n})} ^\mathrm{T} \textbf{{K}}_+ A_xz_I^{n} \right)\\ &+ \Delta x_2 \left( S(A_xz_{I'}^{n+1}) - S(A_xz_{I'}^{n}) + {(\delta _x^+ z_{I'}^{n+1})} ^\mathrm{T} \textbf{{K}}_+ A_xz_{I'}^{n+1} - {(\delta _x^+ z_{I'}^{n})} ^\mathrm{T} \textbf{{K}}_+ A_xz_{I'}^{n} \right), \\[0.3cm] \textbf{VI} = & \Delta t \sum \limits_{j \neq {I,I'}} \left( {(\delta _t^+z_{j+1}^n)}^\mathrm{T} \textbf{{K}}_+A_tz_{j+1}^n - {(\delta _t^+z_j^n)}^\mathrm{T} \textbf{{K}}_+A_tz_j^n \right) \\ &+ \Delta t \left( {(\delta _t^+z_{I',-}^n)}^\mathrm{T} \textbf{{K}}_+A_tz_{I',-}^n - {(\delta _t^+z_I^n)}^\mathrm{T}\textbf{{K}}_+A_tz_I^n \right) \\ &+ \Delta t \left( {(\delta _t^+z_{I+1}^n)}^\mathrm{T} \textbf{{K}}_+A_tz_{I+1}^n - {(\delta _t^+z_{I',+}^n)}^\mathrm{T} \textbf{{K}}_+A_tz_{I',+}^n \right) \bigg). \end{array} \end{equation}$

结合跳跃条件 (3.8) 和边界条件 (3.13), $\textbf{VI}$ 式可化简为

$\begin{equation}\notag \begin{array}{rl} & \Delta t \left( {(\delta _t^+z_{I',-}^n)}^\mathrm{T} \textbf{{K}}_+A_tz_{I',-}^n - {(\delta _t^+z_{I',+}^n)}^\mathrm{T} \textbf{{K}}_+A_tz_{I',+}^n \right)\\ = & \Delta t \left( (\delta _t^+ v_{I',-}^n - \delta _t^+ v_{I',+}^n)A_tp_{I',-}^n + (\delta_t^+ w_{I',-}^n - \delta_t^+w_{I',+}^n)A_tq_{I',-}^n \right)\\ =& 0, \end{array} \end{equation}$

其中 $ p $ 和 $ q $ 在界面点处连续. 由此可得定理成立.

3.3.2 局部能量守恒算法 II (LEPS II)

接下来, 对于节点 $ (x_j,t_n) $, 其中 $ j\notin \left\lbrace I,I',I+1 \right\rbrace $, 使用辛 Euler 方法格式和 AVF 方法分别离散方程组 (3.5a) 的空间方向和时间方向, 可以得到

$ \begin{equation}\label{LEPS II a} \left\lbrace \begin{array}{ll} \textbf{{K}}_+ \delta_x^+ A_t z_j^n + \textbf{{K}}_-\delta_x^-A_tz_j^n = A_t m_j^n,\\ \textbf{{M}} \delta_t^+ z_j^n = \displaystyle \int_{0}^{1} \left( \nabla_zS((1-\xi) z_j^n + \xi z_j^{n+1}) - ((1-\xi) m_j^n + \xi m_j^{n+1})\right) \mathrm{d}\xi. \end{array} \right. \end{equation}$

对于节点 $ (x_I,t_n) $, $ (x_{I'},t_n) $ 和 $ (x_{I+1},t_n) $, 分别可以得到

$\begin{equation}\label{LEPS II b} \left\lbrace \begin{array}{ll} \textbf{{K}}_+ \delta_{x_1}^+ A_t z_I^n + \textbf{{K}}_-\delta_x^-A_tz_I^n = A_t m_I^n,\\ \textbf{{M}} \delta_t^+ z_I^n = \displaystyle \int_{0}^{1} \Big( \nabla_zS((1-\xi) z_I^n + \xi z_I^{n+1}) - ((1-\xi) m_I^n + \xi m_I^{n+1})\Big) \mathrm{d}\xi, \end{array} \right. \end{equation}$
$\begin{equation}\label{LEPS II c} \left\lbrace \begin{array}{ll} \textbf{{K}}_+ \delta_{x_2}^+ A_t z_{I'}^n + \textbf{{K}}_-\delta_{x_1}^-A_tz_{I'}^n = A_t m_{I'}^n,\\ \textbf{{M}} \delta_t^+ z_{I'}^n = \displaystyle \int_{0}^{1} \Big( \nabla_zS((1-\xi) z_{I'}^n + \xi z_{I'}^{n+1}) - ((1-\xi) m_{I'}^n + \xi m_{I'}^{n+1})\Big) \mathrm{d}\xi \end{array} \right. \end{equation}$

和

$\begin{equation}\label{LEPS II d} \left\lbrace \begin{array}{ll} \textbf{{K}}_+ \delta_x^+ A_t z_{I+1}^n + \textbf{{K}}_-\delta_{x_2}^-A_tz_{I+1}^n = A_t m_{I+1}^n,\\ \textbf{{M}} \delta_t^+ z_{I+1}^n = \displaystyle \int_{0}^{1} \left( \nabla_zS((1-\xi) z_j^n + \xi z_{I+1}^{n+1}) - ((1-\xi) m_{I+1}^n + \xi m_{I+1}^{n+1})\right) \mathrm{d} \xi. \end{array} \right. \end{equation}$

消除辅助变量 $ m $, 得到

$\begin{equation}\label{LEPS II} \left\lbrace \begin{array}{ll} \textbf{{M}} \delta_t^+ z_j^n + \textbf{{K}}_+\delta_x^+ A_t z_j^n + \textbf{{K}}_-\delta_x^-A_tz_j^n = \displaystyle \int_{0}^{1} \left( \nabla_zS((1-\xi) z_j^n + \xi z_j^{n+1}) \right) \mathrm{d}\xi, \quad j \notin \left\lbrace I,I',I+1 \right\rbrace, \\ \textbf{{M}} \delta_t^+ z_I^n + \textbf{{K}}_+\delta_{x_1}^+ A_t z_I^n + \textbf{{K}}_-\delta_x^-A_tz_I^n = \displaystyle \int_{0}^{1} \Big( \nabla_zS((1-\xi) z_I^n + \xi z_I^{n+1}) \Big) \mathrm{d}\xi,\\ \textbf{{M}} \delta_t^+ z_{I'}^n + \textbf{{K}}_+\delta_{x_2}^+ A_t z_{I'}^n + \textbf{{K}}_-\delta_{x_1}^-A_tz_{I'}^n = \displaystyle \int_{0}^{1} \Big( \nabla_zS((1-\xi) z_{I'}^n + \xi z_{I'}^{n+1}) \Big) \mathrm{d}\xi,\\ \textbf{{M}} \delta_t^+ z_{I+1}^n + \textbf{{K}}_+\delta_x^+ A_t z_{I+1}^n + \textbf{{K}}_-\delta_{x_2}^-A_tz_{I+1}^n = \displaystyle \int_{0}^{1} \Big( \nabla_zS((1-\xi) z_{I+1}^n + \xi z_{I+1}^{n+1}) \Big) \mathrm{d}\xi, \end{array} \right. \end{equation}$

其中界面点 $ x_{I'} $ 处需满足跳跃条件 (3.8). 类似地, 当 $ \Delta x_1 = 0 $ 时, 意味着公式 (3.45b) 与 (3.45c) 自然失效, 仅保留公式 (3.45d).

定理 3.14 算法 (3.46) 满足如下离散弱局部能量守恒律

(a) $ j \notin \left\lbrace I,I',I+1\right\rbrace $,

$\begin{equation} \label{LECL II a} \delta _t^+ \left( S(z_j^n) + {(\delta _x^+ z_{j-1}^n)} ^\mathrm{T} \textbf{K}_+ z_j^n \right) + \delta _x^+\left( -{(\delta _t^+z_{j-1}^n)}^\mathrm{T} \textbf{K}_+A_tz_j^n \right) = 0; \end{equation}$

(b) $ j=I $,

$\begin{equation}\label{LECL II b} \delta _t^+ \left( S(z_I^n) + {(\delta _x^+ z_{I-1}^n)} ^\mathrm{T} \textbf{K}_+ z_I^n \right) + \delta _x^+\left( -{(\delta _t^+z_{I-1}^n)}^\mathrm{T} \textbf{K}_+A_tz_I^n \right) = 0; \end{equation}$

(c) $ j=I' $,

$\begin{equation}\label{LECL II c} \delta _t^+ \left( S(z_{I'}^n) + {(\delta _x^+ z_{I}^n)} ^\mathrm{T} \textbf{K}_+ z_{I'}^n \right) + \delta _x^+\left( -{(\delta _t^+z_{I}^n)}^\mathrm{T} \textbf{K}_+A_tz_{I'}^n \right) = 0; \end{equation}$

(d) $ j=I+1 $,

$\begin{equation}\label{LECL II d} \delta _t^+ \left( S(z_{I+1}^n) + {(\delta _x^+ z_{I'}^n)} ^\mathrm{T} \textbf{K}_+ z_{I+1}^n \right) + \delta _x^+\left( -{(\delta _t^+z_{I'}^n)}^\mathrm{T} \textbf{K}_+A_tz_{I+1}^n \right) = 0. \end{equation} $

证 对于 $ j \notin \left\lbrace I,I',I+1\right\rbrace $ 的情况, 将算法 (3.46) 中的第一个表达式与 $ \delta_t^+ z_j^n $ 做内积, 可以得到

$\begin{equation} \notag {(\delta _t^+ z_j^n)} ^\mathrm{T}\textbf{{K}}_+\delta_x^+A_tz_j^n + {(\delta_t^+z_j^n)} ^\mathrm{T} \textbf{{K}}_- \delta_x^-A_tz_j^n = {(\delta _t^+z_j^n) } ^\mathrm{T} \int_{0}^{1} \nabla_zS((1-\xi) z_j^n + \xi z_j^{n+1})\mathrm{d}\xi, \end{equation}$

其中,

$\begin{equation}\notag \begin{array}{rl} & {(\delta _t^+ z_j^n)} ^\mathrm{T}\textbf{{K}}_+\delta_x^+A_tz_j^n + {(\delta_t^+z_j^n)}^\mathrm{T} \textbf{{K}}_-\delta_x^-A_tz_j^n \\ = & {(\delta _t^+ z_j^n)} ^\mathrm{T}\textbf{{K}}_+\delta_x^+A_tz_j^n - {(\delta_x^+A_tz_{j-1}^n)}^\mathrm{T}\textbf{{K}}_+\delta_t^+z_j^n\\ = & \delta_t^+\left( -{(\delta_x^+z_{j-1}^n)}^\mathrm{T} \textbf{{K}}_+z_j^n \right) + \delta_x^+\left( {(\delta_t^+z_{j-1}^n)}^\mathrm{T}\textbf{{K}}_+(A_tz_j^n) \right). \end{array} \end{equation}$

此外,

$\begin{equation}\notag {(\delta_t^+ z_j^n) } ^\mathrm{T} \int_{0}^{1}\nabla_zS((1-\xi) z_j^n + \xi z_j^{n+1}) \mathrm{d}\xi = \delta_t^+\left( S(z_j^n)\right), \end{equation}$

由此可得等式 (3.47) 成立. 类似地, 只需将上述推导中的空间步长 $ \Delta x $ 分别替换为 $ \Delta x_1 $ 和 $ \Delta x_2 $, 即可得到结论 (b), (c) 和 (d).

注 3.7 由于离散 Leibnitz 准则的特殊性, 在证明结论 (b), (c) 和 (d) 时要求空间步长保持一致, 即需满足 $ \Delta x =\Delta x_1 = \Delta x_2 $ (或 $ \Delta x = \Delta x_1, \Delta x_2 = 0 $, 或 $ \Delta x=\Delta x_2, \Delta x_1=0 $). 此时界面点 $ x=0 $ 恰好与网格点重合.

定理 3.15 在边界条件 (3.13) 下, 算法 (3.46) 保持如下的离散弱全局能量守恒律

$\begin{equation} \Delta x \sum \limits_{j=0}^J \left( S(z_j^{n+1}) + {(\delta _x^+ z_{j-1}^{n+1})} ^\mathrm{T} \textbf{K}_+ z_j^{n+1} \right) = \Delta x \sum \limits_{j=0}^J \left( S(z_j^{n}) + {(\delta _x^+ z_{j-1}^{n})} ^\mathrm{T} \textbf{K}_+ z_j^{n} \right). \end{equation} $

证 首先, 对于 $ j \notin \left\lbrace I,I',I+1\right\rbrace $ 的情况, 将式 (3.47) 展开, 可以得到

$\begin{equation} \notag \begin{array}{rl} 0= & \dfrac{S(z_j^{n+1})-S(z_j^n) + {(\delta _x^+ z_{j-1}^{n+1})} ^\mathrm{T} \textbf{{K}}_+ z_j^{n+1} - {(\delta _x^+ z_{j-1}^n)} ^\mathrm{T} \textbf{{K}}_+ z_j^n}{\Delta t} \\ &+ \dfrac{-{(\delta _t^+z_{j}^n)}^\mathrm{T} \textbf{{K}}_+A_tz_{j+1}^n + {(\delta _t^+z_{j-1}^n)}^\mathrm{T} \textbf{{K}}_+A_tz_j^n}{\Delta x}. \end{array} \end{equation}$

类似地, 对于 $ j = I,I' $ 和 $ I+1 $, 我们可以得到

$\begin{equation} \notag \begin{array}{rl} 0=& \dfrac{S(z_I^{n+1})-S(z_I^n) + {(\delta _x^+ z_{I-1}^{n+1})} ^\mathrm{T} \textbf{{K}}_+ z_I^{n+1} - {(\delta _x^+ z_{I-1}^n)} ^\mathrm{T} \textbf{{K}}_+ z_{I}^n}{\Delta t} \\ &+ \dfrac{-{(\delta _t^+z_{I}^n)}^\mathrm{T} \textbf{{K}}_+A_tz_{I',-}^n + {(\delta _t^+z_{I-1}^n)}^\mathrm{T} \textbf{{K}}_+A_tz_I^n}{\Delta x}, \end{array} \end{equation}$
$\begin{equation} \notag \begin{array}{rl} 0=& \dfrac{S(z_{I'}^{n+1})-S(z_{I'}^n) + {(\delta _x^+ z_{I}^{n+1})} ^\mathrm{T} \textbf{{K}}_+ z_{I'}^{n+1} - {(\delta _x^+ z_{I}^n)} ^\mathrm{T} \textbf{{K}}_+ z_{I'}^n}{\Delta t} \\ &+ \dfrac{-{(\delta _t^+z_{I',+}^n)}^\mathrm{T} \textbf{{K}}_+A_tz_{I+1}^n + {(\delta _t^+z_{I}^n)}^\mathrm{T} \textbf{{K}}_+A_tz_{I',-}^n}{\Delta x} \end{array} \end{equation}$

和

$\begin{equation} \notag \begin{array}{rl} 0=& \dfrac{S(z_{I+1}^{n+1})-S(z_{I+1}^n) + {(\delta _x^+ z_{I'}^{n+1})} ^\mathrm{T} \textbf{{K}}_+ z_{I+1}^{n+1} - {(\delta _x^+ z_{I'}^n)} ^\mathrm{T} \textbf{{K}}_+ z_{I+1}^n}{\Delta t} \\ &+ \dfrac{-{(\delta _t^+z_{I+1}^n)}^\mathrm{T} \textbf{{K}}_+A_tz_{I+2}^n + {(\delta _t^+z_{I',+}^n)}^\mathrm{T} \textbf{{K}}_+A_tz_{I+1}^n}{\Delta x}. \end{array} \end{equation}$

将上述四个方程的等号两侧分别乘以 $ \Delta x \Delta t $ 后, 对指标 $ j $ 进行求和, 可得

$ \mathbf{V I I}+\mathbf{V I I I}=0, $

其中,

$\begin{equation} \notag \begin{array}{rl} \mathbf{V I I} = & \Delta x \sum \limits_{j=0}^J \left( S(z_j^{n+1}) - S(z_j^{n}) + {(\delta _x^+ z_{j-1}^{n+1})} ^\mathrm{T} \textbf{{K}}_+ z_j^{n+1} - {(\delta _x^+ z_{j-1}^{n})} ^\mathrm{T} \textbf{{K}}_+ z_j^{n} \right), \\[0.3cm] \mathbf{V I I I} = & \Delta t \sum \limits_{j \neq {I,I',I+1} } \left( {(\delta _t^+z_{j}^n)}^\mathrm{T} \textbf{{K}}_+A_tz_{j+1}^n - {(\delta _t^+z_{j-1}^n)}^\mathrm{T} \textbf{{K}}_+A_tz_{j}^n \right)\\ &+ \Delta t \left( {(\delta _t^+z_{I}^n)}^\mathrm{T} \textbf{{K}}_+A_tz_{I',-}^n - {(\delta _t^+z_{I-1}^n)}^\mathrm{T} \textbf{{K}}_+A_tz_I^n \right) \\ &+ \Delta t \left( {(\delta _t^+z_{I',+}^n)}^\mathrm{T} \textbf{{K}}_+A_tz_{I+1}^n - {(\delta _t^+z_{I}^n)}^\mathrm{T} \textbf{{K}}_+A_tz_{I',-}^n \right) \\ &+ \Delta t \left( {(\delta _t^+z_{I+1}^n)}^\mathrm{T} \textbf{{K}}_+A_tz_{I+2}^n - {(\delta _t^+z_{I',+}^n)}^\mathrm{T} \textbf{{K}}_+A_tz_{I+1}^n \right). \end{array} \end{equation}$

结合边界条件 (3.13), $\textbf{VIII}$ 式可化简为

$\begin{equation}\notag \begin{array}{rl} & \Delta t \left( {(\delta _t^+z_{I}^n)}^\mathrm{T} \textbf{{K}}_+A_tz_{I',-}^n \right) + \Delta t \left( {(\delta _t^+z_{I',+}^n)}^\mathrm{T} \textbf{{K}}_+A_tz_{I+1}^n \right) \\ &- \Delta t \left( {(\delta _t^+z_{I}^n)}^\mathrm{T} \textbf{{K}}_+A_tz_{I',-}^n \right) + \Delta t \left( {(\delta _t^+z_{I',+}^n)}^\mathrm{T} \textbf{{K}}_+A_tz_{I+1}^n \right)\\ =& 0. \end{array} \end{equation}$

如注 3.5 所述, 此处步长统一用 $ \Delta x $ 表示. 由此可得定理成立.

注 3.8 在本研究中, 我们采用隐式中点格式与辛 Euler 方法为带 Delta 势的 NLS 方程构造了一系列局部保结构算法. 值得注意的是, 隐式中点格式采用隐式方法, 可通过迭代法处理 Delta 势函数, 因此无需 $ x=0 $ 精确位于网格点上. 与之相反, 辛 Euler 方法采用显式方法, 要求 $ x=0 $ 精确处于网格点才能保证辛结构的保持.

注 3.9 本文基于弱多辛形式, 为带 Delta 势的 NLS 方程给出一系列局部保结构算法构造的统一框架. 然而, 这些算法的适用性并不局限于该特定方程. 事实上, 任何能够表示为弱多辛形式的 PDEs, 均适用于该框架. 这意味着, 我们的算法不仅针对特定方程, 更重要的是为所有可弱多辛化的 PDEs 提供了统一的数值求解框架, 从而为该类问题提供了统一有效的数值求解途径.

4 数值实验

为探究所构造算法的性能, 我们将展示齐次 Dirichlet 边界条件下带 Delta 势的 NLS 方程的数值实验. 实验内容主要包括: (1) 数值模拟效果与全局最大误差; (2) 不同算法中的质量守恒与能量守恒效果.

为更有效地展示守恒特性, 现分别定义离散全局质量与能量的相对误差如下

$\begin{equation}\notag GN=|\frac{\mathcal{N}^n-\mathcal{N}^0}{\mathcal{N}^0}|, \quad GE=|\frac{\mathcal{E}^n-\mathcal{E}^0}{\mathcal{E}^0}|,\quad n=0,1,\cdots,N, \end{equation}$

其中 $ \mathcal{N}^n $ 和 $ \mathcal{E}^n $ 分别表示 $ t = t_n $ 时刻的离散全局质量和能量. 此外, 为展示所构造算法对应的波函数数值解的最大误差, 我们采用

$\begin{equation}\notag error=\max \limits_{0\leq j \leq J}\big(\max(|p_j^n-p(x_j,t_n)|,|q_j^n-q(x_j,t_n)|)\big) \end{equation}$

来衡量 $ t=t_n $ 时刻波函数在整个空间域中的最大误差, 其中 $ p(x_j,t_n) $ 和 $ p_j^n $ 分别表示 $ p(x,t) $ 在点 $ (x_j,t_n) $ 处的精确解与数值解.

我们选择包含精确解见文献[38] 的实验方案

$\begin{equation}\label{solution} u(x,t)=\psi(x) \mathrm{exp} (i(2\gamma-1)^2t/8), \end{equation}$

其中

$\begin{equation}\notag \psi(x)= \left\lbrace \begin{array}{lr} k \mathrm{sech}(k(x-x_0)), \quad x> 0,\\ k \mathrm{sech}(k(x+x_0)), \quad x<0, \end{array} \right. \end{equation}$

参数 $ k $ 满足

$\begin{equation}\notag k=1/2-\gamma, \quad \mathrm{tanh}(kx_0)=\gamma/k. \end{equation}$

对于相应的稳态 NLS 方程, 当 $ \gamma<1/4 $ 时, 存在唯一束缚态解. 本文研究的解析解涉及该束缚态随时间演化的过程. 根据束缚态的存在条件, 我们在数值实验中选取 $\gamma = 0.1$. 由于精确解 $u(x, t)$ 在任意给定的 $ t $ 时刻, 随 $ x $ 远离 $ 0 $ 时迅速衰减, 为便于数值实现, 对问题 (3.1) 施加齐次 Dirichlet 边界条件

$\begin{equation} \psi(a)=\psi(b)=0, \end{equation}$

其中, 空间计算域取 $ a=-50 $, $ b=50 $, 空间步长固定为 $ \Delta x=0.1 $ 且确保界面点 $ x=0 $ 位于网格节点上, 时间步长固定为 $ \Delta t=0.01 $. 在时间离散的每一步中, 鉴于产生的非线性代数系统均为隐式形式, 我们采用不动点迭代法求解, 并设定统一的收敛准则: 相邻两次迭代向量的差值数量级小于 $ \mathcal{O}(10^{-15}) $ 时停止迭代.

4.1 数值模拟与全局最大误差

首先, 我们在时间区间 $ t \in [0,3] $ 内展示算法 $\textbf{MS I}$ (3.7), $\textbf{MS II}$ (3.16) 及 $\textbf{LEPS II}$ (3.46) 的数值模拟结果, 并与该方程的精确解 (4.1) 进行对比, 如图 1 所示. 从图形表现来看, 各算法的数值解曲线与精确解曲线高度一致, 波在指定区间内的传播行为与理论预期相符, 且波的形状在传播过程中得到了良好的保持. 这里需要说明的是, 算法 $\textbf{MS III}$ (3.24), $\textbf{MS IV}$ (3.31) 和 $\textbf{LEPS I}$ (3.38) 所得数值解的图形特征与算法 $\textbf{MS II}$ (3.16) 相似, 鉴于篇幅限制, 此处不再重复展示.

图 1

新窗口打开| 下载原图ZIP| 生成PPT

图 1   不同算法的数值解与精确解对比


其次, 以算法 $\textbf{MS I}$, $\textbf{MS III}$, $\textbf{LEPS I}$ 和 $\textbf{LEPS II}$ 为例, 我们在图 2-3 中分别呈现了从初始时刻 $ t=0s $ 至$ t=10s $ 的数值解的最大误差, 时间步长设置为 $ \Delta t=0.01 $. 根据图 2-3 所示结果, 在多辛算法与局部能量守恒算法的长时间演化过程中, 数值解的最大误差始终被控制在 $ \mathcal{O}(10^{-3}) $ 数量级, 充分体现了这些算法在长时间数值模拟中良好的稳定性.

图 2

新窗口打开| 下载原图ZIP| 生成PPT

图 2   带 Delta 势的 NLS 方程的多辛算法的数值解的最大误差


图 3

新窗口打开| 下载原图ZIP| 生成PPT

图 3   带 Delta 势的 NLS 方程的局部能量守恒算法的数值解的最大误差


4.2 不变量的守恒特性

本节对所提出算法在质量和能量守恒特性方面进行了数值验证. 图 4 展示了多辛算法 $\textbf{MS II}$ 和 $\textbf{MS IV}$ 从初始时刻 $ t=0s $ 至 $ t=10s $ 的不变量的离散全局相对误差的演化结果, 取 $ \Delta t=0.01 $ 的时间步长. 从图 4 观察到, 两类多辛算法的离散全局质量保持得非常好, 离散全局质量的相对误差分别稳定在 $ \mathcal{O}(10^{-13}) $ 和 $ \mathcal{O}(10^{-14}) $ 数量级, 这与计算机舍入误差相当, 充分验证了所构造多辛算法的守恒性能. 鉴于其他多辛算法 ($\textbf{MS I}$ 和 $\textbf{MS III}$) 在离散全局质量守恒性质上表现出相似精度, 此处不再重复展示.

图 4

新窗口打开| 下载原图ZIP| 生成PPT

图 4   带 Delta 势的 NLS 方程的离散全局质量的相对误差


接下来, 我们对两类局部能量守恒算法的不变量守恒特性进行检验. 图 5-8 展示了算法 $\textbf{LEPS I}$ 和 $\textbf{LEPS II}$ 从初始时刻 $ t=0s $ 至 $ t=10s $ 的不变量结果, 采用 $ \Delta t=0.01 $ 的时间步长. 从图 5-6 可以看出, 对于两类局部能量守恒算法的不变量的离散全局相对误差变化情况, 算法 $\textbf{LEPS I}$ 的离散全局质量和能量均保持良好, 分别达到 $ \mathcal{O}(10^{-11}) $ 和 $ \mathcal{O}(10^{-9}) $ 数量级. 相比之下, 算法 $\textbf{LEPS II}$ 在离散全局质量和能量守恒方面表现更优, 尤其是离散全局能量相对误差被严格限制在 $ \mathcal{O}(10^{-16}) $ 的舍入误差水平. 但离散全局能量误差呈现近似线性增长的趋势, 我们推测这可能与迭代算法的误差累积有关. 尽管算法 $\textbf{LEPS I}$ 在能量守恒方面略逊于算法 $\textbf{LEPS II}$, 但仍处于可接受范围内.

图 5

新窗口打开| 下载原图ZIP| 生成PPT

图 5   算法 $\textbf{LEPS I}$ 对带 Delta 势的 NLS 方程的模拟结果


图 6

新窗口打开| 下载原图ZIP| 生成PPT

图 6   算法 $\textbf{LEPS II}$ 对带 Delta 势的 NLS 方程的模拟结果


图 7

新窗口打开| 下载原图ZIP| 生成PPT

图 7   在精确解所满足边界下, 带 Delta 势的 NLS 方程的离散局部能量残差


图 8

新窗口打开| 下载原图ZIP| 生成PPT

图 8   在精确解所满足边界下, 带 Delta 势的 NLS 方程的离散全局能量的相对误差


图 7-8 分别给出了不采用齐次 Dirichlet 边界条件 (如: 精确解所满足的边界) 下的局部能量守恒算法 $\textbf{LEPS I}$ 和 $\textbf{LEPS II} $的局部守恒律的残差分布和离散全局能量的相对误差. 可以观察到, 两类算法的离散局部能量都保持得非常好, 特别地, 算法 $\textbf{LEPS II}$ 的局部能量残差控制在 $ \mathcal{O}(10^{-13}) $ 数量级. 相比之下, 算法 $\textbf{LEPS I}$ 和 $\textbf{LEPS II}$ 的离散全局能量的相对误差在精确解所满足的边界条件下会出现显著误差积累, 最终导致算法无法达到预想的精度. 这进一步表明, 局部保结构算法不依赖于边界条件.

综上所述, 本文提出的局部保结构算法在理论推导与数值实效方面均表现出良好的一致性, 所有数值实验结果均验证了其有效性和可靠性.

5 结论

以带 Delta 势的 NLS 方程为例, 本文基于弱多辛形式, 给出了一系列局部保结构算法构造的统一框架. 该框架包含四类多辛算法和两类局部能量守恒算法, 这些算法不依赖于特定的边界条件, 适用于任意时空区域. 我们从理论角度严格证明了相应算法的离散弱局部守恒律, 以及在适当边界条件下的离散弱全局守恒律. 数值实验验证了所提算法在长时间模拟中的有效性与稳定性. 该框架可推广至任何可表示为弱多辛形式的 PDEs, 为此类问题的数值求解提供了普适性的计算途径, 也为计算数学专业的教学研究提供了典型的教学案例.

参考文献

Wang Z, Zhu J, Wang C.

Finite difference alternative unequal-sized weighted essentially non-oscillatory schemes for hyperbolic conservation laws

Physics of Fluids, 2022, 34: 1-29

[本文引用: 1]

Elgindi T, Masmoudi N.

$L_{\infty}$ ill-posedness for a class of equations arising in hydrodynamics

Archive for Rational Mechanics and Analysis, 2020, 235: 1979-2025

DOI:10.1007/s00205-019-01457-7      [本文引用: 1]

We give a new approach to studying norm inflation (in some critical spaces) for a wide class of equations arising in hydrodynamics. As an application, we prove strong ill-posedness of the n-dimensional Euler equations in the class C1 n L2( ) and also in Ckn L2( ) where can be thewhole space, a smooth bounded domain, or the torus. We also apply our method to the Oldroyd B, surface quasi-geostrophic, and Boussinesq systems.

Kundu P, Almusawa H, Fahim M.

Linear and nonlinear effects analysis on wave profiles in optics and quantum physics

Results in Physics, 2021, 23: 1-9

[本文引用: 1]

DeWitt S, Thornton K.

Phase field modeling of microstructural evolution

Computational Materials System Design, 2018, 22: 67-87

[本文引用: 1]

Feng K, Qin M.

The symplectic methods for the computation of Hamiltonian equations

Numerical Methods for Partial Differential Equations, 2006, 1297: 1-37

[本文引用: 1]

Bridges T, Reich S.

Multi-symplectic integrators: numerical schemes for Hamiltonian PDEs that conserve symplecticity

Physics Letters A, 2001, 284: 184-193

DOI:10.1016/S0375-9601(01)00294-8      URL     [本文引用: 1]

Marsden J, Patrick G, Shkoller S. Multisymplectic geometry,

variational integrators, and nonlinear PDEs

Communications in Mathematical Physics, 1998, 199(2): 351-395

DOI:10.1007/s002200050505      URL    

Reich S.

Multi-symplectic Runge-Kutta collocation methods for Hamiltonian wave equations

Journal of Computational Physics, 2000, 157: 473-499

DOI:10.1006/jcph.1999.6372      URL     [本文引用: 2]

Cohen D, Hairer E, Lubich C.

Conservation of energy, momentum and actions in numerical discretizations of non-linear wave equations

Numerische Mathematik, 2008, 110: 113-143

DOI:10.1007/s00211-008-0163-9      URL     [本文引用: 1]

Chen Y, Sun Y, Tang Y.

Energy-preserving numerical methods for Landau-Lifshitz equation

Journal of Physics A, 2011, 44: 307-329

Wang T, Guo B, Xu Q.

Fourth-order compact and energy conservative difference schemes for the nonlinear Schrödinger equation in two dimensions

Journal of Computational Physics, 2013, 243: 382-399

DOI:10.1016/j.jcp.2013.03.007      URL    

Yi N, Liu H.

An energy conserving local discontinuous Galerkin method for a nonlinear variational wave equation

Communications in Computational Physics, 2018, 23: 747-772

DOI:10.4208/cicp.OA-2016-0189      URL    

Cai J, Wang Y, Jiang C.

Local structure-preserving algorithms for general multi-symplectic Hamiltonian PDEs

Communications in Computational Physics, 2019, 235: 210-220

[本文引用: 1]

Yusuf A, Alshomrani A, Sulaiman T.

Extended classical optical solitons to a nonlinear Schrodinger equation expressing the resonant nonlinear light propagation through isolated flaws in optical waveguides

Optical and Quantum Electronics, 2022, 54: 1-13

DOI:10.1007/s11082-021-03373-1      [本文引用: 1]

Peregrine D.

Water wave, nonlinear Schrödinger equations and their solutions

J Austral Math Soc Ser B, 1983, 25(1): 16-43

DOI:10.1017/S0334270000003891      URL     [本文引用: 1]

Du D, Wu Y, Zhang K.

On the blow-up for nonlinear Schrödinger equation with inhomogeneous perturbation

Mathematical Methods in the Applied Sciences, 2024, 47(7): 5664-5676

DOI:10.1002/mma.v47.7      URL    

Li Z, Tian S, Yang J.

Riemann-Hilbert approach and soliton solutions for the higher-order dispersive nonlinear Schrödinger equation with nonzero boundary conditions

East Asian Journal on Applied Mathematics, 2021, 11(2): 369-388

DOI:10.4208/eajam      URL     [本文引用: 1]

Bao W, Cai Y.

Uniform error estimates offinite difference methods for the nonlinear Schrödinger equation with wave operator

SIAM Journal on Numerical Analysis, 2012, 50: 492-521

DOI:10.1137/110830800      URL     [本文引用: 1]

Wang T, Guo B, Xu Q.

Fourth-order compact and energy conservative difference schemes for the nonlinear Schrödinger equation in two dimensions

Journal of Computational Physics, 2013, 243: 382-399

DOI:10.1016/j.jcp.2013.03.007      URL     [本文引用: 1]

Akrivis G, Dougalis V, Karakashian Q.

On fully discrete Galerkin methods of second-order temporal accuracy for the nonlinear Schrödinger equation

Numerische Mathematik, 1991, 59: 31-53

DOI:10.1007/BF01385769      URL     [本文引用: 1]

Chen J, Qin M.

Multi-symplectic Fourier pseudospectral method for the nonlinear Schrödinger equation

Electronic Transactions on Numerical Analysis, 2001, 12: 193-204

[本文引用: 1]

Wang Y, Wang B, Qin M.

Local structure-preserving algorithms for partial differential equations

Science in China Series A: Mathematics, 2008, 51: 2115-2136

[本文引用: 2]

Cai J, Wang Y, Liang H.

Local energy-preserving and momentum-preserving algorithms for coupled nonlinear Schrödinger system

Journal of Computational Physics, 2013, 239: 30-50

DOI:10.1016/j.jcp.2012.12.036      URL     [本文引用: 1]

Wang J, Wang Y, Liang D.

Construction of local structure-preserving algorithms for the general multi-symplectic Hamiltonian system

Communications in Computational Physics, 2020, 27: 828-860

DOI:10.4208/cicp      URL     [本文引用: 2]

Wang J, Zhou Z, Wang Y.

Local structure-preserving algorithms for the Klein-Gordon-Zakharov equation

Acta Mathematica Scientia, 2023, 43: 1211-1238

DOI:10.1007/s10473-023-0313-2      [本文引用: 1]

Bridges T, Reich S.

Multi-symplectic integrators: numerical schemes for Hamiltonian PDEs that conserve symplecticity

Physics Letters A, 2001, 284: 184-193

DOI:10.1016/S0375-9601(01)00294-8      URL     [本文引用: 1]

Dalfovo F, Giorgini S, Pitaevskii L.

Theory of Bose-Einstein condensation in trapped gases

Review of Modern Physics, 1999, 71: 463-512

DOI:10.1103/RevModPhys.71.463      URL     [本文引用: 1]

Tornberg A, Engquist B.

Regularization techniques for numerical approximation of PDEs with singularities

Journal of Scientific Computing, 2003, 19: 527-552

DOI:10.1023/A:1025332815267      [本文引用: 1]

Qian X, Fu H, Song S.

Conservative modified Crank-Nicolson and time-splitting wavelet methods for modeling Bose-Einstein condensates in delta potentials

Applied Mathematics and Computation, 2017, 307: 1-16

DOI:10.1016/j.amc.2017.02.037      URL     [本文引用: 1]

Zhou X, Cai Y, Cai X.

Accurate and efficient numerical methods for the nonlinear Schrödinger equation with Dirac delta potential

Calcolo, 2023, 60(4): 1-28

DOI:10.1007/s10092-022-00496-z      [本文引用: 1]

LeVeque R, Li Z.

The immersed interface method for elliptic equation with discontinuous coefficients and singular sources

SIAM Journal on Numerical Analysis, 1994, 31: 1019-1044

DOI:10.1137/0731054      URL     [本文引用: 1]

Rutka V, Li Z.

An explicit jump immersed interface method for two-phase Navier-Stokes equations with interfaces

Computer Methods in Applied Mechanics and Engineering, 2008, 197: 2317-2328

DOI:10.1016/j.cma.2007.12.016      URL     [本文引用: 1]

Cheng B, Ya M, Chuan F.

Nonlinear Schrödinger equation with a Dirac delta potential: finite difference method

Communications in Theoretical Physics, 2020, 72: 1-6

[本文引用: 1]

Bai J, Wang L.

EJIIM for the stationary Schrödinger equations with delta potential wells

Applied Mathematics and Computation, 2015, 254: 113-124

DOI:10.1016/j.amc.2014.12.095      URL     [本文引用: 1]

Bai J, Li C, Liu X.

Weak multi-symplectic reformulation and geometric numerical integration for the nonlinear Schrödinger equations with delta potentials

IMA Journal of Numerical Analysis, 2018, 38(1): 399-429

DOI:10.1093/imanum/drw062      URL     [本文引用: 4]

Bai J.

Multi-symplectic Runge-Kutta-Nyström methods for nonsmooth nonlinear Schrödinger equations

Journal of Mathematical Analysis and Applications, 2016, 444: 721-736

DOI:10.1016/j.jmaa.2016.06.060      URL     [本文引用: 3]

Bai J, Ullah H, Li C.

Energy-preserving methods for non-smooth nonlinear Schrödinger equations

Applied Numerical Mathematics, 2023, 185: 188-202

DOI:10.1016/j.apnum.2022.11.017      URL     [本文引用: 1]

Witthaut D, Mossmann S, Korsch H.

Bound and resonance states of the nonlinear Schrödinger equation in simple modle systems

Journal of Physics A, 2005, 38: 1777-1792

DOI:10.1088/0305-4470/38/8/013      URL     [本文引用: 1]

/

〈 〉

版权所有 © 《数学物理学报》杂志社
地址：武汉市71010号信箱(430071)　电话: 027-87199206(中文版)　027-87199087(英文版)
E-mail: actams@wipm.ac.cn
本系统由北京玛格泰克科技发展有限公司设计开发