非线性薛定谔方程数值解的matlab仿真.docx

ADMIN非线性薛定谔方程数值解的MATLAB仿真利用分步快速傅里叶变换对光纤中光信号的传输方程进行数值求解1、非线性薛定谔方程非线性薛定谔方程NONLINEARSCHRODINGEREQUATION，NLSE是奥地利物理学家薛定谔于1926年提出的，应用在量子力学系统中。由于量子力学主要研究粒子的动力学运动状态，所以不能运用牛顿力学公式来表示。通常在量子力学中，研究系统的状态一般通过波函数X，T来表示。而对波函数的研究主要是求解非线性薛定谔方程。本文主要研究光脉冲在光纤中传输状态下的演变。一般情况下，光脉冲信号在光纤中传输时，同时受到光纤的色散和非线性效应的影响。通过MAXWELL方程，考虑到光纤的色散和非线性效应，可以推导出光信号在光纤中的传输方程，即非线性薛定谔方程。NLSE是非线性偏微分方程，一般很难直接求出解析解，于是通过数值方法进行求解。具体分为两大类1分布有限差分法SPLITSTEPFINITEDIFFERENCE，SSFD；2分步傅里叶变换法SPLITSTEPFOURIERTRANS，SSFT。一般情况，在达到相同精度，由于分步傅里叶变换法采用运算速度快的快速傅里叶变换，所以相比较有限差分法运算速度快一到两个数量级。于是本文介绍分步傅里叶变换法来对光纤中光信号的传输方程，即非线性薛定谔方程进行数值求解。并通过MATLAB软件对结果数值仿真。非线性薛定谔方程的基本形式为2|TXIUU其中U是未知的复值函数目前，采用分步傅立叶算法SPLITSTEPFOURIER求解非线性薛定谔方程的数值解应用比较多。分步傅立叶方法最早是在1937年开始应用的，这种方法己经被证明是相同精度下数值求解非线性薛定愕方程最快的方法，部分原因是它采用了快速傅立叶变换算法FASTFOURIERTRANSALGORITHM。基于MATLAB科学计算软件以及MATLAB强大的符号计算功能，完全可以实现分步傅立叶数值算法来对脉冲形状和频谱进行仿真。一般情况下，光脉冲沿光纤传播时受到色散和非线性效应的共同作用，假设当传输距离很小的时候，两者相互独立作用，那么，根据这种思想可建立如下分步傅立叶数值算法的数学模型(I)把待求解的非线性薛定谔方程写成以下形式ˆUDNZ其中是线性算符，代表介质的色散和损耗，是非线性算符，它决定了脉冲传输ˆDˆ过程中光纤的非线性效应。一般来讲，沿光纤的长度方向，色散和非线性是同时作用的。分步傅立叶法假设在传输过程中，光场每通过一小段距离H，色散和非线性效应可以分别作用，得到近似结果。也就是说脉冲从Z到ZH的传输过程中分两步进行。第一步，只有非线性作用，方程II式中的0；第二步，再考虑线性作用，方程II式中的0ˆDˆN这样方程2在这两步中可分别简化为ˆUDZ得到了上面两个方程(III)，就可以分别求解非线性作用方程和线性作用方程，然后讨论分步傅立叶法的数值算法。由于方程(III)是一个偏微分方程，需要通过傅立叶变换把偏微分方程转换为代数方程，进行运算。傅立叶变换的定义如下1,,,EXP1,,,2FUZTZUZTIDIT在计算时一般采用快速傅立叶变换(FFT)算。为了保证精度要求，一般,Z还需要反复调整纵向传输步长Z和横向脉冲取样点数T来保证计算精度。2、分步傅立叶数值算法的MATLAB实现现待求解的非线性薛定谔方程如下2204AIIAZT其中，AZ,T是光场慢变复振幅，Z是脉冲沿光纤传播的距离；(II)(III)(IV)(V),,VG是群速度；是色散系数；是非线性系数；1TTZ//PSKM1/WKM是光纤损耗系数，它与用分贝表示的损耗系数的关系为/KMDB43DB首先，可以将方程(V)归一化振幅，是入射脉冲的峰值功0,/UAZTP率，此时方程V可改写为2204IIZT为了使用分步傅立叶法求解方程VI，将方程VI写成以下形式ˆUDNZ进一步，可以得出如下方程(VII)220ˆITNPIU然后，按照步骤1和步骤2，依次计算方程(VII)的线性算符和非线性算符。最后在步骤3中，运行步骤1和步骤2的MATLAB程序，得出线性算符和非线性算符的精确数值解及其仿真曲线。步骤1线性算符方程的求解线性算符的方程如下2IUTZ用傅立叶变换方法，得到一个常微分方程(IX)24IUZ解方程IX得2,0,EXPIZZ式中是初值的傅立叶变换，将进行反傅立叶变换就得到了0,U,T,U(VI)(VII)(VIII)(IX)(X)。方程X的求解公式为,UZT2,{EXP0,}IZUZTFFUT其中和分别表示傅立叶变换和反傅立叶变换运算。F步骤2非线性算符方程的求解非线性部分的方程如下20UPIZ同STEP1的方法，解方程(XII)，得到20,,EXP,ITZ式中是初值的傅立叶变换，将进行反傅立叶变换就得到了0,U0TZ。方程XIII的求解公式为ZT20,{EXP,0,}ZFPIUTZFT其中和分别表示傅立叶变换和反傅立叶变换运算。F步骤3算法在MATLAB中的实现在MATLAB中，设有限时长序列的长度为，它对应于一个频域内的长XN1NN度为N的有限长序列，的角频，其中T是1XKXK21KN序列的采样时间间隔XN这种正反离散傅立叶变换的关系式为12EXP1NJJXKDFTXNJKNXNIKXJNN然后用MATLAB中的离散傅立叶变换DFT函数FFT和离散傅立叶反变换IDFT的函数IFFT来实现方程VIII，XII式中的傅立叶和反傅立叶变换运算。进一步，得到方程XI，XIV的数值解及仿真曲线。最后，通过测试一组参数，得到方程(V)在该算法下的MATLAB运算结果。MATLAB总共用时3426S，求得的的结果曲线如下图所示。结果表明，算法正确而且精度(XI)(XII)(XIII)(XV)(XIV)也比较高，能够在非线性薛定谔方程的求解中广泛应用。附录MATLAB的脚本文件源代码PO200输入光强，单位WALPHA35光纤损耗值，单位DB/KMGAMMA20光纤非线性参数TO1初始脉冲宽度，单位秒C50第一次计算输入的啁啾参数B21000波数的倒数CPUTIME0TICLN1ISQRT1PI31415926535ALPHALPHA/4343LDTO2/ABSB2扩散长度，单位是MAOSQRTPO光振幅TAU4096E121E124095E12DT1E12H1000步长FORII010115不同的光纤长度不同，这个量可变ZIILDUAOEXP1IC/2TAU/TO2FIGURE1PLOTABSU, R TITLE PULSE XLABEL TIME YLABEL AMPLITUDE GRIDONHOLDONLMAXSIZEUFWHM1FINDABSUABSMAXU/2FWHM1LENGTHFWHM1DW1/L/DT2PIW1L/21L/21DWUFFTSHIFTU零延迟对中的谱WFFTSHIFTW零延迟对中的谱SPECTRUMF