计算脉冲在非线性耦合器中演化的Matlab 程序 q_!3<.sf D3eK!'qS % This Matlab script file solves the coupled nonlinear Schrodinger equations of
QeGU]WU{ % soliton in 2 cores coupler. The output pulse evolution plot is shown in Fig.1 of
/)~M cP3 % Youfa Wang and Wenfeng Wang, “A simple and effective numerical method for nonlinear
Z EW`?6 % pulse propagation in N-core optical couplers”, IEEE Photonics Technology lett. Vol.16, No.4, pp1077-1079, 2004
V5=Injs * fYwumx`J %fid=fopen('e21.dat','w');
^VA)vLj@ N = 128; % Number of Fourier modes (Time domain sampling points)
3'8~H]<W M1 =3000; % Total number of space steps
il: ""x7^y J =100; % Steps between output of space
4WLB,<b} T =10; % length of time windows:T*T0
=uHTpHR T0=0.1; % input pulse width
h <?Vzl MN1=0; % initial value for the space output location
_b+3;Dy dt = T/N; % time step
sviGS&J9h n = [-N/2:1:N/2-1]'; % Index
_$r+*nGDz t = n.*dt;
W*P/~U= u10=1.*sech(1*t); % input to waveguide1 amplitude: power=u10*u10
{|qz> u20=u10.*0.0; % input to waveguide 2
k7j;'6 u1=u10; u2=u20;
<3 i!{"} U1 = u1;
)pg?Z M9 U2 = u2; % Compute initial condition; save it in U
nF=h|rN ww = 4*n.*n*pi*pi/T/T; % Square of frequency. Note i^2=-1.
# 6JG#!W w=2*pi*n./T;
zDX-}t_'q g=-i*ww./2; % w=2*pi*f*n./N, f=1/dt=N/T,so w=2*pi*n./T
[xHK^JP 8F L=4; % length of evoluation to compare with S. Trillo's paper
tYnNO K*| dz=L/M1; % space step, make sure nonlinear<0.05
<|v]9`' for m1 = 1:1:M1 % Start space evolution
DwoO([&I u1 = exp(dz*i*(abs(u1).*abs(u1))).*u1; % 1st sSolve nonlinear part of NLS
'C(YUlT2?P u2 = exp(dz*i*(abs(u2).*abs(u2))).*u2;
.2`S07Z ca1 = fftshift(fft(u1)); % Take Fourier transform
Jg@PhN<9 ca2 = fftshift(fft(u2));
<=WQs2 c2=exp(g.*dz).*(ca2+i*1*ca1.*dz); % approximation
7uYJ _R c1=exp(g.*dz).*(ca1+i*1*ca2.*dz); % frequency domain phase shift
Hg<]5 u2 = ifft(fftshift(c2)); % Return to physical space
9T)-|fja_ u1 = ifft(fftshift(c1));
zT.qNtU% if rem(m1,J) == 0 % Save output every J steps.
nP] ~8ViS U1 = [U1 u1]; % put solutions in U array
gVO[R6C 5C U2=[U2 u2];
&P rx=L` MN1=[MN1 m1];
Z O&5C6qa z1=dz*MN1'; % output location
8xLvpgcZ end
r.[9/'> end
XJ.vj+XXb
hg=abs(U1').*abs(U1'); % for data write to excel
G" wy? ha=[z1 hg]; % for data write to excel
%{ax oGd t1=[0 t'];
iJU]|t hh=[t1' ha']; % for data write to excel file
$cnIsyKWY %dlmwrite('aa',hh,'\t'); % save data in the excel format
E NygD figure(1)
m+zzhv1 waterfall(t',z1',abs(U1').*abs(U1')) % t' is 1xn, z' is 1xm, and U1' is mxn
~i(X{ ^,3 figure(2)
5MT$n4zKu waterfall(t',z1',abs(U2').*abs(U2')) % t' is 1xn, z' is 1xm, and U1' is mxn
k\A8Z[ _L9`bzZj
非线性超快脉冲耦合的数值方法的Matlab程序 b3W@{je c{z QX0 在研究脉冲在非线性耦合器中的演变时，我们需要求解非线性偏微分方程组。在如下的
论文中，我们提出了一种简洁的数值方法。 这里我们提供给大家用Matlab编写的计算程序。
.^ soX} 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
NeQ/#[~g G;MmD?VJ g =j6f/8 c5pF?k FaD % This Matlab script file solves the nonlinear Schrodinger equations
6Z0@4_Y@B6 % for 3 cores nonlinear coupler. The output plot is shown in Fig.2 of
Jc/*w % Youfa Wang and Wenfeng Wang, “A simple and effective numerical method for nonlinear
LNtBYdB`pK % pulse propagation in N-core optical couplers”, IEEE Photonics Technology lett. Vol.16, No.4, pp1077-1079, 2004
(]1n ! 4Z,MqG> C=1;
.hXxh)F M1=120, % integer for amplitude
'`I&g8I\ M3=5000; % integer for length of coupler
J;HkR9<C N = 512; % Number of Fourier modes (Time domain sampling points)
6Gwk*%sb dz =3.14159/(sqrt(2.)*C)/M3; % length of coupler is divided into M3 segments, make sure nonlinearity<0.05.
DR;rK[f T =40; % length of time:T*T0.
uL`;KD dt = T/N; % time step
pri=;I(2A n = [-N/2:1:N/2-1]'; % Index
eNR>W>;' t = n.*dt;
*M gl X< ww = 4*n.*n*pi*pi/T/T; % Square of frequency. Note i^2=-1.
7' ]n_-fu w=2*pi*n./T;
/0IvvD!7N g1=-i*ww./2;
z1K@AaRx g2=-i*ww./2; % w=2*pi*f*n./N, f=1/dt=N/T,so w=2*pi*n./TP=0;
2Gd.B/L6 g3=-i*ww./2;
)[i0 ~o[ P1=0;
nDhr;/"i P2=0;
. _Bejh P3=1;
r3*0`Rup P=0;
|wZcVc t~ for m1=1:M1
QX-%<@ p=0.032*m1; %input amplitude
B agO0# s10=p.*sech(p.*t); %input soliton pulse in waveguide 1
8>%:MS" s1=s10;
9Ra*bP ]1 s20=0.*s10; %input in waveguide 2
/_rEI,[k s30=0.*s10; %input in waveguide 3
JH C 6l s2=s20;
g1UP/hNJ\8 s3=s30;
B&3oo p10=dt*(sum(abs(s10').*abs(s10'))-0.5*(abs(s10(N,1)*s10(N,1))+abs(s10(1,1)*s10(1,1))));
`7 [z% cuK %energy in waveguide 1
`fYICp p20=dt*(sum(abs(s20').*abs(s20'))-0.5*(abs(s20(N,1)*s20(N,1))+abs(s20(1,1)*s20(1,1))));
.SzP ig %energy in waveguide 2
pUi|&F K"> p30=dt*(sum(abs(s30').*abs(s30'))-0.5*(abs(s30(N,1)*s30(N,1))+abs(s30(1,1)*s30(1,1))));
ya .!zGH %energy in waveguide 3
M{w[ hV for m3 = 1:1:M3 % Start space evolution
0]p!
Bscaf s1 = exp(dz*i*(abs(s1).*abs(s1))).*s1; % 1st step, Solve nonlinear part of NLS
LQ(z~M0B s2 = exp(dz*i*(abs(s2).*abs(s2))).*s2;
Q8OA{EUtq s3 = exp(dz*i*(abs(s3).*abs(s3))).*s3;
e=e^;K4 sca1 = fftshift(fft(s1)); % Take Fourier transform
/%fBkA#n sca2 = fftshift(fft(s2));
Jr+~' sca3 = fftshift(fft(s3));
1j"_@?H[ sc1=exp(g1.*dz).*(sca1+i*C*sca2.*dz); % 2nd step, frequency domain phase shift
7L)edR [ sc2=exp(g2.*dz).*(sca2+i*C*(sca1+sca3).*dz);
B WRAz*V sc3=exp(g3.*dz).*(sca3+i*C*sca2.*dz);
Q$u&/g3NvL s3 = ifft(fftshift(sc3));
?tx%K U\3 s2 = ifft(fftshift(sc2)); % Return to physical space
_kGJqyYV s1 = ifft(fftshift(sc1));
q% *-4GP end
/Y| y0iK p1=dt*(sum(abs(s1').*abs(s1'))-0.5*(abs(s1(N,1)*s1(N,1))+abs(s1(1,1)*s1(1,1))));
6:_@ ;/03% p2=dt*(sum(abs(s2').*abs(s2'))-0.5*(abs(s2(N,1)*s2(N,1))+abs(s2(1,1)*s2(1,1))));
J8IdQ:4^l p3=dt*(sum(abs(s3').*abs(s3'))-0.5*(abs(s3(N,1)*s3(N,1))+abs(s3(1,1)*s3(1,1))));
>v- -R8I * P1=[P1 p1/p10];
-hL 0}Wy$N P2=[P2 p2/p10];
`Tw DR6& P3=[P3 p3/p10];
?'SHt9b3| P=[P p*p];
RI.6.f1dy end
xgeDfp F' figure(1)
Lxz!>JO> plot(P,P1, P,P2, P,P3);
vz$-KT4e^ TNun)0p 转自：
http://blog.163.com/opto_wang/
