蘇 文， 李盈盈

(1.武漢大學 動力與機械學院， 武漢 430072； 2.國核(北京)科學技術研究院， 北京 100029)

?

第一性原理研究C/Si比對SiC理想剪切強度的影響

蘇 文1*， 李盈盈2

(1.武漢大學 動力與機械學院， 武漢 430072； 2.國核(北京)科學技術研究院， 北京 100029)

使用第一性原理計算方法研究C/Si比對β-SiC的理想剪切強度的影響，建立化學計量比下理想β-SiC分子模型，用C原子替換其中的Si原子，經蒙特卡羅方法優(yōu)化得到非化學計量比下的SiC分子模型，這些晶體的C/Si比的變化范圍是1～1.52.每個模型樣品都沿[100]方向剪切，得到其剪切應力應變關系，最終了解C/Si比對SiC理想剪切強度的影響.計算結果表明，隨著C/Si比增加，SiC的理想強度總體上呈下降趨勢.在剪切載荷逐漸增大情況下，晶體會在C團簇周圍的C-C鍵或者C-Si鍵處斷裂.C原子的聚集，會使SiC的剪切強度減少，因為C團簇內部存在應變，而且C團簇越大，SiC的剪切強度下降得越劇烈.

理想剪切強度； 第一性原理； 蒙特卡羅模擬； SiC； C團簇； 非化學計量比

SiC由于擁有眾多的優(yōu)點，如熱導系數高、熱膨脹系數低、抗熱沖擊能力強、高溫下抗腐蝕性能好、中子吸收截面低以及高溫下能保留裂變氣體等，而被認為是未來理想的核材料[1].然而純SiC也有其致命的缺點，如斷裂韌性差、對微觀結構中的缺陷非常敏感、脆性及可靠性低，這些缺點是其應用在實際工程之中的最大阻礙.為了克服純SiC的缺點，科學家們將純SiC通過連續(xù)SiC纖維的強化，形成了SiC/SiC陶瓷基復合材料.這種SiC/SiC復合材料與純SiC相比，強度、硬度、抗損傷能力都有所增強，并且對于缺陷和刻痕毫不敏感，這使得它能夠適用于惡劣的裂變反應堆堆芯環(huán)境[2].同時SiC/SiC復合材料的誘導放射性低，固有活性低、衰變熱低及氚滲透性低，使得其也能應用于聚變系統中[3].種種的優(yōu)點都使得SiC/SiC復合材料在第四代反應堆和聚變堆設計中都備受矚目.

經過幾十年的發(fā)展，SiC/SiC復合材料已經形成了具有特色的三代，每一代性能的增強都是通過改進SiC纖維的性能來達到的，因為SiC纖維對復合材料的抗輻照性、高溫強度和穩(wěn)定性都至關重要[4]，因此更深入了解SiC的性質對于制造出高性能的SiC/SiC復合材料是必不可缺的.科研人員們?yōu)榱颂岣逽iC的抗輻照性和高溫下的性能，這些年來研發(fā)SiC纖維的趨勢是尋找靠近化學計量比、氧含量低、結晶度高的SiC纖維[5-7]，因為這樣的SiC纖維在中子輻照下其機械性能不會降低[6].但是，有些在化學計量比附近、氧含量低的SiC纖維(如Hi-Nicalon Type-S)較其他類型的SiC纖維(如Hi-Nicalon)，強度卻有所下降[4].換而言之，C/Si比對SiC的機械性質的影響我們還不清楚.

理想脆性斷裂定義為原子鍵的斷裂，因此理想強度是材料在斷裂臨界點處所受的應力值[8].理想強度的計算能夠讓我們從原子成鍵本質上了解材料的機械性能，因此，多年來不少科研工作者用計算的方法研究了SiC的機械性質.Lambrecht 等采用Orowan方程的方法研究了SiC的楊氏模量和剪切模量等性質，這些彈性常量的計算值與實驗值比較接近[9].王自強等用第一性原理計算方法研究了β-SiC的彈性、穩(wěn)定性、理想強度等性質[10].Yoshitaka等人也用第一性原理的方法研究了SiC的理想拉伸和剪切強度[11].雖然世界使用計算方法研究理想SiC性質的科研人員不在少數，但是大部分都是研究化學計量比下SiC的性質，并沒有過多關注非化學計量比的SiC得機械性能，因此對C/Si比對SiC性能的影響也知之甚少.

李盈盈等[12]用第一性原理的方法研究了C/Si比對SiC理想拉伸強度的影響，得到C/Si比逐步增大，SiC的理想拉伸強度總體呈下降趨勢.本文的研究內容，是在李的研究[12]上更進一步，使用第一性原理軟件研究C/Si對SiC[100]方向上的理想剪切強度的影響.通過蒙特卡羅模擬的方法產生不同C/Si的非化學計量的SiC模型，再用第一性原理計算方法計算每個模型的應力應變曲線，最后了討論C/Si比對SiCl理想剪切強度的影響.

1 計算方法

為了研究C/Si比對SiC理想強度的影響，首先需要建立不同C/Si比的SiC模型.相對于形成空位或者間隙原子這樣的缺陷，C原子取代Si原子的形成能低得多，而Si原子取代C原子的形成能遠高于其他缺陷[13]，因此我們的不同C/Si的SiC模型的建立是用C原子取代SiC晶體中的Si原子，經蒙特卡羅模擬優(yōu)化形成的.由于工業(yè)實際用到的SiC中C/Si比均是大于1的[4]，因此我們構建了比例從1到1.52的10種比例的模型.計算采用3×2×2的超晶胞，每個單晶胞有4個C原子和4個Si原子，總共96個原子.用數量為1～10的C原子取代晶胞中的Si原子，構建10組非化學計量比的SiC晶胞，這些晶胞的C/Si比分別是1.04、1.09、1.13、1.18、1.23、1.29、1.34、1.4、1.46、1.52.每組根據取代的C原子所在晶胞中位置的不同，取5個不同的結構作為計算樣品.

C原子取代Si原子的初始結構中，C原子只是簡單地占據了Si原子的位置，此時的C-C鍵和C-Si鍵的長度相等，這并不是一個穩(wěn)定的狀態(tài).因此，我們使用蒙特卡洛方法來優(yōu)化C原子在結構中的位置[14]，使系統處于總能相對降低的亞穩(wěn)定狀態(tài).系統總能是用大規(guī)模分子并行計算模擬器(LAMMPS)來計算的，蒙卡優(yōu)化結構的過程如下：首先，隨機地選擇結構中目標數量Si原子，在選擇的Si原子的位置上放上等量的C原子取代原來的Si原子；然后，隨機選擇相鄰的兩個原子，如果這兩個原子的類型不同，一個C原子，一個Si原子，交換兩個原子的位置后計算系統總能，如果交換位置之后總能下降，保留新結構，如果能量上升，保留原結構；如上重復1 000步以上，收斂的標準是平均每個原子能量小于10 meV，將最后的晶體結構用于剪切應力的計算，每個C/Si比比選擇5個不同結構的計算樣品以上.

圖1所示的是C/Si比為1.09的晶胞結構.這個晶胞中有46個Si原子和50個C原子，黃顏色的大分子為硅原子，灰顏色的小分子為C原子，C原子數量為50,其中，A和B兩個碳原子取代了硅原子的位置,晶體沿[100]方向剪切(箭頭方向).

圖 1 C/Si比為1.09的SiC超晶胞Fig.1 A SiC super cell with 50 C atoms and 46 Si atoms

系統總能的計算使用的是基于密度泛函理論(DFT)計算方法的第一性原理計算軟件VASP[15-16]，核心區(qū)域附近的波函數是用投影增強波法[17]處理的.晶體波函數用平面波基展開，平面波動能的截止能是520 eV.對于所有C/Si比不同的樣品，均選取(1×1×1)均勻網格的Monkhorst-Pack 類型的k點進行布里淵區(qū)的積分.

剪切應力應變曲線的計算過程中，應變過程是連續(xù)的，每一步初始結構均是在上一步弛豫完成后的結構上增加應變后得到的，每步增加的剪切應變?yōu)?.01.弛豫過程中，使用共軛梯度法最小化Hellmann-Feynman力，剪切[100]方向上的應變是固定的，使用準牛頓法使所有原子在其它方向上完全弛豫后，直到作用在所有原子上的力均小于0.05 GPa.從VASP 輸出結果中， 可以直接得到弛豫后的原子位置、總能量、應力、應變等數據，剪切試驗一直進行到有原子鍵斷裂為止.

2 計算結果

2.1 化學計量比下SiC計算結果

我們首先計算了固體3C-SiC的平衡晶格常數、剪切模量、理想剪切強度，并將與實驗數據和計算結果做比較，比較結果見表1.由于我們計算所用的晶體是沒有任何缺陷的，而實際材料總存在缺陷，因此計算所得的強度值與實驗值相比會高出一兩個數量級.文章后面的部分也會介紹到非化學計量比下SiC的理想強度，這個值比沒有缺陷的晶體材料值低得多，這也說明了缺陷確實會降低材料的理想剪切強度.

表1 β-SiC的平衡晶格常數a0，彈性常量C11、C12、C44，剪切模量G100，剪切強度τ100及相應的應變ε100

圖2所示的是理想無缺陷的3C-SiC的應力應變曲線.應力開始隨著應變的增加開始呈線性關系，在應變ε=0.13的時候出現一個極大值，隨后應力隨應變線性增加.忽略極大值出現的微小擾動，在ε<0.34時，應力隨應變是線性變化的，發(fā)生的是彈性形變，彈性形變區(qū)的剪切模量為169 GPa.當ε>0.35時，此時材料出現非彈性形變，應力隨應變增加至屈服強度后急劇下降，此時原子鍵發(fā)生斷裂.β-SiC的理想強度為60 GPa，其發(fā)生原子鍵斷裂相應的應變ε為0.44，在最大應力下，C-Si鍵的鍵長為2.15 ?，而在完全沒有應力的情況下，C-Si鍵的鍵長為1.90 ?.

2.2 非化學計量比下SiC的計算結果

實際工業(yè)制造的SiC中，C濃度是要比Si的濃度高，為了研究C含量對SiC的影響，用C原子取代晶體中的Si原子，構建了10組C/Si比大于1的SiC模型.96個原子的超晶胞中，含有48個C原子和48個Si原子，分別用1～10個C原子取代等數量的Si原子，構建了C/Si比從1.04～1.52的10組SiC模型.每個比例選取5個不同結構作為計算的樣品，每個結構經過蒙羅卡羅方法優(yōu)化使系統能量最低，最后用第一性原理軟件包VASP優(yōu)化，優(yōu)化后的同一C/Si比的計算樣品的能量差別很小，小于10-3eV每個原子.

圖2 理想3C-SiC和C/Si比為1.04的SiC樣品的沿[100]方向上剪切的應力應變曲線Fig.2 The stress-strain curve of ideal 3C-SiC and a sample with the C/Si ratio 1.04 sheared along [100] direction

對于每一個樣品，通過VASP計算得到沿[100]方向剪切的應力應變曲線.計算結果顯示，C/Si比對SiC的剪切模量幾乎沒什么影響，從圖2中也可以看出，不同C/Si比的應力應變曲線的斜率幾乎是一致的.圖2中所示的是C/Si為1.04的應力應變曲線和沒有缺陷的β-SiC晶體的應力應變曲線.從圖中可以比較得出，C/Si比為1.04的SiC的理想剪切強度比C/Si為1的明顯低了許多.除此之外，無缺陷SiC在彈性形變之后，有一段非彈性形變，當應力超過非彈性形變的屈服強度后，原子鍵才發(fā)生斷裂，應力突降，而C/Si比為1.04的SiC應力達到了彈性應變的最大之后，原子鍵就發(fā)生了斷裂.這說明非化學計量比的SiC較化學計量比的SiC不僅強度下降了，而且變得更脆.每個比例的理想強度取五個樣品計算結果的平均值，表2及圖3列出了不同C/Si比SiC理想強度的平均值和標準差.

表2 不同C/Si比的SiC沿[100]方向上的理想剪切強度及相應的最大應變值

表2中列出了不同C/Si比SiC的理想剪切強度及相應的最大應變值.表中的理想強度為同一C/Si比5個不同結構理想剪切強度的平均值，強度的標準差為同一C/Si比5個不同樣品的最大剪切應力的標準差，最大應變?yōu)橥籆/Si比的5個樣品最大剪切應力對應的應變的平均值，最大應變標準差則為這5個應變的標準差.從表中可以看出，β-SiC的理想強度隨著C/Si比的增加總體上是呈現下降趨勢的，但是不是嚴格的遞減規(guī)律，原子鍵斷裂的最大應變隨C/Si比的增加總體上也是減小的.

圖3所示的是SiC理想剪切強度隨C/Si比變化曲線，曲線橫坐標是C/Si比，縱坐標是SiC的理想剪切強度，圖中以理想強度為中點的用黑線段也表示出了理想強度的偏差范圍.從圖中可以再出C/Si比從1增加到1.09的范圍內，SiC的理想強度是急劇下降的，理想強度在為C/Si比1.09～1.40范圍內的變化基本不大，曲線中出現了第1個平臺，之后又下降，在C/Si比1.46后又出現第2個平臺.總體上來看，SiC的理想剪切強度隨C/Si比增加是呈現下降趨勢的.

圖3 SiC沿[100]方向理想剪切強度隨C/Si比的變化曲線Fig.3 The ideal shear strength along [100] direction of SiC as a function of the C/Si ratio

3 結果分析

3.1 鍵長分析

固體SiC中，每個C原子周圍是4個Si原子，每個Si原子周圍是4個C原子.C原子一旦取代Si原子必然會形成至少5個C原子聚集的C團簇(見圖1)，這個C團簇的結構與金剛石的結構相似.而金剛石中的C-C鍵鍵長為1.55?，而SiC晶體中的Si-C鍵鍵長為1.89?，前者比后者要短18%.圖1所示的C/Si比為1.09的結構中，AD(C-C)原子之間的鍵長為1.66?，CD(Si-C)原子之間的鍵長為1.93 ?.C原子取代Si原子之后，取代后C原子周圍的鍵長均會發(fā)生變化，首先取代的C原子與周圍的C原子形成了C-C鍵，C-C鍵比原來的Si-C短，C團簇進一步收縮，與C團簇相連的Si-C鍵在C團簇的收縮應力下被拉長.因此，經C原子取代后Si原子后的SiC中，C-C鍵要比金剛石中的C-C要長，而Si-C鍵也比晶體SiC中的Si-C鍵長，在C團簇周圍形成了內部拉伸應力.這個內部應力會使SiC的機械性能降低，也使得結構在剪切過程中最先發(fā)生原子鍵斷裂的地方出現在C團簇中被拉長的C-C鍵或者Si-C鍵.從能量的角度看，金剛石中C-C鍵的結合能能是8.55 eV[19]，而SiC中Si-C鍵的結合能則是2.35 eV[20]，也就是說，形成金剛石結構較形成SiC結構釋放更多的熱量.總之，SiC中C原子聚集成C團簇，導致SiC產生內應力，進一步導致SiC機械性能變差.

圖4所示的是圖1結構達到剪切強度后，原子鍵斷裂后的結構.圖中可以看出斷裂的鍵均出現在2個取代的碳原子A原子和B原子周圍，如AD鍵，還有未標記出的B原子周圍的C-C鍵和Si-C鍵.這也證實了由于C團簇的出現形成了內應力，致使斷裂的地方都出現在C團簇附近.

圖4 C/Si比為1.09的SiC沿[100]方向上剪切后，原子鍵斷裂后的結構Fig.4 The configuration of broken SiC sample with the C/Si ratio of 1.09

3.2 C團簇大小的影響

當2個以上C原子取代的Si原子比較靠近時，這時SiC晶體中就會形成較大的C團簇，C團簇的大小對理想剪切強度也有很大的影響.圖3中同一C/Si比下不同結樣品的理想強度偏差的出現，由于取代的C原子所在位置不同和C原子團簇的大小不同，導致5個計算樣品結構不同，計算出的理想強度結果出現偏差.關于C團簇大小對理想強度的影響在文獻[12]有更深入的討論，由于理想強度的定義為原子鍵發(fā)生斷裂所受的應力值，因此剪切與拉伸強度的本質一樣，也是C團簇越大，則剪切強度下降的越多.

4 結論

本文用第一性原理的方法研究了C/Si比對SiC[100]方向上理想剪切強度的影響，從計算結果來看，隨著C/Si比的增加，SiC在[100]方向上的理想剪切強度總體上是呈現下降趨勢的.所選取的計算樣品中，由于C原子取代Si原子在SiC內部集聚成C團簇，使得SiC在C團簇附近出現內應力.內應力的出現使SiC的理想剪切強度下降，并且在剪切作用下發(fā)生斷裂的原子鍵均為受內應力影響的C-C鍵或者Si-C鍵.相比于均勻分布的C團簇，大的C團簇能進一步減少SiC的理想剪切強度.

[1] NEWSOME G, SNEAD L L, HINOKI T, et al. Evaluation of neutron irradiated silicon carbide and silicon carbide composites[J]. Journal of Nuclear Materials, 2007, 371(1-3): 76-89．

[2] KATOH Y, SNEAD L L, HENAGER J R, et al. Current status and recent research achievements in SiC/SiC composites[J]. Journal of Nuclear Materials, 2014, 455(1-3): 387-397．

[3] 于海蛟, 周新貴, 張 煒, 等. 聚變堆結構材料SiC/SiC的研究進展[J]. 材料導報, 2009(23): 104-108．

[4] HINOKI T, LARA-CURZIO E, SNEAD L L. Mechanical Properties of High Purity SiC Fiber-Reinforced CVI-SiC Matrix Composites[J]. Fusion Science and Technology, 2003, 44(1): 211-218．

[5] SNEAD L L, KATOH Y, KOHYAMA A, et al. Evaluation of neutron irradiated near-stoichiometric silicon carbide fiber composites[J]. Journal of Nuclear Materials, 2000, 283: 551-555．

[6] HINOKI T, KATOH Y, KOHYAMA A. Effect of Fiber Properties on Neutron Irradiated SiC/SiC Composites[J]. Materials Transactions, 2002, 43(4): 617-621．

[7] OZAWA K, NOZAWA T, KATOH Y, et al. Mechanical properties of advanced SiC/SiC composites after neutron irradiation[J]. Journal of Nuclear Materials, 2007, 367: 713-718．

[8] LI T, MORRIS J W, CHRZAN D C. Ideal tensile strength of B2 transition-metal aluminides[J]. Physical Review B, 2004, 70(5): 054107(1-5)．

[9] LAMBRECHT W R L, SEGALL B, METHFESSEL M, et al. Calculated elastic constants and deformation potentials of cubic SiC[J]. Physical Review B, 1991, 44(8): 3685-3694．

[10] LI W, WANG T. Elasticity, stability, and ideal strength of β-SiC in plane-wave-based ab initio calculations[J]. Physical Review B, 1999, 59(6): 3993-4001．

[11] UMENO Y, KUBO A, NAGAO S. Density functional theory calculation of ideal strength of SiC and GaN: Effect of multi-axial stress[J]. Computational Materials Science, 2015, 109: 105-110．

[12] LI Y, XIAO W. First principles study of the C/Si ratio effect on the ideal tensile strength of β-SiC[J]. Computational Materials Science, 2015, 110: 215-220．

[13] GAO F, BYLASKA E J, WEBER W J et al. Ab initio and empirical-potential studies of defect properties in 3C-SiC[J]. Physical Review B, 2001, 64(24): 245208(1-7)．

[14] ERHART P, ALBE K. Analytical potential for atomistic simulations of silicon, carbon, and silicon carbide[J]. Physical Review B, 2005, 71(3)：035211(1-14)．

[15] KRESSE G, HAFNER J. Ab initio molecular dynamics for liquid metals[J]. Physical Review B, 1993, 47(1): 558-561．

[16] KRESSE G, HAFNER J. Ab initio molecular-dynamics simulation of the liquid-metal-amorphous-semiconductor transition in germanium[J]. Physical Review B, 1994, 49(20): 14251-14269．

[17] KRESSE G, JOUBERT D. From ultrasoft pseudo-potentials to the projector augmented-wave method[J]. Physical Review B, 1999, 59(3): 1758-1775．

[18] GUO H, WANG B, JIA P, et al. In-plane shear behaviours of a 2D-SiC/SiC composite under various loading conditions[J]. Ceramics International, 2015, 41(9): 11562-11569．

[19] ZAVODINSKY V G. Ab initio simulation of diamond epitaxial growth on copper[J]. Computational Materials Science, 2006, 36(1-2): 139-142．

[20] AGRAWAL B K, YADAV P S, YADAV R K, et al. Ab initio theoretical study of SiC microclusters[J]. Progress in Crystal Growth and Characterization of Materials, 2006, 52(1-2): 15-20．

First-principles study of C/Si ratio effect on the ideal shear strength of SiC

SU Wen1， LI Yingying2

(1.School of Power and Mechanical Engineering， Wuhan University， Wuhan 430072；2.State Nuclear Power Research Institute， Beijing 100029)

The C/Si atom ratio effect on the ideal shear strength of SiC is investigated with first-principle calculation. The stoichiometric ideal β-SiC model is established. Then the Si atom of the model is replaced by C atom and Monte Carlo (MC) simulations are used to generate non-stoichiometric SiC samples. The C/Si atom ratio of which ranges from 1 to 1.52. Each non-stoichiometric SiC super cells is sheared along the [100] direction and the stress-strain relationship is obtained. The calculation results show that the ideal shear strength of the SiC decreases generally with increase of the C/Si atom ratio. Under the shear loading, the sample breaks at the Si-C bonds or the C-C bonds around the carbon cluster. The concentration of C atoms reduces the ideal shear strength because of the internal strain around the carbon clusters. Furthermore, the larger the size of carbon cluster is, the more dramatically the ideal shear strength reduces.

ideal shear strength； first-principle calculation； Monte Carlo simulation； SiC； carbon cluster； non-stoichiometric

2015-11-03.

1000-1190(2016)02-0196-06

O482.1

A

*E-mail: suwen@whu.edu.cn.