| tianmen |
2011-06-12 18:33 |
求解光孤子或超短脉冲耦合方程的Matlab程序
计算脉冲在非线性耦合器中演化的Matlab 程序 oC7#6W:@w }?fa+FQGp % This Matlab script file solves the coupled nonlinear Schrodinger equations of 8dA/dMQ % soliton in 2 cores coupler. The output pulse evolution plot is shown in Fig.1 of @tj0Ir v % Youfa Wang and Wenfeng Wang, “A simple and effective numerical method for nonlinear ]#:xl}'LS % pulse propagation in N-core optical couplers”, IEEE Photonics Technology lett. Vol.16, No.4, pp1077-1079, 2004 _-!6@^+ FBY~Z$o0. %fid=fopen('e21.dat','w'); GzB%vsv95 N = 128; % Number of Fourier modes (Time domain sampling points) 6op\g].P M1 =3000; % Total number of space steps YD+C1*c! J =100; % Steps between output of space -+PPz?0 T =10; % length of time windows:T*T0 ,@!d%rL:4] T0=0.1; % input pulse width P)\f\yb MN1=0; % initial value for the space output location ^|K*lI/ dt = T/N; % time step ffB]4 n = [-N/2:1:N/2-1]'; % Index n9J>yud| t = n.*dt; _:K}DU'6 u10=1.*sech(1*t); % input to waveguide1 amplitude: power=u10*u10 (w^&NU'e u20=u10.*0.0; % input to waveguide 2
g8x8u| u1=u10; u2=u20; _$cBI_eA7 U1 = u1; _x'StD U2 = u2; % Compute initial condition; save it in U 8/F2V?iT ww = 4*n.*n*pi*pi/T/T; % Square of frequency. Note i^2=-1. 2 sOc]L:9 w=2*pi*n./T; ^7''x,I g=-i*ww./2; % w=2*pi*f*n./N, f=1/dt=N/T,so w=2*pi*n./T (i?^g & L=4; % length of evoluation to compare with S. Trillo's paper uB^]5sqfk dz=L/M1; % space step, make sure nonlinear<0.05 K5ph x for m1 = 1:1:M1 % Start space evolution N-Z 9
u1 = exp(dz*i*(abs(u1).*abs(u1))).*u1; % 1st sSolve nonlinear part of NLS Grub1=6l u2 = exp(dz*i*(abs(u2).*abs(u2))).*u2; vOj$-A--qU ca1 = fftshift(fft(u1)); % Take Fourier transform tQ|I$5jNJ ca2 = fftshift(fft(u2)); 5;Z~+$1 c2=exp(g.*dz).*(ca2+i*1*ca1.*dz); % approximation *&PgDAQ c1=exp(g.*dz).*(ca1+i*1*ca2.*dz); % frequency domain phase shift Sh*P^i.]+ u2 = ifft(fftshift(c2)); % Return to physical space <|_Ey)1
6 u1 = ifft(fftshift(c1)); A^Zs?<C- if rem(m1,J) == 0 % Save output every J steps. \mDm*UuG
U1 = [U1 u1]; % put solutions in U array {Eb6. U2=[U2 u2]; ie,{C MN1=[MN1 m1]; ;'Q{ ywr z1=dz*MN1'; % output location GkC88l9z end =@z"k'Vl` end C;ye%&g> hg=abs(U1').*abs(U1'); % for data write to excel xV6j6k ha=[z1 hg]; % for data write to excel ={Hbx>p t1=[0 t']; Mzd}9x$'J hh=[t1' ha']; % for data write to excel file *,pqpD> %dlmwrite('aa',hh,'\t'); % save data in the excel format t(yv figure(1) [~o3S$C&7 waterfall(t',z1',abs(U1').*abs(U1')) % t' is 1xn, z' is 1xm, and U1' is mxn xle29:?l figure(2) ?`XKaD!
f waterfall(t',z1',abs(U2').*abs(U2')) % t' is 1xn, z' is 1xm, and U1' is mxn Cnr48ukq ~;W%s 非线性超快脉冲耦合的数值方法的Matlab程序 5eLPn ~h/U ;Da 在研究脉冲在非线性耦合器中的演变时,我们需要求解非线性偏微分方程组。在如下的论文中,我们提出了一种简洁的数值方法。 这里我们提供给大家用Matlab编写的计算程序。 pi{ahuI#_o Youfa Wang and Wenfeng Wang, “A simple and effective numerical method for nonlinear pulse propagation in N-core optical couplers”, IEEE Photonics Technology lett. Vol.16, No.4, pp1077-1079, 2004 3IkG*enI LQjqwsuN{ 2P=;r:cx ;1 fM L,8 % This Matlab script file solves the nonlinear Schrodinger equations 'yNp J' % for 3 cores nonlinear coupler. The output plot is shown in Fig.2 of pDLo`F}A % Youfa Wang and Wenfeng Wang, “A simple and effective numerical method for nonlinear (`p(c;"*C! % pulse propagation in N-core optical couplers”, IEEE Photonics Technology lett. Vol.16, No.4, pp1077-1079, 2004 0d2%CsMS"D }VGiT~2$ C=1; ]VME`]t` M1=120, % integer for amplitude m+=!Z|K M3=5000; % integer for length of coupler D4|_?O3|m N = 512; % Number of Fourier modes (Time domain sampling points) qrkT7f dz =3.14159/(sqrt(2.)*C)/M3; % length of coupler is divided into M3 segments, make sure nonlinearity<0.05. i&$uG[&P T =40; % length of time:T*T0. 8f.La dt = T/N; % time step P33E\O n = [-N/2:1:N/2-1]'; % Index vn<S" t = n.*dt; C0%%@
2+ ww = 4*n.*n*pi*pi/T/T; % Square of frequency. Note i^2=-1. UPYM~c+} w=2*pi*n./T; hA8 zXk/'8 g1=-i*ww./2; X`b5h}c g2=-i*ww./2; % w=2*pi*f*n./N, f=1/dt=N/T,so w=2*pi*n./TP=0; -ZB"Yg$l g3=-i*ww./2; d#\n)eGr P1=0; >`:+d'Jv0 P2=0; v`&Z.9!Tz^ P3=1; FScQS.qF P=0; +0 MKh for m1=1:M1 m
C Ge*V} p=0.032*m1; %input amplitude Y"qY@` s10=p.*sech(p.*t); %input soliton pulse in waveguide 1 UVK"%kW#( s1=s10; g&>Hy!v, s20=0.*s10; %input in waveguide 2 {/<& s30=0.*s10; %input in waveguide 3 J%)2,szn0 s2=s20; !UTJ) & s3=s30; U\ued=H p10=dt*(sum(abs(s10').*abs(s10'))-0.5*(abs(s10(N,1)*s10(N,1))+abs(s10(1,1)*s10(1,1)))); zTLn*? %energy in waveguide 1 eW0=m:6 p20=dt*(sum(abs(s20').*abs(s20'))-0.5*(abs(s20(N,1)*s20(N,1))+abs(s20(1,1)*s20(1,1)))); meR2"JN' %energy in waveguide 2 _ LNPB$P p30=dt*(sum(abs(s30').*abs(s30'))-0.5*(abs(s30(N,1)*s30(N,1))+abs(s30(1,1)*s30(1,1)))); N6;Z\\&0^q %energy in waveguide 3 0, /x# for m3 = 1:1:M3 % Start space evolution .a*$WGb s1 = exp(dz*i*(abs(s1).*abs(s1))).*s1; % 1st step, Solve nonlinear part of NLS +QX>:z s2 = exp(dz*i*(abs(s2).*abs(s2))).*s2; TpgBS4q s3 = exp(dz*i*(abs(s3).*abs(s3))).*s3; ydQS"]\g sca1 = fftshift(fft(s1)); % Take Fourier transform OX4D' sca2 = fftshift(fft(s2)); F]YKYF'1I sca3 = fftshift(fft(s3)); Q\L5ZJ%y/ sc1=exp(g1.*dz).*(sca1+i*C*sca2.*dz); % 2nd step, frequency domain phase shift SK52.xXJ sc2=exp(g2.*dz).*(sca2+i*C*(sca1+sca3).*dz); fP58$pwu sc3=exp(g3.*dz).*(sca3+i*C*sca2.*dz); !\1 W*6U8; s3 = ifft(fftshift(sc3)); d&`j8O s2 = ifft(fftshift(sc2)); % Return to physical space ;L2bC3 s1 = ifft(fftshift(sc1)); tHbPd.^ end )\vHIXnfJ1 p1=dt*(sum(abs(s1').*abs(s1'))-0.5*(abs(s1(N,1)*s1(N,1))+abs(s1(1,1)*s1(1,1)))); or}*tSKX p2=dt*(sum(abs(s2').*abs(s2'))-0.5*(abs(s2(N,1)*s2(N,1))+abs(s2(1,1)*s2(1,1)))); L?x?+HPY. p3=dt*(sum(abs(s3').*abs(s3'))-0.5*(abs(s3(N,1)*s3(N,1))+abs(s3(1,1)*s3(1,1)))); 31& .Lnq P1=[P1 p1/p10]; j_@3a)[NY P2=[P2 p2/p10]; #~;8#!X P3=[P3 p3/p10]; x-&v|w ' P=[P p*p]; Jr)`shJ" end t[hocl/6 figure(1) Aw;vg/#~md plot(P,P1, P,P2, P,P3); `aL|qyrq# |+-i'N9 转自:http://blog.163.com/opto_wang/
|
|