魏煒，余新剛

(中國科學院大學工程科學學院， 北京 100049) (2020年3月18日收稿； 2020年4月20日收修改稿)

近年來，隨著世界范圍內(nèi)的能源危機日趨嚴重，受控核聚變技術(shù)受到越來越多的關(guān)注。經(jīng)過數(shù)十年的發(fā)展，磁約束聚變技術(shù)已經(jīng)取得了突破性的進展，但在科學和技術(shù)上仍然存在許多問題亟待解決，面對等離子體材料的選擇就是其中之一。目前主流的面向等離子體材料包括碳、鈹、鎢等[1]。然而，實驗表明這些材料都在一定程度上存在難以克服的缺陷，包括物理化學濺射、力學性能退化、氚滯留等一系列問題[2-3]。相比固體材料，液態(tài)鋰第一壁現(xiàn)已在很多磁約束核聚變裝置上被證實具有獨特的優(yōu)勢，不僅可以解決大部分固體第一壁所面臨的問題，同時還可以顯著提高芯部等離子體的約束性能[4-5]。

然而，作為具有一定化學活性的堿金屬，液態(tài)鋰對大部分金屬材料都具有強烈的腐蝕作用，其中包括各種不銹鋼材料，而不銹鋼是目前核聚變實驗裝置中最常用的回路結(jié)構(gòu)材料。因此，液態(tài)鋰與鐵基合金的相容性研究就顯得尤為必要。實際上，早在20世紀80年代，人們就已經(jīng)認識到液態(tài)鋰可以作為聚變裝置中的氚增殖材料，因而提出使用液態(tài)鋰作為聚變堆包層中的冷卻流體，并在隨后開展了大量關(guān)于液態(tài)鋰與包層結(jié)構(gòu)材料的腐蝕實驗研究[6-8]。結(jié)果表明液態(tài)鋰對鐵基合金的腐蝕是個復(fù)雜的物理化學過程，其主要形式包括以下3種：1)鋰原子在合金中的滲透擴散；2)合金組分元素在液態(tài)鋰中的溶解；3)液態(tài)鋰以及雜質(zhì)元素與合金組分元素之間發(fā)生化學反應(yīng)生成中間產(chǎn)物。具體腐蝕過程受多種因素的影響，包括合金微觀晶體結(jié)構(gòu)、組分元素，液態(tài)鋰中的雜質(zhì)元素，環(huán)境溫度，以及液態(tài)鋰流速等[6]。顯然，3種形式的腐蝕過程都會改變材料的微觀結(jié)構(gòu)，從而引起材料宏觀力學性能的變化。Popovich等[9]對液態(tài)鋰腐蝕前后的純鐵進行了單向拉伸試驗，結(jié)果表明腐蝕后其拉伸強度和韌性都顯著降低，降低程度與溫度和應(yīng)變率有關(guān)，當溫度為400 ℃，應(yīng)變率為0.002 7/s時，純鐵拉伸強度下降接近15%，斷裂應(yīng)變下降85%，表現(xiàn)出典型的脆化效應(yīng)。Chopra和Smith[10]開展了鐵素體鋼HT-9在液態(tài)鋰中的疲勞實驗，結(jié)果表明腐蝕可弱化材料的晶界，進而導致材料疲勞壽命的顯著降低。Mustari和Takahashi[11]使用由316L不銹鋼制成的波紋管在液態(tài)鋰中進行低周疲勞實驗，在斷裂后的試件內(nèi)同時觀察到晶界裂紋和穿晶裂紋，表明液態(tài)鋰的腐蝕不僅對晶界有削弱作用，同時對晶粒內(nèi)部也有破壞作用。但是，受客觀條件的限制，目前的實驗研究僅能給出實驗現(xiàn)象，尚不能很好地揭示腐蝕和材料損傷演化的動態(tài)過程和微觀機理。

分子動力學(molecular dynamics,MD)方法是從原子尺度上對材料進行動態(tài)數(shù)值模擬的方法，對于揭示材料損傷演化的微觀機理具有顯著優(yōu)勢，目前已被廣泛應(yīng)用于多種材料相互作用的動態(tài)模擬研究。例如Webb等[12]模擬鉛液滴在銅表面的動態(tài)鋪展過程，解釋了不同晶格取向?qū)︿佌顾俣鹊挠绊懸?guī)律。Palafox-Hernandez等[13]通過計算應(yīng)力、勢能、擴散系數(shù)、密度等參數(shù)，證實液態(tài)鉛在不同的銅晶面上具有不同的合金化行為。然而，針對液態(tài)鋰對鐵基合金的腐蝕問題，由于缺乏相應(yīng)的原子間相互作用勢函數(shù)，相關(guān)研究尚不多見。Gan等[14]開發(fā)了基于改進分析型嵌入原子法(modified analytical embedded atom method,MAEAM)理論的鐵鋰勢函數(shù)，使得對液態(tài)鋰腐蝕鐵基合金的大規(guī)模分子動力學模擬成為可能。

基于液態(tài)鋰腐蝕鐵基合金的工程問題，本文借助Gan等[14]發(fā)展的鐵鋰勢函數(shù)，對含鋰的單晶鐵在單向拉伸載荷作用下的動態(tài)響應(yīng)進行大規(guī)模分子動力學模擬，重點考察鋰原子的存在對鐵晶體彈塑性力學行為的影響規(guī)律。

1 模型和方法

藍色和黃色分別表示Fe原子和Li原子

所有的模型都首先在等溫等壓系綜下采用Nose-Hoover熱浴法馳豫50 ps，并確保系統(tǒng)能量已趨于穩(wěn)定。溫度取為室溫300 K，壓力為0，時間步長1 fs。

拉伸載荷的施加通過在每個時間步內(nèi)對所有原子沿拉伸方向(z軸)的坐標根據(jù)設(shè)定的應(yīng)變率進行調(diào)整來實現(xiàn)，并在加載過程中保持其他2個方向的壓力為0，受分子動力學模擬尺度的限制，應(yīng)變率取為5×108/s。

系統(tǒng)的應(yīng)力利用Virial定理[16]進行計算，其具體公式為

在前期工作中，我們完全采用Gan等[14]開發(fā)的鐵鋰勢函數(shù)進行分子動力學模擬，結(jié)果表明該勢函數(shù)并不能很好地描述α-Fe晶體中位錯等缺陷的動態(tài)演化過程。因此，為了獲得更加準確的結(jié)果，改用由Mendelev等[15]開發(fā)的EAM勢函數(shù)來模擬Fe-Fe之間的相互作用，采用Nichol和Ackland[17]開發(fā)的EAM勢函數(shù)模擬Li-Li之間的相互作用，這2個勢函數(shù)已經(jīng)被證實可以很好地處理單晶α-Fe和單晶Li中的缺陷。Fe-Li之間的相互作用則由Gan等[14]的MAEAM勢來描述。

本文使用LAMMPS軟件[18]進行MD模擬，使用OVITO軟件[19]進行原子構(gòu)型可視化，采用公共近鄰分析方法(common neighbor analysis,CNA)[20]識別系統(tǒng)中的缺陷，同時配合位錯提取算法(dislocation extraction algorithm，DXA)[21]進行位錯識別。

2 結(jié)果與討論

2.1 應(yīng)力-應(yīng)變曲線

圖2展示了不同Li原子數(shù)分數(shù)下的典型應(yīng)力-應(yīng)變曲線，可以看到在拉伸的初始階段，應(yīng)力與應(yīng)變呈線性關(guān)系。隨著應(yīng)變的進一步增大，應(yīng)力應(yīng)變關(guān)系表現(xiàn)出明顯的非線性效應(yīng)。這種非線性效應(yīng)是由于在單晶體中不存在位錯、晶界等缺陷，材料的本構(gòu)關(guān)系主要由原子間的相互作用決定，通常都具有非線性特性，類似的現(xiàn)象在實驗[22]和模擬[23]中都有被觀察到。

圖2 不同Li原子數(shù)分數(shù)下的應(yīng)力-應(yīng)變曲線

通過擬合線彈性階段(應(yīng)變范圍[0, 0.03])內(nèi)的應(yīng)力-應(yīng)變曲線，可以得到不同Li原子濃度下的楊氏模量，結(jié)果列于表 1。從表中可以看出，我們計算的純α-Fe的楊氏模量平均值為195.7 GPa，與文獻[24]中的MD模擬結(jié)果187 GPa非常接近，驗證了我們模擬過程的可靠性。至于含Li原子的α-Fe晶體的楊氏模量，目前尚無可靠數(shù)據(jù)可用來比較。我們的計算結(jié)果顯示楊氏模量受Li原子數(shù)分數(shù)的影響并不明顯，最大相差僅為4.1%。從圖2也可以看到，不同Li原子數(shù)分數(shù)下的應(yīng)力-應(yīng)變曲線在彈性階段內(nèi)幾乎重合。然而，屈服應(yīng)力受Li原子數(shù)分數(shù)的影響相對比較明顯，具體分析在2.3節(jié)進行。

表1 不同Li原子數(shù)分數(shù)下<110>晶向α-Fe的楊氏模量

2.2 塑性變形機理

圖3所示為純α-Fe晶體屈服后的原子構(gòu)型，從圖中可以看出，晶體的屈服是由bcc與hcp之間的相變所控制。當應(yīng)力達到臨界值時，球形的hcp團簇開始在晶體內(nèi)成核，如圖3(a)所示。隨著載荷的增加，hcp團簇在晶體內(nèi)部迅速擴展，當hcp團簇的尺寸增大到一定程度，邊界上的應(yīng)力不足以維持hcp結(jié)構(gòu)，隨即hcp團簇發(fā)生破裂，局部的hcp結(jié)構(gòu)重新變回bcc結(jié)構(gòu)，同時伴隨著大量位錯的產(chǎn)生。圖3(b)中顯示了晶體內(nèi)局部的由hcp原子轉(zhuǎn)變而來的二次bcc相。當繼續(xù)增大載荷時，晶體的塑性變形主要由位錯的運動支配，圖3(c)所示為應(yīng)變達到0.15時晶體內(nèi)部的復(fù)雜位錯線分布。這種應(yīng)力驅(qū)動的相變過程與DFT(第一性原理)的結(jié)果[25]保持一致。

(a)和(b)中的原子根據(jù)CNA著色，其中綠色、藍色、天藍色和紅色分別代表bcc、fcc、hcp和未知結(jié)構(gòu)原子；(c)中綠色、紫色和藍色線分別表示1/2<111>、<100>和<110>位錯；為了清晰(a)中的bcc原子被隱藏

Li原子的引入對晶體的塑性變形行為具有顯著的影響。圖4所示為Li原子數(shù)分數(shù)為0.1%的模型在拉伸載荷作用下的變形過程。從圖4(a)中可以看出，置換Li原子使得所在Fe晶格發(fā)生了一定的變形，部分Fe原子在CNA分析中被識別為非bcc原子。圖4(b)和4(c)為塑性變形初期的原子構(gòu)形，可以看出，相變幾乎已經(jīng)被完全壓制，位錯成為塑性變形的主導機制。

(a)和(b)中的黃色原子是放大的Li原子，其他所有原子的著色方案和圖3相同；為了顯示位錯線(c)中原子被縮小，位錯線著色方案和圖3相同

圖5(a)所示為不同Li原子數(shù)分數(shù)模型中的hcp原子數(shù)達到最大占比、位錯首次出現(xiàn)對應(yīng)的應(yīng)變。Hcp原子數(shù)量達到最大占比的時刻，對應(yīng)于hcp團簇破裂的瞬間。從圖中可以看出，隨著Li原子數(shù)分數(shù)的增加，hcp團簇破裂時的應(yīng)變逐漸減小，產(chǎn)生位錯的應(yīng)變也逐漸減小。此外，對于不含Li原子的純α-Fe晶體，位錯的出現(xiàn)與hcp團簇的破裂幾乎同時發(fā)生，而一旦引入Li原子，位錯的出現(xiàn)則明顯先于hcp團簇的破裂。這一現(xiàn)象說明，Li原子的存在促進了位錯的形成，而位錯的滑動一方面釋放了晶體內(nèi)的應(yīng)力，導致應(yīng)力驅(qū)動的bcc到hcp的相變無法持續(xù)發(fā)展；另一方面，位錯與hcp團簇的相互作用，導致團簇的提前破裂。因此，隨著Li原子數(shù)分數(shù)的增加，位錯變得越來越容易產(chǎn)生，hcp團簇的體積也越來越小，如圖5(b)所示。

圖5 不同Li原子數(shù)分數(shù)下hcp相原子數(shù)達到最大占比、位錯首次出現(xiàn)對應(yīng)的應(yīng)變，以及hcp相原子數(shù)最大占比

因此，可以得出結(jié)論，置換Li原子顯著抑制bcc到hcp的相變，隨著Li原子數(shù)分數(shù)的增加，α-Fe晶體的塑性變形逐漸由位錯支配。

2.3 屈服應(yīng)力

2.3.1 Li原子數(shù)分數(shù)的影響

本文中的屈服應(yīng)力取為應(yīng)力-應(yīng)變曲線上彈性階段應(yīng)力的峰值，對應(yīng)著相變或位錯產(chǎn)生時的應(yīng)力。圖6所示為屈服應(yīng)力隨Li原子數(shù)分數(shù)的變化曲線。隨著Li原子數(shù)分數(shù)增加，屈服應(yīng)力呈單調(diào)減小趨勢，即使只存在0.05%濃度的Li原子時，屈服應(yīng)力仍然出現(xiàn)了大幅度的減小，當濃度達到0.3%時減小幅度已經(jīng)達到26.3%。這一結(jié)果與Popovich[9]的實驗結(jié)果一致，即Li原子的存在可以顯著降低α-Fe的強度。

圖中誤差條為標準差，后文相同

Li原子導致屈服應(yīng)力減小的原因可以從兩個方面來解釋。一方面，Li原子的引入增加了缺陷潛在成核位點的數(shù)量，而這被認為是影響屈服應(yīng)力的一個重要因素[26-27]。另一方面，雖然Li的密度比α-Fe小得多，但是實際上Li原子的半徑比Fe原子半徑更大，這種原子尺寸上的不匹配，導致在Li原子周圍產(chǎn)生局部應(yīng)變場，勢能增加，F(xiàn)e-Fe之間的鍵能強度降低，導致位錯更容易發(fā)生，在宏觀上表現(xiàn)為屈服應(yīng)力變小。這種原子尺寸不匹配可以用尺寸因子來衡量，按Hepburn等[28]的定義，溶質(zhì)原子X的尺寸因子為ΩSF=ΔV/Vave，其中ΔV表示用一個溶質(zhì)原子X置換一個溶劑原子所引起的體積變化，Vave表示每個溶劑原子的平均體積。查閱文獻可知Li原子在α-Fe中的尺寸因子ΩSF為7.5%[29]，意味著置換Li原子在α-Fe晶格中實際上起著膨脹作用，這無疑會使Li原子周圍的原子的內(nèi)應(yīng)變能增加，導致只需更少的機械變形就可使它們在熱振動中離開平衡位置，進而引起屈服。

2.3.2 Li原子分布的影響

本文中Li原子在α-Fe晶體中隨機分布，通過仔細觀察原子構(gòu)形可發(fā)現(xiàn)，即使在低Li原子數(shù)分數(shù)的晶體中也存在彼此非常靠近的Li原子。為了考察Li原子的分布規(guī)律對屈服應(yīng)力的影響，我們在引入Li原子時，限制任意2個Li原子的間距不小于3a，這個距離超過了Li原子相互作用的截斷半徑。為敘述簡便，這里將Li原子分布被約束的系統(tǒng)稱為約束隨機分布，而無約束的系統(tǒng)稱為完全隨機分布。

兩種分布的屈服應(yīng)力隨Li原子數(shù)分數(shù)的變化如圖7所示?？梢郧宄乜吹?，當Li原子的分布被約束時，其屈服應(yīng)力相對未被約束時大幅增加，這表明鄰近Li原子比單個Li原子會造成更大的屈服應(yīng)力降低。這其實并不難理解，因為從直觀上說鄰近Li原子會引起更嚴重的局部晶格畸變，而晶格的變形必然引起Fe原子之間鍵能的變化，從而導致位錯更容易產(chǎn)生。此外，約束隨機分布屈服應(yīng)力的Li濃度敏感性也顯著降低，具體表現(xiàn)為隨著Li原子數(shù)分數(shù)增加屈服應(yīng)力減小的速率明顯小于完全隨機分布。對于完全隨機分布，在較低Li濃度下鄰近Li原子出現(xiàn)的概率較低，所以數(shù)量較少。但隨著Li原子數(shù)分數(shù)增加，鄰近Li原子越來越有可能形成，數(shù)量也隨之增加。因此，在完全隨機分布中Li濃度在很大程度上決定了鄰近Li原子的數(shù)量，最終顯著影響屈服應(yīng)力。相反，由于完全不存在彼此非?？拷腖i原子，約束隨機分布的屈服應(yīng)力在任何Li濃度下都完全受單個Li原子控制，因此Li濃度的影響相對較小。

圖7 不同Li原子分布下的屈服應(yīng)力

2.3.3 溫度的影響

已有的研究結(jié)果顯示，許多fcc和bcc金屬晶體的屈服應(yīng)力都表現(xiàn)出顯著的溫度依賴性，一般趨勢為屈服應(yīng)力隨著溫度的升高而減小，但是減小的規(guī)律與材料性質(zhì)有關(guān)，最常見的兩種規(guī)律分別是σy～T1/2和σy～T，其中σy和T分別為屈服應(yīng)力和溫度[30-31]。為了研究溫度對Fe-Li固溶體屈服應(yīng)力的影響，我們在1 ～ 900 K溫度范圍內(nèi)分別對純α-Fe和Li原子數(shù)分數(shù)為0.05%的晶體進行了拉伸模擬。圖8(a)展示了兩種晶體在不同溫度下典型的應(yīng)力-應(yīng)變曲線，可以看到隨著溫度的升高，兩種晶體的彈性模量逐漸降低，這一點與現(xiàn)有的實驗結(jié)果一致[32]。此外，從圖中還可以看出，隨著溫度升高，屈服點明顯前移，晶體屈服發(fā)生得越來越早。不同溫度下的屈服應(yīng)力變化如圖8(b)所示，在1 K時純α-Fe的屈服應(yīng)力平均值為32.2 GPa，與DFT計算得到的結(jié)果32.3 GPa[33]保持一致。隨著溫度的升高，兩種晶體的屈服應(yīng)力都呈單調(diào)減小的趨勢，這可以用簡單的熱力學觀點來解釋：溫度越高，原子動能越大，熱振動幅度越大，這意味著在相同載荷作用下原子更容易離開其平衡位置，產(chǎn)生塑性變形。

然而，兩種晶體屈服應(yīng)力的減小規(guī)律卻表現(xiàn)出了不同之處。如圖8(b)所示，純α-Fe的屈服應(yīng)力擬合曲線為σy=32.9-0.66T1/2，符合前文提到的σy～T1/2規(guī)律，這個結(jié)果也與文獻[34]結(jié)果一致。但是，在含Li晶體中，屈服應(yīng)力隨溫度變化卻符合σy～T規(guī)律，擬合結(jié)果為σy=22.4-0.01T。實際上，材料的屈服同時由原子動能和鍵能決定。當溫度較低時，原子動能較低，屈服應(yīng)力主要受原子間的鍵能影響，Li原子的存在明顯降低了Fe原子之間鍵的強度，因此屈服應(yīng)力下降明顯。隨著溫度的升高，原子的動能逐漸增加，對屈服應(yīng)力的影響也越來越顯著，相對而言，Li原子引起的Fe原子鍵能變化所帶來的影響越來越小。綜合起來，使得α-Fe晶體的屈服應(yīng)力與溫度的關(guān)系轉(zhuǎn)變?yōu)榫€性關(guān)系。此外，從圖8(b)中還可以看出，隨著溫度的升高，兩種晶體的屈服應(yīng)力逐漸趨于一致，證實Li原子的影響越來越小。

圖8 不同溫度下的應(yīng)力-應(yīng)變曲線和屈服應(yīng)力

3 結(jié)論

本文使用分子動力學方法對單晶α-Fe的拉伸力學行為進行了數(shù)值模擬，分析置換Li原子對其彈塑性力學行為的影響，主要結(jié)論如下：

1)塑性變形機理受Li原子的影響顯著，在純α-Fe中塑性變形首先由bcc到hcp的相變引起，隨后由位錯控制，Li原子的存在可以顯著地抑制相變的發(fā)生，塑性變形主要通過位錯的運動實現(xiàn)；

2)屈服應(yīng)力隨Li原子數(shù)分數(shù)增大而單調(diào)減小，并且還會受到Li原子分布的顯著影響，鄰近Li原子相比單個Li原子對屈服應(yīng)力影響更大；

3)純α-Fe和含Li原子的α-Fe晶體的屈服應(yīng)力都隨溫度升高而減小，但是Li原子的存在改變了屈服應(yīng)力與溫度的具體關(guān)系。

在實際的金屬材料中，具有單晶結(jié)構(gòu)的比較少，通常都含有大量的晶界、位錯等缺陷，這些缺陷對于材料的性能具有至關(guān)重要的作用。在未來的工作中，我們將逐步分析金屬材料中的典型缺陷與Li原子的相互作用，以便對液態(tài)Li腐蝕造成的α-Fe力學性能退化有更加全面的了解。