(*Программа «Two-temperature Thermal Conductivity in Solid Metal due to only Electron-Phonon Scattering» -

сокращенное название суть “Electron-Phonon 2T-sigma”*)

 

(*Вычисляет частоту столкновений s-электрона с импульсом p с фононами и затем интегрированием по p коэффициент теплопроводности*)

(*температура электронов TeK, ионов TiK*)

(*алюминий*)

 

 

TeK=800. ;

(*температура электронов в К*)

 

mu=  -0.0000129462;

(*химпотенциал электронов в а.е.*)

 

TiK=800.

(*температура ионов в К*)

 

Ti=TiK/11605./27.2 ;

(*температура электронов в а.е.*)

 

Te=TeK/11605./27.2 ;

(*температура электронов в а.е.*)

 

rho=2.6989;

(*плотность в г/см^3*)

 

A=26.982 ;

(*молярная масса *)

 

z=3.

(*заряд иона в а.е.*)

 

M=1836.*A ;

(*масса атома в а.е.*)

 

na= rho/A*6.022*0.148*0.1;

(*концентрация атомов в а.е.*)

 

re=M*na ;

(*плотность в атомных единицах*)

 

qD=(6*Pi*Pi*na)^(1./3)

(*квазиимпульс Дебая в атомных единицах*)

 

si=6.260 ;

(*скорость продольного звука в км/сек*)

 

vat=3.*(10^5)/137 ;

(*а.е. единица скорости*)

 

s1=si/vat

s=s1/5.1

(*скорость звука в атомных единицах*)

 

ens= -11.1/27.2;

(*ens - дно s-зоны в а.е.*)

 

mu0=0;

(*положение химпотенциала электронов при Te=0*)

 

m=((3.*z*na*Pi*Pi)^(2./3))/(2*(mu0-ens))

(*масса s-электрона в а.е.*)

 

(*вычисление константы томас-фермиевского экранирования кулоновского взаимодействия за счет s-электронов. Это kap*)

Xx=100

 

te=(ens-mu)/Te ;

 

hs[tk_]= Sqrt[tk]*Exp[te-tk]/((Exp[te]+Exp[-tk])^2);

 

UK= NIntegrate[hs[tk],{tk,0.,Infinity}];

 

dnsdm= Sqrt[2.]/(Pi*Pi)*((Sqrt[m*Te])^3)*UK/Te ;

 

kap=Sqrt[4*Pi*dnsdm]

(*обратная длина томас-фермиевского экранирования в а.е.*)

 

(*вычисление производной химпотенциала по электронной температуре dmu/dTe*)

 

mus=mu-ens ;

(*химпотенциал, отсчитанный от дна зоны проводимости*)

 

ef1[x_]=Exp[-x*x/Te]/( Exp[-x*x/Te]+Exp[-mus/Te]) ;

 

is1= NIntegrate[ef1[x],{x,0.,Infinity}] ;

 

ef2[x_]=(3*x*x-mus)*ef1[x] ;

 

is2= NIntegrate[ef2[x],{x,0.,Infinity}] ;

 

mut= -is2/is1/Te

(*производная химпотенциала по электронной температуре dmu/dTe в а.е.*)

 

 

xx=11111

(*вычисление теплоемкости s-электронов на единицу объема cvs*)

 

fcv[x_]=(3*x*x*(mut-mus/Te)+5*x*x*x*x/Te)*ef[x] ;

 

icv= NIntegrate[fcv[x],{x,0.,Infinity}];

 

cvs=Sqrt[2.]/Pi/Pi*((m)^(3/2.))*icv

(*теплоемкость единицы объема s-электронов в а.е.)

 

(*вычисление среднего квадрата скорости s-электронов vs2*)

ezv=Exp[-xv] ;

 

hes=Exp[(ens-mu)/Te] ;

 

hs[xv_]=Sqrt[xv]*ezv/(hes+ezv);

 

ihs=NIntegrate[hs[xv],{xv,0.,Infinity}];

 

hvs[xv_]=xv*hs[xv] ;

 

ihvs=NIntegrate[hvs[xv],{xv,0.,Infinity}];

 

vs21=2*Te/m*ihvs/ihs

(*средний квадрат скорости s-электронов vs2 в а.е.*)

 

pF=Sqrt[2*m*(mu0-ens)]

(*Ферми-импульс в а.е.*)

 

en[p_]=ens+p*p/(2*m) ;

(*закон дисперсии электронов*)

 

(*вычисление частот столкновений электрона с импульсом p с фононами*)

 

cnu=Pi/(re*s*s)*s*((4*Pi*z*na)^2)*2*Pi*m/((2*Pi)^3) ;

(*постоянный множитель для этих частот*)

 

w[q_]=cnu*q*q/((q*q+kap*kap)^2) ;

(*функция импульса фонона q, содержащая квадрат матричного элемента и множитель от фазового объема фонона*)

 

(*вычисляет случай 1 с испусканием фонона*)

 

p1=Sqrt[p*p] ;

 (*модуль импульса электрона после испускания фонона. Здесь пренебрегаем изменением абсолютной величины импульса электрона в результате испускания фонона.*)

(*Выражение для частоты столкновений электрона с импульсом с фононами при их испускании (индекс i означает irradiation). Это выражение представляет собой интеграл. Причем пределы интегрирования зависят от значения импульса p электрона. Как эти пределы зависят от импульса, поясняется в сопроводительном тексте «Пояснения к программе “Electron-Phonon 2T-sigma & kappa”». Кроме того пояснения приведены ниже при описании рис. 1. В программе подынтегральная функция в выражении для  записывается в следующем виде:*)

 

wp[p_,q_]=w[q]/p*(q*q-(p-p1)*(p-p1))/(2*p*p1)*Exp[s*q/Ti] /(Exp[s*q/Ti]-1) ;

 

(*Комбинация из нескольких постоянных сомножителей, входящая в интеграл , была обозначена выше как «cnu» (constanta nu)*)

 

 Рис. 1. На рисунке изображена область интегрирования в интеграле  в плоскости . Это плоскость электронного импульса и импульса фонона q. Метка  на оси q на рисунке соответствует дебаевскому импульсу. Метка фермиевский импульс, где z – число электронов в зоне проводимости в расчете на один атом. Прямая  соответствует точному рассеянию электрона назад. Квазиимпульс фононов ограничен значением . Поэтому область интегрирования ограничивают две прямые:  и .

 

 

 (*     0<p<0.5*qD      *)

nup1[p_]= NIntegrate[wp[p,q],{q,0.0001, 2*p }]  ;

(*В этой строке находится частота столкновений при испускании фононов электрона с импульсом 0<p<0.5*qD. См. рис. 1. Из-за излома на границе области интегрирования мы разбиваем интеграл  на два – один с 0<p<0.5*qD и другой с  0.5*qD <p<Infinity             *)

 

 

(*        0.5*qD <p<Infinity         *)

nup3[p_]= NIntegrate[wp[p,q],{q, 0.0001, qD}] ;

(*частота столкновений при испускании фононов электрона с импульсом 0.5*qD <p<Infinity *)

 

 

xx=123

 

 (*вычисляет случай 2 с поглощением фонона*)

 

pm=Sqrt[p*p] ;

(*импульс электрона после поглощения фонона*)

 

wm[p_,q_]=w[q]/p*(q*q-(p-pm)*(p-pm))/(2*p*pm)/(Exp[s*q/Ti]-1) ;

(*подынтегральная функция для вычисления частот столкновения при поглощении фонона. Эта строка соответствует интегралу *)

 

(*     0<p<0.5*qD    *)

num2[p_]= NIntegrate[wm[p,q],{q, 0.0001, 2*p }]  ;

(*частота столкновений при поглощении фононов электрона с импульсом 0<p<0.5*qD *)

 

 

(*        0.5*qD <p<Infinity         *)

num4[p_]= NIntegrate[wm[p,q],{q, 0.0001, qD}] ;

(*частота столкновений при поглощении фононов электрона с импульсом 0.5*qD<p<Infinity *)

 

xx=12345

 

nuF=nup3[pF]+num4[pF] ;

(*частота столкновений с фононами электрона с импульсом Ферми в а.е.*)

 

nuFa=nuF*41

(*частота столкновений с фононами электрона с импульсом Ферми в единицах 10^15 1/сек*)

 

(*вычисление коэффициента теплопроводности*)

 

faw[p_] = 1./(( Exp[(en[p]-mu)/(2*Te)]+ Exp[-(en[p]-mu)/(2*Te)])^2)*

((en[p]-mu)/Te+mut)/Te*(en[p]-mu)*

p*p/m/m*2*Pi*p*p*2/((2*Pi)^3)*(2/3.);

 (*faw - часть подынтегральной функции для вычисления коэффициента электропроводности. Коэффициент электропроводности представляет собой интеграл

.

faw  равно подынтегральному выражению без множителя с частотой , где *)

 

iaw2=NIntegrate[faw[p]/(nup1[p]+num2[p]),{p, 0, 0.5*qD}]

(*полный вклад в теплопроводность электронов с импульсом 0<p<0.5*qD *)

 

iaw4=NIntegrate[faw[p]/(nup3[p]+num4[p]),{p, 0.5*qD, 50. }]

(*полный вклад в теплопроводность электронов с импульсом 0.5*qD <p<Infinity*)

 

kappa=  iaw2+  iaw4

(*коэффициент теплопроводности в а.е.*)

 

kappaSI=kappa1*1.090523*(10^4)

(*это – коэффициент теплопроводности в Вт/(м сек К) (система СИ) *)

 

 

nuat=(1/3.)*cvs*vs2/iaw ;

(*это – частота s-i столкновений в атомных единицах*)

 

nus=nuat*41

(*это – частота в единицах 10^15   1/сек*)