| 查看: 1197 | 回复: 0 | |||
[求助]
用分步傅里叶法求解光纤中非线性薛定谔方程
|
|
不知这里哪位虫子有用分步傅里叶法求解光纤中非线性薛定谔方程的MATLAB程序啊？我用Agrawal那本非线性光纤光学后面的程序加以改动，但是得出的图效果不理想，有没有哪位同学有程序啊？真心求交流。。。。。 附程序： %传输距离为1到4，其他参数不变，最后得出不同传输距离下频域光脉冲形状的比较 clear all; for m=1:4 %---Input Field profile distance = 10*m; %传输距离 beta2 =-1;%二阶色散参数 chirp = 0; % 初始啁啾参数 peakpower=1;%峰值功率 T0=1;%初始脉宽 r=1 ; %非线性系数 N=sqrt(r*peakpower*T0^2) ; alpha=1;%光纤损耗 %---set simulation parameters nt =1024; Tmax = 10; % FFT points and window size,设置采样点和时间窗口 step_num = round(20*distance*N^2) % No. of z steps to deltaz = distance/step_num; % step size in z dtau = (2*Tmax)/nt; % step size in tau %---tau and omega arrays tau = (-nt/2:nt/2-1)*dtau; % temporal grid omega = (pi/Tmax)*[(-nt/2:nt/2-1)]; % frequency grid %---store dispersive phase shifts to speedup code dispersion = exp((1i*0.5*beta2*omega.^2-alpha/2)*deltaz);% phase factor hhz = 1i*N^2*deltaz; % nonlinear phase factor uu =sqrt(peakpower)*exp(-0.5*tau.^2.*(1+1i*chirp)/T0^2); % gaussian temp = uu.*exp((abs(uu).^2).*hhz/2); % note hhz/2 subplot(2,1,1) if m==1 plot(fftshift(omega)/(2*pi),-log10(abs(temp).^2),'r') hold on else if m==2 plot(fftshift(omega)/(2*pi),-log10(abs(temp).^2),'y') hold on else if m==3 plot(fftshift(omega)/(2*pi),-log10(abs(temp).^2),'b') hold on else plot(fftshift(omega)/(2*pi),-log10(abs(temp).^2),'k') end end end %*********[ Beginning of MAIN Loop]*********** % scheme: 1/2N -> D -> 1/2N; first half step nonlinear for n=1:step_num f_temp = ifft(temp).*dispersion; uu = fft(f_temp); temp = uu.*exp((abs(uu).^2).*hhz/2); end uu = temp.*exp(-(abs(uu).^2).*hhz/2);% Final field temp = fftshift(ifft(uu)).*(nt*dtau)/sqrt(2*pi);%Final spectrum subplot(2,1,2) if m==1 plot(fftshift(omega)/(2*pi),log10(abs(temp)).^2,'r') hold on else if m==2 plot(fftshift(omega)/(2*pi),log10(abs(temp)).^2,'y') hold on else if m==3 plot(fftshift(omega)/(2*pi),log10(abs(temp)).^2,'b') hold on else plot(fftshift(omega)/(2*pi),log10(abs(temp)).^2,'k') end end end %***************[ End of MAIN Loop ]************** %----Plot output pulse shape and spectrum %TEMP(m,  =temp;end %for m=1:4 plot(fftshift(omega)/(2*pi),-20*(log10(abs(TEMP(m, ).^2)));% plot(fftshift(omega)/(2*pi),-20*(log10(abs(TEMP(m, ).^2)));%end |
回复此楼» 猜你喜欢
284求调剂
已经有10人回复
-
一志愿山东大学药学学硕求调剂
已经有4人回复
-
07化学280分求调剂
已经有4人回复
-
298-一志愿中国农业大学-求调剂
已经有12人回复
-
求材料，环境专业调剂
已经有3人回复
-
335求调剂
已经有5人回复
-
求调剂
已经有7人回复
-
一志愿吉大化学322求调剂
已经有4人回复
-
环境学硕288求调剂
已经有8人回复
-
341求调剂(一志愿湖南大学070300)
已经有6人回复
点击这里搜索更多相关资源
