Below are my codes
n = 5*10e12;vf = 1.0*10e6;mo = 9.109*10e-31;m = 0.44*mo;tau = 10e-12;e = 1; %1.6022e-19;
E0 = 10; % Unit in KV/m
s = 3.5*10e3;hBar = 1; %6.52*10e-16;
Ed = 50;Eo = 8.85*10-12;B = [0.1, 0.2, 0.3]; % Unit in Tesla
W = linspace(0,10, 101); % However many you want.
k = W./s;a0 = (2*pi*(hBar^2)*Eo)/m/e^2;s0 = ((e^2)*n*tau)/m;jo = 1./(e*n*vf);[L] = meshgrid(-10000:200:10000); % 1:0.1:3 span for 'q'
legendStrings = cell(length(B), 1);for k1 = 1:length(B) thisN = B(k1); Wc = ((e.*thisN)./m); g = 1 + (1i.*Wc.*tau); gstar = 1 - (1i.*Wc.*tau); BF = repmat(((k.*vf)./Wc),101,1); for i = 1:length(W) R = (L.*(besselj(L,BF).^2))./(1-1i.*(W(i)-(L.*Wc)).*tau); Rx = ((Wc./k).*sum(R(:))); Sx = ((L.*besselj(L,BF)).^2)./(1-1i.*(W(i)-(L.*Wc)).*tau); Sxx = ((2*s0)./(BF.^2)).*sum(Sx(:)); gkw = 1+(1i.*(1./(Eo*(Ed+1))).*(Sxx./(s - Rx))); J = besselj(L,BF)./(1-1i.*(W(i)-(L.*Wc)).*tau).*(L+((k.*a0.*Sxx)./(Wc.*tau)./(Eo.*(s-Rx)))).*((g.*(L+1).*besselj(L+1,BF))+(gstar.*(L-1).*besselj(L-1,BF))); Jx = jo.*((abs((s0*E0)./gkw)).^2).*(1./(BF.^2)./(1+(Wc.*tau).^2)).*sum(J(:)); end legendStrings{k1} = sprintf('B = %.3f', thisN); plot(W, real(Jx), '-', 'LineWidth', 2, 'MarkerSize', 15); hold on; drawnow;endgrid on;fontSize = 15;xlabel('\omega (ns^{-1})', 'FontSize', fontSize)ylabel('\it j_{x} (\mu A/m)', 'FontSize', fontSize)title('\it j_{x} (\mu A/m) vs. \omega (ns^{-1}) ', 'FontSize', fontSize)legend(legendStrings, 'Location', 'Best')
Best Answer