賀朝會(huì), 唐 杜, 李 奎

(西安交通大學(xué) 核科學(xué)與技術(shù)學(xué)院, 西安 710049)

輻射在材料中產(chǎn)生位移損傷是一個(gè)從微觀到宏觀的過(guò)程，一般經(jīng)歷離位級(jí)聯(lián)、熱峰階段、缺陷初步演化和缺陷遷移4個(gè)階段[1]，最終導(dǎo)致材料宏觀性能的改變。實(shí)驗(yàn)中可觀察的是宏觀性能的變化過(guò)程，微觀變化過(guò)程只有通過(guò)數(shù)值模擬才能揭示。位移損傷效應(yīng)的模擬難度很大，因?yàn)樗婕耙粋€(gè)多尺度的變化過(guò)程，材料尺寸從亞納米到米量級(jí)，時(shí)間從亞皮秒到年量級(jí)，如圖1所示，這需要用到多種不同的方法和軟件進(jìn)行多尺度的模擬研究[2]，如圖2所示。

圖1 材料性能多尺度變化示意圖Fig.1 Multi-scale sketch map of material performance change

圖2 多尺度數(shù)值模擬方法示意圖Fig.2 Multi-scale sketch map of numerical simulation methods

1 第一性原理[3-9]

第一性原理方法(first-principle method)也被稱為從頭算方法(ab initio methods)，是指基于量子力學(xué)理論，完全由理論推導(dǎo)而得，不使用除基本物理常數(shù)和原子量以外的實(shí)驗(yàn)數(shù)據(jù)以及經(jīng)驗(yàn)或者半經(jīng)驗(yàn)參數(shù)的求解薛定諤方程的計(jì)算方法?；舅枷胧牵簩⒍嘣訕?gòu)成的體系理解為電子和原子核組成的多粒子系統(tǒng)，根據(jù)原子核和電子互相作用的原理及其基本運(yùn)動(dòng)規(guī)律，運(yùn)用量子力學(xué)基本原理最大限度地對(duì)問(wèn)題進(jìn)行“非經(jīng)驗(yàn)性”處理。

量子力學(xué)描述n粒子體系的波函數(shù)包含3n個(gè)坐標(biāo)，相應(yīng)的薛定諤方程是含3n個(gè)變量的偏微分方程，而針對(duì)絕大部分問(wèn)題涉及的體系，這就會(huì)是一個(gè)非常復(fù)雜甚至無(wú)法解決的問(wèn)題。

密度泛函理論(density-functional theory，DFT)采用電子密度分布n(r)作為基本變量，而非波函數(shù)，研究原子、分子、固體等多電子體系的基態(tài)性質(zhì)。DFT建立在Hohenberg和Kohn提出的非均勻電子氣理論基礎(chǔ)之上，可以歸結(jié)為兩個(gè)基本定理。定理1：不計(jì)自旋的全同費(fèi)米子系統(tǒng)非簡(jiǎn)并基態(tài)的所有性質(zhì)都是粒子密度函數(shù)的唯一泛函。該定理保證了粒子密度作為體系基本物理量的合法性，同時(shí)也是DFT名稱的由來(lái)。定理2：DFT的變分法中，對(duì)于一個(gè)給定的外勢(shì)，真實(shí)電子密度使能量泛函取得最小值。由于DFT中融入了統(tǒng)計(jì)的思想，不必考慮每個(gè)電子的行為，只需計(jì)算總的電子密度，計(jì)算量大為減少，從而大大簡(jiǎn)化了應(yīng)用量子力學(xué)計(jì)算材料性質(zhì)所涉及的復(fù)雜數(shù)學(xué)問(wèn)題。DFT把多電子問(wèn)題成功轉(zhuǎn)化為單電子問(wèn)題，近似求解薛定諤方程，可用于描述原子、分子等微觀粒子體系的量子行為，分析其電子結(jié)構(gòu)，計(jì)算總能量。已廣泛應(yīng)用在材料設(shè)計(jì)、合成、模擬計(jì)算和評(píng)價(jià)等諸多方面，成為計(jì)算凝聚態(tài)物理、計(jì)算材料科學(xué)和量子化學(xué)的重要基礎(chǔ)和核心技術(shù)。

廣義的第一性原理包括兩大類，以Hartree-Fork自洽場(chǎng)計(jì)算為基礎(chǔ)的ab initio從頭算和DFT計(jì)算。也有人主張，ab initio專指從頭算，而第一性原理和所謂量子化學(xué)計(jì)算特指DFT計(jì)算。

第一性原理主要應(yīng)用于1)材料物性計(jì)算，包括晶格常數(shù)、彈性模量、相變、電介常數(shù)、熱膨脹系數(shù)、光學(xué)特性等；2)電子結(jié)構(gòu)計(jì)算、力學(xué)行為預(yù)測(cè)；3)化學(xué)成鍵機(jī)理描述；4)半導(dǎo)體材料能帶計(jì)算，包括缺陷形成能與結(jié)合能、晶格原子的離位閾能、雜質(zhì)或缺陷的遷移能；5)分子、團(tuán)簇、晶體等的結(jié)構(gòu)優(yōu)化及新構(gòu)相預(yù)測(cè)；6)納米器件設(shè)計(jì)；7)極端條件下的材料行為模擬等。例如，CdSxTe1-x三元半導(dǎo)體材料的基態(tài)物性計(jì)算，硅材料中的間隙缺陷結(jié)構(gòu)及擴(kuò)散路徑，硅材料中不同電荷態(tài)空位和間隙原子的形成能，半導(dǎo)體材料AlSb點(diǎn)缺陷遷移率退化等。

2 分子動(dòng)力學(xué)[1-14]

分子動(dòng)力學(xué)(molecular dynamics，MD)是一種確定性方法，按照體系內(nèi)部的內(nèi)稟動(dòng)力學(xué)規(guī)律來(lái)確定位形的轉(zhuǎn)變，跟蹤系統(tǒng)中每個(gè)粒子的個(gè)體運(yùn)動(dòng)。根據(jù)統(tǒng)計(jì)物理規(guī)律，給出分子的坐標(biāo)、速度等微觀量與溫度、壓力、比熱容、彈性模量等宏觀可觀測(cè)量之間的關(guān)系。

對(duì)于多體系問(wèn)題的嚴(yán)格求解，需要建立并求解體系的薛定諤方程，根據(jù)波恩-奧本海默近似，原子核的運(yùn)動(dòng)可以用經(jīng)典動(dòng)力學(xué)方法處理。把薛定諤方程的求解簡(jiǎn)化為牛頓運(yùn)動(dòng)方程的計(jì)算。具體步驟包括：1)建立一個(gè)由n個(gè)粒子(分子)組成的模型體系；2)解n個(gè)粒子(分子)組成的模型體系的牛頓運(yùn)動(dòng)方程直至平衡；3)平衡后，進(jìn)行材料性能的計(jì)算，對(duì)模擬結(jié)果進(jìn)行分析。

MD方法只考慮多體系中原子核的運(yùn)動(dòng)，而不考慮電子的運(yùn)動(dòng)，忽略量子效應(yīng)。適用于經(jīng)典近似，在很多材料體系應(yīng)用中都較精確。不適用于涉及電荷重新分布的化學(xué)反應(yīng)、鍵的形成與斷裂、解離、極化以及金屬離子的化學(xué)鍵等領(lǐng)域。不適用的領(lǐng)域可用量子力學(xué)方法替代，也有結(jié)合經(jīng)典分子動(dòng)力學(xué)和量子力學(xué)的新方法，叫做第一性原理分子動(dòng)力學(xué)。

進(jìn)行MD運(yùn)算有7個(gè)必備步驟：1)建立計(jì)算模型；2)設(shè)定計(jì)算模型的初始坐標(biāo)和初始速度；3)選定合適的時(shí)間步長(zhǎng)；4)選取合適的原子間相互作用勢(shì)函數(shù)，便于進(jìn)行力的計(jì)算；5)選擇合適的算法、邊界條件和外界條件；6)計(jì)算；7)對(duì)計(jì)算數(shù)據(jù)進(jìn)行統(tǒng)計(jì)處理。MD運(yùn)算主要用于研究材料微觀性能。如單晶硅輻照級(jí)聯(lián)碰撞過(guò)程分子動(dòng)力學(xué)模擬。圖3為初級(jí)反沖原子(primary knock-out atom, PKA)產(chǎn)生的缺陷及其隨時(shí)間演化的過(guò)程。圖4為2 keV 硅PKA產(chǎn)生的點(diǎn)缺陷的數(shù)目隨時(shí)間的變化情況。主要經(jīng)歷3個(gè)階段：初始階段——缺陷數(shù)目急劇上升；熱峰階段——缺陷數(shù)目達(dá)到頂峰并開(kāi)始下降；趨衡階段——缺陷數(shù)目和構(gòu)型趨于穩(wěn)定。

3 動(dòng)力學(xué)蒙特卡羅[14-20]

MD方法可以模擬材料中產(chǎn)生位移損傷的離位級(jí)聯(lián)、熱峰階段和缺陷初步演化的前3個(gè)階段，時(shí)間尺度到納秒量級(jí)、尺寸尺度達(dá)亞微米量級(jí)。而對(duì)于第4階段的缺陷遷移的模擬，MD方法則無(wú)能為力，只有依靠動(dòng)力學(xué)蒙特卡羅(dynamics Monte Carlo，KMC)方法。KMC方法是將普通的蒙特卡羅方法中的事件賦予時(shí)間的尺度，模擬系統(tǒng)隨時(shí)間演化的過(guò)程。可以模擬界面、可動(dòng)粒子、擴(kuò)展缺陷、缺陷團(tuán)簇和復(fù)雜團(tuán)簇的遷移、擴(kuò)散、復(fù)合、發(fā)射粒子、捕獲和解離等，時(shí)間尺度達(dá)105s。模擬中所有事件的發(fā)生頻率遵從Arrhenius方程：

vij=v0·exp(-Eij/kBT)

(1)

其中,v0為指前因子；Eij為某粒子從狀態(tài)i轉(zhuǎn)變?yōu)闋顟B(tài)j所需越過(guò)的勢(shì)壘；kB為玻耳茲曼常數(shù)；T為溫度。圖5為KMC模擬241.9 keV硅PKA產(chǎn)生的不同缺陷隨時(shí)間的演化。缺陷退火分3個(gè)階段：Ⅰ.1×10-11s ≤t< 2×10-3s：點(diǎn)缺陷成團(tuán)；Ⅱ.2×10-3s ≤t< 2×102s：缺陷團(tuán)內(nèi)間隙原子和空位發(fā)生復(fù)合反應(yīng)；Ⅲ.t≥ 2×102s：小缺陷團(tuán)緩慢發(fā)射點(diǎn)缺陷。

4 待研究問(wèn)題

雖然已經(jīng)開(kāi)展了一些研究工作，但仍然存在許多亟待研究的問(wèn)題。

1) 基本參數(shù)研究

在應(yīng)用LAMMPS進(jìn)行MD模擬計(jì)算時(shí)發(fā)現(xiàn)一些材料的勢(shì)函數(shù)還存在問(wèn)題，如GaAs，SiGe等；應(yīng)用KMC進(jìn)行缺陷演變規(guī)律研究，模擬的缺陷種類較少，一些類型的缺陷團(tuán)的形成和解離反應(yīng)還無(wú)法模擬，需添加各種缺陷反應(yīng)類型和反應(yīng)參數(shù)。

2) 不同輻射源產(chǎn)生的缺陷演化規(guī)律異同性研究

圖3 PKA產(chǎn)生的缺陷及其隨時(shí)間演化的過(guò)程Fig.3 Defects produced by primary knock-out atom and their evolution with time

圖4 2 keV Si PKA產(chǎn)生點(diǎn)缺陷的數(shù)目隨時(shí)間的變化情況Fig.4 The evolution of the number of point defects produced by 2 keV Si PKA with time

(a) The evolution of interstitials and vacacies with time

(b) The evolution of two,three and four vacacies with time

(c) The evolution of impurity related defects with time 圖5 241.9 keV硅PKA產(chǎn)生的不同缺陷隨時(shí)間的演化Fig.5 The evolution of the different defects produced by 241.9 keV Si PKA with time

針對(duì)不同輻射源產(chǎn)生的位移損傷，已有初步的等效性研究結(jié)果，主要是研究材料的宏觀性能的改變，其深層次的原因有待揭示?；谌毕莸亩喑叨饶M方法，可以研究中子、質(zhì)子、重離子和電子產(chǎn)生的缺陷及其演化規(guī)律的異同和對(duì)Si器件性能影響的差異；從缺陷的產(chǎn)生和演變的角度分析不同輻射源之間的等效性，揭示其內(nèi)在機(jī)理；研究均勻分布缺陷和非均勻分布缺陷以及不同濃度的缺陷演化規(guī)律的異同。

3) 單粒子位移損傷研究

已經(jīng)開(kāi)始了Si器件的單粒子位移損傷研究，可以進(jìn)一步開(kāi)展二元半導(dǎo)體材料，如SiC、GaAs、GaN、SiGe等的單粒子位移損傷研究、光電材料和器件的單粒子位移損傷研究，以及半導(dǎo)體輻射探測(cè)器的位移損傷漏電流研究。

4) 位移損傷缺陷長(zhǎng)時(shí)間演化的溫度效應(yīng)研究

由于不同溫度下缺陷會(huì)有不同的演化行為，例如V+O在室溫下形成VO后是穩(wěn)定的，不會(huì)進(jìn)一步形成V2O和VO2缺陷，但當(dāng)溫度高于300 ℃時(shí)，可能形成V2O和VO2兩種缺陷。因此，可以利用MD模擬PKA在不同溫度下的缺陷初態(tài)分布；并將不同溫度下的缺陷初態(tài)分布引入KMC模擬中，進(jìn)一步研究不同溫度下各種缺陷的長(zhǎng)時(shí)間演化行為。

5) 軟件研究

要開(kāi)展多尺度的模擬研究，不同軟件之間的銜接是個(gè)問(wèn)題，目前僅實(shí)現(xiàn)了MD結(jié)果到KMC輸入的轉(zhuǎn)換，如圖6所示，而其他軟件之間的銜接難度和工作量都很大，且MD結(jié)果到KMC輸入的轉(zhuǎn)換目前是通過(guò)獨(dú)立編程的軟件實(shí)現(xiàn)，還沒(méi)有實(shí)現(xiàn)自動(dòng)化。此外，軟件的國(guó)產(chǎn)化更是一個(gè)長(zhǎng)期、艱巨的任務(wù)，可以從編制用于二元半導(dǎo)體材料如GaAs、GaN、SiC、SiGe等的KMC軟件入手，首先實(shí)現(xiàn)KMC軟件的國(guó)產(chǎn)化。

實(shí)際要研究的問(wèn)題還有很多，需要進(jìn)一步細(xì)化分析。

(a) The introduction of the number and spatial distributionof defects from MD to KMC

(b) The evolution of the number of defects with time 圖6 MD與KMC的銜接Fig.6 The connection between MD and KMC

5 結(jié)論

應(yīng)用第一性原理、分子動(dòng)力學(xué)和動(dòng)力學(xué)蒙特卡羅等方法，可以跨尺度，尺寸從亞納米到米量級(jí)，時(shí)間從亞皮秒到年量級(jí)，模擬材料中的位移損傷效應(yīng)，揭示輻照效應(yīng)機(jī)理。多尺度數(shù)值模擬方法為材料輻照效應(yīng)研究提供了有力的工具，但仍然存在許多問(wèn)題亟待研究。

[1]NORDLUND K, DJURABEKOVA F. Multiscale modelling of irradiation in nanostructures[J]. Journal of Computational Electronics, 2014, 13 (1): 122-141.

[2]DE LA RUBIA T D, CATURLA M J, ALONSO E A, et al. The primary damage state and its evolution over multiple length and time scales: Recent atomic-scale computer simulation studies[J]. Radiation Effects and Defects in Solids, 1999, 148 (1): 95-126.

[3]KRESSE G, HAFNER J. First-principles study of the adsorption of atomic H on Ni (111), (100) and (110)[J]. Surf Sci, 2000, 459(3): 287-302．

[4]JI Y F, WANG B, LUO Y. First principles study of O2adsorption on reduced rutile TiO2-(110) surface under UV illumination and its role on CO oxidation[J]. J Phys Chem C, 2013, 117(2): 956-961.

[5]SEBBARI K, ROQUES J, SIMONI E, et al. First-principles molecular dynamics simulations of uranyl ion interaction at the water/rutile TiO2(110) interface[J]. Surf Sci, 2012, 606 (15/16): 1 135-1 141.

[6]BENDAVID L I, CARTER E A. First-principles predictions of the structure, stability and photocatalytic potential of Cu2O surfaces[J]. J Phys Chem B, 2013, 117 (49): 15 750-15 760.

[7]ISLAM M M, DIAWARA B, MAURICE V, et al. First principles investigation on the stabilization mechanisms of the polar copper terminated Cu2O(111) surface[J]. Surf Sci, 2009, 603 (13): 2 087-2 095.

[8]LI K, ZHAO Y L, DENG J, et al. Adsorption of radioiodine on Cu2O surfaces: A first-principles density functional study[J]. Acta Phys Chim Sin, 2016, 32 (9): 2 264-2 270.

[9] 李奎. 碘分子、碘陰離子在Cu2O和CuFeS2表面吸附的機(jī)理研究[D]. 西安: 西安交通大學(xué), 2016. (LI Kui. Study on adsorption mechanisms of iodine anions on Cu2O and CuFeS2surfaces[D]. Xi’an: Xi’an Jiaotong University, 2016.)

[10] PLIMPTON S. Fast parallel algorithms for short-range molecular-dynamics[J]. Journal of Computational Physics, 1995, 117 (1): 1-19.

[11] BORODIN V A. Molecular dynamics simulation of annealing of post-ballistic cascade remnants in silicon[J]. Nucl Instrum Methods Phys Res B, 2012, 282: 33-37.

[12] CATURLA M J, DE LA RUBIA T D, GILMER G H. Disordering and defect production in silicon by keV ion irradiation studied by molecular dynamics[J]. Nucl Instrum Methods Phys Res B, 1995, 106 (1/4): 1-8.

[13] SANTOS I, MARQUéS LA, PELAZ L, et al. Temperature effect on damage generation mechanisms during ion implantation in Si: A classical molecular dynamics study[J]. AIP Conf Proc, 2012, 1 496 (1): 229-232.

[14] 唐杜. 硅基器件的單粒子翻轉(zhuǎn)及單粒子位移損傷的數(shù)值模擬研究[D]. 西安：西安交通大學(xué), 2016. (TANG Du. Numerical simulations on single event upset and single particle displacement damage in silicon devices[D]. Xi’an: Xi’an Jiaotong University, 2016.)

[15] FICHTHORN K A, WEINBERG W H. Theoretical foundations of dynamical Monte Carlo simulations[J]. J Chem Phys, 1991, 95 (2): 734-750.

[16] MARTIN-BRAGADO I, ZOGRAPHOS N. Indirect boron diffusion in amorphous silicon modeled by kinetic Monte Carlo[J]. Solid-State Electronics, 2011, 55 (1): 25-28.

[17] MARTIN-BRAGADO I, TIAN S, JOHNSON M, et al. Modeling charged defects, dopant diffusion and activation mechanisms for TCAD simulations using kinetic Monte Carlo[J]. Nucl Instrum Methods Phys Res B, 2006, 253 (1/2): 63-67.

[18] FISICARO G, PELAZ L, ABOY M, et al. Dopant dynamics and defects evolution in implanted silicon under laser irradiations: A coupled continuum and kinetic Monte Carlo approach[C]// International Conference on Simulation of Semiconductor Processes and Devices (SISPAD). Glasgow, Scotland, 2013: 33-36.

[19] TANG D, MARTIN-BRAGADO I, HE C H, et al. Time dependent modeling of single particle displacement damage in silicon devices[J]. Microelectronics Reliability, 2016, 60: 25-32.

[20] 唐杜, 賀朝會(huì), 臧航, 等. 硅單粒子位移損傷多尺度模擬研究[J]. 物理學(xué)報(bào), 2016, 65(8): 084209. (TANG Du, HE Chao-hui, ZANG Hang, et al. Multi-scale simulations of single particle displacement damage in silicon[J]. Acta Phys Sin, 2016, 65(8): 084209.)