第四章 固体能带理论

4.3赝势方法

在大多数情况下，芯态与价态的本征谱在能量上可以明显地区分开来。化学环境对芯态波函数一般只有微小的影响，在固体能带中他们构成非常狭窄的、几乎没有色散的能带，它们的能量位置可以因化学环境而有位移。由于这一特点，在芯能级谱中常作为区分原子或化学环境的特征。然而，固体（金属、半导体、绝缘体）的电子结构性质主要是由费米能级附近的电子态决定的。在能带理论研究中，计算位于深能级的被填满的芯态代价是很昂贵的：一方面，大大增加了能带的数量；另一方面，一个全电子的、没有被屏蔽的晶体势以及芯态的波函数是坐标空间定域性极强的，因而在动量空间收敛很慢。此外，由于离子实的总能量基本不随晶体结构变化，因此，在同样的计算精度下，局限于价态、类价态的总能量计算绝对精度要比全电子方法高得多。于是，能带计算中局限于价态、类价态的方法是非常有价值的，也是非常实用的。

1 赝势的导出

赝势的导出不是唯一的。原始的赝势方法是建立于正交化平面波方法上的。对一个由许多原子组成的固体，坐标空间根据波函数的不同特点可分成两部分：(1) 近原子核局域，所谓芯区。波函数由紧束缚的芯电子波函数组成，与近邻的原子的波函数相互作用很小；(2) 其余区域，价电子波函数相互交叠、相互作用。尽管芯区的势很强地吸引价电子，但是正交化平面波方法中对价态和芯态正交的要求而产生的动能，对价态的贡献就如同一个有效的排斥势。两者的和是价态的有效势。于核的库仑势相比，这种有效势较弱。图4.3.1表示晶体中赝势、赝波函数与周期势、布洛赫波函数的关系。下面就按照这种想法来导出赝势。

图4.3.1 晶体中周期势()V a 、布洛赫波函数()b ψ与赝势()ps V c 、赝波函数()ps d ψ比较

如果用V φ和c φ分别表示晶体哈密顿算符H 的精确的价态V E 和芯态c E 的波函数，满足：

V V V H E φφ= (4.3.1)

和

c c c H E φφ= (4.3.2)

用类似正交化平面波方法构造晶体价态波函数V φ：

ps

V V cV c c

φψμφ=+∑ (4.3.3)

与正交化平面波方法不同，这里c φ是真正的晶体芯态波函数。正交化平面波中的平面波现被ps

V ψ取代，

后面就会看到这就是赝波函数。作c V 0φφ=，可得系数

ps

cV c V

μφψ=- (4.3.4) 现将V H E -作用于ps

V ψ上，有

()()()()ps ps ps ps

V V

V V c c V

V c c V V c c V c c c H E H E H E H E ψφφφψφφψφφψ??-=-+=-=- ???

∑∑∑

(4.3.5)

就有

()ps V c c V V c 0H H E E φφψ??+--=????

∑ (4.3.6)

将哈密顿算符写成

H T V =+ (4.3.7)

如果令

()ps V c c c

V V H E φφ=+-∑ (4.3.8)

则形式上就给出

ps ps V V 0T V E ψ??+-=?? (4.3.9)

ps V 就是赝势，式 (4.3.9) 就是赝波函数满足的方程。

赝势是核的库仑吸引势V 加上一个短程的、非厄米的排斥势()V c c c

H E φφ+-∑，两项之和使总的势

减弱，变得比较平坦。对这样的赝系统，用平面波展开赝波函数可以很快收敛。值得指出的是，虽然ps

V

ψ是赝波函数，但由此得到的能量并非“赝能量”，而是相应于真实晶体波函数真实价态的本征能量V E 。

赝势是非局域的，可以表示成局域的()ps L V r 和非局域的()ps

NL ,'V r r 两项之和：

()()()ps ps ps

L NL ,'','V V V δ=-+r r r r r r (4.3.10)

如果考虑原子球对称性，利用球谐函数，赝势的非局域部分可表示成

()()()()()ps NL NL ,,',,;',','Y ,,'Y ','l l lm lm l

l m

V r r r r υθ?θ?θ?υθ?*==∑∑r r (4.3.11)

一般l υ多取成径向为局域的，即

()()(),''l l r r r r r υυδ=- (4.3.12)

角部分为非局域的，这样非局域赝势的径向部分仅与轨道量子数l 有关，

()()()()()ps NL ,,,'Y ,,'Y ','l l lm lm l m

l m

V r r r lm θ?υθ?υ*==∑∑r r (4.3.13)

这种形式的赝势称为是半非局域的。

经验赝势方法

经验赝势方法 (empirical pseudopotential method, EPM) 是用实验数据拟合有限几个()V K 的值。这时，晶体势被假定表示成原子势的线性组合：

()(),a V V =--∑R r r R τ

τ (4.3.14)

在倒格矢空间展开：

()()()i a a V V S e ?=∑K r K

r K K (4.3.15)

其中是结构因子：

()a i S e ?=∑a

a

K K ττ (4.3.16)

这里()a V K 是原子势的形状因子。

对大的K ，()a V K 一般很小，再加上考虑对称性，这种原子势的形状因子只需拟合有限的几个离散的值。例如：对金刚石结构 (如C ，Si ，Ge 和-Sn α等)，原子形状因子()a V K 一般考虑23,8,11=K [单位

()

2

2a π]就已足够；而对闪锌矿结构(如GaAs, GaP 等)，阳离子和阴离子的原子形状因子()a V K 一般考

虑23,4,11=K [单位()2

2a π]。

经验赝势的拟合过程是：选取初始的()a V K ，解薛定谔方程得到()n E k 和(),n ψk r ，与实验数据 (一般是能带、态密度、响应函数等，对金属来说还可以是费米面) 作比较，修改()a V K ，重复上述过程直至得到与实验接近的结果。

经验赝势在60年代和70年代是研究半导体、金属等材料的电学、光学性质的主要的理论工具。但是，经验赝势基本上不能解决不同化学环境中的应用问题，特别是存在电荷转移的情况在经验赝势方法中是很难考虑的。

在密度泛函理论LDA 的框架下，目前经验赝势的首要用途是在现代从头计算原子赝势自洽迭代计算中作初始值使用。因此，它的原子形状因子不再作为有限的几个倒格矢的离散值，而是将它拟合成动量的解析表达式。一个典型的形式是由S. G. Louie (雷干城) 等建议的：

()()

()2122

34exp 1

a

b q b V

q b q b -=

??-+?? (4.3.17)

对一些主要的半导体材料的拟合参数i b 已列成表。这种赝势的最大优点是它的易于处理的简单形式，能以最少的代价用于任何晶体的能带计算。虽然这种经验赝势产生的电荷分布与真实的不完全一致，但它已正确地体现晶体的对称性，在一个耗费很大的第一性原理赝势自洽从头计算中给出了一个相当好的初始值。

2 模型赝势和模守恒赝势

现代能带理论中，自洽求解Kohn-Sham 方程是个在应用中非常有实际意义的基本课题。为此，构造能用于自洽计算和不同化学环境中的原子赝势是势在必行的事。模型赝势和模守恒赝势就是能用于自洽计算的两类常用的原子赝势。前者是半经验的，而后者则是第一性原理从头计算的。

离子赝势 对原子的全电子势问题，总的电子数密度()ρr 可以分成芯电子数密度()c ρr 和价电子数密度

()V ρr 之和，即

()()()c V

2

c V 1

Z Z i i ρρρψ+==+=

∑

r r r (4.3.18)

这里c Z 表示芯电子数，而V Z 则表示价电子数。薛定谔方程则为

()()()()()c V 2Coul XC

2,,i i i Z Z V V r ρρψεψ+??

-?-++=????

r r r r (4.3.19) 尝试将问题限制在原子赝势及原子赝波函数，则赝原子的薛定谔方程为：

()()()()2ps ps ps ps ps coul XC ,,i i i i V V V ρρψεψ??-?+++=??

r r r r (4.3.20)

其中

()()

occ

2

ps

ps 1

i i ρ

ψ==∑r r (4.3.21)

赝波函数()ps i ψr 可用一个么正变换ij T 与全电子波函数()i ψr 联系起来

()()ps i ij j j

T ψψ=∑r r (4.3.22)

就可以从原子的全电子的解形式上得到离子赝势 ()()()()()()()ps ps ps V c Coul c XC c 22,,i j j i ij i j i

Z Z

V T V V r r

ψεεψψρρ?=--

+-++??∑r r r r r r ()()()()()()()}

()ps ps ps XC XC c XC V Coul Coul XC XC ,,,,,,,j V V V V V V V ρρρρρρρψ+--+-+-r r r r r r r r (4.3.23) 上式第一项表示这种赝势一般是与能级有关的并且是非定域的。

模型赝势 可以看出，离子赝势的形式是非常复杂的，而最简单的离子赝势形式是定域的、且与能量无关的：

()()()V c ps V c c i Z r r r V r Z r r r ?->?=?-≤??

(4.3.24)

这里c r 作为一个可调参数用来拟合原子数据。另一相似的形式是

()()()V

c ps c i Z r r r V r A r r ?->?

=?≤??

(4.3.25) 对应一个在芯区为常数A 的赝势。这种离子赝势用动量空间变量q 表示的形式近似为

()()2V

c 24cos i Z V q r q q

π- (4.3.26)

据此，加入一些可调参数可以构造模型赝势，称为半经验赝势。它的一般形式为

()()()ps 41

2342cos exp i a V q a q a a q q

=

+???? (4.3.27) 加入()

44exp a q 因子使我们可以通过适当地选择()40a <以保证这种赝势的傅里叶展开很快收敛。常用的一些半导体材料的离子赝势的参数i a 可以从文献中得到。这种形式的赝势通常称为软芯 (softcore) 势。上述形式的赝势也是定域的。一般地离子赝势应是非定域的，但上述形式已是一个可行的较好的近似。

模守恒赝势 从密度泛函理论的观点，人们致力于确定没有任何附加经验参数的赝势，所谓第一性原理从头算原子赝势。现在在能带理论计算中最常用的是D. R. Hamann 等提出的模守恒赝势 (norm conserving pseudopotential, NCPP) 就是一种产生原子赝势的具体方法。这种赝势所对应的波函数不仅与真实势对应

的波函数具有同样的能量本征值，而且在c r 以外，与真实波函数的形状和幅度都相同 (模守恒)，另外在c r 以内变化缓慢、没有大的动能，如图所示。这种赝势能产生正确的电荷密度，适合作自洽计算。

图4.3.2 从头算原子赝势方法示意图。在c r r <的芯区内， 全电子势V 、波函数ψ被赝势ps V 和赝波函数ps ψ代替。

模守恒赝势的产生步骤如下：

(1) 求解定域密度泛函理论框架下的孤立原子全电子势()V r 的薛定谔方程。这个方程可以是

()()()()()()222

2

1202l l l

l l d u r dV r u r l l du M V u r M

dr

dr

r dr r αε+????

-

-

-

++-=????????

(4.3.28) 其中()210.5l M V αε=--，α是精细结构常数，()l ru r 是径向波函数，l ε是本征值。

(2) 选择芯半径c l r ，它必须位于全电子的价态波函数节点以外。c l r 越大，赝势越平坦；c l r 越小，赝波函数越精确。并选一截断函数f ：

()()c c exp l l f r r r r λ

??=-??

(4.3.29)

这里λ的最佳选择是3.5。

(3) 从全电子势()V r 构造一过渡原子赝势：

()()()()1c c 1l l l l V r f r r V r c f r r =-+???? (4.3.30)

这里l c 是一可调参数，以使全电子势的能量本征值与赝势的一致，并使过渡赝波函数()1l w r 除一幅度因子()l γ外，在c l r r >以外与真实波函数()l u r 一致，即当c l r r >

()()1l l l w r u r γ→ (4.3.31)

(4) 对过渡赝波函数再加上一个短程的模守恒修正项以产生模守恒的赝波函数()2l w r ；

()()()1

21c l l l l l l w r w r r f r r γδ+??=+?? (4.3.32)

其中l δ由()2l w r 的归一化条件

2

20

1l

R l drw =?

(4.3.32)

确定。l R 是求解方程 (4.3.28) 时积分所要取的最大的范围，一般取 c 2.5l l R r 就可以了。

(5) 由()2l w r 得到模守恒的原子赝势()2l V r 为

()()()()()212

c 21122

2c c 21222l l l l l l l l l l l r f r r l r r V r V r V r w r r r r r λ

λ

γδλλλλε+??++??

????

=+

-+-????

????

??????

(4.3.33) (6) 然后扣除上述赝势中的价电子屏蔽部分，得到离子赝势为

()()()()()()

XC ion 2220

4'''4'''r

l l E r V r V r dr r r dr r r r

r δρπ

ρπρδρ∞

????=-

--

?

? (4.3.34)

这里()r ρ定义为

()()2

2l l l

r n w r r ρ=????∑ (4.3.35)

l n 是价态占据数。如果所有的电子全被处理成芯电子，就没有屏蔽，则

()()ion 2l l V r V r = (4.3.36)

由上述产生从头算赝势过程可以看出，本质上从头算原子赝势是核与芯电子联合产生三有效势，是从原子的薛定谔方程从头计算得到的，这种赝势可以给出价电子或类价电子 (包括部分芯电子，如果需要的话) 的正确的电荷分布，因此适合作自洽计算。它具有较好的传递性，可用在不同的化学环境中。但它的定域性较强，使得动量空间的展开收敛较慢。已提出了一些平滑的方案部分解决这一问题。可以在文献中找到。

需要指出的是，从头算赝势的产生不是唯一的。首先，依赖于对c r 的选取。一般来说，较小的c r 定域性较强，平面波收敛慢，但传递性好，可以用于不同化学环境。

模守恒赝势的解析形式 为了应用方便，从头计算得到的模守恒赝势的数值结果还需拟合成解析函数的形式，分定域部分和非定域部分

()()ps ion ion c ,,l l l m

l m

V r V lm V r V lm lm ==+?∑∑ (4.3.37)

这里l 求和我02 或3。其中

(

)core V c 2erf

i

Z V r c

r

=-

∑ (4.3.38)

而

()()2

3

ion

231

2i

r

l l l

i i i V r A r A e α-+=?=+∑ (4.3.39)

值得一提的是，赝势的定域部分和非定域部分的选择不是唯一的。

计算中常需用到动量空间的表示，式 (4.3.38) 的动量空间表示用参数积分的方法容易得到

()22

core V

c 2

core 1

c 8exp 4i

i i Z q V q c

q πΩα=??

=-

- ???

∑ (4.3.40)

而对式 (4.3.39)，利用平面波的展开式

()()()()(),exp i 4i j Y Y l

l lm lm l m

qr π*

?=∑q r r q (4.3.41)

和球谐函数的加法定理

()()()21

Y 'Y P cos 4lm

lm

l

m

l θπ*+=∑q q (4.3.42)

其中θ为'q 与q 之间的夹角。于是可以得到

()()()()3

ion

31

c

2

,'P cos ,,',,'l l

l

l i l i i l i i V q q A Q q q A R q q θααΩ+=???=

+?

?∑ (4.3.43) 其中j l 和P l 分别是球贝塞耳函数和勒让德多项式，而

()()()()220,,'exp j j 'l i i l l Q q q drr r qr q r αα∞

=-? (4.3.44)

()()()()420

,,'exp j j 'l i i l l R q q drr r qr q r αα∞

=-? (4.3.45)

借助于调制贝塞耳函数I v 完成积分，得到

(

)221'',,'exp 24l i l i i qq q q Q q q ααα+????

+=

- ? ?????

(4.3.46)

和

(

)22221212''''1',,'I I 422422l i l l i i i i i q q qq qq q q qq R q q l αααααα-+????????

??++??=

--+-?? ? ? ? ?????????????

(4.3.47)

对0q =或'0q =，上式除0l =外为零。而

()2203',,'exp 416i i i q q Q q q π

ααα??+=- ???

(4.3.48)

和

()222205''3,,'exp 44216i i i i q q q q R q q π

αααα????++=-+ ??

????

? (4.3.49) G . B. Bachelet, D. R. Hamann 和M. Schl üter 已计算了从H (氢) 到Pt (铂) 所有原子的模守恒赝势并列成表

格形式，常称为BHS 赝势。需要指出的是，式 (4.3.43) 中的i A 需从该表格中的i C 转换：

6

1i ij j j C Q A ==-∑ (4.3.50)

其中

()()()1121110 0.5 i ij ij ki k i ij ki kj ij k i j Q S Q i j S Q Q Q i j -=-=?

?>?

????=-=? ???????

?-< ???

??∑∑ (4.3.51) 而

()()20ij i j S drr r r ΦΦ∞

=? (4.3.52)

且

()()()2

2

2 1,2,

3 4,5,6i i r i r e i r r e

i ααΦ--?=?=?=?? (4.3.53)

最后

6

1

1i l il

i A C Q -==-∑ (4.3.54) 保持计算精度的这种转换的Fortran 源码可以在文献中找到。