含自相位调制非线性薛定谔方程的同伦分析

Homotopy Analysis of Nonlinear Schrödinger Equation with Self-Phase Modulation

摘要: 基于同伦分析方法研究了含自相位调制的非线性薛定谔方程。该方程可以用来描述光信号在光纤传输过程中因损耗、色散等导致的体系非线性效应。求出了方程的孤子解和周期解，并讨论了体系的二维和三维演化行为。

Abstract: This paper investigates the nonlinear Schrödinger equation with self-phase modulation based on the homotopy analysis method. This equation can be used to describe the nonlinear effects of the system caused by loss, dispersion, and other factors during optical signal transmission in optical fibers. The soliton solutions and periodic solutions of the equation are obtained, and the two-dimensional and three-dimensional evolution behaviors of the system are presented.

1. 引言

自然界中许多问题都是非线性的，线性方程通常仅用来描述一些理想化的情况，而无论是科学研究还是工程实际中，由于各种因素的影响，无法达到理想状态，因此非线性方程在研究中更具广泛性。例如，在研究水流中的孤立波的运动时，就可以使用KdV方程[1] [2]、KP方程[3] [4]、非线性薛定谔方程[5] [6]等非线性微分方程。除了水波的研究，在其他领域也广泛采用非线性微分方程，比如在光纤[7] [8]、神经网络[9] [10]以及通信工程[11]，甚或经济学[12]等领域。过去很长一段时间内，求解非线性微分方程发展出许多十分成熟的有效方法，一类是数值分析[13] [14]，另一类包括Darboux变换[15] [16]、摄动法[17]和物理信息神经网络法[18]等。但对于更复杂的非线性现象，普通的方法在求解过程中总会遇到计算量大、计算时间长，乃至无法求解等问题。此外，渐近方法在传统上对于强非线性问题都是难以求得解析解的，当非线性变强时非线性问题的求解往往失效，因此，本文采用近来发展的同伦分析法(HAM) [19]-[21]。该方法的核心思想就是通过将原始的非线性问题转化为一系列简单的线性问题来逼近解。研究表明，在逻辑上同伦分析法还内蕴了如Adomian分解法[22]和Lyapunov人工小参数法[23]等一些传统的方法。因此，它可以被视为一种统一的、更具普遍性的一般方法。

本文使用同伦分析法研究含有自相位调制的非线性薛定谔方程[24]。非线性薛定谔方程在物理学中扮演者非常重要的角色。薛定谔方程是描述量子系统动态演化的核心方程。线性薛定谔方程一般适用于描述单一量子态的演化，而在实际情况下，量子系统之间存在的相互作用会导致波函数的非线性行为。因此，使用非线性薛定谔方程可以更好地描述更复杂的物理现象，如孤立子、波动和相互作用等。除此之外，在光纤中，由于光纤中光的传播是一个非线性的过程，光脉冲会因为非线性效应而变形，其中自相位调制(Self-Phase Modulation, SPM)就是一种重要的非线性效应，常见于光纤通信和光脉冲传输中，SPM机制与光信号的传播过程中的非线性效应密切相关，具体来说，自相位调制由非线性介质引起，当光波通过非线性介质(如光纤)传播时，光的强度变化会影响介质的折射率。即当光强增加时，折射率也会随之增加，这会导致光的相位发生变化。这个相位变化与光波的强度(光场)成正比，进而引起相位的非线性调制，导致频谱扩展、脉冲展宽、信号失真和相位噪声等。因此，就需要更复杂的信号处理算法来恢复信号。非线性薛定谔方程能够有效地描述这些脉冲的变化，包括形成孤立子，保持形状在传播过程中不变，从而减少信息传输的失真。由于光信号在光纤中传输演化时受损耗、色散以及非线性等多种物理效应影响，导致光脉冲在传输中会产生各种畸变。为了更准确地描述和解释这些非线性问题，我们通过求解非线性薛定谔方程[25]来进行研究。该方程的表达形式为：

$\frac{\partial E}{\partial z}+\frac{\alpha }{2}E+\frac{i}{2}{\beta }_1\frac{{\partial }^{2}E}{\partial t^2}-\frac{1}{6}{\beta }_2\frac{{\partial }^{3}E}{\partial t^3}-\gamma {\left|E\right|}^{2}E=0$ (1.1)

其中$E\left(z,t\right)$ 为光信号的复包络；$z$ 和$t$ 分别为传输距离和参照系中的时间量度。传播参数$\alpha$ ，${\beta }_1$ ，${\beta }_2$ ，$\gamma$ 分别反映了信号传输过程中的功率损耗、群速度色散(GVD)、三阶色散(TOD)、自相位调制(SPM)等效应。本文的研究分以下几个部分：在第2节进行方法介绍；在第3节用同伦分析法求解非线性薛定谔方程(1.1)的孤子解，随后在第4节求解周期解。第5节是结论。

2. 方法介绍

考虑一般形式表示的非线性演化方程

$N\left[u\left(r,t\right)\right]=0$ (2.1)

其中N为非线性算子，$u\left(r,t\right)$ 为未知函数，$r$ 和t分别代表空间和时间的自变量。方法的主要步骤如下

步骤1：构造简单同伦

$\begin{array}{l}H\left[\Phi \left(r,t;q\right);u_0\left(r,t\right),T\left(r,t\right),h,q\right]\hfill \\ =\left(1-q\right)\left\{L\left[\Phi \left(r,t;q\right)-u_0\left(r,t\right)\right]\right\}\hfill \\ -qhT\left(r,t\right)N\left[\Phi \left(r,t;q\right)\right]\hfill \end{array}$ (2.2)

其中辅助参数$h\ne 0$ ，辅助函数$T\left(r,t\right)\ne 0$ ，$L$ 为辅助线性算子。

步骤2：令同伦(2.3)为零，得到零阶形变方程

$\left(1-q\right)\left\{L\left[\Phi \left(r,t;q\right)-u_0\left(r,t\right)\right]\right\}=qhT\left(r,t\right)N\left[\Phi \left(r,t;q\right)\right]$ (2.3)

其中，$\Phi \left(r,t;q\right)$ 为上述方程的解，其解不仅依赖于初始猜测解$u_0\left(r,t\right)$ 、辅助线性算子$L$ 、辅助参数$h$ 和辅助函数$T\left(r,t\right)$ ，还依赖于嵌入变量$q\in \left[0,1\right]$ 。从(2.3)式中可以看到，当嵌入变量$q$ 从0增大到1时，$\Phi \left(r,t;q\right)$ 从初始猜测解$u_0\left(r,t\right)$ 连续变化到非线性方程的解$u\left(r,t\right)$ ，在同伦理论中，我们把这种连续变化的过程称为形变。

步骤3：定义m阶形变导数

$u_0^{\left(m\right)}\left(r,t\right)={\frac{{\partial }^m\Phi \left(r,t;q\right)}{\partial q^m}|}_{q=0}$ (2.4)

将$\Phi \left(r,t;q\right)$ 进行泰勒展开，并令

$u_m\left(r,t\right)=\frac{u_0^{\left(m\right)}\left(r,t\right)}{m!}=\frac{1}{m!}{\frac{{\partial }^m\Phi \left(r,t;q\right)}{\partial q^m}|}_{q=0}$ (2.5)

联合(2.4)和(2.5)得到

$\Phi \left(r,t;q\right)=u_0\left(r,t\right)+{\displaystyle \sum _{m=1}^{+\infty }{u}_m}\left(r,t\right)q^m$ (2.6)

当选取合适的初始猜测解$u_0\left(r,t\right)$ 、辅助线性算子$L$ 、非零辅助参数$h$ 和辅助函数$T\left(r,t\right)$ 时，就可以得到对应的级数解的表达式

$u\left(r,t\right)=u_0\left(r,t\right)+{\displaystyle \sum _{m=1}^{+\infty }{u}_m}\left(r,t\right)$ (2.7)

步骤4：将零阶形变方程(2.3)对$q$ 求导m次，再除以$m!$ ，再令$q=0$ ，得到m阶形变方程

$L\left[u_m\left(r,t\right)-X_m{u}_{m-1}\left(r,t\right)\right]=hT\left(r,t\right)R_m\left[u_{m-1}\left(r,t\right)\right]$ (2.8)

其中

$X_m=\{\begin{array}{l}0,m=1\hfill \\ 1,m=2,3,4,\cdots ,n\hfill \end{array}$ (2.9)

同时，有

$R_m\left(u_{m-1},r,t\right)=\frac{1}{\left(m-1\right)!}{\frac{{\partial }^{m-1}N\left[\Phi \left(r,t;q\right)\right]}{\partial q^{m-1}}|}_{q=0}$ (2.10)

需要注意的是，所有的高阶形变方程(2.8)都具有相同的线性算子$L$ ，并且$R_m\left(u_{m-1},r,t\right)$ 可通过定义的非线性算子$N$ 给出，同时，我们发现高阶形变方程的右端只依赖于矢量$u_{m-1}$ ，因此，当m从1到n依次取值时，就可以得到n阶近似解。

下面我们具体运用同伦分析法，探讨含自相位调制的光信号在光纤中的传输行为。

3. 孤子解

现在考虑方程(1.1)的孤子解。利用行波变换$E=E\left(z,t\right)=u\left(z,t\right)+iv\left(z,t\right)$ ，$\xi =z-\omega t$ ，其中$\omega$ 为常数。将时空变量$z$ 和$t$ 的导数都变换为新变量$\xi$ 的导数后，可得到一组$u$ 和$v$ 的耦合方程组，具体形式如下：

$\begin{array}{c}\frac{\partial u}{\partial \xi }+\frac{\alpha }{2}u-\frac{{\beta }_1}{2}{\omega }^2\frac{{\partial }^{2}v}{\partial {\xi }^2}+\frac{{\beta }_2}{6}{\omega }^3\frac{{\partial }^{3}u}{\partial {\xi }^3}+\gamma \left(u^2+v^2\right)v=0\\ \frac{\partial v}{\partial \xi }+\frac{\alpha }{2}v+\frac{{\beta }_1}{2}{\omega }^2\frac{{\partial }^{2}u}{\partial {\xi }^2}+\frac{{\beta }_2}{6}{\omega }^3\frac{{\partial }^{3}v}{\partial {\xi }^3}-\gamma \left(u^2+v^2\right)u=0\end{array}$ (3.1)

想要求得方程(1.1)的孤子解和周期解，我们必须得到$R_m\left[u_{m-1}\right]$ 和$R_m\left[v_{m-1}\right]$ 的表达式。首先求方程(1.1)的孤子解。

步骤1：设定初始猜测解为

$u_0\left(\xi \right)=2{\text{e}}^{-\xi }-{\text{e}}^{-2\xi },v_0\left(\xi \right)=2{\text{e}}^{-\xi }-{\text{e}}^{-2\xi }$ (3.2)

步骤2：设定辅助线性算子为

$L_u=\frac{{\partial }^{2}u}{\partial {\xi }^2}+u,L_v=\frac{{\partial }^{2}v}{\partial {\xi }^2}+v$ (3.3)

步骤3：根据上节步骤得到$R_m\left[u_{m-1}\right]$ 和$R_m\left[v_{m-1}\right]$ 的表达式为

$\begin{array}{c}{R}_m\left[u_{m-1}\right]={\displaystyle \sum _{i=0}^{m-1}\frac{\partial u_{m-1-i}}{\partial \xi }}+\frac{\alpha }{2}{\displaystyle \sum _{i=0}^{m-1}{u}_{m-1-i}}-\frac{{\beta }_1}{2}{\displaystyle \sum _{i=0}^{m-1}{\omega }_{m-1-i}}{\displaystyle \sum _{j=0}^i{\omega }_{i-j}}\frac{{\partial }^2{v}_j}{\partial {\xi }^2}\\ +\frac{{\beta }_2}{2}{\displaystyle \sum _{i=0}^{m-1}{\omega }_{m-1-i}}{\displaystyle \sum _{j=0}^i{\omega }_{i-j}}{\displaystyle \sum _{k=0}^j{\omega }_{j-k}}\frac{{\partial }^3{u}_k}{\partial {\xi }^3}+\gamma {\displaystyle \sum _{i=0}^{m-1}{v}_{m-1-i}}\\ \times \left({\displaystyle \sum _{j=0}^i{u}_{i-j}}{u}_j+{\displaystyle \sum _{j=0}^i{v}_{i-j}}{v}_j\right)\end{array}$ (3.4)

$\begin{array}{c}{R}_m\left[v_{m-1}\right]={\displaystyle \sum _{i=0}^{m-1}\frac{\partial v_{m-1-i}}{\partial \xi }}+\frac{\alpha }{2}{\displaystyle \sum _{i=0}^{m-1}{v}_{m-1-i}}+\frac{{\beta }_1}{2}{\displaystyle \sum _{i=0}^{m-1}{\omega }_{m-1-i}}{\displaystyle \sum _{j=0}^i{\omega }_{i-j}}\frac{{\partial }^2{u}_j}{\partial {\xi }^2}\\ +\frac{{\beta }_2}{2}{\displaystyle \sum _{i=0}^{m-1}{\omega }_{m-1-i}}{\displaystyle \sum _{j=0}^i{\omega }_{i-j}}{\displaystyle \sum _{k=0}^j{\omega }_{j-k}}\frac{{\partial }^3{v}_k}{\partial {\xi }^3}-\gamma {\displaystyle \sum _{i=0}^{m-1}{u}_{m-1-i}}\\ \times \left({\displaystyle \sum _{j=0}^i{u}_{i-j}}{u}_j+{\displaystyle \sum _{j=0}^i{v}_{i-j}}{v}_j\right)\end{array}$ (3.5)

步骤4：利用微分方程的解的性质，求出所需的未知参数。

通过计算，可得出参数$\omega$ 随收敛参数$h$ 变化的关系，如图1。

Figure 1. The change curve of parameter ω with convergence parameter h

图1. 收敛参数ω随收敛参数h的变化曲线

从图1中可以看到，当选取$\alpha =2$ ，${\beta }_1=0.1$ ，${\beta }_2=3.5$ ，$\gamma =1$ ，参数$\omega$ 随着$h$ 变化时在区间$\left(-0.1,0.7\right)$ 时出现一个平台，此时$\omega$ 的值恒定为$-1.17$ 。这种情况说明收敛参数$h$ 在区间$\left(-0.1,0.7\right)$ 时，能保证高阶近似解$u_m\left(\xi \right)$ 和$v_m\left(\xi \right)$ 是收敛的，再通过迭代的方法，就可求出高阶近似解$u_1\left(\xi \right)$ 和$v_1\left(\xi \right)$ 的形式

$\begin{array}{c}{u}_1\left(\xi \right)=\frac{2}{37}{\text{e}}^{-6\xi }h-\frac{6}{13}{\text{e}}^{-5\xi }h+\frac{24}{17}{\text{e}}^{-4\xi }h-\frac{8}{5}{\text{e}}^{-3\xi }h\\ -0.611323{\text{e}}^{-2\xi }h+1.20704{\text{e}}^{-\xi }h\\ v_1\left(\xi \right)=-\frac{2}{37}{\text{e}}^{-6\xi }h+\frac{6}{13}{\text{e}}^{-5\xi }h-\frac{24}{17}{\text{e}}^{-4\xi }h+\frac{8}{5}{\text{e}}^{-3\xi }h\\ -0.540589{\text{e}}^{-2\xi }h-0.0551302{\text{e}}^{-\xi }h\end{array}$ (3.6)

根据图1，当令参数$h=-0.01$ ，可得方程(1.1)的孤子图，如图2。

(a) (b)

Figure 2. Soliton graphs of equation (1.1) are solved by homotopy analysis method. (a) Evolution curve of homotopy soliton solution with ξ; (b) A three-dimensional diagram of the relationship between homotopy approximate solutions with z and t

图2. 利用同伦分析法求解方程(1.1)的孤子图。(a) 同伦孤子解随ξ的演化曲线；(b) 同伦近似解随z和t的变化关系三维图

可以发现，体系具有标准的孤子解，其行为符合标准的孤子特征，随横坐标具有一个包络，体现了振幅以及强度，说明我们这里采用的方法非常好得描述了非线性薛定谔方程的解。

4. 周期解

现在，我们考虑方程(1.1)的周期解。求解依照如下步骤。

步骤1：设定初始猜测解为

$u_0\left(\xi \right)=2\sqrt{2}\cos \left(\xi \right)-\sqrt{2}\cos \left(2\xi \right),v_0\left(\xi \right)=2\sqrt{2}\cos \left(\xi \right)-\sqrt{2}\cos \left(2\xi \right)$ (3.7)

步骤2：设定辅助线性算子为

$L_u=\frac{{\partial }^{2}u}{\partial {\xi }^2},L_v=\frac{{\partial }^{2}v}{\partial {\xi }^2}$ (3.8)

步骤3：根据上一小节的步骤，可以得到$R_m\left[u_{m-1}\right]$ 和$R_m\left[v_{m-1}\right]$ 的表达式为

$\begin{array}{c}{R}_m\left[u_{m-1}\right]={\displaystyle \sum _{i=0}^{m-1}\frac{\partial u_{m-1-i}}{\partial \xi }}+\frac{\alpha }{2}{\displaystyle \sum _{i=0}^{m-1}{u}_{m-1-i}}-\frac{{\beta }_1}{2}{\displaystyle \sum _{i=0}^{m-1}{\omega }_{m-1-i}}{\displaystyle \sum _{j=0}^i{\omega }_{i-j}}\frac{{\partial }^2{v}_j}{\partial {\xi }^2}\\ +\frac{{\beta }_2}{2}{\displaystyle \sum _{i=0}^{m-1}{\omega }_{m-1-i}}{\displaystyle \sum _{j=0}^i{\omega }_{i-j}}{\displaystyle \sum _{k=0}^j{\omega }_{j-k}}\frac{{\partial }^3{u}_k}{\partial {\xi }^3}+\gamma {\displaystyle \sum _{i=0}^{m-1}{v}_{m-1-i}}\\ \times \left({\displaystyle \sum _{j=0}^i{u}_{i-j}}{u}_j+{\displaystyle \sum _{j=0}^i{v}_{i-j}}{v}_j\right)\end{array}$ (3.9)

$\begin{array}{c}{R}_m\left[v_{m-1}\right]={\displaystyle \sum _{i=0}^{m-1}\frac{\partial v_{m-1-i}}{\partial \xi }}+\frac{\alpha }{2}{\displaystyle \sum _{i=0}^{m-1}{v}_{m-1-i}}+\frac{{\beta }_1}{2}{\displaystyle \sum _{i=0}^{m-1}{\omega }_{m-1-i}}{\displaystyle \sum _{j=0}^i{\omega }_{i-j}}\frac{{\partial }^2{u}_j}{\partial {\xi }^2}\\ +\frac{{\beta }_2}{2}{\displaystyle \sum _{i=0}^{m-1}{\omega }_{m-1-i}}{\displaystyle \sum _{j=0}^i{\omega }_{i-j}}{\displaystyle \sum _{k=0}^j{\omega }_{j-k}}\frac{{\partial }^3{v}_k}{\partial {\xi }^3}-\gamma {\displaystyle \sum _{i=0}^{m-1}{u}_{m-1-i}}\\ \times \left({\displaystyle \sum _{j=0}^i{u}_{i-j}}{u}_j+{\displaystyle \sum _{j=0}^i{v}_{i-j}}{v}_j\right)\end{array}$ (3.10)

步骤4：利用微分方程的解的性质，求出所需的未知参数。

通过计算，给出参数$\omega$ 随收敛参数$h$ 的变化图3。

Figure 3. The change curve of parameter ω with convergence parameter h

图3. 收敛参数ω随收敛参数h的变化而变化

从图3中可以看到，当选取$\alpha =0.2$ ，${\beta }_1=13.5$ ，${\beta }_2=-1$ ，$\gamma =-0.1$ 时，参数$\omega$ 随着$h$ 的增大先变大，然后进入一个平台，再变小。该平台处于区间$\left(-0.04,0.15\right)$ 时，$\omega$ 有一恒定值$-0.49$ 。即当收敛参数$h$ 在区间$\left(-0.04,0.15\right)$ 时，可以保证解$u_m\left(\xi \right)$ 和$v_m\left(\xi \right)$ 一定是收敛的，接着通过迭代的方法，求出解$u_1\left(\xi \right)$ 和$v_1\left(\xi \right)$ 的形式

$\begin{array}{c}{u}_1\left(\xi \right)=1.55815h\cos \left(\xi \right)-1.41421h\cos \left(2\xi \right)\\ -0.219989h\cos \left(3\xi \right)+0.106066h\cos \left(4\xi \right)\\ -0.339411h\cos \left(5\xi \right)+0.00392837h\cos \left(6\xi \right)\\ v_1\left(\xi \right)=-1.62886h\cos \left(\xi \right)+1.48492h\cos \left(2\xi \right)\\ +0.219989h\cos \left(3\xi \right)-0.106066h\cos \left(4\xi \right)\\ +0.339411h\cos \left(5\xi \right)-0.00392837h\cos \left(6\xi \right)\end{array}$ (3.11)

然后令$h=-0.01$ ，可得方程(1.1)的周期解图4。

(a) (b)

Figure 4. The period graph of equation (1.1) is solved by homotopy analysis method. (a) Evolution curve of homotopy periodic solution with ξ; (b) A three-dimensional diagram of the relationship between homotopy approximate solutions with z and t

图4. 利用同伦分析法求解方程(1.1)的周期图。(a) 同伦周期解随ξ的演化曲线；(b) 同伦近似解随z和t的变化关系三维图

我们可以发现，体系周期解的形式类似于正弦曲线的震荡，波函数幅度恒为正，最小值为零。在两个大的波峰之间还有两个小的波峰，其幅值约为大峰的1/3，但两个小峰间的极小值大于零，表明此时依然有光强。体系的震荡行为在纵方向并非完全对称但周期性出现。

5. 结论

在本研究中，我们采用同伦分析法讨论了含自相位调制的非线性薛定谔方程。作为一种行之有效的解析近似方法，同伦分析法的独特之处在于求解过程不依赖于任何小参数，为我们提供了极大的灵活性来选取初始猜测解和辅助参数，并且求解过程也体现出了一定的科学美感。这种方法的优势在于能够有效处理强非线性问题，特别是在传统微扰法失效的情况下。我们成功获得了非线性薛定谔方程的孤子解和周期解。孤子解描述了在非线性介质中能够保持形状不变的波包，而周期解则揭示了系统中的周期性行为。这些解析解不仅在理论上具有重要意义，还为我们理解复杂非线性系统的动力学行为提供了深刻洞察。

所得结果的意义体现在两个方面：首先，从理论角度来看，它们为非线性偏微分方程的求解提供了一种新的视角；其次，从应用角度来看，这些结果在光纤通信领域具有重要的实际价值。特别是在光纤中的脉冲传输问题上，我们的方法可以精确描述自相位调制等非线性效应对信号传播的影响，为优化光通信系统设计提供了理论基础。当然，该方法也具有它的一定局限性。例如我们在选取初始猜测解和辅助参数时过于自由，无法用一个严密的数学理论指导它们的选取。但相较于解决实际问题时所带来的便利，这有时又可忽略不计。期望本研究中运用的方法能够推广到更广泛的非线性问题领域。此外，研究过程也体现出，适当在数学公共课的教学中引入前沿知识，有助于创新性人才的自主探究及解决超出课本知识范围的研究问题。

基金项目

本文受到国家自然科学基金(项目号11804396)和2024年江苏省高校“高质量公共课教学改革研究”专项课题(项目号2024GZJX143)的资助。

NOTES

^*通讯作者。