固体能带计算方法

第23卷第1期贵州大学学报(自然科学版)V o.l23N o.1 2006年 2月Journa l o f G uizhou U n i ve rsity(N atura l Sc i ences)F eb.2006文章编号 1000-5269(2006)01-0068-05

固体能带计算方法*

闫万,谢 泉,朱林山,金石声

(贵州大学计算机科学与工程学院,贵阳550025)

摘 要 作者综述了固体能带常用的计算方法,并采用第一性原理赝势平面波方法计算了Si和Ge的电子能带,对计算结果进行了分析。

关键词 固体能带;赝势平面波方法;第一性原理

中图分类号 O481 1 文献标识码 A

0 引言

早期,人们对于晶体中电子能带结构的了解只是定性的,随着科技的发展情况就完全不同了。一方面由于纯单晶的制备,低温、微波等技术的发展,人们可以通过回旋共振、磁阻、磁光、光谱等实验来直接测定若干晶体的电子能带结构;另一方面,由于电子计算机的广泛应用,在理论上也可以对能带结构进行更精确的计算。能带理论在材料的研究中已经成为非常有用的工具[1,2]。

对于三维的周期场中的单电子问题,通常采用各种近似方法求解,首先,选取某个具有布洛赫函数形式的完全集,把晶体电子态的波函数用此函数集展开,然后代入薜定谔方程,确定展开式的系数所必须满足的久期方程,据此可求得能量本征值,再依照逐个本征值确定波函数展开的系数。选择不同的函数集合,有不同的计算方法。借助快速大容量的电子计算机,使实际的计算更快捷。本文将对这些常用的计算方法进行论述,并采用第一性原理分子动力学方法对S i和Ge的能带进行了计算。

1 固体能带计算的方法

1 1 原子轨道线性组合法(LCAO)[3]

紧束缚方法(T igh-t B anding)又称为原子轨道的线性组合法,简称LCAO方法。第一次由F.B l o ch在1929年推出的,其中心思想就是用原了轨道的线性组合来作为一组基函数,由此而求解固体的薛定谔方程。这个方法是基于这样的物理图像,即认为固体中的电子态与其组成的自由原子态差别不大。按B loch 观点,当晶体电子离某个原子较近时,原子势场对电子有较强的束缚作用,此时电子又有共有化运动,它有等同的几率出现在每个原子范围,就是说,电子的共有化运动带有强烈的原子态的色彩。因此就用原子轨道的线性组合来构成电子的共有化轨道,即B l o ch函数。

紧束缚方法在绝缘体的能带结构研究中是很成功的。由于原子轨道处在不同的格点上,由它们组成的基函数一般是非正交和。因此,必然会遇到中心积分的计算问题,而且本征方程形式也不简便。为了克服这些困难,人们道过Slater-Koster参量法,W o Ao H arr i s on键轨道近似以及O lc Ao方法解决。

1 2 正交化平面波方法和赝势方法

1 2 1 平面波方法

平面波方法是计算方法中最简洁,物理定义较明显的晶体电子能带计算方法,其方程形式为

K ?(k+K)2-E n(k)」 KK +V(K -K) c n(k,K)=0

*收稿日期:2005-11-11

基金项目:贵州省教育厅重点基金(05JJ002),贵州大学人才引进基金资助项目(04RCJJ001),湖南省自然科学基金(04JJ3030),湖南省教育厅青年基金资助项目(03B006),国家人事部留学归国人员择优资助项目(200403)。

作者简介:闫万(1978-),女,贵州大学计算机科学与工程学院研究生,研究方向:能带计算。

齐次线性方程组有非平庸解的条件是

det

(k +K )2-E n (k ) KK +V(K -K )=0

对角化这个行列式,即可确定本征值,进而得到展开系数。

1 2 2 正交化平面波方法

[4]平面波方法计算固体能带的优点是概念简单明了,缺点是收敛很慢,即使借助现代高速大容量的计算机作计算也不现实,于是C .H erring 提出了正交化平面波方法,其基本思想是展开波函数的基不但含有动量k +K 较小的平面波成分,还有原子核附近具有大动量的孤立原子波函数的成分,并且与孤立原子芯态波函数组成布洛赫正交,这种基函数称为正交化平面波。

定义正交化平面波| k +K >为

k+K >=

k +K >- c c (k ,r)>< c (k,r)k +K >

久期方程为

det < k+K H k+K >-E k =0一个正交化平面波在远离原子核处的行为象一个平面波,而在近核处具有原子波函数的迅速变化的特征,这样,就能用这种基较好的描写价态的特征。OP W 方法已被广泛应用于金刚石结构、立方闪锌矿结构和六角纤维结构半导体的能带计算。

1 2 3 赝势方法

[5]1959年,基于正交化平面波方法,Ph illips 和K lei n m an 提出了赝势的概念。正交化平面波具有在势函数中加入排斥势的效应,并证明在整个金属中有效势显得很小,这就能解释电子在金属中的行为,也使我们有了一个将电子波函数采取微扰展开的模式,因为电子导带与离子核之间的相互作用 赝势是很弱的。总的晶格赝势就是单个离子与导带相互作用的赝势之和。

赝势:V ps =V +

c (E V -E C ) c >< c v >和 c >分别表示晶体哈密顿算符的精确价态E c 和芯态E v 的波函数。

赝势波函数 ps V >所满足的方程

?T +V p s -E V 」| ps V >=0

对势场的近似主要指局域密度交换关联近似和型状近似。在固体中,对原子的坐标空间波函数根据其特点,分为两部分:

1)近原子区域所谓的芯区。波函数由紧束缚的芯电子波函数组成,与近邻的原子的波函数相互作用很小。

2)其余区域,价电子波函数相互交叠,有较强的相互作用。发生相作用电子中主要贡献来自包围芯区的价电子,其余电子的贡献可以忽略。

赝势实际上是一种算子,但可以近似的将它处理成原子间距离的简单函数,叫做局部赝势。这样做的结果是,只要给定周期性排列的原子的位置,系统的能量就可以算出来。

确定赝势的方法不是唯一的,主要有以下几种:

1)经验赝势方法,利用实验数据拟合有限几个V (K )的值,主要用途是在现代从头算原子赝势自洽迭代计算中作初始值使用。

2)模型赝势是半经验的,能用于自洽计算的原子赝势。

3)模守恒赝势[6]

是第一性原理从头算原子赝势,是核与芯电子联合产生的有效势,是从原子的薛定谔方程从头计算得到的,它可以给出价电子或类价电子的正确电荷分布,适于作自洽计算。

4)模非守恒赝势(超软芯赝势),其优点是容易选择芯区的赝势波函数,减少了必须的平面波函数的数目,较大的减轻了计算工作量。

69 第1期闫万等:固体能带计算方法*

1 3 原子球近似

1 3 1 M uffi n -tin 势[7]

主要思想:把原胞分为两个区域,以原子为中心的球形区及球外的区域。在球内,取球对称势;在球外取常数势,通常选取适当的能量零点,使此常数为零(如图1)

。

图1 M u ffi n -tin 球区域Ⅰ及球外区域Ⅱ的外法线方法在一个原胞中势场可以表示为

V( )=V( ) ( < V )0 ( V

)M u ffi n -tin 势的选取的常用方法:对一个原子周围

的势场贡献最大的是中心原子势场,然后还有它的最近

邻原子对这部分空间的势场。次近邻和远近邻原子的

贡献逐渐减弱。外围原子的效应,在很大程度上通过

M u ffi n -tin 球间区域的平缓势场自动引入M uffin -ti n 球

的势函数。

1 3

2 缀加平面波(AP W )[8]基于M uffi n -ti n 势的选取,可以建立起一套缀加平面波。将原胞分为两个区域,写出AP W 的基函数为

(k ,r)=4 e

ik r V l=0 +l m =-l i l j L (k V )Y *l m (k )Y lm ( )R l (E , )/R l (E , V ) ( V )e ik r

( V

)求解能量本征值及本征函数中的展开系数的久期方程为det H ij +S ij -E ij =0.

AP W 矩阵元 M ij =H ij +S ij -E ij = c (k 2

j -E ) ij -4 v 2v e ik ij r v G v ij .

1 3 3 线性化的缀加平面波方法(LAP W )

LAP W 解决了AP W 方法的在概念和计算上的困难:

(1)建立M uffi n -ti n 球的径向解时,隐含了能量E ,使久期方程成为E 的超越方程,对E 的依赖关系复杂;

(2)当径向解节点落在M uffin -ti n 球面上时,久期方程有奇性;

(3)基函数的导数在球面上不连续。

LAP W 基函数:

L (k, )= lm

a l m R L (E )+

b l m R l (E )Y lm ( ) ( < v ) -

12c d ik ( V )LAP W 与AP W 基函数不同之处在于:球内的函数不再是能量本征值的函数,而是某个确定值E ,有待选定。

LAP W 方法的久期方程和矩阵元:det H ij -ES ij

=0其中H ij = d 3r *(k i ,r),H (k j ,r )S ij =

d 2r *(k i ,r) (k i ,r),1 3 4 全电势LAP W 方法(FLAP W )

[9]前述几种能带计算方法,采用的都是MT 势,在原子球内有球对称形式,但这种势对于非密堆积的开结构,如半导体、薄膜等不是很好的近似。FL AP W 的基本思想:一个局域电荷外面的势场仅通过多极矩而依赖于电荷,为了求得晶体中球间区域内的势场,只需要知道缓变的球间区内的电荷密度和各个球内电荷的多极矩(快速收敛的)傅立叶表 70 贵州大学学报(自然科学版)第23卷

示。

多极矩并不和唯一的电荷密度对应,可以把球内的真实电荷用一个赝电荷来代替,使之和真实电荷分布有相同的多极矩,但有快速收敛性的傅立叶表示。这样,将M uffin -ti n 球内(尤其是近核处)快速变化的电荷分布用另一个平滑的赝电荷分布来代替。然后,求解泊松方程,得到球外空间及球面上的势场分布。最后,用球内的真实电荷密度,用上面求得的球边界上的势场作为边界条件,通过格林函数方法可以求得球内势场的简谐展式。

迄今为止FLAP W 方法已经逐渐发展成为基于DFT 基础上的计算最为准确的一种第一性原理计算方法。

2 第一性原理计算硅和锗的能带

基于密度泛函理论(DFT)[10],使用第一性原理赝势平面波方法对S i 和Ge 的能带进行了计算,该方法有机地结合了赝势和平面波基组,是目前计算机模拟实验中最先进、最重要的赝势能带方法之一。

DFT 中单电子运动的薛定谔方程可表示为

- 22

- q Z q r -R q + (r )

r -r d r +V XC i (r)= i i (r)(1) (r)= i n i i (r)

2(2)(2)式中 i (r )表示单电子波函数,n i 表示本征态的电子占据数, (r)表示多电子密度。(1)式中第一项表示体系中有效电子动能;第二项表示体系中各原子核对电子的吸引库仑势,其具体形式采用了规范保守赝势(nor m -conserv i n g pseudopoten tial)

[11];第三项表示电子库仑势;第四项表示交换关联势,其具体形式采用了局域密度近似(LDA)[12]。计算中应用周期性边界条件,将单电子的轨道波函数用平面波展开,

我们测试后表明平面波截断能量E cu t =300eV 使单胞能量收敛到10-5eV ,选取ScreenedH artree -Fock +Lo -

ca lDensity Approx i m ation

[13]来处理交换关联能,交换关联势的公式采用了CA-PZ [14]的形式,布里渊区积分采用M onkhors-t Pack

[15]形式的特殊高对称k 点方法,使用了

10个特殊k 点。经过计算得到S i 和Ge 的能

带如图2所示。图2 Si(a)和G e(b)的能带

比较Ge 和S i 的能带可知,Ge 和S i 的价带是相似的,价带宽带约为12e V,顶点都在G 点,但是它们的 71 第1期闫万等:固体能带计算方法*

导带结构不同:

(1)S i 的导带底在 轴靠近X 点处的 1态,而Ge 的导带底在L 点的L1态;

(2)在G 点靠近能隙的导带G 2 和G 15两个能级的顺序对两种材料来说是不相同的:

对S,i G 2 >G 15,而对于Ge ,G 2

通过对能带计算结果的分析和比较可知,采用第一性原理赝势平面波方法可以得到较为准确的结果。3 结束语

自量子力学发展以来,固体能带理论成为深入研究固体物理的理论基础,对学科和技术的发展起到了很大的推动作用。随着科学的飞速发展,新材料和新现象层出不穷。与此同时,能带理论在方法上,相应地也有不少的进展,再加上计算机技术的日新月异,方法中的创新得以在计算中实现。能带理论的计算结果,不仅能阐明一些新的现象,而且也可以预测新的材料和有关的物理性质,因此,借助大容量高速度的计算机对固体能带进行深入研究已经成为必然趋势。

参考文献

[1]L ars on P ,M ahan ti S D.E lectron i c S truct u re and T ransport of B i 2Te3and BaB i T e3[J].P hy sical R e v i e w B ,2000,61(12):8162~8171.

[2]L ars on P ,M ahan ti S D.S truct u ral S t ab ili ty ofN -i con tai n i ng H al-f heu sl er C o m pounds[J].Phys i cal Revie w B,2000,62(19):12754~12763.

[3]B l och F .Z .P hys .,1928,52:555

[4]H erri ng C .P hys .R e v .,1940,57:119;H erri ng C and H illA G .P hys .R e v .,1940,58;152

[5]J . C.Ph illi ps and L .K l ei nm an,P hys .R e v .116,287(1959).

[6]H a m ann D R ,Sch l t erM and Ch iang C.Phys .R ev .L ett .,1979,B15:2154

[7]S l ater J C.P hy s .Rev .,1937,51:846;1953,92:603

[8]Loucks T L .Aug m en t ed PlaneW aveM et hod[M ].N e w York:W.a .B en j a m i n Inc .,1967

[9]W i m m er E,K rakauer H,W ei nertM and Free m an A ,J .Phys .R e v .,1981,B24:864

[10]H ohenberg ,P .;Kohn,W.Inho m ogeneous el ectron gas .Phys .R ev .B ,136,864-871(1964).

[11]H a m ann ,D .R .;Sch l u ter ,M.;C h i ang ,C .Nor m-Con s ervi ng Pseudopoten tial s[J].Phys .R ev .L ett .,43,1494-1497(1979).

[12]Perde w,J .P.;Zunger , A.Sel-f i n teracti on correcti on to d ens i ty -f unctional approxi m ations f or many -electron s yste m s .Phys .Rev .B ,23,5048

-5079(1981)

[13]Sei d ,l A.;Gorli ng ,A .;V og,l P .;M aj e w sk,i J . A.;Levy ,M. Generalized Kohn -Sh a m sch e m es and t he band-gap p rob le m ,P hy s .Rev .

B ,53,3764-3774(1996).

[14]C eperl ey ,D .M.;A lder ,B .J .Ground State of t he E l ectron Gas by a S t ochasticM et hod[J].P hys .R e v .L ett .,45,566-569(1980).

[15]M onkhorst ,H.J .;Pack ,J . D.Speci al poi n ts for Brill ou i n-zone i n tegrations[J].Phys .R ev .B,13,5188-5192(1976).

TheM et hods of Soli d Band Structure Calculation

YAN W anjun ,X I E Quan ,Z HU Linshan ,JI N Sh isheng

(Coll ege of Co m puter S ci en ce&Engi neering ,Gu iz hou Un i versity ,Gu i yang 550025)

Abstract The calcu lati o n m et h ods of the so li d band struct u re are rev i e w ed,the band str ucture o f Si and G e are calcu lated by t h e first princi p le m ethods based on plane -w ave pseudopotential theo r y ,and the ca-l culati o n resu lts are ana l y zed .

K ey w ords so li d band structure;p lane -w ave pseudopotentialm ethod ;first pri n c i p le 72 贵州大学学报(自然科学版)第23卷