How to plot the graph for the given code ie Vth vs ls??

if true
 if true
clear all;
close all;
clc;
K=8.61*10^-5;
T=300;
q=1.6*10^-19;
Vt=K*T/q;
ni=4.15;
Na=2000;
ls=3;
Nsd_p=1*10^20;
for x=1:70
Nsd(x)=Nsd_p*exp(-x^2/(2*(ls)^2));
end
figure,plot(Nsd);
for x=1:70
ei(x)=0.054*(1-((Nsd(x)/(2.7*10^19))^(1/3)));
end
figure,plot(ei);
del_eg=1.11;
ni_eff=sqrt(ni*ni*exp(del_eg/(K*T)));
%ni_eff=6.05*10^-20;
Eg=1.11;
Eg_eff=1.17-1.11;
for x=1:70
Ed(x)=Eg_eff-ei(x);
end
figure,plot(Ed);
for x=1:70
Ef(x)=(Eg/2)+K*T*log(Nsd(x)/(ni_eff));
end
figure,plot(Ef);
SD=4;
for x=1:70
%Nsd__(x)=4.91*10^23;
Nssd(x)=(Nsd_p*exp((-x^2)/(2*ls^2))+(Nsd_p*exp(-(ls-x)^2)/(2*ls^2))/(1+SD*exp(Ef(x)-Ed(x)/K*T)));
end
figure,plot(Nssd)
tsi=20;
tox=2;
epsilonsi=2;
epsilonox=4;
lambda=((epsilonsi*tsi*tox)/(2*epsilonox))
Vgs=0.1;
Vfb=1;
PI=(Vgs-Vfb)-(lambda^2)*q*(Na-Nssd(x))/(epsilonsi)
Nde=2.7*10^19;
Seff=log(Nde/Nsd_p)*(-2*ls^2)
Lg=70;
Leff=Lg-2*Seff
Deff=Seff+Leff
Vbi=Vt*log(Nde*Na/(ni_eff)^2)
Vds=1;
c2=((Vbi-PI)*(1-(exp(-Deff/lambda))/(exp(-Seff/lambda))+Vds))/(exp(Deff/lambda))-(exp(Seff/lambda))/(exp(-Seff/lambda))*(exp(-Deff/lambda));
c1=(Vbi-c2*exp(Seff/lambda)/exp(-Seff/lambda))*exp(-Deff/lambda);
fi=Vt*log(Na/ni);
epsilon=2;
xmin=((lambda/2)*log(c1/c2));
Vth=Vfb+(2*fi)-(c1*exp(-xmin/lambda))-(c2*exp(xmin/lambda))+((lambda^2)*q*(Na-Nssd(1))/epsilonsi)
end
end

Answer 1

[Walter Roberson]

lsvals = [3 3.2 4.1 4.18];
for J = 1 : length(lsvals)
  ls = lsvals(J);
  ...
  ...
  Vth(J) = ....
end

plot(lsvals, Vth)

[ES]

What is the range of Is? I think lsvals = [3 3.2 4.1 4.18]; was just an example in the answer.

What is the range you want to plot for?

[ARTI]

ls varies from 2nm to 5nm and still after replacing the values the graph is a straight line. please help