| tianmen |
2011-06-12 18:33 |
求解光孤子或超短脉冲耦合方程的Matlab程序
计算脉冲在非线性耦合器中演化的Matlab 程序 !T{+s
T #-0e0 % This Matlab script file solves the coupled nonlinear Schrodinger equations of Bf utmI % soliton in 2 cores coupler. The output pulse evolution plot is shown in Fig.1 of YOl$sgg} % Youfa Wang and Wenfeng Wang, “A simple and effective numerical method for nonlinear @/z\p7e % pulse propagation in N-core optical couplers”, IEEE Photonics Technology lett. Vol.16, No.4, pp1077-1079, 2004 mZ+!8$1X >JpBX+]5m %fid=fopen('e21.dat','w'); R}nvSerVb N = 128; % Number of Fourier modes (Time domain sampling points) D:z'`v0j M1 =3000; % Total number of space steps e\%,\uV} J =100; % Steps between output of space xfYKUOp/ T =10; % length of time windows:T*T0 5\Q Tm; T0=0.1; % input pulse width aAg Qv* MN1=0; % initial value for the space output location Y^fw37b dt = T/N; % time step &jE\D^>ko n = [-N/2:1:N/2-1]'; % Index ,!#Am13 t = n.*dt; 7(Fas(j3 u10=1.*sech(1*t); % input to waveguide1 amplitude: power=u10*u10 6{h\CU}" u20=u10.*0.0; % input to waveguide 2 {wqT$( (< u1=u10; u2=u20; ={g)[:(C. U1 = u1; F&d!fEHU U2 = u2; % Compute initial condition; save it in U s<I)THC ww = 4*n.*n*pi*pi/T/T; % Square of frequency. Note i^2=-1. ;UQGi}?CD w=2*pi*n./T; R"B{IWQi g=-i*ww./2; % w=2*pi*f*n./N, f=1/dt=N/T,so w=2*pi*n./T ;uBGB
h< L=4; % length of evoluation to compare with S. Trillo's paper c4H6I~2Na dz=L/M1; % space step, make sure nonlinear<0.05 n7t}G'*Y!^ for m1 = 1:1:M1 % Start space evolution KF%BX~80C u1 = exp(dz*i*(abs(u1).*abs(u1))).*u1; % 1st sSolve nonlinear part of NLS hb`9Vn\-E u2 = exp(dz*i*(abs(u2).*abs(u2))).*u2; y=Y k$:-y ca1 = fftshift(fft(u1)); % Take Fourier transform )z[C= ca2 = fftshift(fft(u2)); [JOa^U= c2=exp(g.*dz).*(ca2+i*1*ca1.*dz); % approximation @:N8V[*u c1=exp(g.*dz).*(ca1+i*1*ca2.*dz); % frequency domain phase shift 7:4c\C0 u2 = ifft(fftshift(c2)); % Return to physical space )N.3Q1g- u1 = ifft(fftshift(c1)); Zbczbnj if rem(m1,J) == 0 % Save output every J steps. vk7IqlEQ U1 = [U1 u1]; % put solutions in U array Z(MZbzY7Hq U2=[U2 u2]; Rhc:szDU MN1=[MN1 m1]; ,rB(WKU z1=dz*MN1'; % output location !>48`o^ end HPtMp#`T end UC`h o%OBF hg=abs(U1').*abs(U1'); % for data write to excel ;\pr05 ha=[z1 hg]; % for data write to excel ![z2]L+TB t1=[0 t']; Cy-p1s hh=[t1' ha']; % for data write to excel file |lNp0b %dlmwrite('aa',hh,'\t'); % save data in the excel format fI1CT)0<e figure(1) 8m0*89HEu waterfall(t',z1',abs(U1').*abs(U1')) % t' is 1xn, z' is 1xm, and U1' is mxn /pF8S!,z figure(2) gC$_yd6m
L waterfall(t',z1',abs(U2').*abs(U2')) % t' is 1xn, z' is 1xm, and U1' is mxn WJ8i=MO67
q0ktABB 非线性超快脉冲耦合的数值方法的Matlab程序 =z. hJu ?o(284sV3 在研究脉冲在非线性耦合器中的演变时,我们需要求解非线性偏微分方程组。在如下的论文中,我们提出了一种简洁的数值方法。 这里我们提供给大家用Matlab编写的计算程序。 c/Pql!h+ 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 k]ZE j/y~ V Rv4p5 ?s, oH
ZX/FIxpy % This Matlab script file solves the nonlinear Schrodinger equations #M!u';bZ % for 3 cores nonlinear coupler. The output plot is shown in Fig.2 of ^_#wo" % Youfa Wang and Wenfeng Wang, “A simple and effective numerical method for nonlinear $~5H-wJ % pulse propagation in N-core optical couplers”, IEEE Photonics Technology lett. Vol.16, No.4, pp1077-1079, 2004 G@P;#l`(D #`y[75<n C=1; / /NV_^$y M1=120, % integer for amplitude '`^~Zy?c M3=5000; % integer for length of coupler tQ@7cjq8bA N = 512; % Number of Fourier modes (Time domain sampling points) P[fy dz =3.14159/(sqrt(2.)*C)/M3; % length of coupler is divided into M3 segments, make sure nonlinearity<0.05. 4L>8RiiQE; T =40; % length of time:T*T0. 8s22VL dt = T/N; % time step Xc[ym n = [-N/2:1:N/2-1]'; % Index QyCrz{/ t = n.*dt; ~Bl,_?CBr ww = 4*n.*n*pi*pi/T/T; % Square of frequency. Note i^2=-1. r)~?5d w=2*pi*n./T; {ccc[G?>.Q g1=-i*ww./2; ?5't1219 g2=-i*ww./2; % w=2*pi*f*n./N, f=1/dt=N/T,so w=2*pi*n./TP=0; :.=:N%3[ g3=-i*ww./2; XR",.3LD P1=0; 2XL^A[? P2=0; ^+-QY\N
j P3=1; R@grY:h P=0; I]n X6=j5 for m1=1:M1 wmV=GV8 d p=0.032*m1; %input amplitude 0#GnmH s10=p.*sech(p.*t); %input soliton pulse in waveguide 1 <mP_K^9c s1=s10; ,eTdQI; s20=0.*s10; %input in waveguide 2 !.%*Tp#k# s30=0.*s10; %input in waveguide 3 _*=4xmB.= s2=s20; 8\E=p+C s3=s30; 8Y% p10=dt*(sum(abs(s10').*abs(s10'))-0.5*(abs(s10(N,1)*s10(N,1))+abs(s10(1,1)*s10(1,1)))); 4 dHGU^#WZ %energy in waveguide 1 +)h# !/ p20=dt*(sum(abs(s20').*abs(s20'))-0.5*(abs(s20(N,1)*s20(N,1))+abs(s20(1,1)*s20(1,1)))); tYMr %energy in waveguide 2 [Cd#<Te3 p30=dt*(sum(abs(s30').*abs(s30'))-0.5*(abs(s30(N,1)*s30(N,1))+abs(s30(1,1)*s30(1,1)))); Jv
5l %energy in waveguide 3 L$a{%]I for m3 = 1:1:M3 % Start space evolution I;AS.y s1 = exp(dz*i*(abs(s1).*abs(s1))).*s1; % 1st step, Solve nonlinear part of NLS .z$UNB(!M s2 = exp(dz*i*(abs(s2).*abs(s2))).*s2; _`C|K>: s3 = exp(dz*i*(abs(s3).*abs(s3))).*s3; g<~ODMCO?W sca1 = fftshift(fft(s1)); % Take Fourier transform StR)O))I sca2 = fftshift(fft(s2)); wmK;0 )|H sca3 = fftshift(fft(s3)); 2N-p97"g sc1=exp(g1.*dz).*(sca1+i*C*sca2.*dz); % 2nd step, frequency domain phase shift P5dD& sc2=exp(g2.*dz).*(sca2+i*C*(sca1+sca3).*dz); W|;`R{<I% sc3=exp(g3.*dz).*(sca3+i*C*sca2.*dz); 6t<[- s3 = ifft(fftshift(sc3)); IY~I=} s2 = ifft(fftshift(sc2)); % Return to physical space 9JMf
T] s1 = ifft(fftshift(sc1)); V[^AV"V end [vBP,_Tjx p1=dt*(sum(abs(s1').*abs(s1'))-0.5*(abs(s1(N,1)*s1(N,1))+abs(s1(1,1)*s1(1,1)))); 1A(f_ 0,.Q p2=dt*(sum(abs(s2').*abs(s2'))-0.5*(abs(s2(N,1)*s2(N,1))+abs(s2(1,1)*s2(1,1)))); W }Ll)7(|T p3=dt*(sum(abs(s3').*abs(s3'))-0.5*(abs(s3(N,1)*s3(N,1))+abs(s3(1,1)*s3(1,1)))); B}y#AVSA P1=[P1 p1/p10]; nRHlHu P2=[P2 p2/p10]; b!QRD'31'j P3=[P3 p3/p10]; 4*n1Xu7^x P=[P p*p]; A[Ce3m end N1>M<N03 figure(1) fP;I{AiN~ plot(P,P1, P,P2, P,P3); 2nFr?Y3g, n68qxD-X 转自:http://blog.163.com/opto_wang/
|
|