结合能的新势函数对高压稠密惰性元素压缩特性的影响

郑兴荣

引用本文:
Citation:

结合能的新势函数对高压稠密惰性元素压缩特性的影响

    作者简介: 郑兴荣(1986-),男,硕士,讲师,主要从事凝聚态理论物理与材料计算、团簇结构研究. E-mail:zhengxingrong2006@163.com;
  • 中图分类号: O521.2; O561.1

A Novel Expression of Cohesive Energy Contributions to the Highly Compressed Characteristic for Rare-Gas Solids

  • CLC number: O521.2; O561.1

  • 摘要: 基于量子理论和原子团簇理论,运用多体展开方法和第一性原理的从头算方法,提出了一种计算稠密惰性元素(氦、氖、氩和氪)原子结合能的新势函数,运用新公式研究了结合能对高压稠密惰性元素高压压缩特性的影响。此公式引入了一个物理参量$\beta $(其值为0.5),使得势函数的表达形式更加简单、准确。对比结果表明,结合能的新势函数能够准确地描述多体相互作用对结合能的贡献,且平均相对误差在5%以内。结合能的新势函数对压缩特性的影响在当前实验压强范围内(氦60 GPa、氖238 GPa、氩114 GPa、氪128 GPa)做出了令人满意的描述,且与实验值及理论计算结果基本完全吻合,平均相对误差在3%以内。最后,以固氩的压强数据为例,验证了势函数的准确性。该势函数不仅适用于更宽密度和更高压强范围,而且对所有惰性元素原子各种状态的结合能、高压压缩特性、定容比热容、熔化曲线和弹性模量的研究具有重要的指导意义。
  • 图 1  U2(M)、Vn(M)和E(M)三者的关系(曲线OE表示函数f(x),虚线OG代表U2(M)x

    Figure 1.  The correlations among U2(M), Vn(M), and E(M) (The curve line OE represents function f(x), and the dash line OG describes U2(M)x.)

    图 2  稠密惰性元素原子结合能的比较

    Figure 2.  Comparison of the cohesive energy for rare-gas solids

    图 3  稠密惰性元素原子的高压压缩特性的比较

    Figure 3.  Comparison of the compressibility of rare-gas solids

    表 1  固氩的压强分量

    Table 1.  The pressure components of solid argon

    RV/(cm3·mol–1)P/GPaError/%
    Exp.[24]Ab initio[11]Eq.(10)
    2.405.887237.61248.95247.043.97
    2.456.262194.11204.44200.723.41
    2.506.653158.51167.59163.253.00
    2.557.061129.38137.12132.812.65
    2.607.484105.53111.96107.982.32
    2.657.92486.0291.2387.701.95
    2.708.38170.0674.1871.121.51
    2.758.85657.0060.2057.571.00
    2.809.34846.3448.7546.490.32
    2.859.85737.6239.4137.54–0.21
    2.9010.38530.5131.8030.31–0.66
    2.9510.93224.7125.6224.65–0.24
    3.0011.49719.9920.6019.79–1.00
    下载: 导出CSV
  • [1] TAN L, ZHANG W, ZHANG F. Research on energy efficiency system of new energy vehicle electric drive [C]//IOP Conference Series: Earth and Environmental Science. IOP Publishing, 2019, 223(1): 012028.
    [2] NADAÏ A, LABUSSIÈRE O. New energy resources in the making [J]. Energy Transitions, 2018, 256: 49.
    [3] MYATT P T, DHAM A K, CHANDRASEKHAR P, et al. A new empirical potential energy function for Ar2 [J]. Molecular Physics, 2018, 116(12): 1598–1623. doi: 10.1080/00268976.2018.1437932
    [4] GORBENKO I I, TROITSKAYA E P, PILIPENKO E A. Elastic properties of compressed rare-gas crystals in a model of deformable atoms [J]. Physics of the Solid State, 2017, 59(1): 132–140. doi: 10.1134/S1063783417010097
    [5] BONNET P. Equation of state for solid rare gases and alkali metals under pressure [J]. Physica B: Condensed Matter, 2016, 492: 50–54. doi: 10.1016/j.physb.2016.04.006
    [6] 田春玲, 刘福生, 蔡灵仓, 等. 多体相互作用对高压固氦状态方程的影响 [J]. 物理学报, 2006, 55(2): 764–769. doi: 10.3321/j.issn:1000-3290.2006.02.051
    TIAN C L, LIU F S, CAI L C, et al. Many-body contributions to the equation of state for highly compressed solid helium [J]. Acta Physica Sinica, 2006, 55(2): 764–769. doi: 10.3321/j.issn:1000-3290.2006.02.051
    [7] 郑兴荣. 多体相互作用对固氖物态方程的影响 [J]. 陇东学院学报, 2015, 26(5): 17. doi: 10.3969/j.issn.1674-1730.2015.05.005
    ZHENG X R. On the effect of multi-body interactions on the equation of state of solid neon [J]. Journal of Longdong University, 2015, 26(5): 17. doi: 10.3969/j.issn.1674-1730.2015.05.005
    [8] ABBASPOUR M, BORZOUIE Z. Equation of state, elastic constants, and melting curve of solid neon using an effective two-body potential including quantum corrections [J]. Fluid Phase Equilibria, 2014, 379: 167–174. doi: 10.1016/j.fluid.2014.07.026
    [9] SCHMIDT K M, VASQUEZ V R. A generalized method for the inversion of cohesive energy curves from isotropic and anisotropic lattice expansions [J]. Physical Chemistry Chemical Physics, 2015, 17(36): 23423–23437. doi: 10.1039/C5CP03792A
    [10] KOHN W, SHAM L J. Self-consistent equations including exchange and correlation effects [J]. Physical Review, 1965, 140(4A): A1133. doi: 10.1103/PhysRev.140.A1133
    [11] TIAN C L, LIU F S, CAI L C, et al. Ab initio calculations of many-body interactions for compressed solid argon [J]. The Journal of Chemical Physics, 2015, 143(17): 174506. doi: 10.1063/1.4935050
    [12] TIAN C L, WU N, LIU F S, et al. Four-body interaction energy for compressed solid krypton from quantum theory [J]. The Journal of Chemical Physics, 2012, 137(4): 044108. doi: 10.1063/1.4737183
    [13] AZIZ R. Accurate intermolecular potential for neon [J]. High Temperature High Pressures, 1980, 12(5): 565–577.
    [14] AZIZ R A, SLAMAN M J. The Ne-Ne interatomic potential revisited [J]. Chemical Physics, 1989, 130(1/2/3): 187–194.
    [15] AZIZ R A, SLAMAN M J. The argon and krypton interatomic potentials revisited [J]. Molecular Physics, 1986, 58(4): 679–697. doi: 10.1080/00268978600101501
    [16] FREIMAN Y A, TRETYAK S M. Many-body interactions and high-pressure equations of state in rare-gas solids [J]. Low Temperature Physics, 2007, 33(6): 545–552. doi: 10.1063/1.2746249
    [17] LOUBEYRE P. Three-body-exchange interaction in dense rare gases [J]. Physical Review B, 1988, 37(10): 5432. doi: 10.1103/PhysRevB.37.5432
    [18] LOUBEYRE P, LETOULLEC R, PINCEAUX J P, et al. Equation of state and phase diagram of solid 4He from single-crystal X-ray diffraction over a large P-T domain [J]. Physical Review Letters, 1993, 71(14): 2272. doi: 10.1103/PhysRevLett.71.2272
    [19] 经福谦. 实验物态方程导引 [M]. 北京: 科学出版社, 1999: 31.
    JING F Q. Introduction to experimental equation of state [M]. Beijing: Science Press, 1999: 31.
    [20] BERNE B J. Statistical mechanics [M]. New York: Plenum Press, 1977.
    [21] 武娜, 田春玲, 刘福生, 等. 固氖高压物态方程的量子理论计算 [J]. 高压物理学报, 2012, 26(1): 41–47. doi: 10.11858/gywlxb.2012.01.006
    WU N, TIAN C L, LIU F S, et al. Equation of state of solid neon from quantum calculation [J]. Chinese Journal of High Pressure Physics, 2012, 26(1): 41–47. doi: 10.11858/gywlxb.2012.01.006
    [22] 郑兴荣, 陈海军, 高晓红, 等. 高压固氩物态方程的量子理论计算 [J]. 高压物理学报, 2017, 31(4): 396–402. doi: 10.11858/gywlxb.2017.04.007
    ZHENG X R, CHEN H J, GAO X H, et al. Quantum calculation for equation of state of compressed solid argon [J]. Chinese Journal of High Pressure Physics, 2017, 31(4): 396–402. doi: 10.11858/gywlxb.2017.04.007
    [23] TANGE Y, NISHIHARA Y, TSUCHIYA T. Unified analyses for P-V-T equation of state of MgO: a solution for pressure-scale problems in high P-T experiments [J]. Journal of Geophysical Research: Solid Earth, 2009, 114: B03208.
    [24] JEPHCOAT A P. Rare-gas solids in the Earth’s deep interior [J]. Nature, 1998, 393(6683): 355. doi: 10.1038/30712
    [25] SCHWERDTFEGER P, HERMANN A. Equation of state for solid neon from quantum theory [J]. Physical Review B, 2009, 80(6): 064106. doi: 10.1103/PhysRevB.80.064106
    [26] SCHMIDT M W, BALDRIDGE K K, BOATZ J A, et al. General atomic and molecular electronic structure system [J]. Journal of Computational Chemistry, 1993, 14(11): 1347–1363. doi: 10.1002/jcc.540141112
    [27] HEMLEY R J, ZHA C S, JEPHCOAT A P, et al. X-ray diffraction and equation of state of solid neon to 110 GPa [J]. Physical Review B, 1989, 39(16): 11820. doi: 10.1103/PhysRevB.39.11820
    [28] DEWAELE A, DATCHI F, LOUBEYRE P, et al. High pressure-high temperature equations of state of neon and diamond [J]. Physical Review B, 2008, 77(9): 094106. doi: 10.1103/PhysRevB.77.094106
    [29] FINGER L W, HAZEN R M, ZOU G, et al. Structure and compression of crystalline argon and neon at high pressure and room temperature [J]. Applied Physics Letters, 1981, 39(11): 892–894. doi: 10.1063/1.92597
    [30] TAKEMURA K, WATANUKI T, OHWADA K, et al. Powder X-ray diffraction study of Ne up to 240 GPa [C]//Journal of Physics: Conference Series. IOP Publishing, 2010, 215(1): 012017.
    [31] ERRANDONEA D, BOEHLER R, JAPEL S, et al. Structural transformation of compressed solid Ar: an X-ray diffraction study to 114 GPa [J]. Physical Review B, 2006, 73(9): 092106. doi: 10.1103/PhysRevB.73.092106
    [32] POLIAN A, BESSON J M, GRIMSDITCH M, et al. Solid krypton: equation of state and elastic properties [J]. Physical Review B, 1989, 39(2): 1332. doi: 10.1103/PhysRevB.39.1332
    [33] ANDERSON M S, SWENSON C A. Experimental equations of state for the rare gas solids [J]. Journal of Physics and Chemistry of Solids, 1975, 36(3): 145–162. doi: 10.1016/0022-3697(75)90004-9
  • [1] 胡金彪经福谦 . 用冲击压缩数据计算物质结合能的一个简便解析方法. 高压物理学报, doi: 10.11858/gywlxb.1990.03.003
    [2] 冉宪文汤文辉 . 一种计算固体冷能、冷压和结合能的新方法. 高压物理学报, doi: 10.11858/gywlxb.2003.01.008
    [3] 刘福生芶清泉经福谦 . 多体作用对氦原子团簇He9压缩特性的贡献. 高压物理学报, doi: 10.11858/gywlxb.1995.04.005
    [4] 刘福生芶清泉经福谦 . 压缩氦原子团簇Hen(n=3,4,5)排斥势及其多体展开分量计算. 高压物理学报, doi: 10.11858/gywlxb.1995.03.005
    [5] 田春玲刘福生蔡灵仓经福谦 . 采用从头计算方法研究液态氦原子间等效对势. 高压物理学报, doi: 10.11858/gywlxb.2006.01.004
    [6] 郑朝阳赵纪军 . 含能材料的从头算分子动力学模拟. 高压物理学报, doi: 10.11858/gywlxb.2015.02.001
    [7] 赵海波肖波柏劲松段书超王刚华阚明先陈芳 . 拉氏方法模拟二维多介质可压缩流体的运动. 高压物理学报, doi: 10.11858/gywlxb.20170694
    [8] 孟川民施尚春董石孙悦焦荣珍杨向东 . 液氮冲击压缩特性的理论计算. 高压物理学报, doi: 10.11858/gywlxb.2002.03.009
    [9] 谭多望孙承纬 . 大锥角罩聚能装药射流理论计算方法. 高压物理学报, doi: 10.11858/gywlxb.2006.03.008
    [10] 柏劲松陈森华 . 重新初始化的LS方法跟踪二维可压缩多介质流界面运动. 高压物理学报, doi: 10.11858/gywlxb.2003.01.004
    [11] 孟川民焦荣珍施尚春董石孙悦杨向东经福谦 . 液氦冲击压缩特性的理论计算. 高压物理学报, doi: 10.11858/gywlxb.2002.01.012
    [12] 陈芳李平刘坤柏劲松林健宇季路成 . 基于PPM的界面压缩方法研究. 高压物理学报, doi: 10.11858/gywlxb.20180663
    [13] 刘福生洪德贵周雪芬 . 碳水混合物冲击压缩特性的理论研究. 高压物理学报, doi: 10.11858/gywlxb.2001.03.004
    [14] 张学莹赵宁王春武 . 多介质流动数值计算中的界面处理方法. 高压物理学报, doi: 10.11858/gywlxb.2006.03.005
    [15] 傅华刘仓理王文强李涛 . 冲击动力学中离散元与有限元相结合的计算方法研究. 高压物理学报, doi: 10.11858/gywlxb.2006.04.007
    [16] 孙燕云刘福生张明建徐利华 . 高温高密度条件下氮分子三体关联效应与冲击压缩特性. 高压物理学报, doi: 10.11858/gywlxb.2009.02.010
    [17] 刘小虎祝涛 . 颗粒体侵彻数值计算方法研究. 高压物理学报, doi: 10.11858/gywlxb.2010.01.004
    [18] 童石磊钟敏柏劲松陈森华 . 处理多介质界面的改进型Front Tracking方法. 高压物理学报, doi: 10.11858/gywlxb.2016.06.006
    [19] 李巧燕施尚春杨金文孙悦 . 钕铁硼的冲击压缩特性. 高压物理学报, doi: 10.11858/gywlxb.2007.02.016
    [20] 杨向东武保剑谢文胡栋经福谦 . 液氦的冲击压缩理论计算及其等效两体作用势的考察. 高压物理学报, doi: 10.11858/gywlxb.1997.03.003
  • 加载中
图(3)表(1)
计量
  • 文章访问数:  261
  • 阅读全文浏览量:  264
  • PDF下载量:  6
出版历程
  • 收稿日期:  2019-02-26
  • 录用日期:  2019-03-22
  • 网络出版日期:  2019-09-22

结合能的新势函数对高压稠密惰性元素压缩特性的影响

    作者简介:郑兴荣(1986-),男,硕士,讲师,主要从事凝聚态理论物理与材料计算、团簇结构研究. E-mail:zhengxingrong2006@163.com
  • 陇东学院电气工程学院物理系,甘肃 庆阳 745000

摘要: 基于量子理论和原子团簇理论,运用多体展开方法和第一性原理的从头算方法,提出了一种计算稠密惰性元素(氦、氖、氩和氪)原子结合能的新势函数,运用新公式研究了结合能对高压稠密惰性元素高压压缩特性的影响。此公式引入了一个物理参量\begin{document}$\beta $\end{document}(其值为0.5),使得势函数的表达形式更加简单、准确。对比结果表明,结合能的新势函数能够准确地描述多体相互作用对结合能的贡献,且平均相对误差在5%以内。结合能的新势函数对压缩特性的影响在当前实验压强范围内(氦60 GPa、氖238 GPa、氩114 GPa、氪128 GPa)做出了令人满意的描述,且与实验值及理论计算结果基本完全吻合,平均相对误差在3%以内。最后,以固氩的压强数据为例,验证了势函数的准确性。该势函数不仅适用于更宽密度和更高压强范围,而且对所有惰性元素原子各种状态的结合能、高压压缩特性、定容比热容、熔化曲线和弹性模量的研究具有重要的指导意义。

English Abstract

  • 随着石油、煤炭等传统能源的枯竭,以太阳能、地热能、风能、海洋能、生物质能和核聚变能等为主的新能源研究迫在眉睫,“能量”的利用效率作为新能源的衡量指标,已成为最热门、最受关注的研究内容之一[1-3]。结合能作为能量最基本的单位,是单个原子结合成相应材料时所放出的能量,放出的能量越多,则材料的结构越稳定,反之,越不稳定。另外,原子能够结合为晶体的根本原因在于,原子结合后整个系统具有更低的能量。因此,以性质最稳定的稠密惰性元素原子(晶体)作为研究对象,对结合能进行深入研究是一项艰巨而重要的任务。结合能作为稠密惰性元素(Rare-Gas Solid,RGS)原子的一个重要物性,可推导出稠密惰性元素原子的许多特性,如热力学参数、弹性系数、熔化曲线、弹性常数和高压压缩特性等[4-8]。另外,通过对原子结合能的理论研究也可以得到物理学、化学、能源学和晶体材料学的许多物质特性,如密度泛函理论[9]、从头算自洽场[10]、单双(三)重激发耦合簇理论CCSD(T)[11-12]、半经验及经验势函数[13-18]等。因此,能否得到简便、准确的结合能对研究其他物性至关重要。近年来,对结合能的研究已趋于白热化,大量文献报道了有关结合能的研究进展及成果[4-11],如通过结合能的能量曲线反转得到原子对势的相互作用、三体相互作用、弹性常数、熔化曲线以及各种体系的密度相关势能等。截至目前,在众多结合能的研究中,有两种典型的结合能公式被广泛应用:一种是基于从头算法的多体势能的简单叠加[10-12],此方法在求晶格原子结合能时只有当多体展开式满足收敛性并达到一定的截断精度时,所得的结合能才是准确、可靠的,因此仅适用于较低压强的稠密惰性元素;另一种是通过能带理论和密度泛函理论得到[19],虽然可应用于更高密度范围,但是准确度不高,而且所得的结合能考虑到温度的影响,情况更加复杂[20]。通过研究发现,宇宙当中,金星、木星和土星等恒星的内部成分中包含大量的惰性元素,其内部压强高达数百吉帕,甚至上千吉帕,因此研究结合能对稠密惰性元素高压压缩特性的影响具有潜在的应用价值。

    基于此,运用原子簇理论和量子力学理论,结合两体势的从头算方法和薛定谔方程,本研究提出一种计算RGS(RGS=He,Ne,Ar,Kr)原子结合能的新势函数,并运用新表达式研究结合能对高压稠密惰性元素状态方程的影响。其中新表达式只用到两体相互作用势U2和总原子相互作用势Vn,避免了对各多体项的分别计算,大大减少了计算量,并且由于U2Vn的收敛最快,所需考虑的近邻原子数较少,因此计算的繁琐程度减轻,计算时间缩短。这些优点将便于我们处理复杂情况下(包括高温、高压、等离子体状态等)的多体相互作用并进行分子动力学模拟。本工作旨在用最简单的势函数解决较复杂的问题,在此过程中两体和总原子势能都得到了合理的强调,同时它也是准确理解压缩特性、光谱散射和体积数据的关键。希望新方案能适用于所有电中性或离子气体系统、固体和液体等。

    • 基于原子团簇理论,对任意一个由n个原子组成的惰性元素原子团簇(RGS)n (RGS=He,Ne,Ar,Kr),该团簇的总相互作用势能Vn

      ${V_n} = {E_n}({r_1},{r_2}, \cdots ,{r_n}) - n{E_0}$

      式中:E0为稠密惰性元素孤立原子的基态能量,En(r1r2,···,rn)代表团簇中n个原子处于r1r2,···,rn位置时团簇体系的基态总能量,可通过求解团簇体系的薛定谔方程得到,即

      $\hat H\psi = E\psi $

      式中:$\hat H$表示n个惰性元素原子所构成团簇的哈密顿算符,可表示为

      $\hat H = - \sum\limits_I {\frac{1}{{2{M_I}}}} {\nabla _I} - \sum\limits_i {\frac{1}{2}} {\nabla _i} + \sum\limits_{I < J} {\frac{{{Z_I}{Z_J}}}{{{r_{IJ}}}}} + \sum\limits_{i < j} {\frac{1}{{{r_{ij}}}}} - \sum\limits_{I,j} {\frac{{{Z_I}}}{{{r_{Ij}}}}} $

      式中:MIZI分别为第I个惰性元素原子核的质量和核电荷数,ZJ为第J个原子核的核电荷数,rIJ为第I个原子核与第J个原子核之间的距离,rij为第i个电子与第j个电子之间的距离。利用Louwdin法求解团簇的总波函数$\psi \left( r \right)$。若第I个原子基态单电子轨道波函数为${\varphi _1}\left( r \right)$,则团簇的总波函数$\psi \left( r \right)$N个惰性元素原子中2N个单电子自旋轨道波函数所构成的反对称Slater行列式

      $\psi \left( r \right) = \frac{1}{{\sqrt {N!} }}\det |{\varphi _1}\left( 1 \right)\alpha \left( 1 \right){\varphi _2}\left( 2 \right)\beta \left( 2 \right){\varphi _2}\left( 3 \right)\alpha \left( 3 \right){\varphi _2}\left( 4 \right)\beta \left( 4 \right)\cdots{\varphi _N}\left( {2N - 1} \right)\alpha \left( {2N - 1} \right){\varphi _2}\left( {2N} \right)\beta \left( {2N} \right)|$

      (4)式中选取的双Slater基态单电子解析波函数为

      $\psi \left( r \right) = \frac{N}{{4{\text{π}}}}\left( {{{\rm e}^{ - \alpha r}} + k{{\rm e}^{ - \beta r}}} \right)$

    • 基于多体展开方法和第一性原理的从头算法[6-7, 11-12, 21],也可得到 (RGS)n中任意一个中心原子M的多体相互势能Un和总相互作用势能Vn,即

      ${V_n}(M) = \sum\limits_{k = 2}^n {{U_k}(M)} = \sum\limits_{i = 1}^{n - 1} {{u_2}} (i,M) + \sum\limits_{1 = i < j}^{n - 1} {{u_3}} (i,j,M) + \sum\limits_{1 = i < j < k}^{n - 1} {{u_4}} (i,j,k,M) + \sum\limits_{1 = i < j < k < l}^{n - 1} {{u_5}} (i,j,k,l,M) + \cdot \cdot \cdot $

      式中:Uk(M)代表k体项的总和,高阶项很小,在这里予以忽略。u2u3u4分别表示团簇中任意两个原子、3个原子和4个原子间的两体势、三体势和四体势,表达式如下

      ${u_2}(i,M) = E\left( {{r_i},{r_M}} \right) - 2{E_0}$

      $\begin{split} {u_3}\left( {i,j,M} \right) & = E\left( {{r_i},{r_j},{r_M}} \right) - 3{E_0} - [{u_2}\left( {i,M} \right)+{u_2}\left( {j,{\rm{M}}} \right){\rm{ + }}{u_2}\left( {i,j} \right)]\\ & = E\left( {{r_i},{r_j},{r_M}} \right) - 3{E_0} - \sum\limits_{1 = i < j}^{n - 1} {{u_2}(i,j)} = E\left( {{r_i},{r_j},{r_M}} \right) - 3{E_0} - {U_2}(M) \end{split}$

      ${u_4}\left( {i,j,k,M} \right) = E\left( {{r_i},{r_j},{r_k},{r_M}} \right) - 4{E_0} - \sum\limits_{1 = i < j < k}^{n - 1} {{u_3}(i,j,k)} - \sum\limits_{1 = i < j}^{n - 1} {{u_2}(i,j)}= E\left( {{r_i},{r_j},{r_k},{r_M}} \right) - 4{E_0} - {U_3}(M) - {U_2}(M)$

      本研究中,所有多体势u2(i, M)、u3(i, j, M)、u4(i, j, k, M)等的计算是基于第一性原理,运用GAMESS计算程序下的从头算方法得到。本研究只用到收敛最快、计算最省时、计算精度最高的两体相互作用势能U2(M),即

      ${U_2}(M) = \sum\limits_{i = 1}^{n - 1} {{u_2}(i,M)} $

    • 温度为T时,稠密惰性元素原子的压强P的表达式为[22]

      $P(V,T) = {P_n}(V) + {P_{\rm {th}}}(V,T) + {P_{\rm {zp}}}(V) = - \frac{{\partial {E_n}(V)}}{{\partial V}} + 9{k_{\rm B}}T \cdot \frac{\gamma }{V}{\left( {\frac{\varTheta }{T}} \right)^{ - 3}} \cdot \int_0^{\varTheta /T} {\frac{{{x^3}}}{{{{\rm e}^x} - 1}}} {\rm d}x + \frac{\gamma }{V}{E_{\rm {zp}}}$

      式中:Pn代表由多体相互作用引起的压强;Pzp表示由零点振动引起的压强;Pth表示由温度引起的压强,通过Debye-Mie-Grüneisen模型获得[23]$\gamma $为Grüneisen参数;$\varTheta $为德拜温度[16],不同体积下的德拜温度$\varTheta $可根据线性公式[24]和Grüneisen参数$\gamma $得到[22]Ezp为晶格的零点振动能,通过德拜近似模型得到。

    • 一般而言,只有当多体项Uk(M)在一定的截断精度下收敛时,多体展开方法才能准确地得到结合能E(M)。然而,随着近邻原子的增加,结合能E(M)的展开和计算变得越来越复杂和繁琐。为了解决该问题,本研究根据几何关系,建立了一个计算结合能的新势函数。

      稠密惰性元素原子的结合能E(M)可表示为[6]

      $E(M) = \frac{1}{2}{U_2}(M) + \frac{1}{3}{U_3}(M) + \frac{1}{4}{U_4}(M) + \frac{1}{5}{U_5}(M) + \cdots $

      根据(12)式,利用多项式函数的特性,构造一个多项式函数f(x),将总原子相互作用势能Vn(M)与结合能E(M)联系起来,即

      $f(x) = {U_2}(M)x + {U_3}(M){x^2} + {U_4}(M){x^3} + {U_5}(M){x^4} + \cdot \cdot \cdot = \sum\limits_{k = 2}^n {{U_k}(M)} {x^{k - 1}}$

      式中:f(x)是变量x 的多项式函数,各项系数Uk(M)(k=2,3,4,5,···)为各多项式的值。

      通过比较(6)式、(12)式和(13)式,有

      ${V_n}(M) = f(1)$

      ${U_2}(M) = {f'}(0)$

      因此,任意一个中心原子M对团簇体系结合能的贡献E(M)对应着f(x)在[0,1]区间的积分,即

      $\begin{split} E(M) & = \frac{1}{2}{U_2}(M) + \frac{1}{3}{U_3}(M) + \frac{1}{4}{U_4}(M) + \frac{1}{5}{U_5}(M) + \cdots \\ & = \left. {\left[ {\frac{1}{2}{U_2}(M){x^2} + \frac{1}{3}{U_3}(M){x^3} + \frac{1}{4}{U_4}(M){x^4} + \frac{1}{5}{U_5}(M){x^5} + \cdots } \right]} \right|_0^1\\ & = \int_0^1 {[{U_2}(M)x + {U_3}(M){x^2} + {U_4}(M){x^3} + {U_5}(M){x^4} + \cdots ]{\rm{d}}x} = \int_0^1 {f(x){\rm{d}}x} \end{split}$

      联立(12)式~(16)式,利用几何关系和定积分的物理意义,通过数据分析软件Origin 8.0得到U2(M)、Vn(M)和E(M)三者之间的关系,如图1所示。

      图  1  U2(M)、Vn(M)和E(M)三者的关系(曲线OE表示函数f(x),虚线OG代表U2(M)x

      Figure 1.  The correlations among U2(M), Vn(M), and E(M) (The curve line OE represents function f(x), and the dash line OG describes U2(M)x.)

      图1中,过点C做直线OD的垂线交OGA点,使${\text{△}}OAF$${\text{△}}EBF$的面积相等,即${S_{{\text{△}}OAF}} = {S_{{\text{△}}EBF}}$,则阴影面积SOFED可转化为${\text{△}}AOC$和梯形BCDE的面积之和,即

      ${S_{OFED}} = {S_{{\text{△}}AOC}} + {S_B}_{CDE}$

      结合(17)式,通过几何关系,得到一个计算积分(16)式的近似公式。该近似公式通过点O和点C的距离参数$\beta $来表示,即结合能E(M)的表达式为

      $\begin{split} E(M) = {S_{OFED}} & = {S_{{\text{△}}AOC}} + {S_{BCDE}} = \frac{1}{2}{\beta ^2}{U_2}(M) + \frac{1}{2}(1 - \beta )[\beta {V_n}(M) + {V_n}(M)]\\ & = \frac{1}{2}[{\beta ^2}{U_2}(M) + (1 - {\beta ^2}){V_n}(M)] \end{split}$

      (18)式即为稠密惰性元素原子结合能的新势函数。由图1可以看出,当总原子势能Vn(M)已知时,结合能E(M)的值可以得到。

    • 稠密惰性元素原子结合能的新势函数为

      $E(M) = \frac{1}{2}[{\beta ^2}{U_2}(M) + (1 - {\beta ^2}){V_n}(M)]$

      (19)式是关于参数$\beta $的新势函数,来自两部分贡献:第1部分是由从头算方法计算得到的两体势${\beta ^2}{U_2}\left( M \right)$,对应着图1中的${\text{△}}OGD$的面积;第2部分是通过求解薛定谔方程得到的总多体势$\left( {1 - {\beta ^2}} \right){V_n}\left( M \right)$,对应着图1中的${\text{△}}OED$的面积。

      图1和(19)式可以得到:若$\beta $=1,则(19)式变为传统的两体相互势能的一半,此时图1OEOG重合,这显然高于E(M)的实际值;若$\beta $=0,则(19)式变为平均总相互作用势能,此时对应图1中的E点,其值低于E(M)的实际值。在精确度要求较低的情况下,这两种方法可以粗略地计算原子的结合能。本研究通过第一性原理的从头算法计算发现,结合能E(M)的准确值介于以上两种情况之间,对应的$\beta $介于0和1之间,与图1中的阴影区域有关。由于U2(M)和Vn(M)的收敛速度非常快,近邻原子的数量也不大,所以相比传统方法的计算时间很短,以便在任何温度下研究多体相互作用和分子动力学模拟。

    • 由于总相互作用势能Vn(M)是一个交错级数,因此通过比较(12)式、(16)式和(18)式发现,如果只考虑两体势,忽略三体及以上多体势,结合能E(M)是参数$\beta $的独立函数,此时$\beta $=1,这也是求解结合能的传统方法;如果考虑两体势和三体势,忽略四体势及四体以上势,可以推导出$\beta $=(1/3)1/2=0.577;如果考虑两体、三体和四体势,忽略五体及五体以上势,则有0<$\beta $<0.577。随着更高多体势被考虑,$\beta $值就越接近0.5。也就是说,可以得到参数$\beta $的一个近似值($\beta $=0.5),它适用于所有的稠密惰性元素原子的结合能。因此,结合能的表达式可简化为

      $E(M) = \frac{1}{2}\left[ {\frac{1}{4}{U_2}(M) + \frac{3}{4}{V_n}(M)} \right]$

      在这个新的结合能公式中,U2(M)可以通过从头算法或者Aziz对势准确得到,且只占结合能E(M)的25%;Vn(M)可通过求解薛定谔方程得到,且占结合能E(M)的75%。利用新公式,得到了结合能E(M)的数值,并将数值结果与前人的理论计算结果、实验值[6-7, 11-12, 16, 18, 25]进行了比较,如图2所示。

      图  2  稠密惰性元素原子结合能的比较

      Figure 2.  Comparison of the cohesive energy for rare-gas solids

      结果表明,新的结合能E(M)与参考文献中的实验结果吻合较好,平均相对误差在3%以内,个别在3%~5%范围内。综上所述,新结合能的表达式具有如下优点:

      (1)新表达式能够更广泛地应用到基于RHF理论的分子体系密度范围之内;

      (2)由于原子多体相互作用势中的坐标独立性,新表达式对于原子晶体的结构转变、相变、弹性系数和高压熔化线等的理解具有特殊的价值;

      (3)依据新表达式的理论计算方法,一些简单的理论和表达式也可以很容易地扩展到原子多体展开理论或分子势理论,如流体微扰理论[19-20]等;

      (4)对于短程排斥势,一些离子体系和分子体系的情况与原子体系的情况类似,因此新表达式和该方案有望应用于这些体系。

    • 根据新结合能的计算结果,运用(11)式可得到稠密惰性元素原子的新结合能对高压压缩特性的影响,并与前人通过从头算理论得到的计算结果[6-7, 11-12, 26]、实验值[15, 18, 24, 27-33]进行比较,如图3所示,其中RV分别为原子间距和对应的摩尔体积。数值结果表明,本研究的数值结果与从头算法结果相近,在目前的实验范围内与实验结果基本一致(氦1~60 GPa[6],氖1~238 GPa[30],氩1~114 GPa[31],氪1~128 GPa[24]),且平均相对误差都在3%以内,极个别在3%~5%。

      图  3  稠密惰性元素原子的高压压缩特性的比较

      Figure 3.  Comparison of the compressibility of rare-gas solids

      为了进一步验证本研究提出的理论,表1列出了固氩的实验压强值、理论计算值,以及本研究的计算值与实验值的相对误差。可以看出,本研究的计算结果与实验值[24]基本吻合,平均相对误差小于3%,说明新方案是可行的。在此工作中,结合能的新表达式准确地解释了惰性元素物质在高压高密度下的稠密特性。

      RV/(cm3·mol–1)P/GPaError/%
      Exp.[24]Ab initio[11]Eq.(10)
      2.405.887237.61248.95247.043.97
      2.456.262194.11204.44200.723.41
      2.506.653158.51167.59163.253.00
      2.557.061129.38137.12132.812.65
      2.607.484105.53111.96107.982.32
      2.657.92486.0291.2387.701.95
      2.708.38170.0674.1871.121.51
      2.758.85657.0060.2057.571.00
      2.809.34846.3448.7546.490.32
      2.859.85737.6239.4137.54–0.21
      2.9010.38530.5131.8030.31–0.66
      2.9510.93224.7125.6224.65–0.24
      3.0011.49719.9920.6019.79–1.00

      表 1  固氩的压强分量

      Table 1.  The pressure components of solid argon

    • 基于量子理论和原子团簇理论,运用多体展开方法和第一性原理的从头算方法,提出了一种新的计算稠密惰性元素(RGS=氦,氖,氩和氪)原子结合能的经验势函数,运用新经验势研究了结合能对高压稠密惰性元素状态方程的影响。此经验势函数由两部分组成:一部分是由从头算法计算得到的两体相互作用势能,另一部分是通过求解薛定谔方程得到的总相互作用势能。此函数引入了一个新的物理参量$\beta $(其值为0.5),使得结合能的经验表达式更加简单、准确。比较结果表明:结合能的新经验势函数能够准确地描述多体相互作用对结合能的贡献,且相对误差在5%以内,能够更广泛地应用到基于RHF理论的分子体系密度范围之内,对原子晶体的结构转变、相变、弹性系数和高压熔化线的理解具有特殊的价值;运用新经验势的理论计算方法,可以将一些简单的理论、表达式很容易地扩展到原子多体展开理论或分子势理论;此外,此表达式和方案做法有望应用于一些离子体系和分子体系的短程排斥势情况。在高压区域,结合能的新势函数对高压压缩特性的影响在当前实验范围内(氦1~60 GPa,氖1~238 GPa,氩1~114 GPa,氪1~128 GPa)做出了令人满意的描述,数值结果与实验结果、理论计算结果基本吻合,平均相对误差在3%以内。最后,以固氩的压强数据为例,验证了该势函数的准确性。该经验势函数不仅适用于更宽的密度和压强范围,而且对所有惰性元素原子的各种状态的高压压缩特性研究具有重要的指导意义。

参考文献 (33)

目录

    /

    返回文章
    返回