基于Materials Studio4.3对α-Al2O3晶体和其Cr离子替位缺陷结构的计算分析

黄浩，戴聪

摘要：利用第一性原理计算掺杂Cr^3前后的α-Al2O3晶体的态密度及能带结构，并分析其光学性质。计算结果表明：α-Al2O3晶体在掺杂了Cr离子后晶体物理结构、能带结构、电子态密度，以及各种光学性质均发生了改变，而Cr离子属于3d电子组态，说明在α-Al2O3替位缺陷结构中其能级的和跃迁导致了光学性质的改变。

关键词：α-Al2O3，第一性原理，态密度，能带结构，光学性质

0 引言：

  单晶α-Al2O3是一种优良的氧化物光学晶体和优良的激光晶体，它具有较高的导热率和较低的线膨胀系数。其结构中存在着缺陷与空位，易通过掺杂和离子注入等技术改变晶体的结构和某些光学性能。当在晶体中掺杂Ti^3,Cr^3,V^3,Ni^2等金属离子后，这些金属离子会取代晶体中的Al^3离子成为中心离子，从而使晶体的局部结构发生改变，以达到改变晶体的光学和电磁学性质的目的。当在α-Al2O3中掺杂入Cr^3离子后，形成的Cr^3：α-Al2O3即俗称的红宝石。红宝石由于其独特的光学性质在实验及生产中得到了广泛的应用，本文将通过/第一性原理计算α-Al2O3的态密度及能带结构，并通过对掺杂后的Cr^3：α-Al2O3的态密度及能带结构计算分析其光学性质。

1.理论方法（第一性原理）：

  量子力学第一性原理：仅需五个物理基本常数 —— 电子质量、电子电量、普郎克常数、光速和玻耳兹曼常数，通过求薛定谔方程得到材料的电子结构，而不依赖于任何经验常数即可以预测微观体系的状态和性质，预测材料的组分、结构、性能之间的关系，进一步设计具有特定性能的新材料。但是在实际计算过程中不难发现仅仅依靠薛定谔方程和五个基本常量是无法求解的，因此第一性原理允许采用一些数学近似，典型的有基于Hartree-Fock近似的“从头计算”方法（ab-initio）以及P Hohenberg和 W Kohn于19年提出的“密度泛函理论”（density functional theory, DFT）。

1.1 “从头计算”方法（ab-initio）

  从头计算方法仅仅利用 Planck 常量、电子质量、电量三个基本物理常数以及元素的原子序数，不再借助于任何经验参数，计算体系全部电子的分子积分，求解 Schrodinger 方程。一般所说的量子化学从头计算是建立在三个基本近似基础上的计算方法，对于三个基本近似进行各种校正的计算方法是高级从头计算。从头计算法中的三个基本近似为：

（1）非相对论近似：认为电子质量等于其静止质量，即，并认为光速接近无穷大。

（2）Born－Oppenheimer 近似（也叫绝热近似）即将核运动和电子运动分离开来处理。由于原子核质量一般比电子的质量约大 103～105 倍，分子中核的运动比电子的运动要慢近千倍。因此在电子运动时，可以把核近似看作不动。

（3）轨道近似（又叫单电子近似，由 Hartree 提出）轨道一词是从经典力学中借用来的概念，在量子化学中指单电子波函数，原子的单电子波函数称为原子轨道，分子的单电子波函数称为分子轨道。轨道近似是把 N个电子体系的总波函数写成 N 个单电子函数的乘积：                   （1.1）

其中每一个电子波函数只与一个电子坐标有关。

  从 头 计 算 方 法 基 于 上 述 三 个 基 本 近 似 后 ， 不 再 作 任 何 近 似 来 求 解Schrodinger 方程，这三个基本近似即为从头计算方法的头，而 Schrodinger方程在引入三个近似后的具体表达形式为 Hartree－Fock－Roothaan 方程。量子化学从头计算对于分子基态性质的研究一般是可靠的，但对于晶体材料来说计算量过大，必须采用进一步的近似。

1.2密度泛函理论（density functional theory, DFT）

  类似于哈特利-福克方法，密度泛函理论（DFT）也把多点字问题简化为电子方程。同时，DFT的表达形式是准确的，没有任何近似过程而且它的有效势是局域的。DFT的基本思想是可以用粒子密度函数来描述原子，分子和固体的基本物理量。因此，可以把求解问题从空间转换到了3维空间。

  密度泛函理论最普遍的应用是通过Kohn-Sham方法实现的。 在Kohn-Sham DFT的框架中，最难处理的多体问题（由于处在一个外部静电势中的电子相互作用而产生的）被简化成了一个没有相互作用的电子在有效势场中运动的问题。这个有效势场包括了外部势场以及电子间库仑相互作用的影响，例如，交换和相关作用。处理交换相关作用是KS DFT中的难点。目前并没有精确求解交换相关能 EXC 的方法。最简单的近似求解方法为局域密度近似(LDA)。LDA近似使用均匀电子气来计算体系的交换能（均匀电子气的交换能是可以精确求解的），而相关能部分则采用对自由电子气进行拟合的方法来处理。

1.3赝势方法（Pseudo-Potential Method）

在实际的固体中原子核附近库仑吸引作用是周期势场偏离平均值很远。在离子实内部势场对电子波函数的影响是很大的，电子波函数变化非常的厉害。这样，势场就不能被看成是很小的微扰势场。但是在近自由电子模型中假定周期势场的起伏是很小的。这样的矛盾必须用引进赝势这个概念来解决。  

所谓赝势，即表示在离子实的内部用假想的势能来代替真实的势能，在求解薛定谔方程的时候，不改变能量本征值和离子实之间区域的波函数。这个假想的势能就称为赝势。其基本思想是，所研究体系的哈密顿算符仅显含价电子部分，而将原子内层的全部电子连同原子核构成的核实对外层价电子的作用用适当的赝势表示，同时引入表示投影算符的势以便将价电子波函数与内层电子波函数分离开来，然后对价电子进行变分并用自洽迭代处理，计算出价电子轨道波函数和能级值等。所引入的两个势均通过全电子的原子从头计算确定，随后用于分子计算。此时，由赝势求出的波函数称为赝波函数。如图1.1所示：

图1.1 真实势，电子波函数，赝势，赝势波函数

从图中可以看出来在离子实之间的区域，真实的势和赝势给出同样的波函数。并且可以看出赝势波函数的与真实电子波函数相比更加的光滑 。

在实际程序设计和应用中，经常会把赝势方法和基函数相结合，可以得到一些常见的模拟计算方法。例如：在固体材料分析中最常见的把赝势和平面波基函数相结合得到平面波赝势方法。赝势方法也渐渐地得到不同类型的表示。例如:模型赝势，模守恒赝势（NCPP）[54]和建立在第一性原理上的超软赝势等等。其中，我们所用的软件VASP[55]主要运用的是超软赝势。

2.计算过程

2.1 理论模型的构建

  本文采用Materials Studio4.3中的CASTEP作计算，CASTEP是专为固体材料科学设计的当前最高水平的量子力学软件包之一，其采用的方法为密度泛函平面波赝势方法，可以对一些物质，如：半导体、陶瓷、金属和矿石等的晶体及其表面特性作第一性原理计算模拟。

  α-Al2O3属三方晶系R-3C空间群，用六方大晶胞表示，a=0.4759nm,c=1.2991nm,按照Pulay规则构建模型如下：

图2.1 a=0.4759nm,c=1.2991nm

2.2收敛性测试

  下面进行收敛性测试以说明结构的合理性，计算采用基于密度泛函理论的第一性原理赝势平面波法，采用广义梯度近似GGA中的PBE泛函模式来确定交换能和相关势，选取超软赝势ultrasoft在倒易空间进行计算，单点能计算精度精度为，总能量采用自洽计算，自洽计算时应用Pulay密度混合法。分别对Al和O的平面波基失能量进行迭代收敛计算，得到结果为：O（2s2 2p4）：-421.4240ev,  Al (3s2 3p1):-53.6159ev。在倒易k空间中平面波的截止能量为450ev。

2.3计算参数测试

  2.3.1平衡晶格常数

  结构优化过程中采用了BFGS算法。

  计算精度为fine；

  计算中设e-cut为340ev；

  总能计算精度为0.1000E-05 eV；

  femi能收敛精度 为0.3448E-07 eV；

  总能优化收-2.86354816E+003 eV，

  各项数据与查阅到的文献中所显示的结果相符。

  Max. Displacemen:  3.393135E-004  A ，

  Max. Force  :7.846966E-003  eV/A 

  Max. Stress:  3.886208E-002  GPa，

  得到优化晶格常数a=0.4816nm,c=1.3139nm,a/c

  与实验值相比误差为：0.05%，

  其中a的误差约为1.2%，

  c的误差为1.1%。

  2.3.2 弹性模量的测试：

      =====================================

      Elastic Stiffness Constants Cij (GPa)

      =====================================

   485.30080   142.415   127.985     0.00000     0.00000     0.00000

   142.415   485.30080   127.985     0.00000     0.00000     0.00000

   127.985   127.985   493.03615     0.00000     0.00000     0.00000

     0.00000     0.00000     0.00000   140.00585     0.00000     0.00000

     0.00000     0.00000     0.00000     0.00000   140.00585     0.00000

     0.00000     0.00000     0.00000     0.00000     0.00000   171.40583

      ========================================

      Elastic Compliance Constants Sij (1/GPa)

      ========================================

   0.0023492  -0.0005678  -0.0004624   0.0000000   0.0000000   0.0000000

  -0.0005678   0.0023492  -0.0004624   0.0000000   0.0000000   0.0000000

  -0.0004624  -0.0004624   0.0022683   0.0000000   0.0000000   0.0000000

   0.0000000   0.0000000   0.0000000   0.0071426   0.0000000   0.0000000

   0.0000000   0.0000000   0.0000000   0.0000000   0.0071426   0.0000000

   0.0000000   0.0000000   0.0000000   0.0000000   0.0000000   0.0058341

  Bulk modulus    =   251.130 +/-  5.083 (GPa)

2.4产生缺陷结构

X 1.37380497 Y 2.37950  Z 6.73799867

将此位置的铝原子替换为Cr原子，进行结构优化选Medium精度，k-point网格5x5x2得到：

（掺杂后）

与我们建立的模型相比已经发生变化：

（原始晶胞）

其中O： 2s2 ，2p4 收敛于-429.2136 eV

Al ：3s2 ，3p1 收敛于-52.7382 eV

Cr ：3s2， 3p6， 3d5， 4s1 收敛于-2457.3110 eV

3.计算结果及分析

图3.1 α-Al2O3的能带结构图

图3.2 α-Al2O3替位缺陷结构的能带结构图

  图3.1，图3.2分别是a-A1203和，a-A1203替位缺陷的能带结构。由图可知，两者的能带结构形状发生了很大变化，图1价带顶和导带底位于G点，属于直接带隙半导体。G点通常称为r点，是布里渊区的中心，即倒格子的原点，k值为(0，0，0)。a-A1203禁带宽度为8．2 eV，a-A1203替位缺陷的禁带宽度为10.2 eV。很明显a-A1203的度更窄，且费米面向高能级方向移动，进入导带，这是由于a-A1203为立方结构，比a-A1203替位缺陷的对称性好。对称性越好，具有相同位置的原子越多相距为Rs(Rs为近邻矢径)的两格点的波函数重叠的越大，使得在具有相同的k值时，能量越低能带向低能级移动，即费米面向高能级方向移动。我们再来看两种结构的总的态密度：

图3.3 α-Al2O3的电子总态密度图

图3.4 α-Al2O3的电子态密度

图3.5 α-Al2O3替位缺陷结构的电子总态密度图

图3.6 α-Al2O3替位缺陷结构的电子态密度图

  由图中总态密度图可知α-Al2O3有两个价带，上价带在-6.2ev~0.2ev之间，下价带在-19.2ev~-15.3ev，导带在5.2ev~15ev，由图4可以很清楚的知道电子态密度贡献度大的是s,p轨道电子，而在α-Al2O3替位缺陷的电子总态密度中我们发现明显有多个价带，且其导带注意由p,d轨道电子贡献，但与α-Al2O3大体类似只是发生价带左移现象，即即费米面向高能级方向移动。

图3.7 α-Al2O3的反射系数图

图3.8  α-Al2O3替位缺陷结构的反射系数图

  图3.5可知反射系数在0~22ev时随着能量的增大而增大，但在22ev时发生急转。图6则是截然不同的结果，在0~2ev的时候反射系数虽能能量的增加急剧减少，2ev~20ev之间有三个波峰和两个波谷，反射谱是能带体系的最直接的表现形式，反射谱中的各个峰归属于不同的带间跃迁，都是固体电子在光电磁波场微扰作用下发生带间跃迁的宏观表现，显然在产生缺陷时对晶体的光学性质有很大的影响。

图3.9 α-Al2O3的吸收系数图

图3.10 α-Al2O3替位缺陷结构的吸收系数图

  由图3.10可知α-Al2O3的替位缺陷结构有明显选择性的吸收，在低能级有0~19.4ev的时候有4个吸收波峰，在20~50ev高能级也有4个吸收波峰，与能带结构相对应，显然这跟电子的激发态有关。α-Al2O3晶体在掺杂了Cr离子后晶体物理结构、能带结构、电子态密度，以及各种光学性质均发生了改变，而Cr离子属于3d电子组态，说明在α-Al2O3替位缺陷结构中其能级的和跃迁导致了光学性质的改变。

4.结论：

  通过本实验，可以得到以下结论：

  1.本文通过第一性原理，利用Materials Studio4.3分别计算了a-Al2O3晶体和其Cr离子的替位缺陷结构的能带结构，态密度，光学性质，通过比较发现掺杂Cr离子使得费米面向着导带移动，下价带，上价带态密度主要由Al和O的s,p轨道电子贡献形成，而导带主要由Cr离子的d轨道电子贡献。

  2.Cr离子的能级和电子跃迁是缺陷结构各项性质发生改变的主要原因。

  3.在晶体中掺杂某些特定的物质可能很大程度上改变其原本性质，这就为我们开发新物质构成，或者按照人们的意愿构造新的物质结构提供了一个方向。

  4.使用Materials Studio4.3软件进行计算，我们不仅可以获得物质结构的各项参数，还可以对未知物质进行各项性质的预测。

参考文献：

1 γB28光学性质的第一性原理研究（马松山等）

2 萤石结构TiO2的电子结构和光学性质（关丽等，物理学报第57卷第1期2008年1月）

3 白宝石及红宝石晶体的电子结构（林理彬，四川大学学报）

4 a-Al2O3电子结构对其力学性能的贡献（刘东亮等，计算机与应用化学第24卷第9期）

5 A1203的能带结构和光学函数的关系（高丽红等，稀有金属材料与工程第38卷增刊2）