高壓下Nb3Sn 單晶的超導(dǎo)相轉(zhuǎn)變

石震天，楊緒佳，王豪陽，喬 力

（太原理工大學(xué)機(jī)械與運(yùn)載工程學(xué)院應(yīng)用力學(xué)研究所，山西 太原 030024）

A15型Nb3Sn超導(dǎo)材料是制造高場(＞10 T)超導(dǎo)磁體線圈的主要材料，被廣泛應(yīng)用于磁約束可控核聚變和高能物理等強(qiáng)磁體領(lǐng)域[1-5]。力學(xué)變形誘導(dǎo)的 Nb3Sn材料超導(dǎo)電性能弱化是強(qiáng)磁場磁體制造需要解決的基礎(chǔ)問題之一，研究高壓下 Nb3Sn單晶的超導(dǎo)相轉(zhuǎn)變行為對于揭示這一弱化機(jī)理具有重要意義。早期對于 Nb3Sn材料超導(dǎo)相轉(zhuǎn)變行為的研究主要集中于工程用復(fù)合超導(dǎo)股線在力學(xué)載荷作用下超導(dǎo)臨界性能參數(shù)退化的實(shí)驗(yàn)和表征[6-7]。基于此，發(fā)展了可用于磁體設(shè)計(jì)的經(jīng)驗(yàn)標(biāo)度關(guān)系[8-12]。伴隨著強(qiáng)磁場制造水平的提升，需要進(jìn)一步探究 Nb3Sn 材料超導(dǎo)相轉(zhuǎn)變力學(xué)效應(yīng)背后的物理機(jī)制。與 Nb3Sn復(fù)合超導(dǎo)股線相比，單晶體和多晶體的組織結(jié)構(gòu)更簡單，且更方便獲得不同載荷模式下的電阻率變化規(guī)律，因此高壓下 Nb3Sn單晶體的超導(dǎo)相轉(zhuǎn)變受到了學(xué)者們的廣泛關(guān)注。

式中： ρ0為材料的殘余電阻率， ρ1、ρ2和T0均為根據(jù)實(shí)驗(yàn)結(jié)果擬合的參數(shù)。該模型形式簡明，在18～850 K溫度范圍內(nèi)與實(shí)驗(yàn)數(shù)據(jù)的誤差小于1%。在Woodard-Cody 模型的基礎(chǔ)上，Qiao 等[18]建立了 Nb3Sn超導(dǎo)材料常態(tài)電阻率經(jīng)驗(yàn)?zāi)Ｐ?，用以描述常態(tài)電阻率、三維應(yīng)變和溫度的耦合效應(yīng)，模型預(yù)測結(jié)果與實(shí)驗(yàn)數(shù)據(jù)一致。這些經(jīng)驗(yàn)?zāi)Ｐ屯ㄟ^簡潔的數(shù)學(xué)公式很好地刻畫了實(shí)驗(yàn)觀測結(jié)果，但是缺乏對超導(dǎo)相轉(zhuǎn)變附近正常態(tài)電阻率靜水壓效應(yīng)機(jī)理的闡述。Webb 等[19]認(rèn)為超導(dǎo)相轉(zhuǎn)變附近常態(tài)電阻率的T2依賴特性是基于電子-聲子散射產(chǎn)生的，但該理論不足以解釋所有的數(shù)據(jù)，仍然需要繼續(xù)探索載荷作用下電阻率T2依賴性的內(nèi)在機(jī)制[20]。為理解這一效應(yīng)，需要定量描述 Nb3Sn單晶體在高壓和極低溫環(huán)境下的變形特征及其誘導(dǎo)的晶體結(jié)構(gòu)和電子結(jié)構(gòu)的變化，以期對實(shí)驗(yàn)現(xiàn)象進(jìn)行準(zhǔn)確分析。

為此，本研究將借助分子動(dòng)力學(xué)模擬方法研究 Nb3Sn單晶在極低溫環(huán)境、高壓載荷下的晶體結(jié)構(gòu)變化，此外，還將建立高壓下 Nb3Sn單晶超導(dǎo)相轉(zhuǎn)變的解析模型，定量給出應(yīng)變誘導(dǎo)費(fèi)米面上電子態(tài)密度 變化在高壓作用下 Nb3Sn單晶超導(dǎo)相轉(zhuǎn)變中的作用。

1 高壓下Nb3Sn 單晶體的分子動(dòng)力學(xué)模擬

利用L A M M P S 程序進(jìn)行靜水壓作用下Nb3Sn單晶體的分子動(dòng)力學(xué)模擬，計(jì)算模型如圖1所示。在LAMMPS 內(nèi)利用自定義構(gòu)建晶格的命令建立 Nb3Sn單晶體的A15 立方結(jié)構(gòu)，晶格常數(shù)a 為5.29 ?，模擬盒子為邊長50a 的正方體，原子數(shù)為1 06。 Nb3Sn單晶體的[001]、[010]、[100]晶軸分別對應(yīng)模擬盒子的x、y、z 3 個(gè)方向。為了使模擬結(jié)果更接近正常大小的晶體，模擬采用周期性邊界條件。利用高斯分布隨機(jī)給出指定溫度下原子的初始速度。

圖1 N b3Sn單晶體計(jì)算模型（a）和A15 型晶格結(jié)構(gòu)（b）Fig. 1 Calculation model of N b3Sn single crystal (a)and A15 type lattice structure (b)

采用組合勢函數(shù)的方法定義 N b3Sn原子間的相互作用，Nb-Nb 的勢函數(shù)采用Zhang 等[21]開發(fā)的用于 Ni62Nb38的原子嵌入勢，Sn-Sn 的勢函數(shù)采用MEAM 勢函數(shù)[22]，Nb-Sn 間的勢函數(shù)為L-J 勢，參數(shù)參考文獻(xiàn)[23]。計(jì)算 N b3Sn單晶體的彈性常數(shù)和晶格常數(shù)，如表1 所示，其中C 為彈性系數(shù)。通過與實(shí)驗(yàn)結(jié)果和第一性原理模擬[24-26]的比較，證明了該勢函數(shù)對于力學(xué)性能參數(shù)表征的可靠性。

表1 Nb3Sn 單晶的力學(xué)性能參數(shù)Table 1 Elastic constants and lattice constant of Nb3Sn single crystal

加載靜水壓之前，在4.2～44.2 K 溫度區(qū)間內(nèi)均勻取41 個(gè)溫度點(diǎn)，分別在這些溫度點(diǎn)且無外部壓強(qiáng)加載條件下用NPT 系綜弛豫40 000 步(時(shí)間步長為0.001 ps)，使晶體結(jié)構(gòu)達(dá)到平衡。弛豫完成之后，對模型連續(xù)施加0～10 GPa 的靜水壓，采用OVITO 軟件進(jìn)行模擬結(jié)果的后處理。

2 Nb3Sn 單晶超導(dǎo)轉(zhuǎn)變隨壓力和溫度的演化關(guān)系

研究高壓下Nb3Sn 超導(dǎo)相的轉(zhuǎn)變行為，除了需要刻畫高壓下Nb3Sn 單晶體的變形外，還需要對超導(dǎo)轉(zhuǎn)變溫度附近的電阻率變化行為進(jìn)行定量描述。Wiesmann 等[27]建立了一個(gè)類似并聯(lián)電阻形式的N b3Sn 超導(dǎo)材料常態(tài)電阻率模型

式中： ρideal為 理想常態(tài)電阻率， ρm為飽和電阻率。受此啟發(fā)，假定超導(dǎo)相轉(zhuǎn)變的電阻率由正常態(tài)電阻率ρn和超導(dǎo)態(tài)電阻率 ρs相互競爭決定，其中 ρs起到“分流” ρn的作用，基于Yang 等[28]的并聯(lián)電阻率模型，N b3Sn 單晶超導(dǎo)相轉(zhuǎn)變的電阻率

式 中： Π =ρn/ρs， 為無量綱電阻率函數(shù)； ε為壓力作用下的應(yīng)變張量。

2.1 壓力-溫度耦合作用下的常態(tài)電阻率ρn(ε,T)

Nb3Sn超導(dǎo)體在低于臨界溫度Tc（約40 K）范圍內(nèi)，常態(tài)電阻率與溫度的關(guān)系可表示為[29-30]

式中： ρ0為不依賴于力學(xué)變形的殘余電阻率；A 為平方項(xiàng)系數(shù)，即正常態(tài)電阻率的T2依賴關(guān)系系數(shù)。實(shí)驗(yàn)揭示了A15 型超導(dǎo)體常態(tài)電阻率的T2依賴性是一種內(nèi)在性質(zhì)。不論是電子-電子散射，還是具有非德拜相似聲子結(jié)構(gòu)的聲子-電子帶間散射，都不能完美地解釋這種性質(zhì)[30]。在高壓作用下，正常態(tài)電阻率的變化依然遵循這一規(guī)律。然而，不同的是，壓強(qiáng)（應(yīng)變）會(huì)導(dǎo)致系數(shù)A 發(fā)生變化，所以該關(guān)系式可以拓 展為

式 中： Aε為與應(yīng)變相關(guān)的平方項(xiàng)系數(shù)，與態(tài)密度的平方成正比，所以這里 Aε可以表示為

式 中： χi、κi和 βij均為根據(jù)第一性原理計(jì)算結(jié)果擬合得到的參數(shù)。將式(7)代入式(6)得到

2.2 無量綱電阻率函數(shù)Π

在高壓-低溫聯(lián)合作用下，當(dāng) Nb3Sn從超導(dǎo)態(tài)向正常態(tài)轉(zhuǎn)變時(shí)，伴隨著超導(dǎo)相消失，材料的電阻率逐漸增加。根據(jù)電阻率-溫度轉(zhuǎn)變曲線變化趨勢可以看出，中點(diǎn)轉(zhuǎn)變溫度 T1/2(ε)和 無量綱電阻率函數(shù) Π具有兩個(gè)特點(diǎn)[28]：（1）材料從正常態(tài)向超導(dǎo)態(tài)轉(zhuǎn)變的過程中，即無量綱溫度 TN=T/T1/2(ε)從1 減小到零的過程中， Π將從1 增加到無窮大；（2）當(dāng) TN從1 開始增加時(shí)， Π從1 到零的轉(zhuǎn)變過程是在1/2 轉(zhuǎn)變溫度范圍內(nèi)逐漸完成的。在超導(dǎo)相轉(zhuǎn)變的區(qū)間內(nèi)，電阻率函數(shù)（兩種極限情形下的電阻率比值）連續(xù)變化。基于 Kim 等[32]和Godeke 等[33]的工作，借助修正的Arrhenius 關(guān)系給出無量綱電阻率函數(shù)的表達(dá)式

式中：無量綱參數(shù) κ ≡ΔTW(ε)/T1/2(ε) 描 述了應(yīng)變作用下超導(dǎo)轉(zhuǎn)變區(qū)域的溫度寬度， Δ TW(ε)表 示 N b3Sn 單晶超導(dǎo)轉(zhuǎn)變過程中完成電阻率躍變所經(jīng)歷溫度區(qū)間的寬度。

為了確定電阻率函數(shù) Π的形式，需要確定壓力載荷下超導(dǎo)臨界溫度 T1/2(ε)和 無量綱參數(shù) κ的表達(dá)形

式。其中，壓力誘導(dǎo)的臨界溫度變化由Maki De Gennes(MDG)關(guān)系式[34]結(jié)合臨界磁場對溫度的偏導(dǎo)數(shù)在臨界溫度處的取值導(dǎo)出

s(ε)為臨界溫度退化的標(biāo)定函數(shù)，可以表示為

超導(dǎo)轉(zhuǎn)變寬度 ΔTW是 描述電阻率轉(zhuǎn)變的另一個(gè)重要參數(shù)， ΔTW=T90%－T10%， T90%、T10%分別是電阻率躍變90%和10%對應(yīng)的溫度[35]。基于多物理場環(huán)境下 Nb3Sn超導(dǎo)臨界曲面表達(dá)式的構(gòu)造方式，通過對 臨界溫度和臨界磁場修正引入應(yīng)變，根據(jù)量綱分析結(jié)果和分離變量方法可以給出超導(dǎo)轉(zhuǎn)變寬度[28]

式中：K 為無量綱比例常數(shù)， g (H,H1/2(0,ε)) 為關(guān)于磁場和上臨界磁場的無量綱函數(shù)。 N b3Sn單晶在一般加載模式下的上臨界磁場可以采用 H1/2(ε)=H1/2(0)s(ε)描述。借鑒Josephson 等[32]和Tinkham 等[35]關(guān)于 高溫超導(dǎo)體電阻率轉(zhuǎn)變的研究，得出轉(zhuǎn)變寬度和磁場有關(guān)，即

進(jìn)而得出無量綱函數(shù) g 與 s(ε)α(1－h(huán))α成比例。通過上述分析可以得出無磁場作用下無量綱參數(shù) κ在一般加 載模式下的形式

3 結(jié)果分析

3.1 高壓下晶體結(jié)構(gòu)的變化

圖2 N b3Sn單晶體積隨溫度的變化Fig. 2 Volume of N b3Sn single crystal varies with temperature

分子動(dòng)力學(xué)模擬可以給出壓力和低溫環(huán)境下Nb3Sn單晶原子尺度的變形和晶格結(jié)構(gòu)變化，對這些細(xì)節(jié)的分析可以彌補(bǔ)實(shí)驗(yàn)觀測的不足。圖2 給出了4.2～44.0 K 溫度區(qū)間內(nèi)模擬單元的體積變化。在無靜水壓作用（0 GPa）時(shí)，隨著溫度的增加，Nb3Sn單晶體的體積增大。當(dāng)環(huán)境溫度從4.2 K 增加到44.0 K 時(shí)，體積增大了0.161%。加載0.2、2.6、5.7 和9.2 GPa 靜水壓時(shí)，隨著環(huán)境溫度從4.2 K 增加到44.0 K 時(shí)，模擬單元體積分別增大了0.151%、0.138%、0.121%和0.106%，體積變化率依賴于單晶所承受的靜水壓力載荷。

在4.2 K 低溫環(huán)境下，當(dāng)加載的靜水壓力從0 GPa 提高到10 GPa 時(shí)，Nb3Sn 單晶體3 個(gè)方向的主應(yīng)變和體積變化分別如圖3(a)和圖3(b)所示。圖3(a)顯示了該靜水壓加載下的3 個(gè)主應(yīng)變，可以看出3 個(gè)主應(yīng)變均隨著靜水壓的增大而增大，且近似為線性關(guān)系。3 個(gè)主應(yīng)變重合，這與經(jīng)典彈性力學(xué)理論吻合，即 ε1=ε2=ε3。隨著加載壓力增大，到達(dá)高壓區(qū)時(shí)與彈性力學(xué)分析得到的線性結(jié)果不同，主應(yīng)變與靜水壓強(qiáng)之間呈弱非線性關(guān)系，這起源于高壓下原子間的非諧相互作用。一般來說，在彈性階段加載應(yīng)力和變形的關(guān)系可以用彈性常數(shù)描述，彈性常數(shù) Cij為 單位體積內(nèi)原子間作用總勢能 U對應(yīng)變分量 的二階偏導(dǎo)數(shù)，即

式中： V0為 模擬體系的體積，U 為模擬體系的勢能， εi和 εj為應(yīng)變分量。當(dāng)施加在Nb3Sn 單晶的靜水壓強(qiáng)較小時(shí)，原子偏離平衡狀態(tài)的位移較小，考慮勢能中原子間相對位移（應(yīng)變）的二次方項(xiàng)即可準(zhǔn)確描述體系的變形，根據(jù)式（16）得到彈性模量為常數(shù)，即加載應(yīng)力和變形之間呈線性關(guān)系。在高壓區(qū)，原子偏離平衡狀態(tài)的位移增大，勢能中的原子間相對位移（應(yīng)變）的高次項(xiàng)（非諧項(xiàng)）不能忽略，從而使得加載應(yīng)力與變形之間呈非線性變化趨勢。由圖3(b)可見，單晶體體積隨靜水壓強(qiáng)的增大而減小，靜水壓達(dá)到10 GPa 時(shí)， Nb3Sn單晶體體積減小了約4%。

圖3 4.2 K 時(shí)Nb3Sn 單晶體在靜水壓作用下的變形情況Fig. 3 Deformation of Nb3Sn single crystal under hydrostatic pressure at 4.2 K

如圖4 所示，借助于可視化軟件OVITO，對靜水壓加載下Nb3Sn 單晶體的變形進(jìn)行可視化處理，給出了 Nb3Sn單晶體靜水壓加載下平面的平均應(yīng)力和整體及局部原子的Mises 應(yīng)力分布。圖4(a)給出了Nb3Sn 單晶體平面的平均應(yīng)力 σxx、 σyy和 σzz，可以看出平均應(yīng)力在10 GPa 左右。圖4(b)為 Nb3Sn單晶體的Mises 應(yīng)力分布云圖，可以看出：在10 GPa 靜水壓作用下， Nb3Sn單晶體發(fā)生了明顯的晶格畸變；但晶體結(jié)構(gòu)完整，沒有發(fā)生損傷；Mises 應(yīng)力主要集中在Nb 原子上。圖4(c)為 Nb3Sn單晶體原子等效應(yīng)力的分布直方圖，可以看出應(yīng)力35 GPa 以上的原子數(shù)是2 GPa 左右的3 倍。從應(yīng)力分布上看，Sn 原子的應(yīng)力明顯小于Nb 原子。模擬結(jié)果表明，Sn 原子所受局部Mises 應(yīng)力接近2 GPa，而Nb 原子的Mises 應(yīng)力在37 GPa 左右。

圖4 靜水壓加載下N b3Sn單晶體的應(yīng)力分布Fig. 4 Stress distribution of N b3Sn single crystal under hydrostatic pressure

圖5 給出了 Nb3Sn 單晶體主應(yīng)力在xz 面上整體及局部原子的主應(yīng)力云圖，其中 σ1、 σ2和 σ3分別為沿著x、y 和z 方向的主應(yīng)力?？梢钥闯?，各方向Sn 原子均受到較大的壓應(yīng)力，且該壓應(yīng)力大于Nb 原子受到的應(yīng)力，約為57 GPa；當(dāng)Nb 原子所在Nb 鏈的方向與主應(yīng)力方向平行時(shí)，會(huì)受到一個(gè)壓應(yīng)力，而其他方向Nb 鏈原子在該方向的應(yīng)力表現(xiàn)為拉應(yīng)力，但圖5 顯示所有Nb 鏈原子所受的應(yīng)力大小大致相同，均在20 GPa 左右。在靜水壓作用下，Sn 原子承受的主應(yīng)力在各方向均為壓應(yīng)力，Nb 原子在所在Nb 鏈方向承受的主應(yīng)力為壓應(yīng)力，其他方向承受的主應(yīng)力為拉應(yīng)力，這也是導(dǎo)致等效應(yīng)力主要集中在Nb 鏈原子的原因。

圖5 N b3Sn 單晶體xz 面的主應(yīng)力分布Fig. 5 Principal stresses distribution of N b3Sn single crystal on the right side under hydrostatic pressure

3.2 高壓下的超導(dǎo)相變

為了研究 N b3Sn 超導(dǎo)材料的壓力敏感性，Ren 等[16]通過實(shí)驗(yàn)研究了高壓下 N b3Sn 單晶體的電阻轉(zhuǎn)變行為?；谏鲜龈邏合翹b3Sn 單晶體的變形分析，借助本研究所建立的超導(dǎo)相轉(zhuǎn)變模型，對實(shí)驗(yàn)給出的經(jīng)驗(yàn)關(guān)系進(jìn)行解釋。模型參數(shù)C 和通過 Ren 等[16]的實(shí)驗(yàn)結(jié)果進(jìn)行最優(yōu)化擬合，其他參數(shù)由實(shí)驗(yàn)結(jié)果直接給出，如表2 所示。計(jì)算中的態(tài)密度模型參數(shù)來源于第一性原理模擬結(jié)果[31]，如表3 所示。

表2 電阻率模型相關(guān)參數(shù)Table 2 Parameters of resistivity model

表3 態(tài)密度函數(shù)相關(guān)參數(shù)Table 3 Parameters of density of state function

圖6 中實(shí)線部分給出了0～44 K 溫度范圍內(nèi) Nb3Sn單晶體在靜水壓作用下的電阻率預(yù)測值?？梢钥闯觯B(tài)電阻率部分與實(shí)驗(yàn)數(shù)值[16]基本吻合，所建立的計(jì)算模型能夠準(zhǔn)確刻畫出實(shí)驗(yàn)觀測到的超導(dǎo)轉(zhuǎn)變溫度附近電阻率的T2依賴性。對于常態(tài)電阻率部分，其斜率隨著靜水壓的增加而減小。在超導(dǎo)轉(zhuǎn)變區(qū)域，臨界溫度隨著靜水壓的增大而減小，由于分子動(dòng)力學(xué)模擬中選取的溫度間隔較大，超導(dǎo)轉(zhuǎn)變過程不夠平緩，并沒有完整的呈現(xiàn)，但是仍然可以看出隨著靜水壓強(qiáng)增大，超導(dǎo)轉(zhuǎn)變過程逐漸向后推移，與實(shí)驗(yàn)數(shù)據(jù)定性吻合。

圖6 靜水壓作用下Nb3Sn 單晶體在低溫區(qū)的電阻率隨溫度的變化Fig. 6 Change of resistivity with temperature in low temperature area for Nb3Sn single crystal under hydrostatic pressure

圖7 靜水壓作用下Nb3Sn 單晶體臨界溫度變化預(yù)測結(jié)果與實(shí)驗(yàn)結(jié)果對比Fig. 7 Comparison between the predicted and measured critical temperature of Nb3Sn single crystal under hydrostatic pressure

4 結(jié) 論

基于分子動(dòng)力學(xué)模擬方法研究了高壓下 Nb3Sn單晶體的晶體結(jié)構(gòu)和局部原子應(yīng)力，在此基礎(chǔ)上，受Wiesmann 并聯(lián)電阻率模型的啟發(fā)，基于 Nb3Sn單晶體常態(tài)電阻率的T2依賴性，建立了高壓下Nb3Sn 單晶體的超導(dǎo)相轉(zhuǎn)變模型。結(jié)果表明：高靜水壓作用下 Nb3Sn單晶體發(fā)生了明顯的晶格畸變，但晶體結(jié)構(gòu)完整，沒有發(fā)生損傷；高靜水壓作用下電阻率轉(zhuǎn)變的模型預(yù)測結(jié)果與實(shí)驗(yàn)觀測結(jié)果吻合較好，證明了壓力誘導(dǎo)費(fèi)米面上電子態(tài)密度的變化在高壓下 Nb3Sn單晶體超導(dǎo)相轉(zhuǎn)變中的作用。相關(guān)研究結(jié)果對于理解高壓下 Nb3Sn 單晶體的超導(dǎo)相變行為具有重要意義，同時(shí)為后續(xù)研究高壓下 Nb3Sn多晶體以及復(fù)合多晶體的相變行為奠定了基礎(chǔ)。