%Fig6_5a.m %Essential Electron Transport for Device Physics %Inelastic LO-phonon scattering rate at finite temperature in GaAs clear all; clf; FS = 12; %label fontsize 18 FSN = 12; %number fontsize 16 LW = 2; %linewidth % Change default axes fonts. set(0,'DefaultAxesFontName', 'Times'); set(0,'DefaultAxesFontSize', FSN); % Change default text fonts. set(0,'DefaultTextFontname', 'Times'); set(0,'DefaultTextFontSize', FSN); colormap(jet); hbar = 1.05457159e-34; % Planck's constant (Js) m2 = 0.07; % Effective electron mass coefficient m0 = 9.109382e-31; % Bare electron mass (kg) wLO=36.3; % longitudinal optic phonon energy (meV) wTO=33.3; % transverse optic phonon energy (meV) hbar1=['\fontname{MT Extra} h \fontname{Times}']; a0=0.529177e-10; % Bohr radius (m) einf =11.1; % high frequency dielectric constant meff = m2*m0; % Effective electron mass (kg) echarge = 1.6021764e-19; % Electron charge (C) kscale=3.0937e+08; %wave vector scale (m^-1) Escale = (hbar^2)*(kscale^2)/(2*meff);% energy scale (J) Escale1 = (Escale*1e3)/echarge; % energy scale in meV wLO2=(wLO/Escale1)^2; wTO2=(wTO/Escale1)^2; gamma=0.005; %Energy broadening (GAMMA / Escale1) 0.01 0.005 eye=complex(0,1); %square root of minus one T = 300; % Absolute temperature 300, 10, 4.2 (K) kB = 1.3806505e-23/echarge; % Boltzmann constant (eV/K) beta = 1/(kB*T); % Inverse thermal energy (eV^-1) Energymax=100; % Electron energy maximum (meV) 100 300 Energystep=Energymax/100; % Electron energy step (meV) /100 ymax=Energymax/Escale1; % max value of energy loss, omega/Ef 18 %max value of scattered wave vector q/kscale set at ymax value xmax=(sqrt(ymax))*(1+sqrt(2));%approximately 15, 6 for n=1e17, 1e18 x1=[gamma/4:gamma/4:xmax];%[0.01:.005:xmax] scattered wave vector step at gamma/4 y=[-ymax:gamma/4:ymax]+eye*gamma;% energy loss step is gamma/4 f1 = zeros(length(y),length(x1)); %************************************************************************* % Calculate the response with phonons %************************************************************************* for ii=1:length(x1) x = x1(ii); epsi = einf*(1.0+((wTO2-wLO2)./(y.*(y-eye*gamma)-wTO2))); f1(:,ii)=-imag(1./epsi)./x; %loss function in the integrand weighted by 1/x end %************************************************************************* % Bose distribution %************************************************************************* Bose_mask = 1.0+(1./(exp(beta*real(y)*Escale1/1000)-1.0)); Bose_indicator = Bose_mask'*ones(1,length(x1)); Emax_tmp = 0.1:Energystep:Energymax; spec_e = zeros(length(Emax_tmp),1); %************************************************************************* % Integration for LO-phonon scattering rate %************************************************************************* % f2 = f1.*Bose_indicator; %************************************************************************* f2 = f1.*Bose_indicator; f2(isnan(f2))=0;%Array elements that are NaN [X1,Y1]=meshgrid(x1,y); for ii = 1:length(Emax_tmp) Emax1 = Emax_tmp(ii)/Escale1; yplot =Emax1-(x1-sqrt(Emax1)).^2; %parabola of integration yplot2=Emax1-(x1+sqrt(Emax1)).^2; %lower boundary of integration %quasi-particle scattering sumindicator=Y1<=ones(length(y),1)*yplot&Y1>=ones(length(y),1)*yplot2; %find values of sum indicator for each column Integral=sum(sum(f2(sumindicator)))/(sqrt(Emax1)); spec_e(ii)=Integral; end spec = spec_e; %************************************************************************* figure(1); surf(x1*kscale,real(y)*Escale1,log10(abs(f1))); az = 0; %set up the angle of view el = 90; %looking from directly above view(az, el); %set the view shading interp; %do color interpolation colorbar; axis([0 xmax*kscale -Energymax Energymax]); hold on; plot(x1*kscale, yplot*Escale1,'r','LineWidth',LW); plot(x1*kscale, yplot2*Escale1,'r','LineWidth',LW); hold off; xlabel('Wave vector, \itq\rm (m^{-1})'); ylabel(['Energy transfer, ',hbar1,'\omega (meV)']); ttl1=['\rmGaAs log_{10}(\itS\rm(\itq\rm, \omega )/\itq\rm), \itk\rm_{scale}=',num2str(kscale,... '%4.1e'),' m^{-1}, \itm\rm^*_e=',num2str(m2,'%5.2f'),... '\times\itm\rm_0, \gamma=',num2str(gamma*Escale1,'%5.2f'),... ' meV, \itT\rm=',num2str(T,'%5.1f'),' K']; title(ttl1); %************************************************************************* figure(2); surf(x1*kscale,real(y)*Escale1,log10(abs(f2))); az = 0; %set up the angle of view el = 90; %looking from directly above view(az, el); %set the view shading interp; %do color interpolation colorbar; axis([0 xmax*kscale -Energymax Energymax]); hold on; plot(x1*kscale, yplot*Escale1,'w','LineWidth',2); plot(x1*kscale, yplot2*Escale1,'w','LineWidth',2); hold off; xlabel('Wave vector, \itq\rm (m^{-1})'); ylabel(['Energy transfer, ',hbar1,'\omega (meV)']); ttl2=['\rmGaAs log_{10}(S(\itq\rm, \omega )\itg\rm/\itq\rm), \itk\rm_{scale}=',num2str(kscale,... '%4.1e'),' m^{-1}, \itm\rm^*_e=',num2str(m2,'%5.2f'),... '\times\itm\rm_0, \gamma=',num2str(gamma*Escale1,'%5.2f'),... ' meV, \itT\rm=',num2str(T,'%5.1f'),' K']; title(ttl2); %************************************************************************* dx1 = x1(2)-x1(1); dy1 = y(2)-y(1); prefactor=(2.0*m2*dx1*dy1*Escale1*0.001*echarge)/(pi*a0*kscale*hbar); %************************************************************************* figure(3); plot(Emax_tmp,spec*prefactor/1e12); ttl3=['\rmFig6.5a, GaAs, \itk\rm_{scale}=',num2str(kscale,... '%4.1e'),' m^{-1}, \itm\rm^*_e=',num2str(m2,'%5.2f'),... '\times\itm\rm_0, \gamma=',num2str(gamma*Escale1,'%5.2f'),... ' meV, \itT\rm=',num2str(T,'%5.1f'),' K']; title(ttl3); ylabel('Inelastic electron scattering rate, 1/\tau_{in} (ps^{-1})'); xlabel('Initial electron energy, \itE_k\rm (meV)'); grid on; axis([0,Energymax,0,10]); max(spec*prefactor/1e12) %*************************************************************************