计算脉冲在非线性耦合器中演化的Matlab 程序 UOY1^wY ?- VN+
d7 % This Matlab script file solves the coupled nonlinear Schrodinger equations of
zEGwQp< % soliton in 2 cores coupler. The output pulse evolution plot is shown in Fig.1 of
07#!b~N % Youfa Wang and Wenfeng Wang, “A simple and effective numerical method for nonlinear
E|TzrH % pulse propagation in N-core optical couplers”, IEEE Photonics Technology lett. Vol.16, No.4, pp1077-1079, 2004
@!S$gTz
y5#_@ %fid=fopen('e21.dat','w');
A g Pg0(G N = 128; % Number of Fourier modes (Time domain sampling points)
c;e ,)$)-| M1 =3000; % Total number of space steps
Z\ Q7#dl J =100; % Steps between output of space
I|M*yObl6 T =10; % length of time windows:T*T0
W) _B(;$] T0=0.1; % input pulse width
"gO5dZ\0 MN1=0; % initial value for the space output location
]+Vcu zq/ dt = T/N; % time step
`7j,njCX. n = [-N/2:1:N/2-1]'; % Index
'71btd1 t = n.*dt;
83 h3C EQ u10=1.*sech(1*t); % input to waveguide1 amplitude: power=u10*u10
$@x kKe" u20=u10.*0.0; % input to waveguide 2
pxF!<nN1, u1=u10; u2=u20;
yx-"YV}5 U1 = u1;
3k/Mig T U2 = u2; % Compute initial condition; save it in U
#7>CLjI ww = 4*n.*n*pi*pi/T/T; % Square of frequency. Note i^2=-1.
!d^` YEfE w=2*pi*n./T;
P TP2QAt g=-i*ww./2; % w=2*pi*f*n./N, f=1/dt=N/T,so w=2*pi*n./T
>"[u.1J_'I L=4; % length of evoluation to compare with S. Trillo's paper
+~@Y#>+./l dz=L/M1; % space step, make sure nonlinear<0.05
7[)(;- for m1 = 1:1:M1 % Start space evolution
9~_6mR< u1 = exp(dz*i*(abs(u1).*abs(u1))).*u1; % 1st sSolve nonlinear part of NLS
W1s|7 u2 = exp(dz*i*(abs(u2).*abs(u2))).*u2;
t0q@]
0B5 ca1 = fftshift(fft(u1)); % Take Fourier transform
RoTT%c P_ ca2 = fftshift(fft(u2));
Px8E~X <@ c2=exp(g.*dz).*(ca2+i*1*ca1.*dz); % approximation
PEWzqZ|!; c1=exp(g.*dz).*(ca1+i*1*ca2.*dz); % frequency domain phase shift
`zL9d lZ u2 = ifft(fftshift(c2)); % Return to physical space
(07d0 <<[ u1 = ifft(fftshift(c1));
*G^]j
)/ if rem(m1,J) == 0 % Save output every J steps.
^#o.WL%4/B U1 = [U1 u1]; % put solutions in U array
O rBFe *2y U2=[U2 u2];
GZ={G2@=I MN1=[MN1 m1];
l0_V-|x z1=dz*MN1'; % output location
j;3o9!.s: end
YV
msWuF end
|2# Ro* hg=abs(U1').*abs(U1'); % for data write to excel
e#('`vGB ha=[z1 hg]; % for data write to excel
v^ tKT & t1=[0 t'];
|`nVr>QF& hh=[t1' ha']; % for data write to excel file
D 4\I;M^ %dlmwrite('aa',hh,'\t'); % save data in the excel format
R<=t{vTJ5 figure(1)
[wR8q,2
waterfall(t',z1',abs(U1').*abs(U1')) % t' is 1xn, z' is 1xm, and U1' is mxn
eBB
D9 SI figure(2)
0TpA3K waterfall(t',z1',abs(U2').*abs(U2')) % t' is 1xn, z' is 1xm, and U1' is mxn
2XtQ"` ) iCS/~[ 非线性超快脉冲耦合的数值方法的Matlab程序 m+g>s&1H
=FnZk J 在研究脉冲在非线性耦合器中的演变时，我们需要求解非线性偏微分方程组。在如下的
论文中，我们提出了一种简洁的数值方法。 这里我们提供给大家用Matlab编写的计算程序。
[x PE?OD 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
f"Iyo:Wt cF2/}m] .tNB07=7 <Va>5R_d< % This Matlab script file solves the nonlinear Schrodinger equations
\K6J{;# L % for 3 cores nonlinear coupler. The output plot is shown in Fig.2 of
S\A[Z&k 0
% Youfa Wang and Wenfeng Wang, “A simple and effective numerical method for nonlinear
wu')Q/v % pulse propagation in N-core optical couplers”, IEEE Photonics Technology lett. Vol.16, No.4, pp1077-1079, 2004
Z ux2VepT s<b7/ ;w' C=1;
wLbngO=VG M1=120, % integer for amplitude
oB9m\o7$ M3=5000; % integer for length of coupler
Q 1Ao65 N = 512; % Number of Fourier modes (Time domain sampling points)
X\%3u PQ dz =3.14159/(sqrt(2.)*C)/M3; % length of coupler is divided into M3 segments, make sure nonlinearity<0.05.
yH^*Fp8V
T =40; % length of time:T*T0.
@Xmk Im dt = T/N; % time step
_H s vF[\[ n = [-N/2:1:N/2-1]'; % Index
bed+Ur& t = n.*dt;
'_)t R;s ww = 4*n.*n*pi*pi/T/T; % Square of frequency. Note i^2=-1.
`vw.~OBl w=2*pi*n./T;
V*}zwm s6 g1=-i*ww./2;
%a `dO EO g2=-i*ww./2; % w=2*pi*f*n./N, f=1/dt=N/T,so w=2*pi*n./TP=0;
w3>|mDA}I g3=-i*ww./2;
AHGcWS\,X P1=0;
iE(grI3 P2=0;
rRYf.~UH@P P3=1;
V{{x~Q9 P=0;
(#]KjpIK
for m1=1:M1
pZxL?N! p=0.032*m1; %input amplitude
$ *A3p s10=p.*sech(p.*t); %input soliton pulse in waveguide 1
d}_c ( s1=s10;
@_3$(*n$~ s20=0.*s10; %input in waveguide 2
lQ"i]};<D s30=0.*s10; %input in waveguide 3
Dl I5} Jh s2=s20;
?W_U{=anl s3=s30;
7g9 ^Jn p10=dt*(sum(abs(s10').*abs(s10'))-0.5*(abs(s10(N,1)*s10(N,1))+abs(s10(1,1)*s10(1,1))));
RZA\-?cO) %energy in waveguide 1
`@7tWX0 p20=dt*(sum(abs(s20').*abs(s20'))-0.5*(abs(s20(N,1)*s20(N,1))+abs(s20(1,1)*s20(1,1))));
'Aj> +H<B %energy in waveguide 2
|T *qAJ8c p30=dt*(sum(abs(s30').*abs(s30'))-0.5*(abs(s30(N,1)*s30(N,1))+abs(s30(1,1)*s30(1,1))));
<J-Z;r(gQN %energy in waveguide 3
ISew]R2 for m3 = 1:1:M3 % Start space evolution
`x)bw s1 = exp(dz*i*(abs(s1).*abs(s1))).*s1; % 1st step, Solve nonlinear part of NLS
HU9y{H s2 = exp(dz*i*(abs(s2).*abs(s2))).*s2;
6l'y s3 = exp(dz*i*(abs(s3).*abs(s3))).*s3;
U I C? S sca1 = fftshift(fft(s1)); % Take Fourier transform
8
-A7 sca2 = fftshift(fft(s2));
$:!T/ *p* sca3 = fftshift(fft(s3));
bl_WN|SQ sc1=exp(g1.*dz).*(sca1+i*C*sca2.*dz); % 2nd step, frequency domain phase shift
zi
.,?Q sc2=exp(g2.*dz).*(sca2+i*C*(sca1+sca3).*dz);
\DK*>
k sc3=exp(g3.*dz).*(sca3+i*C*sca2.*dz);
()?co<@(l s3 = ifft(fftshift(sc3));
Xkom@F~] s2 = ifft(fftshift(sc2)); % Return to physical space
`g N68:B s1 = ifft(fftshift(sc1));
3:lp"C51 end
nD\os[ 3 p1=dt*(sum(abs(s1').*abs(s1'))-0.5*(abs(s1(N,1)*s1(N,1))+abs(s1(1,1)*s1(1,1))));
u^%')Nc p p2=dt*(sum(abs(s2').*abs(s2'))-0.5*(abs(s2(N,1)*s2(N,1))+abs(s2(1,1)*s2(1,1))));
]bb}[#AY p3=dt*(sum(abs(s3').*abs(s3'))-0.5*(abs(s3(N,1)*s3(N,1))+abs(s3(1,1)*s3(1,1))));
3ohcHQ/a P1=[P1 p1/p10];
yuEOQ\!(u P2=[P2 p2/p10];
+Q31K7G r P3=[P3 p3/p10];
J5_Y\@ P=[P p*p];
F
t/
x 5 end
"B3: m-' figure(1)
W y*7j B plot(P,P1, P,P2, P,P3);
:<k|u!b}y % T \N@ 转自：
http://blog.163.com/opto_wang/
