缺陷形成能的定义及相关讨论
(2010-02-12 01:06)
• 标签： - 分类： Theoretical Basis
[转]缺陷形成能的定义及相关讨论
第一性原理计算方法 2009-08-21 03:59:54 阅读124 评论0 字号：大中小
原帖地址：http://www.materialssimulation.com/node/258
satchel1979
笔名: 陈征征
个人主页:
所在院校或者研究机构: California State University Northridge
研究方向:
主要研究方向为缺陷对于结构材料力学性质的影响，并发展起跨尺度相关算法。具体研究体心立方材料中线缺陷及空位、杂质等的相互作用。此外也承担一部分第一原理表面科学计算，包括表面态对于吸附原子/空位的扩散势垒的影响。对磁性问题也比较感兴趣。最常用的软件为VASP、DMOL、LAMMPS等。正在积极学习各类算法，以期实现TB、
MD、MC及至平面波第一原理计算软件包。
实际工作中，各类缺陷的形成能EF是一个非常重要的概念。它的数值可以直接反应特定缺陷形成的难易程度，材料合成环境对于缺陷形成的影响，以及复合缺陷体系的稳定性等等。因此，如何尽量准确的确定EF对于计算材料学而言是比较核心的任务之一。本文试图讨论利用第一原理方法计算EF时需要注意的几个方面，以期帮助读者理解文献，为进一步的学习和研究奠定基础。
一般而言，缺陷体系相较于完美体系，也即参考体系，经常会有粒子数目上的变化（这里粒子是指构成材料的原子以及电子）。因此可以将其视为巨正则系综。按照最为直接的定义，其缺陷形成能EF可表示为[1]
EF = G - \\sum_i Ni μi - ne μe
其中G是缺陷体系的Gibbs自由能。求和遍历该体系中包含的所有元素种类。Ni和μi分别为第i种元素的原子个数以及当前体系下该类原子的化学势。上式将电子单列出来，ne和μe分别为变化的电子数（相对于中性参考体系而言）以及电子的化学势。
这个定义虽然普适，但是在实际应用中，特别是以第一原理为基础的研究中并不常见，目前的文献都利用参考体系对上式进行了改写。首先是用体系总能E代替了 G。因此熵的贡献被忽略了。这个忽略在低温下是可以接受的。其次，变化的电子数应为该缺陷体系所处价态q的负数。很显然用价态表示EF更为直接和易懂。最后，各类元素的化学势按体系的化学式求和，应等于该种物质最小构成单元的总能。例如Y2O3，氧和钇的化学势满足2μY + 3μO = E[Y2O3]。参考体系通常是由若干最小构成单元组成的超单胞，而缺陷体系可视为参考体系增减若干原子所得。综上所述，EF可重新写为
EF(q) = E[defect,q]-E[ref]- \\sum_i ni μi + q*μe
其中ni是缺陷体系第i种原子减去参考体系中的同类原子所得。这个定义式有两个难点需要进一步讨论：其一是各类原子的化学势如何界定，我们将在后面涉及；其二是电子的化学势μe如何决定。我们需要假定：电子的得失由体系价带顶能级Ev’和能量为Eres的费米海交换电子完成。因此可以定义
μe = Ef = Eres - Ev’，
其中Ef是费米能级，其参考点是Ev’，为了保证和体系总能的能量零点一致，Ev’也必须计入EF的表达式。因此
EF(q) = E[defect,q]-E[ref] - \\sum_i ni μi + q(Ev’ + Ef)
至此，我们得到了文献上经常看到的带电缺陷的形成能表达式。该式也非常清楚的表明EF取决于缺陷的价态（q）、载流子类型及密度（Ef的位置）以及材料合成环境（μi）这三个条件。据此可以画出EF随各条件变化的关系曲线图，进而可以展开非常细致和详尽的讨论。具体实例请参考文献[2]。
特别需要注意的是q不为零时的一些修正方法。首先是如何确定Ev’。缺陷的存在可以严重的改变其附近的能带结构，使得Ev’无法直接应用[3]。一种可行的解决方法是用参考体系的价带顶Ev代替Ev’。但是这种方法需要额外考虑如何使二者对齐。常用的周期性边界条件PBC使得带电缺陷和自己的映像产生静电相互作用。这个额外的能量项使得整个实空间中的位势产生了变化，从而导致缺陷体系的能级相对于参考体系有一个平移。
这个平移量无法直接得出。因此 Laks和Van de Walle首先提出对齐二者能级的方法[3]：分别计算两个体系的静电势分布，然后进行二维的平均（也即对yz平面进行平均，得到V(x)）。取所有平面中距离缺陷最远的那个平面x0，定义“位势矫正”ΔV如下：
ΔV = V(x0)defect – V(x0)ref
因此，EF可重新写为
EF(q) = E[defect,q] - E[ref] - \\sum_i n_i μi + q(Ev+Ef+ΔV)
第二个修正也与PBC的使用有关。带电缺陷及其映像的相互作用会产生专属于缺陷轨道/悬挂键交叠的能带Ed。Ed的存在导致电子的得失不再与Ev’而转而与Ed相关。这与孤立带电缺陷的情况有所偏差。因为在后一种情况下，Ev’对应于Γ点处的最高占据态。因此Wei建议加入“色散修正”[4]：
Ecorr = q*(Ed [Γ]-Ed [MP k-mesh])
其中前者是单Γ点的计算结果，后者是标准的Monkhorst-Pack多k点计算结果。
最后一个常用的修正就是Makov-Payne修正，即电偶极矩修正:
EMP = q2* α *(1/εL)
其中α是Madelung常数，ε是介电常数，L是超单胞的边长。Makov-Payne指出，这个修正项可以有效的提高能量计算相对于超单胞大小的收敛速度[5]。
上述的三个修正都与有限大小的超单胞以及PBC有关。也就是说，如果体系足够大，这几种修正都可以不要。Castleton、H?glund和 Mirbt针对这种论点作了非常细致的研究[6]。现行计算能力下不可能采用“足够大”的体系，因此他们利用多项式拟和计算无限大体系下的缺陷形成能 EF0:
EF(L) = EF0+a/L+b/L3
EF(L)是采用不同大小的超单胞得出的形成能。而所有这些超单胞必须保证形状完全相同（也即三个方向重复次数相同）。通过对InP中11种不同缺陷不同价态的研究，他们的结论是：(1) 多项式拟和是最为可靠和准确的方法，当然计算量也极为庞大；(2) 对于单个超单胞计算，位势矫正在大多数情况下可以给出合理的答案；(3) 色散修正并不可靠，有时候甚至修正到相反的方向[7]；(4) Makov-Payne修正很多情况下并不能给出比未修正的结果更好的答案，它最能发挥效用的情况是不考虑原子弛豫的体系。
最后回过头来讨论一下原子的化学势问题。单质比较简单，μ的取值为该元素稳定状态下单个原子的能量。但是对于化合物而言，μ可以在一个范围内变化。其上限与该种元素的单质情况相同，超出这个限值，原子将会形成单质沉积下来，而不会形成化合物。此外，如上所述，各类元素的化学势按体系的化学式求和，应等于该种物质最小构成单元的总能，因此其下限也可以计算出来。同样以二元体系Y2O3为例：
μYmin = (E[Y2O3] – 3μOmax)/2
三元或更高阶的化合物情况更为复杂。因为每两种元素都可能形成沉积物。因此条件更为严格。对于杂质来讲，其下限为负无穷，对应于完全没有杂质的情况（因为此时杂质形成能为无穷大）。而其上限一般决定于与化合物中某类元素最易形成沉积物的条件。
确定方法同上。
至此，我们系统的讨论了缺陷形成能的定义式及各种修正方法。虽然大的原则相同，但是不同课题组在实际工作中往往有不同的修正方式。这使得这个问题更为复杂化。本文较倾向于文献[6]中的结论。而进一步的验证或者新修正方法的提出更应引起我们的注意。最后，我希望对本文所阐述的各个方面的理解都没有问题。如果某些地方出现了偏差或者误解。希望大家指出。
参考文献
[1] S.B. Zhang, J.E. Northrup, Phys. Rev. Lett. 67, 2339 (1991).
[2] C.G. Van de Walle, J. Neugebauer, J. App. Phys. 95, 3851 (2004).
[3] D.B. Laks, C.G. Van de Walle, P.E. Blochl, S.T. Pantelides, Phys. Rev. B45, 10965 (1992).
[4] S.H. Wei, Comput. Mater. Sci. 30, 337 (2004).
[5] G. Makov, M.C. Payne, Phys. Rev. B51, 4014 (1995).
[6] C. Castleton, A. H?glund, Mirbt, Phys. Rev. B73, 035215 (2006).
[7] 文献[2]中指出，色散修正只应该从产生浅能级的缺陷形成能，且q=0的价态上减去（与文献[5]中的修正恰好相反）。关于这个分歧我的理解是在有限。望大家指正