低压缩硬质材料氮碳酰亚胺的第一性原理计算 本文关键词：亚胺，压缩，原理，材料，计算
低压缩硬质材料氮碳酰亚胺的第一性原理计算 本文简介：硬质材料碳氮酰亚胺的电子结构特征摘要：采用第一性原理计算方法研究了氮碳酰亚胺C2N2(NH)的晶体结构﹑力学和电子特性。计算得到的晶格参数、实验数据均与先前的理论结果很好吻合。且计算表明常压下正交结构的C2N2(NH)力学性质稳定。计算得到的高体弹性模量和剪切模量表明C2N2(NH)可做为超强抗压缩
低压缩硬质材料氮碳酰亚胺的第一性原理计算 本文内容：
硬质材料碳氮酰亚胺的电子结构特征
摘要：采用第一性原理计算方法研究了氮碳酰亚胺C2N2(NH)
的晶体结构﹑力学和电子特性。计算得到的晶格参数、实验数据均与先前的理论结果很好吻合。且计算表明常压下正交结构的C2N2(NH)力学性质稳定。计算得到的高体弹性模量和剪切模量表明C2N2(NH)可做为超强抗压缩硬质材料的潜在候选者。此外，文章还讨论了其各向异性和德拜温度的大小。另外电子能态密度和电子局域函数分析表明，存在于CN4四面体中的C-N强共价键是其具有高体弹模量和剪切模量以及小泊松比的主要原因。
关键词：第一性原理计算
力学和电子性质
弹性各向异性
电子局域函数
Abstract
First
principles
calculations
are
performed
to
investigate
the
structural,mechanical,and
electronic
properties
of
C2N2(NH).
Our
calculated
lattice
parameters
are
in
good
agreement
with
the
experimental
data
and
previous
theoretical
values.
Orthorhombic
C2N2(NH)
phase
is
found
to
be
mechanically
stable
at
an
ambient
pressure.
Based
on
the
calculated
bulk
modulus
and
shear
modulus
of
polycrystalline
aggregate,C2N2(NH)
can
be
regarded
as
a
potential
candidate
of
ultra-incompressible
and
hard
material.
Furthermore,the
elastic
anisotropy
and
Debye
temperatures
are
also
discussed
by
investigat-ing
the
elastic
constants
and
moduli.
Density
of
states
and
electronic
localization
function
analysis
show
that
the
strong
C–N
covalent
bond
in
CN4
tetrahedron
is
the
main
driving
force
for
the
hing
bulk
and
shear
moduli
as
well
as
small
Poisson’s
ratio
of
C2N2(NH)
Keywords:
First
principles
calculations
Mechanical
and
electronic
Elastic
anisotropy
Electronic
localization
function
1.
引言
探索和设计超强抗压缩超硬材料是一个具有挑战性的长期课题，在基础科学和技术应用方面有着非常重要的意义[1]。
超硬材料中除了已知的金刚石和立方氮化硼（c-BN），大量的实验和理论工作一直致力于新型硬质材料设计与合成。众所周知，由轻元素B,C,N和O所组成的强三维共价化合物是超硬材料的强有力候选者，例如BC2N[2]﹑BC5[3]﹑B6O[4]等。在这些共价化合物中，CN化合物（C3N4，CN，C3N等）因拥有与金刚石可相比的高体弹性模量和硬度值而备受关注。尤其是β-C3N4
[5]，其理论体弹模量值高达~
430
GPa，因此实验上一直在致力合成具有这种结构的碳氮相。虽然许多研究人员试图通过使用各种合成技术，包括高压和高温（HP-HT）法去合成CxNy晶体，但是目前还没有可靠的证据来证实这种碳氮相的成功合成。化学计量N:C比值为1.3-1.5的CxNy晶体在以往的工作中报道过[6,7]。然而，这些材料均是是非晶或微晶材料。它们的结构和化学成分没有得到很好的表征。此外，在高温高下合成的材料一般很难在常压条件在保存。
最近，E.
Horvath-Bordon等人[8]首次在高温高压条件下，利用激光加热金刚石对顶砧技术，合成出N:C比为3:2的单晶碳氮酰亚胺C2N2(NH)。这个新的氮碳单晶致密相可以恢复保留至常压环境。使用TEM﹑EELS﹑SIMS技术，人们发现C2N2(NH)具有类似于
Si2N2(NH)的缺陷纤锌矿结构。此外，E.
Horvath-Bordon等人计算得到C2N2(NH)的体弹性模量为277
GPa。最近，这一具有缺陷纤锌矿结构的碳氮酰亚胺C2N2(NH)
在采用双氰胺为原材料使用激光加热的高压实验中，被同步辐射X射线衍射和拉曼散射再次确定[9]。在卸压时，样品中没有新相变的发生，C2N2(NH)
稳定存在至常压。然而在理论上关于这种新的单晶碳氮相的物理性质和化学性质很少研究，如弹性常数﹑热力学性质﹑电子特性等等。因此，作为一种新的氮碳相，人们可能会期望其具有优异的机械和其它新颖的物理性质和化学性质。在本文中，我们采用第一性原理总能量计算发现C2N2(NH)是一个宽带隙绝缘体。此外，计算表明C2N2(NH)在常压下具有高体弹性模量和较大的剪切模量。C2N2(NH)的这些优异的物理性质归因于其内部有趣的化学成键。此外，文章还还研究了具有类似结构的Si2N2(NH)来进行比较。
2.
计算方法
密度泛函理论（DFT）[10，11]计算中电子相互关联势采用局域密度近似（LDA）[12]和广义梯度近似（GGA）[13]，结构局域优化以及电子性质的计算由VASP软件包完成[14]。电子和离子之间的作用采用投影缀加波（PAW）“冷冻核”全势[15]，其中H：1s1,C：2s22p2,N:
2s22p3﹑Si：3s23p2视为价电子。对于总能量计算，平面波展开截断能为520
eV，布里渊区特殊k点取样采取的是Monkhorst-Pack（MP）方法[16]，C2N2(NH)
和Si2N2(NH)
的k点网格为10×
10
×
10。采用“应变-应力”方法应变应力张量，计算晶格弹性常数。相应的多晶体弹模量、剪切模量、杨氏模量以及泊松比等通过Voigt-Reuss-Hill方法
[17]。原子受力收敛到0.001
eV/?，原子受力的最大应变为
0.2％。
3.
结果与讨论
3.1
晶体结构
实验发现单晶碳氮酰亚胺C2N2(NH)
在P
>
27
GPa和T
>
2000
K时，晶体结构所属的空间群为Cmc21。单胞中具有四个分子式（见图1），其中所有的C原子形成CN4四面体结构，两个不等价的N1和N2原子与三个近邻原子以三个共价键相连接。我们优化C2N2(NH)和Si2N2(NH)
晶格的几何形状和离子的位置。表1中给出了采用GGA和LDA近似计算的平衡晶格参数与键长结果，其他理论值计算值以及其相应的实验数据[8，9，18]也在表中给出。显然，采用LDA近似得到的结构参数比GGA近似小，非常合理。对
C2N2(NH)
而言，计算的晶格常数a，b，和c与实验值吻合非常好，偏离误差为0.5％，0.2％和0.6％，同时其与以前的采用LDA近似的理论值也符合的很好。此外我们计算的键长结果也与实验结果（表1）相一致。对于Si2N2(NH)，计算的晶格常数与实验数据的最大误差约1.2％（LDA）和0.9％（GGA）。
表1.
C2N2(NH)和Si2N2(NH)的平衡晶格参数﹑体弹模量及偏导﹑键长.
Structure
Method
a
b
c
B0
B0
dC/Si-N1
dC/Si-N2
dN2-H
C2N2(NH)
GGA
7.6827
4.5186
4.0546
256
3.50
1.47,1.48
1.46
1.06
LDA
7.5738
4.4429
4.0038
281
3.49
1.45,1.46
1.44
1.05
Exp.a
7.5362
4.4348
4.0298
1.45,1.46
1.43
1.05
Exp.b
7.618
4.483
4.038
258
6.3
Theoreticala
7.5726
4.4425
4.0036
Theoreticalb
271
3.97
Theoreticalb
288
3.94
Si2N2(NH)
GGA
9.2720
5.4640
4.8550
1.74,1.75
1.74
1.03
LDA
9.1371
5.3437
4.7852
1.72,1.73
1.72
1.04
Exp.c
9.1930
5.4096
4.8190
依据文献9，为了进一步的探索C2N2(NH)
的高压行为，我们计算了不同压力点下的C2N2(NH)
的总能。然后将计算得到的E-V数据进行三阶Birch-Murnaghan
物态方程拟合，得到的体弹性模量B0和压力导数B0
在GGA/LDA近似下的结果分别是256/281
GPa
和
3.50/3.49，与表1中的实验数据及其它理论结果相一致。此外，图2
给出了压力作用下C2N2(NH)
单胞体积随压力的变化情况。引人注目的是，实验测量数据完全位于我们采用GGA和LDA近似计算的两条曲线之间。更有趣的是，我们采取GGA和LDA近似的结果与之前理论工作所采取的DFT-
PBE0和DFT-
B3LYP近似的结果在图2中非常一致。因此，上述的计算结果与实验数据以及其他理论结果的一致性极大的验证了我们计算的准确性和可靠性。
图1.
C2N2(NH)的晶体结构示意图.
图2.
压力下C2N2(NH)的体积变化.
3.2
力学性质
3.2.1弹性性质
固体的弹性性质描述了固体经历应变，变形，然后恢复并返回到原来的形状的行为。据我们所知，实验和理论上迄今都没有关于C2N2(NH)
弹性常数的报道。我们希望本文的工作能为今后的实验测量提供一定的参考。弹性常数计算采用应力应变法[20]，对晶体施加一个下的应变，然后充分优化晶体结构和原子的位置，然后有计算得到的弹性常数。表2列出了弹性常数的计算结果。对于稳定的正交晶体，9个独立的Cij弹性常数应满足波恩力学稳定标准[21]，即C11
>
0,C22
>
0,C33
>
0,C44
>
0,C55
>
0,C66
>
0,[C11
+
C22
+
C33
+
2(C12
+
C13
+
C23)]
>
0,(C11
+
C22
–
2C12)
>
0,(C11
+
C33
–
2C13)
>
0,and
(C22
+
C33
–
2C23)
>
0。
显然，计算的关于C2N2(NH)和Si2N2(NH)
的弹性常数Cij满足上述的稳定标准。因此，具有正交结构的C2N2(NH)和Si2N2(NH)
在常压下力学性质稳定。
基于计算的弹性常数Cij，相应的多晶的体弹性模量，剪切模量可以有由Voigt-Reuss-Hill近似确定[17]。此外，杨氏模量EH和泊松比VH可以有下面的公式得到：,(1),(2)
下标H代表Hill近似。计算的体积弹性模量，剪切模量，杨氏模量和泊松比均列于表3。GGA/
LDA近似下，C2N2(NH)
的体弹性模量是276/282
GPa，稍高于超硬材料B6O（270
GPa）[22]的实验数据，表明C2N2(NH)具有超高的不可压缩性。此外，这些计算的体弹模量与直接拟合三阶Birch-Murnaghan
物态方程得到的结果一致，进一步说明我们关于其弹性性质计算的可靠性。作为表征材料硬度的一个更有效的指标，剪切模量是材料抵抗剪切变形的能力。值得注意的是，C2N2(NH)
的剪切模量是非常大，在GGA/
LDA近似下为265/270
GPa。因此，其具有很强的抗剪应变能力。此外，材料中的原子定向成键大小对其硬度具有非常重要的作用，可以通过G/B比率确定。C2N2(NH)
的G/B计算值是0.960，这比金刚石（1.10）稍小，但是和理论预言的超硬材料α-C3N2（0.96）和β-C3N2（1.07）接近[23]。表明在C2N2(NH)
中的C-N键具有非常强的方向性。所有的这些优异的力学性质表明C2N2(NH)
具有极强的抗压缩性，是潜在的超硬材料。与C2N2(NH)
相比，Si2N2(NH)
弹性模量明显偏低，这可能是由于较弱的Si-N成键。
表2.
常压下C2N2(NH)和Si2N2(NH)的弹性常数
structure
method
C11
C22
C33
C12
C23
C13
C44
C55
C66
C2N2(NH)
GGA
597
567
804
89
79
107
335
221
222
LDA
616
576
871
80
71
103
347
221
218
Si2N2(NH)
GGA
285
195
313
38
26
43
125
62
75
LDA
302
255
375
83
47
79
147
86
109
3.2.1弹性各向异性
晶体的弹性各向异性对材料的物理机制，如各向异性的塑性变形，裂纹行为和弹性不稳定可以产生很大的影响。因此，计算材料弹性的各向异性对提高其机械耐久性是非常重要的
[24]。剪切各向异性因子表征不同晶面原子之间的成键的强度。位于
和
方向的{100}剪切面的各向异性因子
.
(3)
位于
和
方向的{010}剪切面的各向异性因子,(4)
位于
和方向的{001}剪切面的各向异性因子
.
(5)
对于各向同性晶体，剪切各向异性因子A1，A2和A3是1.0，而任何值小于或大于1.0则体现了晶体各向异性程度。此外，正交晶体除了所谓的剪切各向异性外，还需考虑线性的体积弹性模量的各向异性。相对于b轴，晶体沿a轴和c轴的体积弹性模量的各向异性由公式
和确定，其中Ba,Bb,和
Bc
是晶体体弹模量沿着各个晶轴方向的分量,定义为
i
=
a,b,和
c.
体弹模量各向异性因子取值为1.0表示其弹性的各向同性，任何偏离1.0表示其具有弹性的各向异性。此外，我们还计算了在多晶材料中剪切模量和体积弹性模量各向异性百分比AG和AB，定义如下：
和，其中下标V和R代表Voigt
和Reuss
近似。AG和AB
的取值为零时对应弹性各向同性，
100％时对应于最大的各向异性。
表3.
常压下C2N2(NH)和Si2N2(NH)的体弹模量﹑剪切模量﹑杨氏模量﹑泊松比及德拜温度.
structure
method
BV
BR
BH
GV
GR
GH
EH
v
GH/
BH
TD
C2N2(NH)
GGA
280
272
276
269
260
265
600
0.137
0.960
1658
LDA
286
277
282
275
265
270
616
0.135
0.957
1662
Theoreticala
277
Si2N2(NH)
GGA
112
107
110
98
90
94
220
0.166
0.855
894
LDA
150
144
147
117
111
114
272
0.192
0.776
978
根据上述提到的关系，我们计算了C2N2(NH)
的各向异性因子，结果如表4所示。很显然，C2N2(NH)
是弹性各向异性的。剪切各向异性结果表明，各向异性沿和方向之间{010}剪切面大于和方向之间的{100}剪切面和沿和方向之间{001}剪切面。此外，有趣的是C2N2(NH)
沿c轴方向的体积弹性模量是大于沿a轴和b轴，这与预测的弹性常数（见表2）是一致的。此外，我们还注意到，剪切模量的百分比AG（1.89％/2.39％，GGA/
LDA）比体积弹性模量百分比AB（1.45％/1.78％，GGA/
LDA）大，表明C2N2(NH)
剪切各向异性强于压缩各向异性。
作为重要的基本参数，Debye温度与许多物理和化学性质密切相关如热传导，动力学特性，融化温度[24]。在低温下其可以从平均声速vm及弹性常数计算,公式如下:,(6)
其中h是普朗克常数，k是波尔兹曼常数，NA是阿伏伽德罗数，n是每个分子式中的原子数，M是每个分子式的质量，ρ为密度。平均声速vm是由,(7)
来决定。vt和vl是多晶材料横向和纵向的弹性波速度,由Navier方程求解[25]，德拜温度的计算值列于表3。计算结果发现，C2N2(NH)
的德拜温度高于Si2N2(NH),表明C2N2(NH)
与Si2N2(NH)
相比具有更高的硬度值。
表4.
常压下C2N2(NH)的弹性各向异性因子.
structure
method
A1
A2
A3
ABa
ABc
Ba
Bb
Bc
AG
AB
C2N2(NH)
GGA
1.1289
0.7288
0.9052
1.060
1.476
760.1
716.8
1058.3
1.89%
1.45%
LDA
1.0835
0.6774
0.8450
1.124
1.651
776.1
690.1
1139.4
2.39%
1.78%
3.3
电子性质
3.3.1
能带结构与电子能态密度
为了深入理解这两种化合物的弹性性质，我们采用GGA近似计算了它们在常压下的能带结构和电子能态密度态密度（DOS），结果如图3所示。C2N2(NH)和Si2N2(NH)
均具有明显的绝缘体特点，能隙分别为?4.36
eV和?5.01
eV（见图3b和3d）。众所周知，密度泛函对带隙计算一般低估至少30％左右。因此，它们精确的带隙值应该分别为5.67
eV和6.51
eV，具有潜在的光学应用价值。
C2N2(NH)
和Si2N2(NH)的原子分立态密度PDOS图如
3a和3c
所示，其中Si2N2(NH)分立态密度的主要特点可以概括如下：（一）在能量较低的范围内，DOS曲线的主要贡献来自于C和N的s电子;（二）费米能级附近的成键态部分主要来自于C-s，C-p,N-s和N-p
电子;（三）费米能级以上主要是反键态部分。能量范围在-23
eV到-16
eV
内，DOS曲线的主要贡献来自于C和N的s电子。能量范围在-15eV到费米能级部分，C-s和C-p轨道电子与N-s和N-p轨道电子之间存在着非常强的杂化作用。表明在CN4的四面体内，C-N共价键作用非常强，这也被后面的电子局域函数结果进一步证实。此外，在DOS导带部分主要是由C-p和N-p电子的贡献。Si2N2(NH)具有与C2N2(NH)相似的电子能态密度特征（参见图3c）。
两种化合物的化学成键机制，可以进一步采用Mulliken原子布局电荷转移来分析，这是判断化合物的内成键一种常有用方法。Mulliken原子布局的绝对值的物理意义不是很大，但是其相对值仍然可以为我们提供一些有用的信息。C2N2(NH)和Si2N2(NH)内电荷转移值在表5中列出，其中在N1和N2是两个不等价的原子如图1所示。显然两种化合物中C和Si的总价电子数分别是是3.57和2.81，小于4;
此外硅原子的价电子总数为6.30
-
6.37个电子，大于C2N2(NH)中C原子的价电子数。这些值表明，这两种化合物的形成也具有离子键的成分。我们也注意到，两种化合物中从H到N的电荷转移为0.4
e。
与C2N2(NH)中从C原子到N原子的电荷转移（0.86
e）相比，Si2N2(NH)中Si原子中更多的电荷转移到N（3.64
e），其中Si所转移的电子主要来自于Si
-2s轨道。因此，Si-
N键的离子性要强于C-N键，这可能也是Si2N2(NH)
具有较低弹性模量的主要原因。因此，这两种化合物中化学键成分是同时具有共价键和离子键。
图2.
常压下C2N2(NH)和Si2N2(NH)的能带结构与电子能态密度.
3.3.2
电子局域函数分析
为了定量确定两种化合物的化学键特征，我们计算了电子局域函数ELF，它表示发现电子的附近的另一个具有相同的自旋电子的概率。
ELF是在正空间中范围从0.0到1.0的等高线图。在ELF是1.0的地区，表示没有找到两个相同的自旋电子的机会。这通常发生在成键分子轨道或孤对原子轨道驻留的地方。ELF
=
0.0对应的地方表示没有电子密度，ELF是0.5表示类似金属中均匀电子气。应该指出，ELF是泡利原理的一种体现，而不是电子密度[28]。如图4为C2N2(NH)
和
Si2N2(NH)
在各自的晶面的ELF。对于C2N2(NH)，我们注意到（001）面中强C-N共价键与C-N共价“点吸引子”的ELF
=
0.9
相对应，同时（100）面中N原子周围的ELF值在H处达到最，体现另一种N和H原子之间的共价相互作用。在Si2N2(NH)中ELF值的分布与C2N2(NH)
相似。然而，在C2N2(NH)
中C和N之间的共价相互作用强于比Si和N之间的相互作用，这是其具有较高体积和剪切模量的主要原因。
表5.
常压下C2N2(NH)和Si2N2(NH)的电荷转移.
structure
s
p
Total
Charge
C2N2(NH)
H
0.6
0.0
0.6
0.4
C
(×2)
1.01
2.56
3.57
0.43
N1
(×2)
1.40
3.95
5.36
-0.36
N2
1.47
4.08
5.55
-0.55
Si2N2(NH)
H
0.6
0.0
0.6
0.4
Si
(×2)
0.79
1.39
2.18
1.82
N1
(×2)
1.66
4.71
6.37
-1.37
N2
1.65
4.66
6.30
-1.30
图4.
常压下C2N2(NH)
和Si2N2(NH)不同晶面的ELF:(a)
C2N2(NH)
(001)面，(b)
C2N2(NH)
(100)面,(c)
Si2N2(NH)
(001)面，(d)
Si2N2(NH)
(100)面.:
4.
结论
基于第一性原理的计算，我们研究C2N2(NH)
和
Si2N2(NH)的力学和电子特性。我们的计算晶格参数与实验数据和以往的理论值吻合非常好。高体积弹性模量和大剪切模量显示C2N2(NH)是一个潜在的超强抗压缩硬质材料。同时，C2N2(NH)化合物显示了不同程度的各向异性。此外，发现C2N2(NH)和Si2N2(NH)是带隙分别为?4.36
eV和5.01
eV的绝缘功能。此外，我们的电子能态密度密度和电子局域函数计算发现CN4的四面体中强的共价C-N键是C2N2(NH)不可压缩性和高硬度的主要原因。我们希望这些计算将有助于今后实验上对这种重要工业材料的性能测试。
References
[1]
R.B.
Kaner,J.
J.
Gilman,S.
H.
Tolbert,Science
308
(2005)
1268.
[2]
V.
L.
Solozhenko,D.
Andrault,G.
Fiquet,M.
Mezouar,and
D.
C.
Rubie,Appl.
Phys.
Lett.
78
(2001)
1385.
[3]
V.
L.
Solozhenko,,O.
O.
Kurakevych,D.
Andrault,Y.
L.
Godec,M.
Mezouar,Phys.
Rev.
Lett.
102
(2009)
015506.
[4]
D.
W.
He,Y.
Zhao,L.
Daemen,J.
Qian,T.
D.
Shen,T.
W.
Zerda,Appl.
Phys.
Lett.
81
(2002)
643.
[5]
A.
Y.
Liu,M.
L.
Cohen,Science
245
(1989)
841.
[6]
E.
Horvath-Bordon,R.
Riedel,A.
Zerr,P.
F.
McMillan,G.
Auffermann,Y.
Prots,W.
Bronger,R.
Kniep,P.
Kroll,Chem.
Soc.
Rev.
35
(2006)
987.
[7]
E.
Kroke,M.
Schwarz,Coord.
Chem.
Rev.
248
(2004)
493.
[8]
E.
Horvath-Bordon,R.
Riedel,P.
F.
McMillan,P.
Kroll,G.
M.
A.v.
Aken,A.
Zerr,P.
Hoppe,O.
Shebanova,l.
McLaren,S.
Lauterbach,E.
Kroke,R.
Boehler,Angew.
Chem.
Int.
Ed.,46
(2007)
1476-1480.
[9]
A.
Salamat,K.
Woodhead,P.
F.
McMillan,R.
Q.
Cabrera,A.
Rahman,D.
Adriaens,F.
Cora,and
J.-P.
Perrillat,Phys.
Rev.
B
80
(2009)
104106
[10]
P.
Hohenberg,W.
Kohn,Phys.
Rev.
136
(1964)
B864.
[11]
W.
Kohn,L.
J.
Sham,Phys.
Rev.
140
(1965)
A1133.
[12]
D.
M.
Ceperley,B.
J.
Alder,Phys.
Rev.
Lett.
45
(1980)
566.
[13]
J.
P.
Perdew,K.
Burke,M.
Ernzerhof,Phys.
Rev.
Lett.
77
(1996)
3865.
[14]
G.
Kresse,D.
Joubert,Phys.
Rev.
B
59
(1999)
1758.
[15]
P.
E.
Blchl,Phys.
Rev.
B,50
(1994)
17953.
[16]
J.
D.
Pack,H.
J.
Monkhorst,Phys.
Rev.
B
16
(1977)
1748.
[17]
R.
Hill,Proc.
Phys.
Soc.
London
65
(1952)
349.
[18]
D.
Peters,H.
Jacobs,J.
Less-Common
Met.
146
(1989)
241.
[19]
F.
Birch,Phys.
Rev.
71
(1947)
809.
[20]
V.
Milman,M.
C.
Warren,J.
Phys.:
Condens.
Matter
13
(2001)
241.
[21]
M.
Born,Proc.
Cambridge
Philos.
Soc.
36
(1940)
160.
[22]
D.
He,S.
R.
Shieh,T.
S.
Duffy,Phys.
Rev.
B
70
(2004)
184121.
[23]
F.
Tian,J.
Wang,Z.
He,Y.
Ma,L.
Wang,T.
Cui,C.
Chen,B.
Liu,G.
Zou,Phys.
Rev.
B
78
(2008)
235431.
[24]
P.
Ravindran,L.
Fast,P.
A.
Korzhavyi,B.
Johansson,J.
Wills,O.
Eriksson,J.
Appl.
Phys.
84
(1998)
4891.
[25]
E.
Schreiber,O.
L.
Anderson,N.
Soga,Elastic
Constants
and
their
Measurements
McGraw-Hill,New
York,1973.
[26]
A.
D.
Becke,K.
E.
Edgecombe,J.
Chem.
Phys.
92
(1990)
5397.
[27]
A.
Savin,B.
Silvi,Nature
371
(1994)
683.
[28]
N.
Ooi,A.
Rairkar,L.
Lindsley,J.
B.
Adams,J.
Phys.:
Condens.
Matter,18
(2006)
97.