i have tried to solve coupled wave equation using split step fourier transform , but getting some issue with the rank
4 次查看(过去 30 天)
显示 更早的评论
given below is the code for pulse compression and decompression for non linear crystal , here i am getting the warning for rank deficiency , and rank =0 , can someone help me out to resolve this issue
clc; clear all; close all; clf; cputime=0; tic; ln=1; i=sqrt(-1); Po=.00064; %input pwr in watts % alpha=0; % Fiber loss value in dB/km % alph=alpha/(4.343); %Ref page#55 eqn 2.5.3 Fiber optic Comm by GP Agrawal % gamma=0.003; %fiber non linearity in /W/m to=125e-12; %initial pulse width in second C=-2; %Input chirp parameter for first calculation b2=-20e-27; %2nd order disp. (s2/m) Ld=(to^2)/(abs(b2)); %dispersion length in meter pi=3.1415926535; Ao=sqrt(Po); %Amplitude %%
% k
n2=5*10^-18;
deltak=(200*pi)/17*10^-3 %delkat k
Lnl=10^-3;
Lgvm=0.6*10^-3;
g=Lnl/Lgvm;
lambda= 800*10^-9 ; % k
n2=6.6508*10^-14;
% deltak=-6;
c=3*1e8;
w1=(2*pi*c)/lambda;
c=3*10^8;
deff=2.0657e-12;
I0=10^9;
%z=1;
tau =- 4096e-12:1e-12: 4095e-12;% dt=t/to
dt=1e-12;
rel_error=1e-5;
h=1000;% step size
for ii=0.1:0.1:1.5 %the various fiber lengths can be varied and this vector can be changed
z=ii*Ld;
u=Ao*exp(-((1+i*(-C))/2)*(tau/to).^2);%page#47 G.P.AGrawal
figure(1)
plot(abs(u),'r');
title('Input Pulse'); xlabel('Time'); ylabel('Amplitude');
grid on;
hold on;
l=max(size(u));
%%%%%%%%%%%%%%%%%%%%%%%
fwhm1=find(abs(u)>abs(max(u)/2));
fwhm1=length(fwhm1);
dw=1/l/dt*2*pi;
w=(-1*l/2:1:l/2-1)*dw;
u=fftshift(u);
w=fftshift(w);
spectrum=fft(fftshift(u)); %Pulse spectrum
for jj=h:h:z
%spectrum=spectrum.*exp(-alph*(h/2)+i*b2/2*w.^2*(h/2)) ;
spectrum=spectrum.*exp(g*w*(h/2)+((n2*deltak)/(w1.*deff))*(h/2)) ;
f=ifft(spectrum);
%f=f.*exp(i*gamma*((abs(f)).^2)*(h));
f=f.*exp((-8*pi*n2^2*I0*Lnl*deltak*exp(1i*deltak*z))/(lambda*w1*deff)*(h/2)+((abs(f)).^2*(h/2)));
spectrum=fft(f);
%spectrum=spectrum.*exp(-alph*(h/2)+i*b2/2*w.^2*(h/2)) ;
spectrum=spectrum.*exp(g*w*(h/2)+((n2*deltak)/(w1.*deff))*(h/2)) ;
% spectrum=spectrum.*exp(-alph*(h/2)+i*b2/2*w.^2*(h/2)) ; % f=ifft(spectrum); % f=f.*exp(i*gamma*((abs(f)).^2)*(h)); % spectrum=fft(f); % spectrum=spectrum.*exp(-alph*(h/2)+i*b2/2*w.^2*(h/2)) ; end f=ifft(spectrum); op_pulse(ln,:)=abs(f);%saving output pulse at all intervals fwhm=find(abs(f)>abs(max(f)/2)); fwhm=length(fwhm); ratio=fwhm/fwhm1; %PBR at every value pbratio(ln)=ratio;%saving PBR at every step size dd=atand((abs(imag(f)))/(abs(real(f)))); phadisp(ln)=dd;%saving pulse phase ln=ln+1; end toc; cputime=toc; figure(2); mesh(op_pulse(1:1:ln-1,:)); title('Pulse Evolution'); xlabel('Time'); ylabel('distance'); zlabel('amplitude'); figure(3) plot(pbratio(1:1:ln-1),'k'); xlabel('Number of steps'); ylabel('Pulse broadening ratio'); grid on; hold on; figure(5) plot(phadisp(1:1:ln-1),'k'); xlabel('distance travelled'); ylabel('phase change'); grid on; hold on; disp('CPU time:'), disp(cputime);
0 个评论
回答(0 个)
另请参阅
Community Treasure Hunt
Find the treasures in MATLAB Central and discover how the community can help you!
Start Hunting!