How to skip NaN arrays outputs of a for loop?
1 visualización (últimos 30 días)
Mostrar comentarios más antiguos
Athira T Das
el 18 de Jul. de 2022
Respondida: Jan
el 18 de Jul. de 2022
clc; clear all; close all;
syms m1 m2 s1 s2 j1 j2
lambda = 1060*10^-9;
M=1
z=linspace(0.00001,8000);
wo = 0.02;
C = 10^(-7);
k=2*pi/lambda;
po=(0.545*C^2*k^2*z).^(-(3/5));
b=0.1;
T1=0;
T2=0;
T3=0;
T4=0;
T5=0;
T6=0;
delta= ((1i*k)./(2*z))+(1./(wo.^2))+(1./(po.^2));
delta_star= subs(delta, 1i, -1i);
etta = delta_star - (1./(delta.*po.^4));
X=((k./(2.*z)).^2).*((1./(2.*1i.*sqrt(delta))).^M).*(1./(16.*delta.*etta));
alpha = (1i.*k./(2.*z)).*(1./(delta.*po.^2)-1);
beta_p = (b./(2.*wo)).*(1./(delta.*po.^2)+1);
beta_n = (b./(2.*wo)).*(1./(delta.*po.^2)-1);
A = (alpha.^2./etta)-((k.^2)-(4.*z.^2.*delta));
B1_p = ((1i.*k.*b)./(2.*z.*wo.*delta))+((2.*alpha.*beta_p)./etta);
B1_n = ((1i.*k.*b)./(2.*z.*wo.*delta))+((2.*alpha.*beta_n)./etta);
B2_p = -(((1i.*k.*b)./(2.*z.*wo.*delta))+((2.*alpha.*beta_p)./etta));
B2_n = -(((1i.*k.*b)./(2.*z.*wo.*delta))+((2.*alpha.*beta_n)./etta));
C1_p = ((b.^2)./(4.*wo.^2.*delta))+(((beta_p).^2)./etta);
C1_n = ((b.^2)./(4.*wo.^2.*delta))+(((beta_n).^2)./etta);
C2_p = ((b.^2)./(4.*wo.^2.*delta))+(((beta_n).^2)./etta);
C2_n = ((b.^2)./(4.*wo.^2.*delta))+(((beta_p).^2)./etta);
m1 = 0:M;
m2=0:M;
s1 = m1./2;
s2=(M-m1)./2;
%%%% j1=m1-2.*s1 %%%%
allVecs_1 = {-2*s1,m1};
sub = cell(1,numel(allVecs_1));
[sub{:}] = ndgrid(allVecs_1{:});
sub = cellfun(@(x)x(:),sub,'UniformOutput', false);
j= cell2mat(sub);
j11 = [ sum(j(:,[1 2]),2) ];
j1=j11(j11>=0);
%%%% j2=(M-m1-(2.*s2)) %%%%
allVecs_1 = {-1*m1,-2*s2};
sub = cell(1,numel(allVecs_1));
[sub{:}] = ndgrid(allVecs_1{:});
sub = cellfun(@(x)x(:),sub,'UniformOutput', false);
j21= cell2mat(sub);
j2212 =M + [ sum(j21(:,[1 2]),2) ];
j2=j2212(j2212>=0);
%%%%%%%%%%%%%%%%%%%%%%%%%%%
[rm1,cm1] = size(m1);
[rm2,cm2] = size(m2);
[rs1,cs1] = size(s1);
[rs2,cs2] = size(s2);
[rj1,cj1] = size(j1);
[rj2,cj2] = size(j2);
for Nm1 = 1:cm1
for Mm1 = 1:rm1
M1= m1(Mm1,Nm1);
T2=0;
for Nm2 = 1:cm2
for Mm2 = 1:rm2
M2= m2(Mm2,Nm2);
T3=0;
for Ns1 = 1:cs1
for Ms1 = 1:rs1
S1= s1(Ms1,Ns1);
T4=0;
for Nj1 = 1:cj1
for Mj1 = 1:rj1
J1= j1(Mj1,Nj1);
T5=0;
for Ns2 = 1:cs1
for Ms2 = 1:rs2
S2= s2(Ms2,Ns2);
T6=0;
for Nj2 = 1:cj2
for Mj2 = 1:rj2
J2= j2(Mj2,Nj2);
e3=(exp(C1_p).*hermiteH(M-M2+J2,((1i.*beta_p)./sqrt(etta))))+(exp(C2_n).*hermiteH(M-M2+J2,((-1i.*beta_n)./sqrt(etta))));
e4=(exp(C1_n).*hermiteH(M-M2+J2,((1i.*beta_n)./sqrt(etta))))+(exp(C2_p).*hermiteH(M-M2+J2,((-1i.*beta_p)./sqrt(etta))));
T6=T6+(gamma(1+M-M1-(2.*S2))/(gamma(J2+1)*gamma(M-M1-(2.*S2)-J2+1))).*((1./(2.*1i.*sqrt(etta))).^(M-M2+J2)).*((1./(po.^2)).^J2).*((b./2.*wo).^(M-M1-2*S2-J2)).*(e4+e3);
end
end
e1=(exp(C1_n).*hermiteH(M2+J1,((1i.*beta_n)./sqrt(etta))))+(exp(C2_p).*hermiteH(M2+J1,((-1i.*beta_p)./sqrt(etta))));
e2=(exp(C1_p).*hermiteH(M2+J1,((1i.*beta_p)./sqrt(etta))))+(exp(C2_n).*hermiteH(M2+J1,((-1i.*beta_n)./sqrt(etta))));
T5=T5+T6.*(((gamma(M-M1+1).*(-1).^S2)./(gamma(S2+1).*gamma(M-M1-2.*S2+1)))).*(((2.*1i)./(delta.^(0.5))).^(M-M1-2.*S2)).*(e1+e2);
T5 = rmmissing(vpa(T5))
end
end
T4=T4+T5.*((gamma(M1-(2.*S2)+1))/(gamma(J1+1)*gamma(1+M1-(2.*S2)-J1))).*((1./(2.*1i.*sqrt(etta))).^(M2+J2)).*((1./(po.^2)).^J1).*((b./2.*wo).^(M1-2.*S1-J1));
T4 = rmmissing(vpa(T4))
end
end
end
T3=T3+T4.*(((gamma(M1+1).*(-1).^S1)./(gamma(S1+1).*gamma(1+M1-2.*S1))).*(((2.*1i)./(delta.^(0.5))).^(M1-2.*S1)));
T3 = rmmissing(vpa(T3))
end
T2=T2+T3*nchoosek(M,M2)*((-1i)^(M-M2));
T2 = rmmissing(vpa(T2))
end
T1=T1+T2*nchoosek(M,M1)*((1i)^(M-M1));
T1 = rmmissing(vpa(T1))
end
end
end
I0 = abs(real(X.*T1))
2 comentarios
Matt J
el 18 de Jul. de 2022
What is your question? The title of your post doesn't convey what problems you're having with the posted code.
Respuesta aceptada
Jan
el 18 de Jul. de 2022
doBreak = false;
for Nj2 = 1:cj2
for Mj2 = 1:rj2
J2= j2(Mj2,Nj2);
e3=(exp(C1_p).*hermiteH(M-M2+J2,((1i.*beta_p)./sqrt(etta))))+(exp(C2_n).*hermiteH(M-M2+J2,((-1i.*beta_n)./sqrt(etta))));
e4=(exp(C1_n).*hermiteH(M-M2+J2,((1i.*beta_n)./sqrt(etta))))+(exp(C2_p).*hermiteH(M-M2+J2,((-1i.*beta_p)./sqrt(etta))));
tmp = (gamma(1+M-M1-(2.*S2))/(gamma(J2+1)*gamma(M-M1-(2.*S2)-J2+1))).*((1./(2.*1i.*sqrt(etta))).^(M-M2+J2)).*((1./(po.^2)).^J2).*((b./2.*wo).^(M-M1-2*S2-J2)).*(e4+e3);
if isnan(tmp)
doBreak = true;
break;
end
T6 = T6 + tmp;
end
if doBreak % The loop over Nj2 also?
break;
end
end
0 comentarios
Más respuestas (0)
Ver también
Categorías
Más información sobre Loops and Conditional Statements en Help Center y File Exchange.
Community Treasure Hunt
Find the treasures in MATLAB Central and discover how the community can help you!
Start Hunting!