计算脉冲在非线性耦合器中演化的Matlab 程序 K~Nx;{{d OIJNOu I % This Matlab script file solves the coupled nonlinear Schrodinger equations of
K G<. s< % soliton in 2 cores coupler. The output pulse evolution plot is shown in Fig.1 of
sB` .G % Youfa Wang and Wenfeng Wang, “A simple and effective numerical method for nonlinear
o1lhVM`15 % pulse propagation in N-core optical couplers”, IEEE Photonics Technology lett. Vol.16, No.4, pp1077-1079, 2004
3N c#6VI w_qX~ d/ %fid=fopen('e21.dat','w');
0"}qND N = 128; % Number of Fourier modes (Time domain sampling points)
#0$fZ M1 =3000; % Total number of space steps
*ThP->&:( J =100; % Steps between output of space
/M!b3bmA T =10; % length of time windows:T*T0
XX&4OV,^%D T0=0.1; % input pulse width
eFKF9m MN1=0; % initial value for the space output location
8! eYax dt = T/N; % time step
RGEgYOO n = [-N/2:1:N/2-1]'; % Index
F3nYMf t = n.*dt;
MTXh-9DA u10=1.*sech(1*t); % input to waveguide1 amplitude: power=u10*u10
8k +^jj u20=u10.*0.0; % input to waveguide 2
!aQb
Kp u1=u10; u2=u20;
Rax]svc U1 = u1;
>|zMN$: U2 = u2; % Compute initial condition; save it in U
(;VlK#rnC ww = 4*n.*n*pi*pi/T/T; % Square of frequency. Note i^2=-1.
sbv2*fno5 w=2*pi*n./T;
| KtI:n4d g=-i*ww./2; % w=2*pi*f*n./N, f=1/dt=N/T,so w=2*pi*n./T
XM1;
>#kz L=4; % length of evoluation to compare with S. Trillo's paper
%9v l dz=L/M1; % space step, make sure nonlinear<0.05
J lp nR#@ for m1 = 1:1:M1 % Start space evolution
IC"Z.'Ph u1 = exp(dz*i*(abs(u1).*abs(u1))).*u1; % 1st sSolve nonlinear part of NLS
q"(b}3 u2 = exp(dz*i*(abs(u2).*abs(u2))).*u2;
lT^/ 8Z<g ca1 = fftshift(fft(u1)); % Take Fourier transform
/U26IbJ ca2 = fftshift(fft(u2));
cl04fqX c2=exp(g.*dz).*(ca2+i*1*ca1.*dz); % approximation
ibH!b S{ c1=exp(g.*dz).*(ca1+i*1*ca2.*dz); % frequency domain phase shift
KE[!{O^(a u2 = ifft(fftshift(c2)); % Return to physical space
"hi d3"G u1 = ifft(fftshift(c1));
*'w?j)}A9g if rem(m1,J) == 0 % Save output every J steps.
_=Z?5{7S > U1 = [U1 u1]; % put solutions in U array
*Xcqnu(' U2=[U2 u2];
&cGa~#-u MN1=[MN1 m1];
y>^ FKN/ z1=dz*MN1'; % output location
2nf<RE> end
m^%@bu, end
;
DXsPp ZC hg=abs(U1').*abs(U1'); % for data write to excel
j+9;Rv t 2 ha=[z1 hg]; % for data write to excel
&&% oazR= t1=[0 t'];
ig x~6G* hh=[t1' ha']; % for data write to excel file
=U7P\s w2 %dlmwrite('aa',hh,'\t'); % save data in the excel format
) >te|@}o figure(1)
"7q!u,u waterfall(t',z1',abs(U1').*abs(U1')) % t' is 1xn, z' is 1xm, and U1' is mxn
}1
,\ *)5 figure(2)
Upa F>,kM waterfall(t',z1',abs(U2').*abs(U2')) % t' is 1xn, z' is 1xm, and U1' is mxn
?wP /l `^ZhxFX 非线性超快脉冲耦合的数值方法的Matlab程序 {8I,uQO Ptm=c6H(' 在研究脉冲在非线性耦合器中的演变时，我们需要求解非线性偏微分方程组。在如下的
论文中，我们提出了一种简洁的数值方法。 这里我们提供给大家用Matlab编写的计算程序。
'8Cg2v5&w 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
{o SdVRI dBw 7l} 3{)!T;W d
2##; [ % This Matlab script file solves the nonlinear Schrodinger equations
GQ(*k)'a % for 3 cores nonlinear coupler. The output plot is shown in Fig.2 of
H +' 6*akV % Youfa Wang and Wenfeng Wang, “A simple and effective numerical method for nonlinear
&@K6;T % pulse propagation in N-core optical couplers”, IEEE Photonics Technology lett. Vol.16, No.4, pp1077-1079, 2004
FI,K 0sO/| e%s1D C=1;
_h+7 K K M1=120, % integer for amplitude
GKc? M3=5000; % integer for length of coupler
D V\7KKJE N = 512; % Number of Fourier modes (Time domain sampling points)
Fr~\ZL dz =3.14159/(sqrt(2.)*C)/M3; % length of coupler is divided into M3 segments, make sure nonlinearity<0.05.
|LW5dtQ T =40; % length of time:T*T0.
x<h|$$4S dt = T/N; % time step
oam$9 q n = [-N/2:1:N/2-1]'; % Index
~x7 CI t = n.*dt;
)T6:@n^]h ww = 4*n.*n*pi*pi/T/T; % Square of frequency. Note i^2=-1.
N a$.VT w=2*pi*n./T;
5vFM0 g1=-i*ww./2;
+ -uQ] ^n g2=-i*ww./2; % w=2*pi*f*n./N, f=1/dt=N/T,so w=2*pi*n./TP=0;
-T}r$A g3=-i*ww./2;
/qKA1-R}4
P1=0;
Wv|CJN;4 P2=0;
mqHcD8X P3=1;
{#st>%i P=0;
Atb`Q'Yrw for m1=1:M1
xax[# Vl4 p=0.032*m1; %input amplitude
SwsJ<Dq^z s10=p.*sech(p.*t); %input soliton pulse in waveguide 1
~s-bA#0S s1=s10;
^&D5J\][ s20=0.*s10; %input in waveguide 2
A!,c@Kv
3 s30=0.*s10; %input in waveguide 3
0BNH~,0u s2=s20;
x <a}*8" s3=s30;
,4S[<(T" p10=dt*(sum(abs(s10').*abs(s10'))-0.5*(abs(s10(N,1)*s10(N,1))+abs(s10(1,1)*s10(1,1))));
h/oun2C %energy in waveguide 1
j,Mb l"P p20=dt*(sum(abs(s20').*abs(s20'))-0.5*(abs(s20(N,1)*s20(N,1))+abs(s20(1,1)*s20(1,1))));
k-H6 c %energy in waveguide 2
*^%+PQ p30=dt*(sum(abs(s30').*abs(s30'))-0.5*(abs(s30(N,1)*s30(N,1))+abs(s30(1,1)*s30(1,1))));
(/2rj[F& %energy in waveguide 3
cRH(@b
Xr for m3 = 1:1:M3 % Start space evolution
B `.aQ s1 = exp(dz*i*(abs(s1).*abs(s1))).*s1; % 1st step, Solve nonlinear part of NLS
DXG`% <ZMn s2 = exp(dz*i*(abs(s2).*abs(s2))).*s2;
X{Fr s3 = exp(dz*i*(abs(s3).*abs(s3))).*s3;
~n8U N< sca1 = fftshift(fft(s1)); % Take Fourier transform
whYk"N sca2 = fftshift(fft(s2));
xT+# K5 sca3 = fftshift(fft(s3));
v-N4&9)%9 sc1=exp(g1.*dz).*(sca1+i*C*sca2.*dz); % 2nd step, frequency domain phase shift
/lbj!\~ sc2=exp(g2.*dz).*(sca2+i*C*(sca1+sca3).*dz);
e`co:HO `# sc3=exp(g3.*dz).*(sca3+i*C*sca2.*dz);
8o[gzW:Q)U s3 = ifft(fftshift(sc3));
HU'w[r 6a s2 = ifft(fftshift(sc2)); % Return to physical space
'j*Q s1 = ifft(fftshift(sc1));
>-\^ )z end
etT9}RbQ p1=dt*(sum(abs(s1').*abs(s1'))-0.5*(abs(s1(N,1)*s1(N,1))+abs(s1(1,1)*s1(1,1))));
cpl Ny?UIC p2=dt*(sum(abs(s2').*abs(s2'))-0.5*(abs(s2(N,1)*s2(N,1))+abs(s2(1,1)*s2(1,1))));
k>F!S`a&m p3=dt*(sum(abs(s3').*abs(s3'))-0.5*(abs(s3(N,1)*s3(N,1))+abs(s3(1,1)*s3(1,1))));
w> 8HS+ P1=[P1 p1/p10];
sZ~03QvkT P2=[P2 p2/p10];
+_ /ys! P3=[P3 p3/p10];
w,X)g{^T P=[P p*p];
)Nqx=ms[(! end
@`)>- k figure(1)
iZ> P>x\ plot(P,P1, P,P2, P,P3);
n-2!<`UFX !@])Ut@tN 转自：
http://blog.163.com/opto_wang/
