计算脉冲在非线性耦合器中演化的Matlab 程序 4E>/*F! J! TK*\a2 % This Matlab script file solves the coupled nonlinear Schrodinger equations of
,P; a/{U % soliton in 2 cores coupler. The output pulse evolution plot is shown in Fig.1 of
Z%HEn$t % Youfa Wang and Wenfeng Wang, “A simple and effective numerical method for nonlinear
Fh!!T%5>C % pulse propagation in N-core optical couplers”, IEEE Photonics Technology lett. Vol.16, No.4, pp1077-1079, 2004
2*7s 9g ym,Ot 1 %fid=fopen('e21.dat','w');
]9yA0,z/ N = 128; % Number of Fourier modes (Time domain sampling points)
<DlanczziF M1 =3000; % Total number of space steps
Zy+QA>d| J =100; % Steps between output of space
2I(@aB+ T =10; % length of time windows:T*T0
v BeU T0=0.1; % input pulse width
/x8C70W^ MN1=0; % initial value for the space output location
C[<\ufclD dt = T/N; % time step
j}?ZsnqV n = [-N/2:1:N/2-1]'; % Index
h
C`p<jp/ t = n.*dt;
XL&eJ u10=1.*sech(1*t); % input to waveguide1 amplitude: power=u10*u10
\$ \(9!= u20=u10.*0.0; % input to waveguide 2
t;f
p<z7N. u1=u10; u2=u20;
\g6 # MNW U1 = u1;
JjO/u>A3;7 U2 = u2; % Compute initial condition; save it in U
Ud(d Wj-/ ww = 4*n.*n*pi*pi/T/T; % Square of frequency. Note i^2=-1.
wqoN@d w=2*pi*n./T;
UmI@":|- g=-i*ww./2; % w=2*pi*f*n./N, f=1/dt=N/T,so w=2*pi*n./T
uKAHJ$% L=4; % length of evoluation to compare with S. Trillo's paper
}m
l bN0v dz=L/M1; % space step, make sure nonlinear<0.05
&b ]KMAo3 for m1 = 1:1:M1 % Start space evolution
f'yd {ihFp u1 = exp(dz*i*(abs(u1).*abs(u1))).*u1; % 1st sSolve nonlinear part of NLS
( L u2 = exp(dz*i*(abs(u2).*abs(u2))).*u2;
JI; i1@| b ca1 = fftshift(fft(u1)); % Take Fourier transform
J-{E`ibGN ca2 = fftshift(fft(u2));
GKDG 5u; c2=exp(g.*dz).*(ca2+i*1*ca1.*dz); % approximation
]mU*Y:< c1=exp(g.*dz).*(ca1+i*1*ca2.*dz); % frequency domain phase shift
^?VT y5yp u2 = ifft(fftshift(c2)); % Return to physical space
~,4Z nuin u1 = ifft(fftshift(c1));
"@|V.d@ if rem(m1,J) == 0 % Save output every J steps.
p q5H{ U1 = [U1 u1]; % put solutions in U array
NOr*+N\ U2=[U2 u2];
p2?+[d MN1=[MN1 m1];
L}pFb@ z1=dz*MN1'; % output location
vK >^#b3 end
I:7,CV end
q q{N; C hg=abs(U1').*abs(U1'); % for data write to excel
~
a&j4E ha=[z1 hg]; % for data write to excel
'bO? =+c t1=[0 t'];
;lt;]7 hh=[t1' ha']; % for data write to excel file
*zht(~% %dlmwrite('aa',hh,'\t'); % save data in the excel format
Z'kYf figure(1)
vw
2@}#\: waterfall(t',z1',abs(U1').*abs(U1')) % t' is 1xn, z' is 1xm, and U1' is mxn
h --!pE+ figure(2)
w`_9 *AF9 waterfall(t',z1',abs(U2').*abs(U2')) % t' is 1xn, z' is 1xm, and U1' is mxn
Jxp'.oo[ m 1;jS| 非线性超快脉冲耦合的数值方法的Matlab程序 p7tC~]r:L 5ZxBmQ 在研究脉冲在非线性耦合器中的演变时，我们需要求解非线性偏微分方程组。在如下的
论文中，我们提出了一种简洁的数值方法。 这里我们提供给大家用Matlab编写的计算程序。
FeMu`| 2 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
R>q'Y mu ~ \ 2b9A' d> )Me&xQTn > `M\xt % This Matlab script file solves the nonlinear Schrodinger equations
9] \vw % for 3 cores nonlinear coupler. The output plot is shown in Fig.2 of
wH< * % Youfa Wang and Wenfeng Wang, “A simple and effective numerical method for nonlinear
ZQ/5]]}3y % pulse propagation in N-core optical couplers”, IEEE Photonics Technology lett. Vol.16, No.4, pp1077-1079, 2004
,{a t?y* n]!H,Q1,T C=1;
jnY4(B
M1=120, % integer for amplitude
35T7g65; M3=5000; % integer for length of coupler
yhmW-#+^e N = 512; % Number of Fourier modes (Time domain sampling points)
1[ Pbs b dz =3.14159/(sqrt(2.)*C)/M3; % length of coupler is divided into M3 segments, make sure nonlinearity<0.05.
z_TK
(;j T =40; % length of time:T*T0.
)M~5F,) dt = T/N; % time step
F Te# @ \I n = [-N/2:1:N/2-1]'; % Index
7Jk.U=vY t = n.*dt;
|99eD gK, ww = 4*n.*n*pi*pi/T/T; % Square of frequency. Note i^2=-1.
|z E7 W w=2*pi*n./T;
P+a&R<Dj4 g1=-i*ww./2;
,*30Q g2=-i*ww./2; % w=2*pi*f*n./N, f=1/dt=N/T,so w=2*pi*n./TP=0;
AXFVsZH"zi g3=-i*ww./2;
s0SB!-Vjm P1=0;
L<n_}ucA P2=0;
YeVhWPn@ P3=1;
r|+Zni] P=0;
RA}PM?D/ for m1=1:M1
9z#IdY$a p=0.032*m1; %input amplitude
)? xg =o/? s10=p.*sech(p.*t); %input soliton pulse in waveguide 1
23PSv8;EM s1=s10;
EifY K s20=0.*s10; %input in waveguide 2
|j;`;"+B s30=0.*s10; %input in waveguide 3
_B2t|uQ s2=s20;
+x` tvo s3=s30;
XB?!V|bno p10=dt*(sum(abs(s10').*abs(s10'))-0.5*(abs(s10(N,1)*s10(N,1))+abs(s10(1,1)*s10(1,1))));
Y+E@afsKs %energy in waveguide 1
|kn}iA@72p p20=dt*(sum(abs(s20').*abs(s20'))-0.5*(abs(s20(N,1)*s20(N,1))+abs(s20(1,1)*s20(1,1))));
]TQjk{X< %energy in waveguide 2
^U1;5+2G+~ p30=dt*(sum(abs(s30').*abs(s30'))-0.5*(abs(s30(N,1)*s30(N,1))+abs(s30(1,1)*s30(1,1))));
*UTk. :G5 %energy in waveguide 3
S9.jc@#.` for m3 = 1:1:M3 % Start space evolution
}v: h EMO s1 = exp(dz*i*(abs(s1).*abs(s1))).*s1; % 1st step, Solve nonlinear part of NLS
oq|K:<l s2 = exp(dz*i*(abs(s2).*abs(s2))).*s2;
Y9Pb s3 = exp(dz*i*(abs(s3).*abs(s3))).*s3;
E'C[+iK6, sca1 = fftshift(fft(s1)); % Take Fourier transform
/w| YNDA]j sca2 = fftshift(fft(s2));
,yC~{ H sca3 = fftshift(fft(s3));
te` 4*t sc1=exp(g1.*dz).*(sca1+i*C*sca2.*dz); % 2nd step, frequency domain phase shift
h0GXN\xI sc2=exp(g2.*dz).*(sca2+i*C*(sca1+sca3).*dz);
tIg_cY_y sc3=exp(g3.*dz).*(sca3+i*C*sca2.*dz);
|i| O9^*% s3 = ifft(fftshift(sc3));
@? t) UE s2 = ifft(fftshift(sc2)); % Return to physical space
Q5Wb ) s1 = ifft(fftshift(sc1));
@E}4LTB end
Mqna0"IYx* p1=dt*(sum(abs(s1').*abs(s1'))-0.5*(abs(s1(N,1)*s1(N,1))+abs(s1(1,1)*s1(1,1))));
ej0q* TH. p2=dt*(sum(abs(s2').*abs(s2'))-0.5*(abs(s2(N,1)*s2(N,1))+abs(s2(1,1)*s2(1,1))));
H.YntFtD' p3=dt*(sum(abs(s3').*abs(s3'))-0.5*(abs(s3(N,1)*s3(N,1))+abs(s3(1,1)*s3(1,1))));
uG/Zpi P1=[P1 p1/p10];
c c@y P2=[P2 p2/p10];
^mH^cP?/ P3=[P3 p3/p10];
G=wJz P=[P p*p];
1v`*%95 end
T8v>J4@t figure(1)
gg<lWeS/3 plot(P,P1, P,P2, P,P3);
>2%!=q3) LnLuWr<;} 转自：
http://blog.163.com/opto_wang/
