计算脉冲在非线性耦合器中演化的Matlab 程序 |8^53*f ? & -L$B
% This Matlab script file solves the coupled nonlinear Schrodinger equations of
dC+WII`V % soliton in 2 cores coupler. The output pulse evolution plot is shown in Fig.1 of
r Q)?Bhf % Youfa Wang and Wenfeng Wang, “A simple and effective numerical method for nonlinear
ramYSX@ % pulse propagation in N-core optical couplers”, IEEE Photonics Technology lett. Vol.16, No.4, pp1077-1079, 2004
QS(aA*D *|W S, %fid=fopen('e21.dat','w');
DmzK* O{ N = 128; % Number of Fourier modes (Time domain sampling points)
lz1RAp0R " M1 =3000; % Total number of space steps
v$~1{}iI5 J =100; % Steps between output of space
! R r k T =10; % length of time windows:T*T0
} )D E T0=0.1; % input pulse width
I)7STzlMj. MN1=0; % initial value for the space output location
{jdtNtw dt = T/N; % time step
rywui10x* n = [-N/2:1:N/2-1]'; % Index
Q8-;w{% t = n.*dt;
%-9?rOr u10=1.*sech(1*t); % input to waveguide1 amplitude: power=u10*u10
RE)! b
u20=u10.*0.0; % input to waveguide 2
E%Tpby}^' u1=u10; u2=u20;
Z[9)
hGh U1 = u1;
(j<FS>## U2 = u2; % Compute initial condition; save it in U
xib?XzxGo ww = 4*n.*n*pi*pi/T/T; % Square of frequency. Note i^2=-1.
Aw?i6d w=2*pi*n./T;
Yf1&"WW4 g=-i*ww./2; % w=2*pi*f*n./N, f=1/dt=N/T,so w=2*pi*n./T
U-^qVlw L=4; % length of evoluation to compare with S. Trillo's paper
|w; hu] dz=L/M1; % space step, make sure nonlinear<0.05
X=C*PWa7 for m1 = 1:1:M1 % Start space evolution
l$[7 pM[ u1 = exp(dz*i*(abs(u1).*abs(u1))).*u1; % 1st sSolve nonlinear part of NLS
;I V u2 = exp(dz*i*(abs(u2).*abs(u2))).*u2;
/Z3 Mlm{ ca1 = fftshift(fft(u1)); % Take Fourier transform
QjT$.pU d ca2 = fftshift(fft(u2));
c_V^~hq c2=exp(g.*dz).*(ca2+i*1*ca1.*dz); % approximation
P" -*'q,9 c1=exp(g.*dz).*(ca1+i*1*ca2.*dz); % frequency domain phase shift
Ygeg[S!7 u2 = ifft(fftshift(c2)); % Return to physical space
|h^[/ u1 = ifft(fftshift(c1));
#3VOC#. if rem(m1,J) == 0 % Save output every J steps.
'%Fg+cZN\ U1 = [U1 u1]; % put solutions in U array
\NZ(Xk U2=[U2 u2];
# <?igtUO MN1=[MN1 m1];
Fw{:fFZC[ z1=dz*MN1'; % output location
&,DZ0xA end
;*{"|l qe end
nm#ISueh hg=abs(U1').*abs(U1'); % for data write to excel
) wZ;}O ha=[z1 hg]; % for data write to excel
]u5B]ZQnA t1=[0 t'];
D;1?IeS hh=[t1' ha']; % for data write to excel file
SI)QX\is8 %dlmwrite('aa',hh,'\t'); % save data in the excel format
pseN!7+or figure(1)
I8x,8}o>V waterfall(t',z1',abs(U1').*abs(U1')) % t' is 1xn, z' is 1xm, and U1' is mxn
T 6ihEb$C figure(2)
g49G 7sk waterfall(t',z1',abs(U2').*abs(U2')) % t' is 1xn, z' is 1xm, and U1' is mxn
MiK
-W '@0Z#A 非线性超快脉冲耦合的数值方法的Matlab程序 %3%bRP }yzCq+ 在研究脉冲在非线性耦合器中的演变时，我们需要求解非线性偏微分方程组。在如下的
论文中，我们提出了一种简洁的数值方法。 这里我们提供给大家用Matlab编写的计算程序。
]3D>ai? 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
N4HIQ\p Wg5<@=x!G ' ]bw37_U, 0#G@F5; < % This Matlab script file solves the nonlinear Schrodinger equations
ayGcc` % for 3 cores nonlinear coupler. The output plot is shown in Fig.2 of
'^|u\$&U % Youfa Wang and Wenfeng Wang, “A simple and effective numerical method for nonlinear
!eu\ShI % pulse propagation in N-core optical couplers”, IEEE Photonics Technology lett. Vol.16, No.4, pp1077-1079, 2004
&xWej2a! v ZiuElxKi C=1;
2RbK##`vC M1=120, % integer for amplitude
C ^IPddw> M3=5000; % integer for length of coupler
}/bxe0px N = 512; % Number of Fourier modes (Time domain sampling points)
=?3b3PZn dz =3.14159/(sqrt(2.)*C)/M3; % length of coupler is divided into M3 segments, make sure nonlinearity<0.05.
T)Y{>wT T =40; % length of time:T*T0.
oBS m>V dt = T/N; % time step
]qd$rX n = [-N/2:1:N/2-1]'; % Index
A+=K<e t = n.*dt;
?S<`*O
+ ww = 4*n.*n*pi*pi/T/T; % Square of frequency. Note i^2=-1.
h}y] Pt? w=2*pi*n./T;
Q]{ `m g1=-i*ww./2;
wi/qI(O! g2=-i*ww./2; % w=2*pi*f*n./N, f=1/dt=N/T,so w=2*pi*n./TP=0;
3<x1s2U g3=-i*ww./2;
;7>k[?'e P1=0;
Yy c fb P2=0;
<*Gd0 v% P3=1;
v]GQb P=0;
\1He 9~6 for m1=1:M1
V8hmfV~=]P p=0.032*m1; %input amplitude
9u;/l#?@T s10=p.*sech(p.*t); %input soliton pulse in waveguide 1
[.Rdq]w6 s1=s10;
^ 'ws/( s20=0.*s10; %input in waveguide 2
rT|wZz9$@ s30=0.*s10; %input in waveguide 3
\
z3>kvk s2=s20;
8w$q4fg0 s3=s30;
J# DN 2y < p10=dt*(sum(abs(s10').*abs(s10'))-0.5*(abs(s10(N,1)*s10(N,1))+abs(s10(1,1)*s10(1,1))));
&J\<"3 %energy in waveguide 1
4 KX\'K p20=dt*(sum(abs(s20').*abs(s20'))-0.5*(abs(s20(N,1)*s20(N,1))+abs(s20(1,1)*s20(1,1))));
(zX75QSKV %energy in waveguide 2
% M*2 j%6 p30=dt*(sum(abs(s30').*abs(s30'))-0.5*(abs(s30(N,1)*s30(N,1))+abs(s30(1,1)*s30(1,1))));
b%QcB[k[WB %energy in waveguide 3
Ya &\ b 6 for m3 = 1:1:M3 % Start space evolution
Z8ds`KZM s1 = exp(dz*i*(abs(s1).*abs(s1))).*s1; % 1st step, Solve nonlinear part of NLS
*.6m,QqJ( s2 = exp(dz*i*(abs(s2).*abs(s2))).*s2;
+-!2nk `"a s3 = exp(dz*i*(abs(s3).*abs(s3))).*s3;
`F$lO2 #k sca1 = fftshift(fft(s1)); % Take Fourier transform
]]NTvr sca2 = fftshift(fft(s2));
!%'"l{R sca3 = fftshift(fft(s3));
P~*'/!@ sc1=exp(g1.*dz).*(sca1+i*C*sca2.*dz); % 2nd step, frequency domain phase shift
e -Z ul.m sc2=exp(g2.*dz).*(sca2+i*C*(sca1+sca3).*dz);
[X 9zrGHt sc3=exp(g3.*dz).*(sca3+i*C*sca2.*dz);
5uX-onP\[ s3 = ifft(fftshift(sc3));
af'gk&% s2 = ifft(fftshift(sc2)); % Return to physical space
7NRm\%^q s1 = ifft(fftshift(sc1));
mndKUI}d end
$H`{wJ?2( p1=dt*(sum(abs(s1').*abs(s1'))-0.5*(abs(s1(N,1)*s1(N,1))+abs(s1(1,1)*s1(1,1))));
N;v]ypak p2=dt*(sum(abs(s2').*abs(s2'))-0.5*(abs(s2(N,1)*s2(N,1))+abs(s2(1,1)*s2(1,1))));
{kghZ ur p3=dt*(sum(abs(s3').*abs(s3'))-0.5*(abs(s3(N,1)*s3(N,1))+abs(s3(1,1)*s3(1,1))));
|=:<[FU P1=[P1 p1/p10];
u}.mJDL P2=[P2 p2/p10];
?IG[W+ M8 P3=[P3 p3/p10];
,u=+%6b)A P=[P p*p];
q?qH7={,eu end
"QvTn= figure(1)
:O7n*lwx plot(P,P1, P,P2, P,P3);
OtbPr F5 [:zP]l.| 转自：
http://blog.163.com/opto_wang/
