Unable to plot array inside a running while loop

9 ビュー (過去 30 日間)
Agnivo Gosai
Agnivo Gosai 2020 年 4 月 21 日
編集済み: Agnivo Gosai 2020 年 4 月 21 日
This is the code
% snyder paper 1981
% parallel waveguides Morishita K.,Optics Letters,Vol 13. No.2
% 15 th April 2020
clearvars
nc1 = 1.449734; % core
nclad1 = 1.44479; % cladding - surrounding
delta1 = 1 - (nclad1/nc1)^2;
a1 = 4.75;% core radius
nc2 = 1.449686; % core
nclad2 = 1.4428; % cladding - surrounding
delta2 = 1 - (nclad2/nc2)^2;
a2 = 4.75; %core radius
d =11.5;% core to core center-distance
lambda = 1.25; % start wavelength in microns
inc = 0.005; % increment for wavelength
%% code to store data
kk = (1.605-lambda)/inc + 1;
kk = round(kk,0);
pp1 = zeros(1,kk);
pp2 = zeros(1,kk);
wave = zeros(1,kk);
coup = zeros(1,kk);
length = zeros(1,kk);
betas = zeros(1,kk);
i=0;
for z = 500:250:1250
z % check value of z
while lambda<=1.605
k0 = 2*3.14/lambda;
V1 = sqrt((a1^2)*(k0^2)*(nc1^2-nclad1^2));
U1 = -0.0024*(V1^4)+0.0493*(V1^3)-0.3757*(V1^2)+1.3444*V1-0.0095;
beta1 = sqrt((1/(a1^2))*((V1^2)/delta1)*(1-delta1*((U1/V1)^2)));
W1 = sqrt((a1^2)*(beta1^2-(k0^2)*(nclad1^2)));
V2 = sqrt((a2^2)*(k0^2)*(nc2^2-nclad2^2));
U2 = -0.0024*(V2^4)+0.0493*(V2^3)-0.3757*(V2^2)+1.3444*V2-0.0095;
beta2 = sqrt((1/(a2^2))*((V2^2)/delta2)*(1-delta2*((U2/V2)^2)));
W2 = sqrt((a2^2)*(beta2^2-(k0^2)*(nclad2^2)));
deltabeta = beta1-beta2;
A1 = ((delta1*delta2)^(1/4))/(sqrt(a1*a2));
A2 = (U1*U2)/((V1*V2)^(3/2));
A3 = W1*d/a1;
A4 = besselk(0,A3,1);
A5 = (besselk(1,W1,1))*(besselk(1,W2,1));
C = A1*A2*A4/A5;
L = 3.14/(2*C);
betad = deltabeta/2;
betab = sqrt(betad^2+C^2);
F = 1/(1+(betad/C)^2);
P1= 1-F*((sin(betab*z))^2);
P2 = F*((sin(betab*z))^2);
i = i +1;
pp1(:,i) = P1*100;
pp2(:,i) = P2*100;
wave(:,i) = lambda*1000;
coup(:,i) = C;
length(:,i) = L;
betas(:,i) = deltabeta*10^3;
figure(1)
plot(wave(:,i),pp1(:,i),'r--') % plot the current value of pp1 against wave
hold on
lambda = lambda + inc;
end
i = 0 % reset i for next iteration
lambda = 1.25 % reset lambda for next iteration
end
hold off

採用された回答

Sriram Tadavarty
Sriram Tadavarty 2020 年 4 月 21 日
編集済み: Sriram Tadavarty 2020 年 4 月 21 日
Hi Agnivo,
A minor update to your code in placing the plot after the while loops solves this. I commented out the lines that are not required.
clearvars
nc1 = 1.449734; % core
nclad1 = 1.44479; % cladding - surrounding
delta1 = 1 - (nclad1/nc1)^2;
a1 = 4.75;% core radius
nc2 = 1.449686; % core
nclad2 = 1.4428; % cladding - surrounding
delta2 = 1 - (nclad2/nc2)^2;
a2 = 4.75; %core radius
d =11.5;% core to core center-distance
lambda = 1.25; % start wavelength in microns
inc = 0.005; % increment for wavelength
%% code to store data
kk = (1.605-lambda)/inc + 1;
kk = round(kk,0);
pp1 = zeros(1,kk);
pp2 = zeros(1,kk);
wave = zeros(1,kk);
coup = zeros(1,kk);
length = zeros(1,kk);
betas = zeros(1,kk);
i=0;
for z = 500:250:1250
z % check value of z
while lambda<=1.605
k0 = 2*3.14/lambda;
V1 = sqrt((a1^2)*(k0^2)*(nc1^2-nclad1^2));
U1 = -0.0024*(V1^4)+0.0493*(V1^3)-0.3757*(V1^2)+1.3444*V1-0.0095;
beta1 = sqrt((1/(a1^2))*((V1^2)/delta1)*(1-delta1*((U1/V1)^2)));
W1 = sqrt((a1^2)*(beta1^2-(k0^2)*(nclad1^2)));
V2 = sqrt((a2^2)*(k0^2)*(nc2^2-nclad2^2));
U2 = -0.0024*(V2^4)+0.0493*(V2^3)-0.3757*(V2^2)+1.3444*V2-0.0095;
beta2 = sqrt((1/(a2^2))*((V2^2)/delta2)*(1-delta2*((U2/V2)^2)));
W2 = sqrt((a2^2)*(beta2^2-(k0^2)*(nclad2^2)));
deltabeta = beta1-beta2;
A1 = ((delta1*delta2)^(1/4))/(sqrt(a1*a2));
A2 = (U1*U2)/((V1*V2)^(3/2));
A3 = W1*d/a1;
A4 = besselk(0,A3,1);
A5 = (besselk(1,W1,1))*(besselk(1,W2,1));
C = A1*A2*A4/A5;
L = 3.14/(2*C);
betad = deltabeta/2;
betab = sqrt(betad^2+C^2);
F = 1/(1+(betad/C)^2);
P1= 1-F*((sin(betab*z))^2);
P2 = F*((sin(betab*z))^2);
i = i +1;
pp1(:,i) = P1*100;
pp2(:,i) = P2*100;
wave(:,i) = lambda*1000;
coup(:,i) = C;
length(:,i) = L;
betas(:,i) = deltabeta*10^3;
% figure(1)
% plot(wave(:,i),pp1(:,i),'r--') % plot the current value of pp1 against wave
% hold on
lambda = lambda + inc;
end
figure
plot(wave,pp1,'r--') % plot the complete wave for the iteration
i = 0 % reset i for next iteration
lambda = 1.25 % reset lambda for next iteration
end
Hope this helps.
Regards,
Sriram

その他の回答 (1 件)

David Hill
David Hill 2020 年 4 月 21 日
If you want everything stored:
nc1 = 1.449734; % core
nclad1 = 1.44479; % cladding - surrounding
delta1 = 1 - (nclad1/nc1)^2;
a1 = 4.75;% core radius
nc2 = 1.449686; % core
nclad2 = 1.4428; % cladding - surrounding
delta2 = 1 - (nclad2/nc2)^2;
a2 = 4.75; %core radius
d =11.5;% core to core center-distance
lambda = 1.25; % start wavelength in microns
inc = 0.005; % increment for wavelength
%% code to store data
kk = (1.605-lambda)/inc + 1;
kk = round(kk,0);
pp1 = zeros(1,kk);
pp2 = zeros(1,kk);
wave = zeros(1,kk);
coup = zeros(1,kk);
length = zeros(1,kk);
betas = zeros(1,kk);
i=0;
count=0;
for z = 500:250:1250
%z % check value of z
count=count+1;
while lambda<=1.605
k0 = 2*3.14/lambda;
V1 = sqrt((a1^2)*(k0^2)*(nc1^2-nclad1^2));
U1 = -0.0024*(V1^4)+0.0493*(V1^3)-0.3757*(V1^2)+1.3444*V1-0.0095;
beta1 = sqrt((1/(a1^2))*((V1^2)/delta1)*(1-delta1*((U1/V1)^2)));
W1 = sqrt((a1^2)*(beta1^2-(k0^2)*(nclad1^2)));
V2 = sqrt((a2^2)*(k0^2)*(nc2^2-nclad2^2));
U2 = -0.0024*(V2^4)+0.0493*(V2^3)-0.3757*(V2^2)+1.3444*V2-0.0095;
beta2 = sqrt((1/(a2^2))*((V2^2)/delta2)*(1-delta2*((U2/V2)^2)));
W2 = sqrt((a2^2)*(beta2^2-(k0^2)*(nclad2^2)));
deltabeta = beta1-beta2;
A1 = ((delta1*delta2)^(1/4))/(sqrt(a1*a2));
A2 = (U1*U2)/((V1*V2)^(3/2));
A3 = W1*d/a1;
A4 = besselk(0,A3,1);
A5 = (besselk(1,W1,1))*(besselk(1,W2,1));
C = A1*A2*A4/A5;
L = 3.14/(2*C);
betad = deltabeta/2;
betab = sqrt(betad^2+C^2);
F = 1/(1+(betad/C)^2);
P1= 1-F*((sin(betab*z))^2);
P2 = F*((sin(betab*z))^2);
i = i +1;
pp1(count,i) = P1*100;
pp2(count,i) = P2*100;
wave(count,i) = lambda*1000;
coup(count,i) = C;
length(count,i) = L;
betas(count,i) = deltabeta*10^3;
lambda = lambda + inc;
end
hold on;
plot(wave(count,:),pp1(count,:),'r--')
i = 0; % reset i for next iteration
lambda = 1.25; % reset lambda for next iteration
end
  1 件のコメント
Agnivo Gosai
Agnivo Gosai 2020 年 4 月 21 日
編集済み: Agnivo Gosai 2020 年 4 月 21 日
Thanks, this is very helpful

サインインしてコメントする。

カテゴリ

Help Center および File ExchangeDeep Learning Toolbox についてさらに検索

タグ

Community Treasure Hunt

Find the treasures in MATLAB Central and discover how the community can help you!

Start Hunting!

Translated by