搜索

x

留言板

尊敬的读者、作者、审稿人, 关于本刊的投稿、审稿、编辑和出版的任何问题, 您可以本页添加留言。我们将尽快给您答复。谢谢您的支持!

姓名
邮箱
手机号码
标题
留言内容
验证码

热峰作用下单斜ZrO2相变过程的分子动力学模拟

赵中华 渠广昊 姚佳池 闵道敏 翟鹏飞 刘杰 李盛涛

引用本文:
Citation:

热峰作用下单斜ZrO2相变过程的分子动力学模拟

赵中华, 渠广昊, 姚佳池, 闵道敏, 翟鹏飞, 刘杰, 李盛涛

Molecular dynamics simulation of phase transition by thermal spikes in monoclinic ZrO2

Zhao Zhong-Hua, Qu Guang-Hao, Yao Jia-Chi, Min Dao-Min, Zhai Peng-Fei, Liu Jie, Li Sheng-Tao
PDF
HTML
导出引用
  • ZrO2陶瓷耐高温、耐腐蚀、抗辐照性能强, 是极具前景的反应堆惰性基质燃料和锕系元素固化材料. 本文联合使用热峰模型和分子动力学方法, 模拟了核辐射环境下ZrO2的相变过程: 基于热峰模型, 从快速重离子注入后能量沉积和传导的多物理过程出发, 建立热扩散方程, 求得ZrO2晶格温度时空演变特性; 然后运用分子动力学方法模拟了该热峰作用下, 单斜ZrO2相变的微观物理过程. 研究发现, 电子能损为30 keV·nm–1的单一快速重离子注入后, ZrO2中心产生一个半径为7 nm的柱形径迹, 径迹中心晶格迅速熔融, Zr原子配位数由7降至4—6, 2 ps时开始结晶并形成空洞, 空洞周围为非晶区, 非晶区外Zr原子配位数变为8, 同时X射线衍射(X-ray diffraction, XRD)计算和分析结果确认发生了单斜相向四方相的转变. 随着热峰能量向周围传递, 相变区逐渐扩大. 经热峰计算和分子动力学模拟, 辐照诱导ZrO2由单斜相转为四方相的快速重离子的电子能损阈值为21 keV·nm–1.
    Owing to its excellent corrosion, radiation and high temperature resistance, ZrO2 has been considered as a strong candidate material for inert fuel for the incineration of actinides. In this paper, a combination of thermal spike model and molecular dynamics is used to simulate the phase transition process of ZrO2 in the nuclear radiation environment. Based on the thermal spike model, two coupled diffusion equations are established with considering the multiple physical process of energy deposition and transmission after the implantation of swift heavy ions into target material. The space-time evolution characteristics of ZrO2 lattice temperature are obtained by solving the coupled diffusion equations numerically. Then the phase transformation of ZrO2 form monoclinic phase to tetragonal phase under the thermal spike is investigated on an atomic scale by means of molecular dynamics. It is found that a cylindrical track with a radius of 7 nm is generated in the center of ZrO2 after the implantation of swift heavy ion with an electronic energy loss of 30 keV·nm–1. The lattice melts immediately in the center of track, accompanied with the coordination number of Zr decreasing from 7 to 4–6. Then at about 2 ps, the melting zone gradually turns cool and recrystallized. And in the center of the melting zone, voids begin to form and are surrounded by a highly disordered amorphous region. Meanwhile, tetragonal phase of ZrO2, whose coordination number of Zr is 8, is formed at the periphery of the amorphous region, which is also confirmed by the XRD calculation results. As energy transfers from track center to the surround, the tetragonal region gradually develops into the whole system, accompanied with the increase of voids size. The simulation results indicate that the irradiation of ZrO2 with swift heavy ions can lead to a transformation from the monoclinic to the tetragonal phase when the deposited electronic energy loss exceeds an effective threshold ~21 keV·nm–1, greater than the experimental value (12 keV·nm–1), which was mainly due to the large difference between the simulated and measured incident ion fluences and the accuracy of the force field used in the molecular dynamics.
      通信作者: 李盛涛, sli@mail.xjtu.edu.cn
    • 基金项目: 国家自然科学基金(批准号: 11690040, 11690041)和强脉冲辐射环境模拟与效应国家重点实验室(西北核技术研究院)(批准号: SKLIPR1709)资助的课题
      Corresponding author: Li Sheng-Tao, sli@mail.xjtu.edu.cn
    • Funds: Project supported by the National Natural Science Foundation of China (Grant Nos. 11690040, 11690041), and the State Key Laboratory of Intense Pulsed Radiation Simulation and Effect, Northwest Institute of Nuclear Technology, China (Grant No. SKLIPR1709)
    [1]

    Lee Y W, Joung C Y, Kim S H, Lee S C 2001 Met. Mater. Int. 7 159Google Scholar

    [2]

    Aoki T, Sagara H, Han C Y 2019 Ann. Nucl. Energy 126 427Google Scholar

    [3]

    Guo T, Wang C, Lü J, Liang T X 2016 J. Nucl. Mater. 481 66Google Scholar

    [4]

    Nandi C, Grover V, Bhandari K, Bhattacharya S, Mishra S, Banerjee J, Prakash A, Tyagi A K 2018 J. Nucl. Mater. 508 82Google Scholar

    [5]

    Kelly P M, Rose L R F 2002 Prog. Mater. Sci. 47 463Google Scholar

    [6]

    Gong X, Ding S R, Zhao Y M, Huo Y Z, Zhang L, LI Y M 2013 Mech. Mater. 65 110Google Scholar

    [7]

    Benyagoub A, Levesque F, Couvreur F, Gibert-Mougel C, Dufour C, Paumier E 2000 Appl. Phys. Lett. 77 3197Google Scholar

    [8]

    Benyagoub A 2006 Nucl. Instrum. Methods Phys. Res., Sect. B 245 225Google Scholar

    [9]

    Benyagoub A 2010 Nucl. Instrum. Methods Phys. Res., Sect. B 268 2968Google Scholar

    [10]

    Schuster B, Fujara F, Merk B, Neumann R, Seidl T, Trautmann C 2012 Nucl. Instrum. Methods Phys. Res., Sect. B 277 45Google Scholar

    [11]

    Manzini A M, Alurralde M A, Giménez G, Luca V 2016 J. Nucl. Mater. 482 175Google Scholar

    [12]

    O'Connell J H, Lee M E, Skuratov V A, Rymzhanov R A 2020 Nucl. Instrum. Methods Phys. Res., Sect. B 473 1Google Scholar

    [13]

    Benyagoub A 2005 Phys. Rev. B 72 094114Google Scholar

    [14]

    Lu F Y, Wang J W, Lang M, Toulemonde M, Namavar F, Trautmann C, Zhang J M, Ewing R C, Lian J 2012 Phys. Chem. Chem. Phys. 14 12295Google Scholar

    [15]

    Sharma A, Varshney M, Shin H J, Kumar Y, Gautam S, Chae K H 2014 Chem. Phys. Lett. 592 85Google Scholar

    [16]

    Wang J W, Lang M, Ewing R C, Becker U 2013 J. Phys. Condens. Matter 25 135001Google Scholar

    [17]

    王成龙, 王庆宇, 张跃, 李忠宇, 洪兵, 苏折, 董良 2014 物理学报 63 153402Google Scholar

    Wang C L, Wang Q Y, Zhang Y, Li Z Y, Hong B, Su Z, Dong L 2014 Acta Phys. Sin. 63 153402Google Scholar

    [18]

    袁伟, 彭海波, 杜鑫, 律鹏, 沈扬皓, 赵彦, 陈亮, 王铁山 2017 物理学报 66 106102Google Scholar

    Yuan W, Peng H B, Du X, Lv P, Shen Y H, Zhao Y, Chen L, Wang T S 2017 Acta Phys. Sin. 66 106102Google Scholar

    [19]

    Toulemonde M, Trautmann C, Balanzat E, Hjort K, Weidinger A 2004 Nucl. Instrum. Methods Phys. Res., Sect. B 216 1Google Scholar

    [20]

    Meftah A, Costantini J M, Khalfaoui N, Boudjadar S, Stoquert J P, Studer F, Toulemonde M 2005 Nucl. Instrum. Methods Phys. Res., Sect. B 237 563Google Scholar

    [21]

    Stodel C, Toulemonde M, Fransen C, Jacquot B, Clément C, Frémont G, Michel M, Dufour C 2018 29th International Conference of the International Nuclear Target Development Society East Lansing, Michigan, United States of America, October 7–12, 2018 p05001

    [22]

    Toulemonde M, Dufour Ch, Meftah A, Paumier E 2000 Nucl. Instrum. Methods Phys. Res., Sect. B 166 903Google Scholar

    [23]

    Coughlin J P, King E G 1950 J. Am. Chem. Soc. 72 2262Google Scholar

    [24]

    Kim W K, Shim J H, Kaviany M 2016 J. Nucl. Mater. 491 126Google Scholar

    [25]

    Sobolev V, Lemehov S 2006 J. Nucl. Mater. 352 300Google Scholar

    [26]

    Mévrel R, Laizet J C, Azzopardi A, Leclercq B, Poulain M, Lavigne O, Demange D 2004 J. Eur. Ceram. Soc. 24 3081Google Scholar

    [27]

    Ziegler J F, Ziegler M D, Biersack J P 2008 Nucl. Instrum. Methods Phys. Res., Sect. B 268 1818Google Scholar

    [28]

    Plimpton S 1995 J. Comput. Phys. 117 1Google Scholar

    [29]

    Buckingham R A 1938 Proc. R. Soc. London, Ser. A 168 264Google Scholar

    [30]

    Houska J 2016 Comput. Mater. Sci. 111 209Google Scholar

    [31]

    López P, Pelaz L, Santos I, Marqués L A, Aboy M 2012 J. Appl. Phys. 111 033519Google Scholar

    [32]

    Nordlund K 1995 Comput. Mater. Sci. 3 448Google Scholar

    [33]

    Coleman S P, Sichani M M, Spearot D E 2014 JOM 66 408Google Scholar

    [34]

    杨辉, 冯泽华, 王贺然, 张云鹏, 陈铮, 信天缘, 宋小蓉, 吴璐, 张静 2021 物理学报 70 054601Google Scholar

    Yang H, Feng Z H, Wang H R, Zhang Y P, Chen Z, Xin T Y, Song X R, Wu L, Zhang J 2021 Acta Phys. Sin. 70 054601Google Scholar

    [35]

    Szenes G 1995 Phys. Rev. B 51 8026Google Scholar

  • 图 1  单斜ZrO2原胞结构模型

    Fig. 1.  The crystal structure of monoclinic ZrO2.

    图 2  热峰注入靶材料示意图(Se = 30 keV·nm–1)

    Fig. 2.  Diagram of thermal spike injection into target material (Se = 30 keV·nm–1).

    图 3  晶格温度分布曲线 (a) 不同位置晶格温度随时间演化; (b) 不同时刻晶格温度的径向分布

    Fig. 3.  Distributions of lattice temperature: (a) Evolution of lattice temperature at different radial distances; (b) radial distribution of lattice temperature at different times.

    图 4  不同时刻ZrO2原子结构变化. 所有原子沿热峰注入方向投影在(010)平面上

    Fig. 4.  Atomic structure of ZrO2 at different times. All the atoms are projected on (010) plane along the injection direction of thermal spike.

    图 5  辐照前后以及标准ZrO2 XRD图谱

    Fig. 5.  XRD patterns of pristine, ion-irradiated and standard ZrO2.

    图 6  不同时刻不同配位数Zr原子空间分布. 所有原子沿离子入射方向投影在(010)平面上, 颜色不同代表配位数不同

    Fig. 6.  Spatial distribution of Zr atoms with different coordination number at various times. All the atoms are projected on (010) plane along the injection direction of thermal spike. Different coordination numbers are distinguished by the atomic color.

    图 7  配位数随电子能损值变化

    Fig. 7.  The change of coordination number with electronic energy loss.

    图 8  不同配位数Zr原子空间分布 (a) 电子能损值为20 keV·nm–1; (b)电子能损值为21 keV·nm–1. 所有原子沿离子入射方向投影在(010)平面上, 颜色不同代表配位数不同

    Fig. 8.  Spatial distribution of Zr atoms with different coordination number: (a) Se = 20 keV nm–1; (b) Se = 21 keV nm–1. All the atoms are projected on (010) plane along the injection direction of thermal spike. Different coordination numbers are distinguished by the atomic color.

  • [1]

    Lee Y W, Joung C Y, Kim S H, Lee S C 2001 Met. Mater. Int. 7 159Google Scholar

    [2]

    Aoki T, Sagara H, Han C Y 2019 Ann. Nucl. Energy 126 427Google Scholar

    [3]

    Guo T, Wang C, Lü J, Liang T X 2016 J. Nucl. Mater. 481 66Google Scholar

    [4]

    Nandi C, Grover V, Bhandari K, Bhattacharya S, Mishra S, Banerjee J, Prakash A, Tyagi A K 2018 J. Nucl. Mater. 508 82Google Scholar

    [5]

    Kelly P M, Rose L R F 2002 Prog. Mater. Sci. 47 463Google Scholar

    [6]

    Gong X, Ding S R, Zhao Y M, Huo Y Z, Zhang L, LI Y M 2013 Mech. Mater. 65 110Google Scholar

    [7]

    Benyagoub A, Levesque F, Couvreur F, Gibert-Mougel C, Dufour C, Paumier E 2000 Appl. Phys. Lett. 77 3197Google Scholar

    [8]

    Benyagoub A 2006 Nucl. Instrum. Methods Phys. Res., Sect. B 245 225Google Scholar

    [9]

    Benyagoub A 2010 Nucl. Instrum. Methods Phys. Res., Sect. B 268 2968Google Scholar

    [10]

    Schuster B, Fujara F, Merk B, Neumann R, Seidl T, Trautmann C 2012 Nucl. Instrum. Methods Phys. Res., Sect. B 277 45Google Scholar

    [11]

    Manzini A M, Alurralde M A, Giménez G, Luca V 2016 J. Nucl. Mater. 482 175Google Scholar

    [12]

    O'Connell J H, Lee M E, Skuratov V A, Rymzhanov R A 2020 Nucl. Instrum. Methods Phys. Res., Sect. B 473 1Google Scholar

    [13]

    Benyagoub A 2005 Phys. Rev. B 72 094114Google Scholar

    [14]

    Lu F Y, Wang J W, Lang M, Toulemonde M, Namavar F, Trautmann C, Zhang J M, Ewing R C, Lian J 2012 Phys. Chem. Chem. Phys. 14 12295Google Scholar

    [15]

    Sharma A, Varshney M, Shin H J, Kumar Y, Gautam S, Chae K H 2014 Chem. Phys. Lett. 592 85Google Scholar

    [16]

    Wang J W, Lang M, Ewing R C, Becker U 2013 J. Phys. Condens. Matter 25 135001Google Scholar

    [17]

    王成龙, 王庆宇, 张跃, 李忠宇, 洪兵, 苏折, 董良 2014 物理学报 63 153402Google Scholar

    Wang C L, Wang Q Y, Zhang Y, Li Z Y, Hong B, Su Z, Dong L 2014 Acta Phys. Sin. 63 153402Google Scholar

    [18]

    袁伟, 彭海波, 杜鑫, 律鹏, 沈扬皓, 赵彦, 陈亮, 王铁山 2017 物理学报 66 106102Google Scholar

    Yuan W, Peng H B, Du X, Lv P, Shen Y H, Zhao Y, Chen L, Wang T S 2017 Acta Phys. Sin. 66 106102Google Scholar

    [19]

    Toulemonde M, Trautmann C, Balanzat E, Hjort K, Weidinger A 2004 Nucl. Instrum. Methods Phys. Res., Sect. B 216 1Google Scholar

    [20]

    Meftah A, Costantini J M, Khalfaoui N, Boudjadar S, Stoquert J P, Studer F, Toulemonde M 2005 Nucl. Instrum. Methods Phys. Res., Sect. B 237 563Google Scholar

    [21]

    Stodel C, Toulemonde M, Fransen C, Jacquot B, Clément C, Frémont G, Michel M, Dufour C 2018 29th International Conference of the International Nuclear Target Development Society East Lansing, Michigan, United States of America, October 7–12, 2018 p05001

    [22]

    Toulemonde M, Dufour Ch, Meftah A, Paumier E 2000 Nucl. Instrum. Methods Phys. Res., Sect. B 166 903Google Scholar

    [23]

    Coughlin J P, King E G 1950 J. Am. Chem. Soc. 72 2262Google Scholar

    [24]

    Kim W K, Shim J H, Kaviany M 2016 J. Nucl. Mater. 491 126Google Scholar

    [25]

    Sobolev V, Lemehov S 2006 J. Nucl. Mater. 352 300Google Scholar

    [26]

    Mévrel R, Laizet J C, Azzopardi A, Leclercq B, Poulain M, Lavigne O, Demange D 2004 J. Eur. Ceram. Soc. 24 3081Google Scholar

    [27]

    Ziegler J F, Ziegler M D, Biersack J P 2008 Nucl. Instrum. Methods Phys. Res., Sect. B 268 1818Google Scholar

    [28]

    Plimpton S 1995 J. Comput. Phys. 117 1Google Scholar

    [29]

    Buckingham R A 1938 Proc. R. Soc. London, Ser. A 168 264Google Scholar

    [30]

    Houska J 2016 Comput. Mater. Sci. 111 209Google Scholar

    [31]

    López P, Pelaz L, Santos I, Marqués L A, Aboy M 2012 J. Appl. Phys. 111 033519Google Scholar

    [32]

    Nordlund K 1995 Comput. Mater. Sci. 3 448Google Scholar

    [33]

    Coleman S P, Sichani M M, Spearot D E 2014 JOM 66 408Google Scholar

    [34]

    杨辉, 冯泽华, 王贺然, 张云鹏, 陈铮, 信天缘, 宋小蓉, 吴璐, 张静 2021 物理学报 70 054601Google Scholar

    Yang H, Feng Z H, Wang H R, Zhang Y P, Chen Z, Xin T Y, Song X R, Wu L, Zhang J 2021 Acta Phys. Sin. 70 054601Google Scholar

    [35]

    Szenes G 1995 Phys. Rev. B 51 8026Google Scholar

  • [1] 冯妍卉, 冯黛丽, 褚福强, 邱琳, 孙方远, 林林, 张欣欣. 纳米组装相变储热材料的热设计前沿. 物理学报, 2022, 71(1): 016501. doi: 10.7498/aps.71.20211776
    [2] 梅涛, 陈占秀, 杨历, 朱洪漫, 苗瑞灿. 非对称纳米通道内界面热阻的分子动力学研究. 物理学报, 2020, 69(22): 224701. doi: 10.7498/aps.69.20200491
    [3] 第伍旻杰, 胡晓棉. 单晶Ce冲击相变的分子动力学模拟. 物理学报, 2020, 69(11): 116202. doi: 10.7498/aps.69.20200323
    [4] 种涛, 王桂吉, 谭福利, 赵剑衡, 唐志平. 窗口声阻抗对锆相变动力学的影响. 物理学报, 2018, 67(7): 070204. doi: 10.7498/aps.67.20172198
    [5] 林长鹏, 刘新健, 饶中浩. 铝纳米颗粒的热物性及相变行为的分子动力学模拟. 物理学报, 2015, 64(8): 083601. doi: 10.7498/aps.64.083601
    [6] 张宝玲, 宋小勇, 侯氢, 汪俊. 高密度氦相变的分子动力学研究. 物理学报, 2015, 64(1): 016202. doi: 10.7498/aps.64.016202
    [7] 张金平, 张洋洋, 李慧, 高景霞, 程新路. 纳米铝热剂Al/SiO2层状结构铝热反应的分子动力学模拟. 物理学报, 2014, 63(8): 086401. doi: 10.7498/aps.63.086401
    [8] 肖红星, 龙冲生. UO2 晶体中低密勒指数晶面表面能的分子动力学模拟. 物理学报, 2013, 62(10): 103104. doi: 10.7498/aps.62.103104
    [9] 饶中浩, 汪双凤, 张艳来, 彭飞飞, 蔡颂恒. 相变材料热物理性质的分子动力学模拟. 物理学报, 2013, 62(5): 056601. doi: 10.7498/aps.62.056601
    [10] 唐翠明, 赵锋, 陈晓旭, 陈华君, 程新路. Al与α-Fe2O3纳米界面铝热反应的从头计算分子动力学研究. 物理学报, 2013, 62(24): 247101. doi: 10.7498/aps.62.247101
    [11] 周婷婷, 黄风雷. HMX不同晶型热膨胀特性及相变的ReaxFF分子动力学模拟. 物理学报, 2012, 61(24): 246501. doi: 10.7498/aps.61.246501
    [12] 杨平, 吴勇胜, 许海锋, 许鲜欣, 张立强, 李培. TiO2/ZnO纳米薄膜界面热导率的分子动力学模拟. 物理学报, 2011, 60(6): 066601. doi: 10.7498/aps.60.066601
    [13] 汪志刚, 吴亮, 张杨, 文玉华. 面心立方铁纳米粒子的相变与并合行为的分子动力学研究. 物理学报, 2011, 60(9): 096105. doi: 10.7498/aps.60.096105
    [14] 陈斌, 彭向和, 范镜泓, 孙士涛, 罗吉. 考虑相变的热弹塑性本构方程及其应用. 物理学报, 2009, 58(13): 29-S34. doi: 10.7498/aps.58.29
    [15] 王海燕, 崔红保, 历长云, 李旭升, 王狂飞. AlAs相变及热动力学性质的第一性原理研究. 物理学报, 2009, 58(8): 5598-5603. doi: 10.7498/aps.58.5598
    [16] 陈贺胜. 带有2+1味道Wilson费米子的格点量子色动力学在有限温度、有限密度下的相变. 物理学报, 2009, 58(10): 6791-6797. doi: 10.7498/aps.58.6791
    [17] 邵建立, 秦承森, 王裴. 动态压缩下马氏体相变力学性质的微观研究. 物理学报, 2009, 58(3): 1936-1941. doi: 10.7498/aps.58.1936
    [18] 邵建立, 王 裴, 秦承森, 周洪强. 冲击加载下孔洞诱导相变形核分析. 物理学报, 2008, 57(2): 1254-1258. doi: 10.7498/aps.57.1254
    [19] 邵建立, 王 裴, 秦承森, 周洪强. 铁冲击相变的分子动力学研究. 物理学报, 2007, 56(9): 5389-5393. doi: 10.7498/aps.56.5389
    [20] 崔新林, 祝文军, 邓小良, 李英骏, 贺红亮. 冲击波压缩下含纳米孔洞单晶铁的结构相变研究. 物理学报, 2006, 55(10): 5545-5550. doi: 10.7498/aps.55.5545
计量
  • 文章访问数:  2235
  • PDF下载量:  131
  • 被引次数: 0
出版历程
  • 收稿日期:  2020-11-06
  • 修回日期:  2021-02-23
  • 上网日期:  2021-06-25
  • 刊出日期:  2021-07-05

热峰作用下单斜ZrO2相变过程的分子动力学模拟

  • 1. 西安交通大学电气工程学院, 电力设备电气绝缘国家重点实验室, 西安 710049
  • 2. 中国科学院近代物理研究所, 兰州 730000
  • 通信作者: 李盛涛, sli@mail.xjtu.edu.cn
    基金项目: 国家自然科学基金(批准号: 11690040, 11690041)和强脉冲辐射环境模拟与效应国家重点实验室(西北核技术研究院)(批准号: SKLIPR1709)资助的课题

摘要: ZrO2陶瓷耐高温、耐腐蚀、抗辐照性能强, 是极具前景的反应堆惰性基质燃料和锕系元素固化材料. 本文联合使用热峰模型和分子动力学方法, 模拟了核辐射环境下ZrO2的相变过程: 基于热峰模型, 从快速重离子注入后能量沉积和传导的多物理过程出发, 建立热扩散方程, 求得ZrO2晶格温度时空演变特性; 然后运用分子动力学方法模拟了该热峰作用下, 单斜ZrO2相变的微观物理过程. 研究发现, 电子能损为30 keV·nm–1的单一快速重离子注入后, ZrO2中心产生一个半径为7 nm的柱形径迹, 径迹中心晶格迅速熔融, Zr原子配位数由7降至4—6, 2 ps时开始结晶并形成空洞, 空洞周围为非晶区, 非晶区外Zr原子配位数变为8, 同时X射线衍射(X-ray diffraction, XRD)计算和分析结果确认发生了单斜相向四方相的转变. 随着热峰能量向周围传递, 相变区逐渐扩大. 经热峰计算和分子动力学模拟, 辐照诱导ZrO2由单斜相转为四方相的快速重离子的电子能损阈值为21 keV·nm–1.

English Abstract

    • 随着核能的不断开发利用, 高放射性核废料持续增多, 如何有效地提高核燃料利用率并减少放射性核废料, 已成为核能可持续发展所必须解决的关键问题. 近年来研究发现, 在反应堆以及加速器驱动的次临界系统中, 采用惰性基质燃料能够有效地促进钚的回收利用和次锕系元素的嬗变, 从源头上降低核废物中长寿命放射性元素的含量[1,2]. 氧化锆(ZrO2)陶瓷耐高温、耐腐蚀、抗辐照性能强, 是惰性基质燃料的候选材料之一[3,4]. 然而, ZrO2在不同温度下呈现三种不同的晶相, 1170 ℃以下为单斜相, 1170—2370 ℃之间为四方相, 2370—2715 ℃之间为立方相. ZrO2在加热冷却过程中, 可发生单斜至四方相变(1100—1250 ℃), 以及四方至单斜相变(1050—700 ℃), 后者的剪切应变为16%—18%, 导致体积增加约5%[5]. 核反应堆的强辐射环境可诱导ZrO2发生相变, 这将对ZrO2在核工业方面的应用产生极大影响, 例如惰性基质燃料相变引起的体积和应力变化将可能损坏燃料棒结构, 损伤核反应装置, 带来核泄漏风险[1,2,6]. 因此, 研究ZrO2在高能粒子辐照下的相变行为具有重要工程意义.

      近20年来, 利用加速器技术对ZrO2在不同种类和能量的快速重离子作用下的抗辐照性能进行了广泛研究[7-12]. 研究发现, 单斜ZrO2晶体在快速重离子辐照下可以转化为四方相, 并且能够在室温下稳定存在. 当前主要采用热峰模型对快速重离子诱导材料相变现象进行机理解释[13-15], 即快速重离子辐照促进材料局部温度升高, 超过材料的相变温度, 从而诱导相变. 但热峰模型只能定量计算辐照后材料内部的温度分布, 不能直观地反映材料微观结构变化[16]. 分子动力学能够模拟几百皮秒和几十纳米时空范围内材料微观结构演化, 已广泛应用于材料的低能辐照损伤研究中[17,18]. 然而快速重离子辐照材料时, 能量损失过程由电子非弹性碰撞主导. 分子动力学无法描述电子与电子以及电子与声子之间的相互作用, 因此不能模拟高能粒子对材料的辐照损伤作用.

      本文结合热峰模型和分子动力学两种方法的优势, 对快速重离子辐照下, ZrO2由单斜相转化为四方相的微观过程进行模拟, 并计算诱导ZrO2相变的快速重离子的电子能损阈值. 首先, 利用热峰模型描述快速重离子入射后, ZrO2内部电子激发、电子-声子能量耦合、电子-电子以及声子-声子能量传递过程, 计算得出晶格温度; 然后, 根据晶格温度标定体系原子动能, 将其作为分子动力学模拟初始条件, 研究ZrO2微观结构演化, 并通过X射线衍射(X-ray diffraction, XRD)计算确定ZrO2相变类型; 随后, 根据Zr原子配位数变化具体分析ZrO2相变过程; 最后, 改变入射离子的电子能损, 进行多次热峰计算和分子动力学模拟, 确定诱导ZrO2发生相变的快速重离子的电子能损阈值.

    • 热峰模型认为快速重离子穿过靶材料时, 首先以离子运动路径为轴线形成圆柱形的电离核心区域, 称为径迹. 在大约10–16至10–15 s时间内将能量传递给靶材电子系统; 随后, 径迹中心电子通过热传导将部分能量传递给径迹周围电子; 同时, 通过电子-声子耦合作用将部分能量传递给原子; 径迹中心原子获得能量后, 通过晶格振动将能量传递给径迹周围原子. 热峰模型的热过程可以通过以下两个耦合微分方程描述[19-21], 分别对应电子和原子系统.

      ${C_{\rm{e}}}({T_{\rm{e}}})\frac{{\partial {T_{\rm{e}}}}}{{\partial t}} \!=\! \frac{1}{r}\frac{{\partial \left[ {r{K_{\rm{e}}}\left( {{T_{\rm{e}}}} \right)\dfrac{{\partial {T_{\rm{e}}}}}{{\partial r}}} \right]}}{{\partial r}} - g\left( {{T_{\rm{e}}} \!-\! {T_{\rm{a}}}} \right) + A\left( {r,t} \right),$

      ${C_{\rm{a}}}({T_{\rm{a}}})\rho \frac{{\partial {T_{\rm{a}}}}}{{\partial t}} \!=\! \frac{1}{r}\frac{{\partial \left[ {r{K_{\rm{a}}}\left( {{T_{\rm{a}}}} \right)\dfrac{{\partial {T_{\rm{a}}}}}{{\partial r}}} \right]}}{{\partial r}} \!+\! g\left( {{T_{\rm{e}}} - {T_{\rm{a}}}} \right).$

      式中, TeTa, CeCa, KeKa分别为电子和原子的温度、比热容和热导率, Ce取1 J·cm–3·K–1, Ke取2 J·cm–1·s–1·K–1[22], Ca[23,24]Ka[24-26]为随温度变化的量, 其表达式分别为

      ${C_{\rm{a}}} \!=\! \left\{ {\begin{aligned} & {\rm{565}}.{{4272 + 61}}.{\rm{164}} \times {\rm{1}}{{\rm{0}}^{ - 3}}{T_{\rm{a}}} \\ & - 114.1728 \times {{10}^5}{T_{\rm{a}}}^{ - 2},&{{T_{\rm{a}}} < 2715{{,}}} \\ &{800}, &{{T_{\rm{a}}} \geqslant 2715{{,}}} \end{aligned}} \right.$

      ${K_{\rm{a}}} = \frac{{2400}}{{{T_{\rm{a}}}}}\left( {\frac{2}{3}\sqrt {\frac{{420}}{{{T_{\rm{a}}}}}} + \frac{{{T_{\rm{a}}}}}{{3 \times 420}}} \right).$

      ρ为晶体质量密度, 取5.68 g·cm–3; g为电子-声子耦合系数, 由电子平均自由程λ决定, g = Ke/λ, 对于ZrO2, λ = 4.5 nm[14]; A(r, t)为入射离子在半径r, 时间t时沉积在电子系统的能量密度, 其表达式为

      $A\left( {r,t} \right) = AD\left( r \right)\alpha {{\rm{e}}^{ - \alpha t}}.$

      式中, D(r)为初始沉积能量空间分布, 指数因子α = 1/τ, τ为电子能量沉积时间, 取10–15 s, 可调参数A表达式为

      $\int_{t = 0}^\infty {\int_{r = 0}^{r = {r_{\rm{m}}}} {A\left( {r,t} \right)2{\rm{\pi}} r{\rm{d}}r{\rm{d}}t = {S_{\rm{e}}}} } .$

      式中, rm为电离径迹半径, Se为入射离子的电子能损值, 对于给定能量和种类的入射离子, 其电子能损值可由SRIM软件包[27]求得.

      数值求解热峰方程, 获取ZrO2原子温度的时空分布.

    • 建立如图1所示的单斜ZrO2原胞模型, 将其扩展为40a0 × 10b0 × 40c0的超晶胞模型, 模型尺寸为20.5816 nm × 5.2075 nm × 21.2428 nm, 共包含1.92 × 106个原子. 其中a0 = 5.1454 Å, b0 = 5.2075 Å, c0 = 5.3107 Å, 为单斜ZrO2的晶格常数.

      图  1  单斜ZrO2原胞结构模型

      Figure 1.  The crystal structure of monoclinic ZrO2.

      计算采用LAMMPS软件[28], ZrO2原子间相互作用力采用Buckingham势函数[29,30]描述, 势能表达式为

      ${U_{ij}} = \frac{{{q_i}{q_j}}}{{{r_{ij}}}} + A\exp \left( {\frac{{ - {r_{ij}}}}{\rho }} \right) - \frac{C}{{r_{ij}^6}}.$

      式中, 第一项为长程库伦势, 后两项为短程Buckingham势; qiqj分别为原子ij的电荷, rij为原子间距, A, ρ, C为可调参数. 在快速重离子辐照下, 热峰中心区域原子处于非平衡态, 原子间距小于平衡距离, 相互作用力为较强的斥力. Buckingham势能够准确模拟平衡态下ZrO2原子间相互作用, 然而, 由其描述的原子间近核作用力为相互吸引力, 与实际不符. 因此, 采用普适的Ziegler-Biersack-Littmark (ZBL)库伦屏蔽势[27]对Buckingham势的近核力进行修正, 并用Fermi函数将两种势函数平滑连接[31]:

      ${U_{ij}} = {f_{\rm{F}}}\left( {{r_{ij}}} \right) \cdot {U_{{\rm{Buckingham}}}} + \left( {1 - {f_{\rm{F}}}\left( {{r_{ij}}} \right)} \right) \cdot {U_{{\rm{ZBL}}}},$

      ${f_{\rm{F}}}\left( {{r_{ij}}} \right){{ = }}{1}/ [{{1{{ + }}{{\rm{e}}^{ - {A_{\rm{F}}}\left( {{r_{ij}} - {r_{\rm{c}}}} \right)}}}}].$

      式中, rij为原子间距; rc, AF为可调参数, rc决定两种势函数间连接点的选取, AF控制连接处的平滑程度.

      模拟过程中采用三维周期性边界条件, 截断半径取1.6 nm. 采用自适应时间步长[32], 即根据原子速度和受力情况, 实时更新积分步长, 限制每个原子在单个步长内最大移动距离为0.5 Å, 最大动能变化为100 eV, 并设置最大允许步长为1 fs.

      热峰注入前, 采用NVT系综进行10 ps弛豫, 使用Nose-Hoover恒温器进行300 K热浴. 随后, 如图2所示, 将热峰沿[010]方向注入, 根据热峰模型计算得出的晶格温度标定原子动能. 采用NVE系综对弛豫阶段进行模拟, 并将XZ方向3 Å厚边界层与300 K Berendsen恒温器连接, 模拟能量通过热传导耗散到环境的过程.

      图  2  热峰注入靶材料示意图(Se = 30 keV·nm–1)

      Figure 2.  Diagram of thermal spike injection into target material (Se = 30 keV·nm–1).

    • 图3展示了不同位置和不同时刻的晶格温度变化趋势, 入射离子的电子能损为30 keV·nm–1, 为实验值(10—50 keV·nm–1)的中位数.

      图  3  晶格温度分布曲线 (a) 不同位置晶格温度随时间演化; (b) 不同时刻晶格温度的径向分布

      Figure 3.  Distributions of lattice temperature: (a) Evolution of lattice temperature at different radial distances; (b) radial distribution of lattice temperature at different times.

      图3可知, 在电子能损为30 keV·nm–1的快速重离子辐照下, ZrO2径迹周围晶格温度迅速升高, 60 fs时各处温度基本已达最大值, 径迹中心温度最高接近12000 K, 远超ZrO2熔点(2715 ℃). 径迹中心温度达到峰值后迅速降低, 距径迹中心较远的晶格温度达到峰值后, 先维持一段时间再缓慢下降, 径迹内各处温差随时间减小, 60 ps时各处温度已接近环境温度300 K. 热峰作用下, 距径迹中心4 nm内晶格温度会超过熔点(2715 ℃), 6 nm内晶格温度会超过单斜相至四方相的相变温度(1170 ℃).

    • 由热峰模型得出的晶格温度径向分布, 反映了快速重离子辐照后, 沉积在不同位置处的原子动能. 根据计算结果, 60 fs时径迹内各处温度基本已达最大值, 故取此时刻晶格温度, 标定不同位置处原子动能. 经分子动力学模拟, 得出不同时刻ZrO2微观结构变化如图4所示.

      图  4  不同时刻ZrO2原子结构变化. 所有原子沿热峰注入方向投影在(010)平面上

      Figure 4.  Atomic structure of ZrO2 at different times. All the atoms are projected on (010) plane along the injection direction of thermal spike.

      图4可知, 0 ps时体系原子处于热平衡态. 热峰注入后, 中心区域原子能量较高, 受热激发离开初始晶格位置, 并与周围原子发生剧烈碰撞. 0.15 ps时, 产生一个熔融的柱形径迹, 距离热峰中心3 nm内原子无序程度较高, 3至7 nm内原子无序程度较低. 2 ps时径迹中心开始形成空洞, 空洞外围为高度无序的非晶区, 非晶区周边原子熔融重排, 形成了两种不同的晶区(晶区一和晶区二). 在(010)投影平面上, 晶区二内Zr, O原子空间分布均匀, 且Zr原子位于临近4个O原子的对称中心, O原子位于临近4个Zr原子的对称中心; 晶区一内Zr与Zr原子间距较远, O与O原子间距较近, 且O原子不处于临近Zr原子的对称中心位置. 由4, 6和100 ps的微观结构可知, 随着能量由径迹中心向四周传递, 重结晶过程逐渐发展至整个体系, 空洞、晶区一和晶区二的尺寸逐渐增大, 100 ps时体系已达平衡态.

      XRD图谱能够反映材料成分、形态以及内部原子排列结构等信息, 广泛应用于物质的晶体结构表征. 取100 ps时晶区二内部分原子坐标, 进行分子动力学计算, 获得辐照后ZrO2晶体XRD图谱[33], 并与辐照前XRD图谱以及单斜(m-ZrO2)、四方(t-ZrO2)标准谱进行对比, 结果如图5所示.

      图  5  辐照前后以及标准ZrO2 XRD图谱

      Figure 5.  XRD patterns of pristine, ion-irradiated and standard ZrO2.

      图5可知, 辐照前XRD图谱在28.2°, 31.5°, 34.2°, 35.3°, 49.3°和50.1°出现衍射峰, 分别对应m-ZrO2的($ {\overline{1}} 11$), (111), (002), (200), (022), (220)晶面衍射(JCPDS. No 83-0936), 表明辐照前ZrO2为单斜相. 经快速重离子辐照后, 单斜相特有的特征峰消失, 并且在30.3°, 34.6°, 50.3°, 50.8°, 59.4°, 60.3°出现新的特征峰, 与t-ZrO2的(101), (002), (112), (200), (103), (211)衍射晶面良好对应(JCPDS. No 88-1007), 证明在快速重离子辐照作用下, ZrO2由单斜相转为了四方相.

      ZrO2晶体中, 单斜相的Zr原子配位数为7, 四方相的Zr原子配位数为8, 因此可根据Zr原子配位数的变化情况研究辐照条件下ZrO2相变过程. 图6展示了不同配位数的Zr原子空间分布随时间演化过程.

      图  6  不同时刻不同配位数Zr原子空间分布. 所有原子沿离子入射方向投影在(010)平面上, 颜色不同代表配位数不同

      Figure 6.  Spatial distribution of Zr atoms with different coordination number at various times. All the atoms are projected on (010) plane along the injection direction of thermal spike. Different coordination numbers are distinguished by the atomic color.

      图6所示, 初始时刻所有Zr原子配位数均为7, 表明ZrO2为单斜相. 0.15 ps时, 距热峰中心7 nm范围内形成一个熔融的柱形径迹, 其中28%的Zr原子配位数降为6, 4.4%的Zr原子配位数降为5, 距径迹中心3 nm范围内, 47.8%的Zr原子配位数降为6, 14.8%的Zr原子配位数降为5, 0.7%的Zr原子配位数降为4, 表明越靠近径迹中心, 原子无序程度越高. 随着能量向四周传递, 径迹中心温度逐渐降低, 原子发生结晶重排, 2 ps时, 径迹中心出现空洞, 其周边由配位数为3或4的Zr原子组成非晶区, 非晶区两侧生长出两种不同的晶相, 分别为配位数为8的四方相以及配位数为6的低配位相. 由4, 6和100 ps微观结构可知, 随着能量继续由中心向四周传递, 四方相以及低配位相逐渐发展至整个体系, 同时空洞范围不断扩大. 4和6 ps时, 四方相区域内仍有部分Zr原子配位数为7, 100 ps时, 四方相区域内所有Zr原子配位数均为8, 表明随着时间推移, 相变区域变得更加致密有序. 100 ps时, 空洞周边部分低配位相逐渐转变为四方相, 将两侧相变区连接起来. 本文模拟所得到的空洞, 在实验中尚未得到验证. 然而, 单斜ZrO2密度为5.65 g/cm3, 四方ZrO2密度为6.10 g/cm3, 一定质量的ZrO2由单斜相转为四方相后, 体积减小约7.4%. 本文分子动力学模拟中采用的系综为NVE系综, 模拟过程中体系体积固定, 因此, 在NVE系综下进行ZrO2由单斜相至四方相的相变模拟, 必然导致空洞的产生. 此外, 在其他材料(如Ge[31], Fe-Cr[34]合金等)的辐照损伤研究中, 曾观察到空洞的产生.

    • 为获得辐照诱导ZrO2由单斜相转为四方相的快速重离子的电子能损阈值, 在电子能损为15—30 keV·nm–1范围内进行了多次热峰计算和分子动力学模拟. 图7展示了Zr原子配位数随电子能损值变化趋势.

      图  7  配位数随电子能损值变化

      Figure 7.  The change of coordination number with electronic energy loss.

      图7所示, 电子能损值为15—20 keV·nm–1时, 配位数为7的Zr原子个数随电子能损值增加而减少, 未出现配位数为8的Zr原子, 表明ZrO2仍保持单斜相. 当电子能损值为21 keV·nm–1时, 72.4% Zr原子配位数由7变为8, 表明ZrO2已经由单斜相转为四方相. 由此得出辐照诱导单斜ZrO2相变的快速重离子的电子能损阈值为21 keV·nm–1. 电子能损值为21—30 keV·nm-1时, 随着电子能损值增加, 配位数为8的Zr原子个数先增加后减小, 占Zr原子总数比例在71.5%—76%之间.

      当入射离子电子能损分别为20和21 keV·nm–1时, 辐照后不同配位数的Zr原子空间分布如图8所示.

      图  8  不同配位数Zr原子空间分布 (a) 电子能损值为20 keV·nm–1; (b)电子能损值为21 keV·nm–1. 所有原子沿离子入射方向投影在(010)平面上, 颜色不同代表配位数不同

      Figure 8.  Spatial distribution of Zr atoms with different coordination number: (a) Se = 20 keV nm–1; (b) Se = 21 keV nm–1. All the atoms are projected on (010) plane along the injection direction of thermal spike. Different coordination numbers are distinguished by the atomic color.

      图8所示, 电子能损值为20 keV·nm–1时, 辐照未能造成ZrO2由单斜相至四方相的转变. 在热峰作用下, 部分原子离开其正常晶格位置, 成为点缺陷, 点缺陷的引入使得部分Zr原子配位数降为6或5, 点缺陷在空间上均匀分布. 电子能损值低于20 keV·nm–1时, 不同配位数Zr原子的空间分布图与图8(a)相似. 当快速重离子电子能损值为21 keV·nm–1时, 72.4% Zr原子配位数由7变为8, 并且在空间上形成连续致密的排列结构, 表明ZrO2已经由单斜相转为四方相. 当电子能损值低于27 keV·nm–1时, 低配位相和四方相中心未形成连续空洞, 与图8(b)相似, 这是由于热峰温度较低, 中心原子动能较小, 不足以引起足够剧烈的碰撞使得热峰中心所有原子运动到足够远的区域, 随着系统快速淬火, 原子迅速冷却, 形成非连续的空洞. 当电子能损值高于27 keV·nm–1时, 低配位相和四方相中心形成连续空洞, 不同配位数Zr原子空间分布与图6(f)相似.

      根据热峰模型[13,35], 快速重离子辐照后, 材料相变区域面积σ与入射离子的电子能损Se存在如下关系:

      $ {\sigma \left( {{S_{\rm{e}}}} \right) = {\rm{\pi }}a_0^2\ln \left( {{S_{\rm{e}}}/{S_{{\rm{et}}}}} \right)},\;\;\;{{S_{{\rm{et}}}} \leqslant {S_{\rm{e}}} \leqslant 2.7{S_{{\rm{et}}}}} ,$

      $ {\sigma \left( {{S_{\rm{e}}}} \right) = {\rm{\pi }}a_0^2{S_{\rm{e}}}/2.7{S_{{\rm{et}}}}},\;\;\;{{S_{\rm{e}}} > 2.7{S_{{\rm{et}}}}} .$

      式中, a0为热峰温度的初始分布半径; Set为入射离子的电子能损阈值.

      经热峰模型计算和分子动力学模拟得出, 辐照致单斜ZrO2相变的快速重离子的电子能损阈值为21 keV·nm–1, 与实验值12 keV·nm–1存在一定误差[7-10], 这主要与辐照剂量以及势函数的准确性有关. 模拟条件为单离子入射, 对应的辐照剂量为2.37 × 1011 ions·cm–2, 而实验中通常辐照剂量较大, 为1012—1015 ions·cm–2, 可能存在同一区域被多个离子辐照的情况. 此外, 分子动力学模拟中采用普适的Ziegler-Biersack-Littmark (ZBL)库伦屏蔽势修正的Buckingham势描述ZrO2原子间相互作用力, 其与真实的原子间作用力相比存在一定差异, 这将导致模拟结果与实验结果存在一定偏差. 未来可开发更为准确的势函数, 开展多个离子同时入射以及多个离子先后入射的研究, 模拟高剂量快速重离子辐照对材料的损伤作用.

    • 本文结合了热峰模型和分子动力学的优势, 模拟了快速重离子辐照下, ZrO2由单斜相转为四方相的微观物理过程, 并计算了辐照诱导ZrO2相变的快速重离子的电子能损阈值. 首先, 基于热峰模型, 考虑快速重离子辐照单斜ZrO2时, 电子激发、电子-声子能量耦合以及电子-电子、声子-声子能量传导多物理过程, 建立热扩散方程, 求得ZrO2晶格温度时空演变特性. 然后, 运用分子动力学方法模拟了热峰作用下ZrO2相变的微观物理过程. 热峰注入后, ZrO2中心形成一个半径为7 nm的径迹, 径迹中心晶格迅速熔融; 2 ps时熔融区逐渐冷却重结晶, 其中心开始形成空洞, 空洞外围为高度无序的非晶区, 非晶区外围生成四方相和低配位相; 随着能量由径迹中心向四周传递, 四方相和低配位相逐渐发展至整个体系, 同时空洞尺寸随之增大. 经热峰计算和分子动力学模拟, 辐照诱导单斜ZrO2相变的快速重离子的电子能损阈值为21 keV·nm–1, 高于实验值12 keV·nm-1, 主要原因为: 1) 入射离子剂量的模拟值小于实验值; 2) 分子动力学中采用势函数描述的原子间相互作用力与真实情况存在一定误差, 未来可开发更精确的势函数描述ZrO2原子间相互作用力, 开展高剂量快速重离子辐照下, 材料辐照损伤的模拟研究. 本文建立的研究方法可有效地模拟高能离子对材料辐照损伤的微观物理过程, 将来可用于其他惰性基质燃料(如YSZ/ZrO2-PuO2等)的高能辐照损伤研究.

参考文献 (35)

目录

    /

    返回文章
    返回