由熱學性質獲取氬晶體原子間各階力常數

黃建平，胡詩一

(湖南師范大學計算機部，量子結構與調控教育部重點實驗室，長沙410081)

1 引 言

根據晶格動力學理論［1］，原子間相互作用決定了晶體的許多物理性質，因此原子間相互作用的研究是一項非常重要的基礎研究，長期以來一直是研究熱點. 例如，最近J?ger 等人［2］通過ab initio 方法計算得到了氬原子對相互作用勢. 然而，由于氬晶體內原子間相互作用不能等同于氬原子對內原子間相互作用，因此J?ger 等人的計算結果不能準確的反映氬晶體內氬原子之間的原子相互作用的實際情況. 雖然通過擬合熱膨脹和比熱的實驗數據，可以獲取晶體中原子間相互作用力常數，但由于現有的比熱公式只考慮到和諧勢能的貢獻［3］，熱膨脹系數公式只考慮到三階非和諧勢能的貢獻［4］，因此我們根據該思路只能對低溫段的熱學性質參數進行擬合而得到二階和三階力常數. 然而在溫度較高情況下，更高階的力常數對熱膨脹和比熱的貢獻不能忽略，因此有必要首先推導得到包含各階力常數的熱膨脹和熱膨脹系數公式，以獲取更為全面、準確的晶體內部原子間相互作用信息.

盡管氬晶體缺乏實際應用價值，但由于其晶體結構簡單，常被用于驗證各種固體理論［5］，因此本文也將以氬晶體作為研究對象. 本文將運用晶格動力學和量子力學定態微擾理論［6］，推導比熱和熱膨脹系數與氬晶體原子間各階力常數的關系公式，再根據這些公式對熱膨脹系數和比熱數據［7］進行擬合，計算出氬晶體原子間各階力常數，并根據這些力常數還原出氬晶體原子間相互作用勢能曲線，并與Morse 勢能［5］曲線進行比較.

2 氬晶體的晶格動力學基礎

其中，σ、ρ、α 和β 可為x、y 或z，且α、β 和σ互不相同. m 為原子質量，k = ( kx，ky，kz)為格波波矢. 當σ = ρ 時δσρ= 1 ，= 0 ;而當σ ≠ρ時δσρ= 0 ，= 1 .

求解晶格動力學方程［1］，得本征值和對應的本征矢e(kj)，其中，j = 1，2，3，e(kj)滿足正交歸一化條件. 原子l 圍繞平衡位置R( )l 進行簡諧振動的瞬時振動位移矢量記為ul( )t . 和諧晶體的晶格原子位移和晶格振動哈密頓可分別表示為其中，N 為原子數，Akj為akj+，和akj是聲子的產生與湮滅算符，?ωkj為聲子能量，nkj是聲子數算符，平均聲子數符合玻色統計. 對(3)式求熱力學平均，可得和諧晶體的晶格振動內能E2，據此可得比熱C2.

勢能中還包含3 階及3 階以上非和諧勢能項，總的非和諧勢能項可表示為

其中，原子ζ 和ζ'互為最近鄰原子，ζ ＜ζ'表示對原子ζ 和ζ' 求和時，避免對勢能的重復計算. 方括號內為原子ζ 和ζ' 之間距離的變化量.

其中，eζζ'為R ζ( )' - R( )ζ 的單位矢量.

3 氬晶體的比熱及熱膨脹系數公式

非和諧勢能較小，在量子力學中可作微擾處理. 在一級近似下，根據量子力學微擾理論，氬晶體晶格振動內能E 可表示為

根據(6)式計算可知，在一級近似下，奇次非和諧勢能對晶格振動內能沒有貢獻. 計算到2n 階非和諧勢能，得晶格振動內能

其中，O 為座標原點. 上式可簡化為

其中，ε2= E2/N 為平均到單個原子的和諧晶體晶格振動內能.

計及至2n 階非和諧勢能對晶格振動內能影響，得到非和諧氬晶體的比熱.

其中，| n ＞或| n' ＞分別是由和諧晶體各種kj 的聲子態| nkj＞或| n'kj＞的直積構成的未微擾態.

先計算3 階非和諧勢能項H3對晶格常數熱膨脹的貢獻Δ3a. 根據(2)式，及聲子產生和消滅算符的性質可知， (10)式中求和項非零的條件是| n ＞和| n' ＞狀態相差一個聲子，設其模為k'j' ，在| n ＞態中該聲子的數目為nk'j'，則| n' ＞態中該聲子的數目為n'k'j'= nk'j'-1 或n'k'j'= nk'j'+1 ，而所有其它模kj 的聲子數在| n ＞和| n' ＞態中是相同的. 根據聲子升降算符的性質，以及以上對| n ＞態和| n' ＞態的分析，可知只有當k1j1、k2j2和k3j3中，有一個為k'j'或- k'j' ，另外的為kj 和- kj 時，＜n'| H3| n ＞才不為零.根據以上分析，計算可得

運用以上方法，我們還計算了其它各階非和諧勢能項對熱膨脹的貢獻，發現只有奇次非和諧勢能項對熱膨脹才有貢獻. 至2n + 1 階的非和諧勢能引起的晶格常數熱膨脹為

根據(9)式和(15)式可知，比熱和熱膨脹系數與原子間各階力常數有關，因而可以利用Peterson 等人［7］的氬晶體比熱和熱膨脹系數的實驗數據計算原子間各階力常數.

4 數值計算結果與討論

在低溫度區，非和諧勢能對氬晶體比熱的貢獻可以忽略，因此我們首先根據和諧氬晶體比熱公式按最小二乘法擬合低溫氬晶體比熱數據，得二階力常數，在此基礎上計算得到全溫區的氬晶體比熱與溫度的關系曲線，如圖1 中的n =2 對應曲線所示.可知，在高溫情形，用和諧晶體晶格模型來描述氬晶體的比熱性質是有較大誤差的，因此必須考慮高偶數階非和諧勢能對比熱的貢獻. 設n = 12，根據(9)式按最小二乘法擬合全溫段的氬晶體比熱數據［7］，得到至12 階的各偶數階力常數.

圖1 氬晶體比熱與溫度關系Fig.1 The heat capacity of argon crystal vs temperature

在低溫度區，首先只考慮到3 階非和諧勢能，根據(15)式按最小二乘法擬合低溫氬晶體熱膨脹系數數據，得3 階力常數，在此基礎上計算得到全溫區的氬晶體熱膨脹系數與溫度的關系曲線，如圖2 中的n = 3 對應的曲線所示. 可知，在高溫區，只計及3 階非和諧勢能對熱膨脹性質的貢獻也會產生較大誤差，必須考慮高奇數階非和諧勢能對熱膨脹的貢獻. 考慮至11 階非和諧勢能，根據(15)式按最小二乘法計擬合全溫段的氬晶體熱膨脹系數數據［7］，由此計算出直至11 階的奇數階力常數.

圖2 氬晶體熱膨脹系數與溫度Fig.2 The thermal expansion coefficient of argon crystal vs temperature

圖3 氬晶體中原子間勢能與距離變化量關系Fig.3 The interatomic potential of argon crystal vs the change of interatomic space

根據各階力常數，得到原子間相互作用勢能δE 與原子間距離變化量δr 之間的關系曲線如圖3所示. 將受Morse 勢［5］作用的晶格原子處于平衡位置時的δr 和δE 分別規定為零，繪制得到Morse勢能曲線如圖3 所示，這與前面得到的勢能曲線較好地吻合，說明本文提供的研究氬晶體中原子間相互作用勢能的方法和結果是正確和可靠的.

［1］ Bottger H. Principles of the theory of lattice dynamics［M］. Weinheim:Physik-Verlag，1983:15.

［2］ J?ger B，Hellmann R，Bich E，et al. Ab intio pair potential energy curve for the argon atom pair and thermophysical properties of the dilute argon gas.I. Argon -argon interatomic potential and rovibrational spectra［J］. Mol. Phys.，2009，107:2181.

［3］ Mohazzabiy P，Behrooziz F. Thermal expansion of solids:a simple classical model［J］. Eur. J. Phys.，1997，18:237.

［4］ Kittel C. Introduction to solid state physics［M］. 8th Ed. New York:John Wiley ＆ Sons，2005:105.

［5］ Jelinek G E. Properties of crystalline argon，krypton，and xenon based upon the Born - Huang method of homogeneous deformations. III. The low-temperature limit［J］. Phys. Rev. B，1972，5:3210.

［6］ Sakurai J J，Napolitano J. Modern quantum mechanics［M］. 2nd ed. San Francisco:Addison - Wesley，2011:303.

［7］ Peterson O G，Batchelder D N，Simmons R O. Measurements of x-ray lattice constant，thermal expansivity，and isothermal compressibility of argon crystals［J］. Phys. Rev.，1966，150:703.