program going beyond matrix dimension
1 Ansicht (letzte 30 Tage)
Ältere Kommentare anzeigen
%first electric field envelop in fourier space
tic;
clc
clear all
close all
Po=9;
Ao=sqrt(Po);
C=0;
%pi=3.1415926535;
x=-60*10^-15:1*10^-15:59*10^-15;
x = 10*x;
i=sqrt(-1);
Wo=120*10^-15;%initial pulse width in second
u=Ao*exp(-(x./Wo).^2)
figure(1)
plot(abs(u),'b');
title('Input Pulse'); xlabel('Time'); ylabel('Amplitude');
hold on
grid on
E1_omega=fft(fftshift(u));
disp(E1_omega)
%second electric field envelop in fourier space
P1=9;
A1=sqrt(P1);
C=0;
%pi=3.1415926535;
% x =-60*10^-15:1*10^-15:59*10^-15
i=sqrt(-1);
Wo=120*10^-15;%initial pulse width in second
u=A1*exp(-(x./Wo).^2)
figure(2)
plot(abs(u),'b');
title('Input Pulse'); xlabel('Time'); ylabel('Amplitude');
hold on
grid on
E2_omega=fft(fftshift(u));
disp(E2_omega)
% patching of the E1 and E2 electric fields with:
E_omega=[E1_omega,E2_omega]
disp(E_omega)
figure(3)
plot(abs(E_omega),'r');
zspan=[0 30*10^-3]
E_omega_0=[0 0.0625]
[z, E_omega] = ode45(@odefun,zspan,E_omega);
E1_omega = E_omega(:, 1:end/2);
E2_omega = E_omega(:, end/2+1:end);
E1_t = ifft(E1_omega, [], 2);
E2_t = ifft(E2_omega, [], 2);
figure(101)
subplot(2, 2, 1)
[X, Y] = meshgrid(1:size(E1_omega, 2), z);
mesh(X, Y, abs(fftshift(E1_omega, 2)).^2);
subplot(2, 2, 2)
[X, Y] = meshgrid(1:size(E1_omega, 2), z);
mesh(X, Y,abs(fftshift(E2_omega, 2)).^2);
subplot(2, 2, 3)
[X, Y] = meshgrid(x, z);
mesh(X, Y,abs(fftshift(E1_t, 2)).^2);
subplot(2, 2, 4)
[X, Y] = meshgrid(x, z);
mesh(X, Y,abs(fftshift(E2_t, 2)).^2);
% g=angle(E_omega)
% K=unwrap(g)
figure(4)
plot(z,E_omega,'linewidth',4)
xlabel('z')
ylabel('phase (radian)')
toc
function dE_omega_dz = odefun(z, E_omega,~,~)
dE_omega_dz=zeros(length(E_omega),1);
Dk=20.94;
l=800; % lambda 800nm
c=3*10^8;
pi=3.1415926535;
omega=2*pi*c/l;
n_2=5*10^-16; %5*10^-16cm^2/W
L_NL=3.7245*10^-3;
LGVM=0.6*10^-3;
I0=50*10^9;
to=125e-12; % initial pulse widthin second
dt=1e-12/to;
dw=1/l/dt*2*pi;
w=(-1*l/2:1:l/2-1)*dw;
% you would have to split the fields back in two:
E1_omega = E_omega(1:end/2);
E2_omega = E_omega(end/2+1:end);
% go back to time space to calculate the nonlinear part:
E1_t = ifft(E1_omega);
E2_t = ifft(E2_omega);
% and calculate the derivatives:
dE_omega_dz(1:length(E_omega)/2) = fft(1i*conj(E1_t).*E2_t.*exp(1i*Dk*z) ...
+ 1i*2*pi*n_2*I0*L_NL/l*(abs(E1_t.^2 + ...
2*abs(E2_t.^2)).*E1_t));
dE_omega_dz(length(E_omega)/2+1:length(E_omega)) = 1i*omega*L_NL/LGVM * E2_t + fft(1i*E1_t.*E1_t.*exp(-1i*Dk*z)....
+ 1i*4*pi*n_2*I0*L_NL/l*(2*abs(E1_t.^2 + ...
abs(E2_t.^2)).*E1_t));
end
% l=max(size(u))
% fwhm1=find(abs(u)>abs(max(u)/2));
% fwhm1=length(fwhm1)
% WWo=(1.699*fwhm1)/2
% max(u)
% 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))
%
% spectrum=spectrum.*exp(-i*w.^2*1.064*10^-6*(h/2*pi))
% f=(ifft((spectrum)));
%
% [t,y]=ode45(@amplitude1,[0 1],[1 1 0]);
% k=1
% y4=y(k)
% k=k+1
% f=f*(-2.02771279053773 + 6.08080900963716i);
% spectrum=fftshift((fft(fftshift((f)))));
% spectrum=spectrum.*exp(-i*w.^2*1.064*10^-6*(h/2*pi));
% f=ifftshift(ifft(ifftshift((((spectrum))))));
% J=max(abs(f))
% j=j+1
% Po(j)=J
%
%
%
%
% plot(abs(f),'r')
% %
%
% fwhm=find(abs(f)>abs(max(f)/2))
% fwhm=length(fwhm)
% s=s+1
% Wo(s)=(1.699*fwhm)/2
% ratio=fwhm/fwhm1
%
% clc
% clear all
% close all
% i=sqrt(-1);
% alph=0
% ln=1
% Po=64;
% 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=20; %initial pulse width in second
% pi=3.1415926535;
% Ao=sqrt(Po); %Amplitude
% %----------------------------------------------------------
% x =- 30:1:29;% dt=t/to
% dt=1/20;
% k=1
% h=0.0077519379844961;% step size
% %the various fiber lengths can be varied and this vector can be changed
%
% u=Ao*exp(-(x/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;
% w1=(-1*l/2:1:l/2-1)*dw
% u=fftshift(u);
% w=fftshift(w1)
% spectrum=fft(fftshift(u)); %Pulse spectrum
% q=h:h:1
% spectrum=spectrum.*exp(-i*1.064*w.^2*(h/4*2*pi));
% f=ifft(spectrum);
%
% [t,y]=ode45(@amplitude1,[0 1],[1 1 0]);
%
% y4=y(k)
% k=k+1
% f=f.*exp(0.335545781012448 + 0.321505706629338i);
% spectrum=fft(f);
% spectrum=spectrum.*exp(-i*1.064*w.^2*(h/4*2*pi)) ;
%
% 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;
% %
%
% plot(abs(f),'b')
3 Kommentare
Adam
am 1 Okt. 2018
Well, unless there is a mistake in your code, that is making your data bigger than it should be, you either need to use less data or get more memory. I'm not familiar with whatever it is your code does so am not going to spot any instructions that accidentally make your arrays bigger than they should be, if there are any such instructions.
Antworten (0)
Siehe auch
Community Treasure Hunt
Find the treasures in MATLAB Central and discover how the community can help you!
Start Hunting!