超硬材料在基础科学和工业应用中具有重要的意义，其设计与合成受到人们的广泛关注^[1-2]。由轻元素（B、C、N和O）形成的强共价键化合物是潜在超硬材料的典型来源，近年来取得了丰硕的研究成果^[3]，并成功合成出B-C^[4-5]、B-O^[6-7]、B-C-N^[8-10]等超硬材料。超硬材料的理论设计和预测也取得了巨大成就，提出了由一元单质^[11]、二元^[12]及三元^[13]化合物组成的新型超硬材料，标志性工作是能够匹敌金刚石硬度的超硬材料$\,\beta$-C₃N₄的预测^[14]。

自1997年采用高压技术通过沃克型多砧装置合成B-C-O化合物（B₆C_1.1O_0.33和B₆C_1.28O_0.31）以来^[15]，三元化合物B-C-O进入了人们的研究视野。2001年，Bolotina等^[16]在5.5 GPa、1400 K的条件下使1∶1混合的B₄C和B₂O₃发生反应，成功制备出B(C,O)_0.1555。实验合成的三元B-C-O化合物证实了其存在的物质基础，自此B-C-O体系化合物作为超硬材料候选物质来源引起了广泛的关注。

以B₂CO（具有与金刚石等电子体的最简配比B-C-O化合物）为例，Li等^[17]开展了对B-C-O化合物的开创性理论研究，预测出两种四方晶系超硬结构，分别为tP4-B₂CO和tI16-B₂CO。实际上，在B-C-O三元体系中存在很多金刚石等电子体化合物，如B₂C_xO（x = 2, 3, 4, ···），随着C含量的增加，B-C-O三元化合物结构中存在C的堆积，进而形成更短更强的sp³杂化C―C共价键。因此，Li等^[17]认为B₂C_xO（x = 2, 3, 4, ···）化合物中极有可能存在比B₂CO更硬的结构。随后，Zhang等^[18]针对B₂C_xO化合物开展了研究，通过系统的粒子群优化算法结构搜索，提出了3种四方晶系的超硬B₂C_xO（x ≥ 2）相，分别为B₂C₂O、B₂C₃O和B₂C₅O，并发现随着碳含量的增加，sp³杂化C―C共价键有利于弹性模量和硬度的提高^[18]，与Li等的观点吻合。2017年，Liu等^[19]从理论上预测出一种类蓝丝黛尔石lonsdalite的斜方晶系B₂CO结构（oP8-B₂CO），这是首次提出非四方结构的B₂CO化合物。基于键阻模型的微观硬度理论分析，预测oP8-B₂CO具有高达47.70 GPa的硬度。随后Qiao等^[20]系统研究了高压下oP8-、tI16-和tP4-B₂CO化合物的结构和电子性质，基于杂化泛函的研究表明，oP8-、tI16-和tP4-B₂CO都是宽的间接带隙半导体材料，其带隙均随压力的增大而增大。

受到近年来所预测的B-C-O结构特征的启发（如tP4-B₂CO和tI16-B₂CO均具有类金刚石结构^[17]，oP8-B₂CO具有类蓝丝黛尔石结构^[19]，B₂C₂O、B₂C₃O和B₂C₅O均具有类金刚石结构），可由不同的tP4-B₂CO超胞通过部分B和O被C取代产生^[18]。基于结构对性质的影响，Liu等^[21]以具有超硬特性的碳同素异形体为模型，通过多步骤手动构建结合第一性原理研究，预测出具有超硬性能的两种正交晶系B₂CO相oP16-B₂CO和oC16-B₂CO，它们分别与Cco-C8和Bct-C4具有相似的结构。

对于B₂CO三元化合物，虽然已经成功预测出一些结构，但是考虑到轻元素的灵活性和报道的宽稳定压力范围，寻找新型未知的热力学稳态相仍然值得关注^[22]。为此，Yan等^[22]采用粒子群优化晶体结构分析软件CALYPSO，探索了B₂CO在环境条件下的热力学稳定结构，发现了一种sp²-sp³杂化共存的新型高硬度正交晶系oI16-B₂CO，比所有先前提出的B₂CO结构在能量上更稳定。受到Yan等研究成果的鼓舞，Chen等^[23]提出了一种四方晶系sp²-sp³ B―C键和B―O键共存的亚稳态B₂CO化合物tP16-B₂CO。研究表明，随着压力的增加，tP16-B₂CO的带隙先增大后减小。通过计算tP16-B₂CO的硬度、应力-应变关系等力学性能，发现其为硬质材料，硬度为23.19 GPa，且最大应力沿[001]方向的各向异性远高于沿[100]方向^[23]。

此外，通过长期探索，非金刚石等电子体的B-C-O化合物研究取得了较大进展。2016年，Wang等^[24]通过从头算变组成演化模拟，探究了0～50 GPa压力范围内的B-C-O体系，发现了一种新型四方相B₄CO₄。B₄CO₄在高于23 GPa的压力下是热力学稳定的，但在室温条件下保持亚稳态，计算显示其硬度为38～41 GPa，表明B₄CO₄具有超硬的潜力。Nuruzzaman等^[25]研究发现，B₄CO₄具有较大的体积模量、剪切模量、杨氏模量、维氏硬度以及较小的泊松比，表明B₄CO₄可能是超硬材料的候选材料。鉴于直观上B₄CO₄结构中心存在较大的空心空间，B₄CO₄是否是一种超硬材料还需要进一步澄清^[26]。Zheng等^[26]采用第一性原理计算方法，详细推导了B₄CO₄的杨氏模量和剪切模量的方向依赖关系，结果表明B₄CO₄具有高度的各向异性。基于应变-应力法，计算了B₄CO₄沿主晶体方向的理想拉伸强度和剪切强度，研究表明，B₄CO₄的塑性变形主要是沿(001)[100]滑移体系的剪切模式，这与B―O键的断裂有关。最弱的理想剪切强度为27.5 GPa，证明了B₄CO₄的高硬度属性。

Wang等^[24]还报道了B-C-O体系中两个低焓亚稳相B₆C₂O₅和B₂CO₂。基于计算的声子色散和弹性常数，B₆C₂O₅和B₂CO₂在环境条件下都满足动力学和弹性力学稳定性。通过计算证实了B₆C₂O₅和B₂CO₂具有高硬属性，硬度值分别为36和40 GPa^[24]。Qiao等^[27]计算了B₂CO₂和B₆C₂O₅的能带结构，发现B₂CO₂和B₆C₂O₅是宽带隙半导体材料，其间接能隙分别为5.66和5.24 eV。弹性各向异性结果表明，B₆C₂O₅的剪切模量具有较大的各向异性，而B₂CO₂的杨氏模量具有较大的各向异性。此外，Lyakhov-Oganov模型的硬度计算表明B₂CO₂是一种潜在的超硬材料。随后，Liu等^[28-29]用第一性原理计算预测了一种四方晶系B-C-O化合物tI12-B₆C₄O₂。tI12-B₆C₄O₂具有较大的力学模量和较高的硬度，是一种硬度为21.9 GPa的高硬材料。各向异性研究表明，沿特定方向的杨氏模量存在E_[110] < E_[100] < E_[001] < E_[111] < E_[011]的关系。

纵观三元B-C-O化合物体系研究发现，虽然取得了丰富的成果，但是缺乏组分对性能调控机理的研究，尤其对于非金刚石等电子体系列的B-C-O化合物，其组分多样，相关的结构报道不多，研究不够系统，尚未揭示组分对其性能的作用规律。为此，本研究开展了非金刚石等电子体系列的B-C-O化合物结构预测，提出一种新型B₄C₆O₄化合物，并基于B₄C₆O₄、B₂CO₂和B₄CO₄，系统研究BC_xO化合物的性质与组分之间的关系。

1.   理论研究细节

通过基于粒子群优化算法的晶体结构预测软件CALYPSO^[30-32]，在常压下通过变组分方式生成非金刚石等电子体B-C-O化合物的潜在结构，然后在CASTEP模块^[33]中采用广义梯度近似的Perdewey-Burke-Ernzerhof（PBE）泛函^[34]进行全弛豫结构优化及性质研究。其中超软赝势^[35]的截止能量设为380 eV，并采用Koelling-Harmon方式进行相对论处理。通过Monkhorst-Pack网格划分产生K点，其中划分密度为2$\text{π}$ × 0.04 Å⁻¹。在几何优化过程中，电子步自洽精度为5 × 10⁻⁴ meV/atom，离子步收敛精度为每原子能量变化小于5 × 10⁻³ meV，指定最大应力的收敛阈值为0.02 GPa，指定最大位移小于5 × 10⁻⁴ Å。计算电子能带结构和态密度时涉及到的能带能量容差低至1 × 10⁻⁵ eV。声子计算采用有限位移法^[36]，计算弹性常数时每个应变产生9个模式，最大应变幅值为0.009，收敛精度为每原子能量变化不超过1 × 10⁻³ meV，最大力为0.002 eV/Å，最大位移不超过1 × 10⁻⁴ Å。

此外，给出了计算过程中涉及的模型类型及其布里渊区高对称点信息。B₄C₆O₄单胞：G(0, 0, 0)→Z(0, 0, 0.5)→T(−0.5, 0, 0.5)→Y(−0.5, 0, 0)→S(−0.5, 0.5, 0)→X(0, 0.5, 0)→U(0, 0.5, 0.5)→R(−0.5, 0.5, 0.5)；B₄C₆O₄原胞：G(0, 0, 0)→Z(0, 0, 0.5)→T(0.5, 0.5, 0.5)→Y(0.5, 0.5, 0)→G(0, 0, 0)→S(0, 0.5, 0)→R(0, 0.5, 0.5)→Z(0, 0, 0.5)；B₄CO₄原胞：Z(−0.5, 0.5, 0.5)→G(0, 0, 0)→X(0, 0, 0.5)→P(0.25, 0.25, 0.25)→N(0, 0.5, 0)→G(0, 0, 0)；B₂CO₂原胞：L(−0.5, 0, 0.5)→M(−0.5, −0.5, 0.5)→A(−0.5, 0, 0)→G(0, 0, 0)→Z(0, −0.5, 0.5)→V(0, 0, 0.5)。

2.   结果与讨论

2.1   优化结构

通过CALYPSO软件预测一种新型非金刚石等电子体的B-C-O化合物，常压下优化得到的结构模型如图1所示。该B-C-O化合物具有(0.5, 0.5, 0)中心对称的正交晶系结构，空间群Cmm2，晶胞参数a、b和c分别为4.587、5.441和3.812 Å。晶胞中B、C、O的原子数分别为4、6、4，这里记为B₄C₆O₄。常压下B₄C₆O₄的密度为3.130 g/cm³，其中C原子有C1和C2两类，结构中各类原子占位情况如表1所示。

图  1  常压下B₄C₆O₄结构（2 × 2 × 2）的超胞模型

Figure  1.  Structure graphs for B₄C₆O₄ (2 × 2 × 2) supercell at ambient pressure

下载: 全尺寸图片 幻灯片

表  1  常压下B₄C₆O₄的原子坐标

Table  1.  Atomic Wyckoff positions of B₄C₆O₄ at ambient pressure

Atom	Wyckoff site	x	y	z
B	4d	0.299	0.500	0.267
C1	2b	0	0.500	0.991
C2	4c	0.250	0.250	0.491
O	4e	0	0.279	0.746

下载: 导出CSV 

| 显示表格

观察发现：结构中存在B原子之间的成键，B―B键长为1.842 Å；O原子与C1原子形成的O―C1键的键长为1.522 Å，该类键与x轴呈空间垂直关系；O原子与C2原子形成的O―C2键的键长为1.511 Å，该类键沿x轴相互串联。C2与B原子形成的键具有两类键长：一类在空间上呈现沿y轴无限连接，键长为1.622 Å；另一类键则在空间上垂直于y轴，键长1.728 Å，并与B―B键以等键数相间的方式连接，沿x轴延伸。整个结构中C原子以4配位方式与周边2个O原子和2个B原子成键，B原子与周边2个C2原子、1个C1原子和1个近邻B原子连接成键，而O原子只以3配位形式与周边3个C原子成键。

2.2   稳定性分析

B₄C₆O₄的弹性稳定性通过计算其独立的弹性常数C_ij，并经波恩弹性稳定性判据验证。B₄C₆O₄的晶体结构属于正交晶系，其独立的弹性常数有9个，波恩弹性稳定性判据^[37]为

计算得到常压和100 GPa高压下B₄C₆O₄的弹性力学常数见表2。经分析发现，独立的弹性常数C_ij满足稳定性判据，证明B₄C₆O₄在常压和100 GPa高压条件下满足弹性稳定性。

表  2  常压和100 GPa下B₄C₆O₄的弹性常数

Table  2.  Elastic parameters of B₄C₆O₄ at ambient pressure and 100 GPa 　GPa

Pressure	C11	C12	C13	C22	C23
Ambient	693.28	15.78	77.13	542.70	167.35
100	1394.79	146.17	335.20	974.72	497.69
Pressure	C33	C44	C55	C66
Ambient	516.07	255.45	197.31	222.21
100	1000.58	517.68	363.80	292.39

下载: 导出CSV 

| 显示表格

虚频会导致晶格畸变，预示晶体的动力学不稳定。图2为常压和100 GPa高压下B₄C₆O₄单胞整个布里渊区的声子色散谱和相关声子态密度（Density of states，DOS）。在其单胞的整个布里渊区没有负声子模式，表明常压和100 GPa高压下B₄C₆O₄动力学稳定。由此可知，预测的B₄C₆O₄同时满足弹性稳定性和动力学稳定性。

图  2  常压和100 GPa下B₄C₆O₄的声子散射谱和声子态密度

Figure  2.  Phonon dispersion spectrum and phonon density of states of B₄C₆O₄ at ambient pressure and 100 GPa

下载: 全尺寸图片 幻灯片

随后基于实际环境开展对B₄C₆O₄的分子动力学模拟：在恒定273 K和标准大气压的系综下进行总时间为3 ps、步长6 fs的模拟研究。模拟结果如图3所示。整个模拟过程中，结构没有明显变化。此外，为了阐明分子动力学模拟过程中结构的微调，这里一并选取了B―B键的键长以及C2与C2、O与O、B与B的原子间距（如图3中红色虚线所示），将B₄C₆O₄的密度$\,\rho$作为观察对象，发现在整个分子动力学模拟过程中B₄C₆O₄结构的局部化学键和原子间距（见表3）存在一定程度的拉伸或压缩，但总体变化微小，模拟初态和末态的体系密度几乎无变化。综上所述，即便在考虑温度和压力的真实环境下，B₄C₆O₄依旧保持其结构完整性和稳定性。

图  3  B₄C₆O₄在恒温273 K和常压下的分子动力学模拟

Figure  3.  Molecular dynamics simulation for B₄C₆O₄ at 273 K and ambient pressure

下载: 全尺寸图片 幻灯片

表  3  恒温常压条件下B₄C₆O₄在分子动力学过程中的动态结构信息

Table  3.  Structural information of B₄C₆O₄ during molecular dynamics under constant temperature and ambient pressure

Time/ps	$\,\rho$/(g·cm−3)	B―B bond length/Å	Atom distance/Å
			C2―C2	O―O	B―B
0	3.130	1.842	2.720	3.038	2.757
0.6	3.175	1.776	2.677	3.013	2.775
1.2	3.137	1.869	2.697	3.026	2.722
1.8	3.066	1.827	2.741	3.010	2.810
2.4	3.123	1.846	2.711	3.063	2.806
3.0	3.128	1.884	2.738	3.033	2.798

下载: 导出CSV 

| 显示表格

2.3   常压性质分析

基于密度泛函理论，B₄C₆O₄在常压下的物理性质如电学、力学等得以研究。鉴于传统PBE泛函存在低估半导体带隙的不足^[38]，为了获取B₄C₆O₄精确的电学性质，基于杂化泛函HSE06^[39]，分别在B₄C₆O₄单胞和原胞模型上进行电学性质研究，结果如图4所示。虽然B₄C₆O₄的单胞和原胞具有不同的倒易空间和布里渊区及高对称点路径，但是研究发现基于原胞和单胞计算所得的价带最高点和导带最低点均位于同一高对称点K上，最高价带（VBM）和最低导带（CBM）之间存在带隙，分别为2.248和2.247 eV。不同结构模型的研究均表明，B₄C₆O₄具有带隙宽度约2.25 eV的直接带隙半导体属性。

图  4  基于HSE06计算得到的常压下B₄C₆O₄的电子能带结构和态密度：(a)原胞，(b)单胞（水平红线、暗青色曲线和金色曲线分别代表费米能级、VBM和CBM）

Figure  4.  Calculated electronic band structures and density of states of B₄C₆O₄ phases via HSE06 with primitive cell (a) and unit cell (b) at ambient pressure (The horizontal red line, dark cyan curve and gold curve represent the Fermi energy level, VBM and CBM, respectively.)

下载: 全尺寸图片 幻灯片

纵观研究提出的B₄C₆O₄以及Wang等^[24]报道的B₂CO₂和B₄CO₄，三者均具有B、O含量相等的特点，即化学式可统一表达为BC_xO（x = 3/2, 1/2, 1/4）。三者的结构中O呈现3配位成键，B和C均呈现4配位成键的诸多结构共性，为此通过理论研究探究三者的电学性质随组分的变化关系。鉴于杂化泛函HSE06对计算设备的硬件要求严苛、计算时间冗长等不足，而PBE泛函可以快速得到电学性质以及外界条件对电学性质的调控趋势，为此这里采用基于原胞模型和PBE交换关联泛函开展对BC_xO（x = 3/2, 1/2, 1/4）的电学性质计算以及压力对电学性质的调控规律研究。

首先研究常压下BC_xO（x = 3/2, 1/2, 1/4）的电学性质，如图5所示。BC_xO（x = 3/2, 1/2, 1/4）均属于半导体材料，其中B₄C₆O₄为直接带隙，而B₂CO₂的价带最高点位于高对称点V，导带最低点位于L与M之间的非高对称点，B₄CO₄的价带最高点位于高对称点G，导带最低点位于X，表明B₂CO₂和B₄CO₄均具有间接带隙。对比研究发现，随着C含量的降低，体系的带隙增大。

图  5  基于PBE计算常压下BC_xO原胞结构的电子能带结构和态密度（水平红线、暗青色曲线和金色曲线分别代表费米能级E_F、VBM和CBM）

Figure  5.  Calculated electronic band structures and density of states of three BC_xO phases with primitive cell via PBE at ambient pressure (The horizontal red line, dark cyan curve and gold curve represent the Fermi energy level E_F, VBM and CBM, respectively.)

下载: 全尺寸图片 幻灯片

一般而言，压力作用下原子间距减小，体积减小，致密度增大。对于3种BC_xO（x = 3/2, 1/2, 1/4）化合物，归一化体积随压力的变化如图6(a)所示。随着x的降低，归一化体积减小；B₄C₆O₄、B₂CO₂和B₄CO₄三者的归一化体积均随压力的增大而持续降低。在所研究的压力范围内，B₄C₆O₄、B₂CO₂、B₄CO₄的归一化体积分别减小了5.221、3.311、3.486 ų，体系分别被压缩了21.955%、20.197%和22.764%。

图  6  基于PBE泛函和原胞模型计算的BC_xO（x = 3/2, 1/2, 1/4）的归一化体积(a)和带隙(b)与压力的关系

Figure  6.  Relationship of the normalized volume (a) and the calculated electronic band gaps (b) of BC_xO (x = 3/2, 1/2, 1/4) with pressure based on PBE funtional and primitive cell

下载: 全尺寸图片 幻灯片

凝聚态材料的电学性质会受压力的影响^[40]。通常而言，高压下原子间距变小，电子重叠度增大，电子不局域于单个原子或键，即形成离域电子，导致带隙减小乃至金属化，典型代表就是金属氢^[41]。然而高压下，电子出现局域化的同时也可能伴随成键态和反键态的形成，导致价带能量更低而导带能量更高，即出现带隙变大现象，典型代表就是钠绝缘化^[42]。很多时候两种作用相互竞争，共同影响材料的电学性质。

为此，研究了BC_xO（x = 3/2, 1/2, 1/4）带隙随压力的变化。基于原胞和PEB算法的研究结果如图6(b)所示。可见，B₂CO₂的带隙有较明显的受压变化趋势，而B₄CO₄和B₄C₆O₄的带隙受压力的影响较小。在所研究的压力范围内，B₂CO₂的带隙随压力的升高呈现持续下降现象，100 GPa内带隙减小425 meV，即在高压的作用下B₂CO₂的晶格参数变小，布里渊区变大，能带展宽，进而导致带隙减小。B₄C₆O₄的带隙也呈现随压力的升高而降低的变化规律，在低于100 GPa的压力范围内带隙降低了168 meV，降幅明显小于B₂CO₂。对于B₄CO₄，带隙在常压到40 GPa区间随压力的升高而持续增大，增幅达62 meV，随后保持带隙最大值5.495 eV直至50 GPa，50 GPa之后带隙随压力的升高而缓慢下降，100 GPa时带隙为5.470 eV。纵观整个升压过程，B₄CO₄的带隙变化微小，即两种作用相互竞争，近乎相持。

通过模拟特定方向的结构应变获得凝聚态材料的内在应力是理解结构变形和强度的常用方法。材料的理想强度（一个完美晶体在弹性力学上由稳定变为不稳定时的临界应力）代表了材料强度的上限。通过模拟拉伸过程的应变-应力关系和键断裂过程，可以从微观原子层面深入理解结构变形乃至失效的理论机制。为此，基于密度泛函理论模拟研究了BC_xO（x = 3/2, 1/2, 1/4）在常见晶向上的拉伸过程，所得应力-应变关系如图7所示。

图  7  模拟得到的3种BC_xO沿特定方向的拉伸应力-应变关系

Figure  7.  Simulated tensile stress-strain relationship along the specific direction for three BC_xO phases

下载: 全尺寸图片 幻灯片

B₄C₆O₄为正交晶系结构，3个常见方向的应力-应变关系如图7(a)所示。[100]方向上应力先随应变的增大而增大，在应变为0.1125时达到峰值42.1 GPa，随后应力略微减小，呈现快速失效特征。[010]和[001]方向上，在低于0.05的小应变区域内应力相近，应变继续增大时应力差异显现，其中[001]方向的应力在应变为0.10时达到最大值40.4 GPa，此后继续增大拉伸应力会导致应变超过极限变形量而造成结构失效，而[010]方向上的应力随应变的增大持续增大到49.5 GPa（应变0.1625），随后应力随应变的增大而略微下降，结构快速失效。

B₂CO₂为单斜晶系结构，常见方向的拉伸模拟研究结果如图7(b)所示。[100]和[001]方向的应力-应变曲线具有相似性，二者的应力先随着应变的增大而增大，均在0.056处达到最大值，分别为21.0和23.5 GPa，随后应力随应变的增大而减小，最大应变均为0.069，应变继续增大会导致结构失效。在[010]方向，低应变区域内应力随着应变的增大而快速增大，高应变区域内应力随应变的增大而缓慢增大，最大应力为66.7 GPa，此时应变也达到极限值0.225，随后若应力继续增大，则结构失效。

B₄CO₄为四方晶系，同样选取了[100]、[110]和[001] 3个方向，其应力-应变曲线如图7(c)所示。在[100]和[001]方向上应力随应变的增大而增大，其中：[100]方向的最大应力为42.3 GPa，此时应变达到极限0.138；[001]方向的最大应力31.0 GPa出现在极限应变0.100处，继续增大应力会导致应变超过极限应变，进而使结构失效塌陷；在[110]方向，应力先随应变的增大而显著增大，当应变为0.100时应力达到最大值55.9 GPa，随后应变继续增大，结构中部分B―O键断裂（图7(c)中S₁到S₂），应力突降至11.2 GPa，此时结构未完全失效，当应变继续增大时应力持续增大，并在0.2875应变处应力达到38.0 GPa，当应变继续增大，结构完全失效。二次变形的存在导致结构在[110]方向的应变高于[100]和[001]方向。

基于模拟的拉伸过程，在结构未失效下开展了3种BC_xO（x = 3/2, 1/2, 1/4）的带隙-应变关系研究，结果见图8。如图8(a)所示，在[001]方向的拉伸过程中，B₄C₆O₄的带隙随应变的增大单调增大；在[100]方向的拉伸过程中，带隙随应变的增大持续下降；而[010]方向的拉伸则引起带隙先缓慢上升，在应变0.063处达到最大值1.068 eV，随后带隙先均匀降低，在临近失效阶段快速降低。如图8(b)所示，B₂CO₂在3个方向的拉伸应力作用下，带隙随结构的改变先快速下降，在达到一定的应变后，[001]方向的拉伸应变引起带隙持续下降，而[100]和[010]方向的应变则均先使体系带隙缓慢上升，达到峰值后，随应变的增大而降低，其中[100]方向的带隙峰值4.155 eV出现在应变0.069处，[010]方向的带隙峰值4.173 eV出现在应变0.050处。B₄CO₄的带隙随[100]、[110]和[001] 3个方向的拉伸应变的变化如图8(c)所示。[001]方向的应变使带隙先升后降，峰值6.218 eV出现在应变0.075处；在[100]方向上带隙随应变先微弱上升后持续降低；而[110]方向上带隙的变化与应力的变化一样，出现两个阶段，即应变在0.100以内带隙快速下降，应变由0.100增大到0.113的过程中带隙出现突降，随后带隙缓慢下降。

图  8  模拟得到的3种BC_xO在沿特定方向拉伸过程中的电阻-应变关系

Figure  8.  Simulated gap-strain relationship along the specific direction for three BC_xO phases during the tensile process

下载: 全尺寸图片 幻灯片

基于密度泛函理论，系统研究了3种BC_xO（x = 3/2, 1/2, 1/4）的力学性质。首先，基于PBE泛函计算了BC_xO（x = 3/2, 1/2, 1/4）的弹性常数C_ij，通过C_ij计算Voigt和Reuss形式的体积模量B和剪切模量G^[43]，并经改进型硬度计算公式计算维氏硬度H_V^[44]。改进型硬度计算公式为

3种BC_xO（x = 3/2, 1/2, 1/4）的体积模量B和剪切模量G与压力的关系如图9(a)和图9(b)所示，其中Voigt和Reuss形式数值分别代表模量的上限和下限。对于体积模量而言，B₄CO₄的体积模量的Voigt和Reuss形式数值非常接近，相差不到0.8 GPa，B₂CO₂和B₄C₆O₄的体积模量的Voigt和Reuss形式数值差值均小于8 GPa。结合图9(a)可知，3种BC_xO（x = 3/2, 1/2, 1/4）均具有较大的体积模量。此外，3种BC_xO（x = 3/2, 1/2, 1/4）的体积模量均呈现与压力正相关的线性关系。如图9(b)所示，对于剪切模量而言，3种BC_xO（x = 3/2, 1/2, 1/4）的Voigt和Reuss形式数值均存在较大差值，且该差值随压力的升高而增大。根据Voigt和Reuss形式的体积模量和剪切模量，计算得到两种形式的维氏硬度H_V及其与压力的关系，如图9(c)所示。常压下3种BC_xO（x = 3/2, 1/2, 1/4）均具有较高的硬度，尤其是新预测的B₄C₆O₄具有40.3～37.3 GPa的常压硬度，即B₄C₆O₄是潜在的常压准“超硬”材料；B₂CO₂和B₄CO₄也分别具有39.2～36.7 GPa、38.9～33.7 GPa的常压硬度，表明二者具有常压高硬性质。除B₄C₆O₄的硬度在常压到10 GPa区间呈现小幅上升外，3种BC_xO（x = 3/2, 1/2, 1/4）的硬度随压力的升高呈现下降趋势，在100 GPa下依旧保持高硬特性。

图  9  3种BC_xO（x = 3/2, 1/2, 1/4）的体积模量(a)、剪切模量(b)、硬度(c)与压力的关系

Figure  9.  Relationship of bulk modulus (a), shear modulus (b), and hardness (c) with pressure for BC_xO (x = 3/2, 1/2, 1/4)

下载: 全尺寸图片 幻灯片

晶体的声子振动频率受其结构中化学键强度的影响，即最大声子振动频率在一定程度上反映了化学键的强弱。基于不同压力下BC_xO（x = 3/2, 1/2, 1/4）的结构，计算得到其声子色散曲线和声子态密度，分析发现在常压到100 GPa的压力范围内BC_xO（x = 3/2, 1/2, 1/4）均满足动力学稳定性。基于声子色散曲线，详细探讨了3种BC_xO（x = 3/2, 1/2, 1/4）的最大振动频率随压力变化的关系。如图10(a)所示，常压下BC_xO的最高声子振动频率均高于30 THz，且呈B₄CO₄ > B₂CO₂ > B₄C₆O₄的关系。此外，随着压力的增大，原子间距减小，成键原子的相互作用力增强，即化学键增强，导致体系的最高声子振动频率随压力同步递增。

图  10  BC_xO（x = 3/2, 1/2, 1/4）的最大声子振动频率(a)和零点振动能(b)与压力的关系

Figure  10.  Relationship of maximum vibrational frenquency (a), and zero-point vibration energy (b) with pressure for BC_xO (x = 3/2, 1/2, 1/4)

下载: 全尺寸图片 幻灯片

零点振动能E_ZP是量子在绝对温度0 K下的振动能量。零点振动的振幅随温度的升高而增大，原子越轻，零点振动越明显。根据声子散射谱和声子态密度，可以计算出零点振动能E_ZP

式中：ћ为普朗克常数，F(ω)为声子态密度。基于相关的声子散射谱，进一步研究了归一化零点振动能与压力的关系。

如图10(b)所示，常压下B₄C₆O₄的归一化零点振动能最大，其次是B₂CO₂，B₄CO₄最小，与三者的归一化相对质量大小密切相关。此外，随着压力的增大，3种BC_xO的零点振动能均增大，且增幅也随x的增大而递减，E_ZP的增幅分别为134.8、87.8和80.1 meV。以上说明C含量对BC_xO的最大振动频率和零点振动能具有显著的影响。

3.   结　论

通过粒子群优化算法生成候选结构，研究了具有正交晶系结构的新型B-C-O化合物B₄C₆O₄。结合第一性原理的声子散射谱和独立弹性常数，证明了B₄C₆O₄结构的常压和高压稳定性。B₄C₆O₄具有带隙宽度约2.25 eV的直接带隙半导体属性。研究同属BC_xO系列且结构具有相似性的B₄C₆O₄、B₂CO₂和B₄CO₄发现，C含量的降低导致带隙增大，三者的分子式体积随C含量的降低而降低，且100 GPa的高压对三者体积均形成高达20%的压缩。高压作用导致B₂CO₂和B₄C₆O₄的带隙持续降低，而B₄CO₄的带隙先升后降，整体变化微弱。应力-应变模拟结果表明，B₄C₆O₄的最大应力49.5 GPa出现在[010]方向的0.1625应变处，[100]和[001]方向能承受的最大应力及应变均低于[010]方向；同样，B₂CO₂在[010]方向能承受的最大拉伸应力66.7 GPa和应变0.225也远大于[100]和[001]方向的应力和应变。B₄CO₄在最大应力和最大应变的[110]方向上出现二次拉伸特质，应变为0.100时应力达到最大值55.9 GPa，增大应变至0.113时结构中出现部分断键但结构未塌陷，继续增大应变到超过0.288时结构的完整性被破坏。应力引起的应变会影响3种BC_xO（x = 3/2, 1/2, 1/4）的带隙。3种BC_xO化合物的力学性能研究表明，它们都具有高弹性模量和高硬度，尤其是新预测的B₄C₆O₄具有40.3～37.3 GPa的常压硬度，是潜在的准“超硬”材料。三者在100 GPa下依旧保持高硬特性。声子散射谱研究表明，常压下BC_xO的最高声子振动频率均高于30 THz，且B₄CO₄的最高声子振动频率最大，B₂CO₂次之，B₄C₆O₄最小，高压作用会引起体系键能的持续增强。