孤立波的Matlab的数值计算和模拟

当V > 0, (N − 1)N < V < N (N + 1) 时，初始条件在 KdV
方程下将给出 N 个孤立子，但也有色散波，限于篇幅， 在此我们就不做演示，读者可自行演示。
由图形我们得到的这些结论，与由反散射方法理论 推导的结果是一致的。
2.5 孤立子在 KdV 方程下的相互作用 KdV 方程只能实现单向运动，速度不同且互相追赶的多孤
y(3:(n-2),3)=-c*(y(5:n,2)-2*y(4:(n-1),2)+… 2*y(2:(n-3),2)-y(1:(n-4),2))+y(3:(n-2),1)... -6*y(3:(n-2),2).*d.*(y(4:(n-1),2)-y(2:(n-3),2));
y(3:(n-2),1)=y(3:(n-2),2); y(3:(n-2),2)=y(3:(n-2),3); if mod(j,40)==0; set(h,'XData',x,'YData',y(:,2)); drawnow; end end 其运行结果如图一；
宽反比于 c 。由计算模拟得到的结果我们不难发现这
些结论的正确性。 下面我们进一步讨论孤立波形状保持不变的物理
内涵和物理意义。 2.2 波形的色散方程和色散关系
色散方程的形式各样，为了对应 KdV 方程，本文取 对应色散方程为：
ut + uxxx = 0;
我们同样利用 matlab 进行模拟（程序附后），为 了和 KdV 的方程结果有所比较，我们仍取初始条件为
− k2 )2 + k2 )2
− (λ1 − (λ1
− λ2 )2 − λ2 )2
；
利用孤波的解析解，我们进行模拟得到图十：
利用双孤立子的解析解，我们进行模拟得到图十一：
由图不难发现，KPⅡ不能产生波包，但是它有着与 KdV 很接近
的性质，对横向扰动的稳定性*。
图九 图 十一
图十
%线性色散图的差分解法 clear clf dt=0.0001;dx=0.1;c=dt/dx^3;x=-100:0.1:100; n=length(x); y(:,1)=2*(sech(x')).^2; plot([-20,10],[0,0],'r'); text(0,0,'水平面'); hold on h=plot(x,y(:,1),'linewidth',2); axis([-20,6,-1,3]) set(h,'EraseMode','xor','MarkerSize',18); y(3:(n-2),2)=-c/2*(y(5:n,1)-2*y(4:(n-1),1) +2*y(2:(n-3),1)-y(1:(n-4),1))+y(3:(n-2),1); for j=3:30000;
包的形状，且碰撞时也是非线性叠加。但和 KdV 方程 不同，它对横向的扰动不具有稳定性，具体可见*。 3.2 KPⅡ方程的模拟
KPⅡ与 KdV 方程一样，它拥有孤波解：
U (x, y,t) =
{ [ ( ) 2k 2 sec h2 k x + λy − 4k 2 + 3λ2 t + δ 0 ]}
双孤立子的解为：
有所不同，一般情况下是采用反散色方法进行解决*，但本文 将采用数值模拟的方式直观的得到一些结论。
w + k 3 = 0;
图三
图四
在程序 1 中，
取 u(x,0) = V sec h2 (x) = 2sec h2 (x) 时，得
到了图一，得到的是 1 个稳定的孤子。
现取 u(x,0) = V sec h2 (x) = 6sec h2 (x) 时，
( ) 改变初始值，令 u(x,0) = 6sec h2 3x ；
运行结果如图二；
图一
图二
两图比较可发现，v1>v2,即 c 越大，移动速度越快，
波高越大，但宽度越窄。这和行波解表达的物理意义是
相符的： c 和 x0 是由初始条件决定的常数, x0 为孤波 的初始位置， c 为孤波解的传播速度，波高为 c /2，波
孤立波的 matlab 的数值计算和模拟
Made by
1 前言 孤立波简称孤波，是一种特殊的波，又是不难在自
然界中出现的波，具有保持其波形和速度不变的特点。 能发生强烈的相互作用，但相互作用后仍能保持其各自 特点、形状、速度不变或只有一些位相改变的那些孤波 称为孤子。孤波被称为自然界的相干结构，混沌运动所 呈现的是非线性系统中奇妙的无序状态，相干结构则反 映了非线性系统中的惊人有序性，因而孤立子理论的产 生与发展是非线性偏微分方程研究中的一个重要组成 部分。虽然随着人们研究的深入，在通信领域，采用光 孤子技术大大提高了通信容量和通信距离，但作为非线 性的两个重要的不同部分，混沌确为更多人所熟知。这 可能是因为孤立子的研究和学习需要更高的数学和物 理知识，对数学物理方法的要求也比较高。但本文将以 matlab 数值模拟为主要方法，简单形象得讲述 KdV 方 程的性质和特点，并以此论证理论推导的结果。 2 KdV 方程的数值模拟 2.1 孤波运动的模拟
h
2
⎡ ⎢
2
⎢⎣
c 2
(x − ct − x0 )⎥⎤ ；
⎥⎦
因为 KdV 方程的行波解都为同向运动的，则在进行 单个孤子子运动的模拟时，我们把 KdV 方程的初始条
件设为： u(x,0) = 2sec h2 (x)；
利用差分法：
03 级 饶泽浪 and 吴辛烨
f ' = f1 + f −1 ; 2Δh
Baidu Nhomakorabea
f ''
=
f1 − 2 f + Δh 2
f −1 ;
f ''' = f 2 − 2 f1 + 2 f −1 − f −2 ; 2Δh 3
在 matlab6.5 下编程模拟，程序 1 如下： dt=0.0001;dx=0.1;c=dt/dx^3;d=dt/dx; x=-6:0.1:25;n=length(x); y(:,1)=2*(sech(x')).^2;%初始条件 h=plot(x(1:end),y(1:end,1),'linewidth',2); axis([-6,25,0,3]) set(h,'EraseMode','xor','MarkerSize',18); y(3:(n-2),2)=-c/2*(y(5:n,1)-2*y(4:(n-1),1)+… 2*y(2:(n-3),1)-y(1:(n-4),1))+y(3:(n-2),1)... -6*y(3:(n-2),1).*d.*(y(4:(n-1),1)-y(2:(n-3),1)); for j=1:100000;
观察汇聚效应和色散效应的模拟图，不难发现，色散效应 和汇聚效应的传播速度相反。这就是孤立波能稳定传播的原 因，一个线性的波动在传播时由于存在色散效应，波动不能持 续。只有当一个线性波动中存在非线性的汇聚效应时，且只有 这两种作用达到某种平衡时，才能出现波形稳定的孤立波。
2.4 利用数值模拟验证 KdV 的初值问题 KdV 方程在不同的初始条件下，孤子的数目和运动情况都
y(3:(n-2),3)=-c*(y(5:n,2)-2*y(4:(n-1),2) +2*y(2:(n-3),2)-y(1:(n-4),2))+y(3:(n-2),1);
y(3:(n-2),1)=y(3:(n-2),2); y(3:(n-2),2)=y(3:(n-2),3); if mod(j,40)==0; set(h,'XData',x,'YData',y(:,2)); drawnow; pause end end %非线性汇聚的差分解法 clear clf lbd=-6;dx=0.1;dt=0.0001;x_begin=-4;x_end=4; x=[x_begin:dx:x_end];xlength=length(x); u(:,1)=2*(sech(x')).^2; u(2:(xlength-1),2)=lbd.*(dt/(2*dx)).* u(2:(xlength-1),1).*(u(3:xlength,1)u(1:(xlength-2),1))+u(2:(xlength-1),1); h1=plot(x,u(:,1),'r','linewidth',2); axis([x_begin x_end -0.2 2.2]) set(h1,'ydata',u(:,2)); for j=3:850
子碰撞。为了达到这个条件，我们依据孤波解的物理意义，取 方程一的初始条件为：
( ) u(x,0) = 6sech2 3(x + 12) + 2sech2 ((x + 4))；
则可得到两孤立子的碰撞模拟，如图七； 再取初始条件为：
( ) u(x,0) = 6sech2 3(x + 21) +
( ) 4sech2 2(x + 14) + 2sech2 ((x + 8))
U (x,
y,t)
=
2
∂2 ∂x 2
ln
F ( x,
y,t)
；
F (x, y,t) = 1+ exp(η1) +
exp(η2) + exp(η1 +η2 + A12) ；
[ ( ) ] ηi = 2ki x + λi y − 4ki 2 + 3λi2 t + δi ；
exp( A12)
=
4(k1 4(k1
度 v 不等。则初始时刻由不同分波合成的波形，因分波的传播 速度不同，波形发生变化，发生弥散，最后消失。 2.3 波动中的非线性汇聚 我们取方程,
ut + 6uut = 0;
推导可见*。 利用 matlab 进行模拟（程序在后），仍取初始条件为
u(x,0) = 2sec h2 (x)；可得到图四：
由图不难发现波速随着振幅的增加而增大，所以随着波的 传播，波包的上端会越来越陡，最后坍塌，这就意味着某种汇 聚效应。
至此我们把孤立子在一维 KdV 方程下的形成过程，及其 在 KdV 方程下的一些最重要的性质，完全以 matlab 数值模拟 的方式做了形象且正确的介绍。
但是到此为止我们都是在一维的 KdV 方程下对孤立波进 行模拟。下面我们将分别模拟在二维的 KdV 方程（即 K-P 方 程）及 Boussinesq 方程下的孤波运动。
(Ut + 6UUx ) +Uxxx x − 3U yy = 0 ；
我们取它的初始条件*为：
u(x, y,0)
= 16
− 4(x −
x0 )2
+
16k
2 i
(
y
−
y0 )2
+ 1/ ki2
[4( x
−
x0 )2
+
16k
2 i
(
y
−
y0 )2
+1/
k
2 i
]2
运行程序，得到相应的图九： 由图我们不难发现，孤立子在 KPI 方程可以出现波
图四 （1,2,3）
图五 （1,2,3）
图六 （1,2,3）
图七 （1,2,3）
3 K-P 方程的数值模拟 作为二维的 KdV 方程，和 KdV 方程一样，K-P 方 程有着多种不同形式，但以两种形式为主。
3.1 KPI 方程的数值计算及模拟 为了和本文的 KdV 方程保持一致。我们取 KPI 方程为：
u(x,0) = 2sec h2 (x)；
运行结果如图三。 由图我们不难发现 初始时刻出现的波形会随着时间的推移发生变化， 发生弥散。其体现的物理内涵为： 色散方程的通解*为:
u(x,t) = ∑ Aexp[i(kx − wt)];
色散方程对应的色散关系*为：
则 v = w = −k 2 ; 我们不难发现对不同的 k 值，其传播速 k
，得到图五，不难发现这猜想的正确性，限于篇幅，读 者可以自行带入其它数值，即可发现结果确是如此。
现在取V = −1, 即 u(x,0) = −sec h2 (x) ，得到图
六。不难发现未产生孤波，只有色散波。当V < 0, 读者 取其它值仍会得到同样的结果，则当V < 0, 不能产生稳
定的孤立子。 用同样的方法，读者代入不同的初始值，则会发现
得到图八，其对应的就是三孤立子碰撞的情况。 由图七，图八，我们不难发现：孤立子在碰撞的时候不满
足一般线性波动的叠加原理。一般线性波动叠加遵循振幅叠加 原理，但孤立子的碰撞，由图可发现振幅非但没有‘加’，反 而‘减’了。碰撞过程就像大的孤子把小的孤子吞掉后，然后 又把小的孤子吐了出来，并且各自都毫发无伤。这种现象很显 然的是一种非线性的叠加，这也正是孤立波最重要的性质之 一。
相应的得到图五 可以发现，在此初始条件下，给出了 2 个孤子，
而 2 = 1× 2;6 = 2 × 3;
则 我 们 可 猜 想 当 V = N (N + 1) ， 即
u(x,0) = N ∗ (N + 1)sec h2 (x)， N 为正常数。在
此初始条件下，将得到 N 个孤立子。为了验证，我们取
u(x,0) = N ∗ (N + 1)sec h2 (x) = 3* 4sech2 (x)
1895 年荷兰数学家科特韦格(D.J.Korteweg)和德 弗里斯(G.de.Vries)在浅水和小振幅假定下得到了浅水 波单向运动方程。KdV 方程经过变换，可得到各中不同 的形式。为了方便，本文取以下形式：
ut + 6uut + uxxx = 0;
此方程的孤波解为*：
u(x,t) =
c
sec