金屬單原子模型催化劑熱穩定性的反應力場(ReaxFF)分子動力學研究

楊文琦 汪 杰 喬園園 王貴昌

(南開大學化學學院，先進能源材料化學教育部重點實驗室，化學科學與工程協同創新中心(天津)，天津 300071)

單原子催化劑(SAC)在催化領域是相對較新的概念，在過去幾年中受到人們的廣泛關注[1-8]。由于它具有很高的化學勢，金屬單原子很容易遷移并聚集形成大的納米顆粒。如果一個金屬單原子想在高溫時穩定存在，那么就要求單原子和金屬基底之間具有強相互作用。較強的相互作用可以阻礙單個原子凝聚燒結，從而保持穩定。如果這個相互作用過于強，催化劑會失去它的催化活性。因此平衡催化劑活性和穩定性這兩者的關系是很有必要的。單原子催化劑燒結過程包括以下幾步：首先原子從基底表面脫離成為小的納米顆粒，之后發生遷移從而形成較大的納米顆粒，也就會發生凝聚燒結現象[9]。一般而言，SAC催化劑往往是過渡金屬負載在過渡金屬氧化物上，金屬-載體的強相互作用有利于SAC的穩定，如Datye課題組[6]發現氧化鈰與鉑的相互作用可以穩定氧化鈰部分晶面上孤立的Pt原子，可以在高溫時保持分散的Pt原子不聚集。另一方面，對于高度配位不飽和的過渡金屬單原子催化劑(如adatom)則研究的還不夠。

Kentaro等通過場離子顯微鏡研究了Pt(111)上Pt團簇的形成以及遷移能壘與溫度的關系，發現在350 K的低溫下，包括五聚體在內的Pt團簇擴散和解離現象都是比較顯著的[10]。Kokalj及其合作者[11]發現Rh1/Rh(111)-ad-tom催化劑可以增強甲烷第一次脫氫形成甲基的反應，同時可以抑制甲基的進一步脫氫反應，從而既保證了甲烷的充分活化，同時也避免了積碳現象的發生。那么，該類ad-atom催化劑在反應溫度下是否可以穩定存在呢?為探討該問題，在本工作中，我們通過反應力場(ReaxFF)結合LAMMPS(large-scale atomic/molecular massively parallel simulator)軟件包進行了大尺度分子動力學模擬，研究了一系列ad-atom過渡金屬催化劑(Cu，Ag，Au，Fe，Ni，Pd，Pt，Ru)在高溫的穩定性。 同時，通過與密度泛函(DFT)計算結果進行對比，驗證Ni/C/H/O力場參數，并探究了金屬單原子催化劑模型Ni1/Ni(111)在H2和O2氣氛下的穩定性行為。該研究工作為人們制備高活性和選擇性ad-atom催化劑提供了一定的理論指導依據。

1 計算方法

采用基于密度泛函理論的從頭量子力學的程序包VASP(Vienna ab-initio simulation package)[12-13]進行計算。基組采用平面波贗勢方法[14-15]，利用PBE關聯梯度修正泛函對模型進行充分優化[16]，采用廣義梯度近似(GGA)描述體系的交換關聯能。平面波的截斷能設置為 400 eV，k 點[17]設置為 3×3×1。 對于含有Fe和Ni的體系，加入了自旋極化。晶胞厚度為4層，擴胞3×3，單原子催化劑模型(M1/M(111))為添加一個額外的原子在fcc或hcp位點。在優化過程中，上3層馳豫，下2層原子保持固定。

大尺度分子動力學模擬軟件包LAMMPS[18]和反應性力場(ReaxFF)用于分子動力學的計算。動力學模擬所用模型為Cu單胞切(111)表面，厚度為4層，擴9×9超胞的周期性slab模型。這4層中下2層固定，上面2層放開，在表層上方有1.5 nm的真空層。時間步長設為0.25 fs，模擬時間設為100 ps，系綜選擇NVT系綜。首先設置初始溫度100 K模擬10 ps穩定，然后用50 ps將系統升溫至所需溫度(最終溫度為 300，400，500和 600 K)， 溫度阻尼為 25.0，最后模擬100 ps來穩定此狀態。SAC體系(Ag，Au，Fe，Ni，Pd，Pt，Ru 和 Cu) 反 應 性 力 場 參 數 來 源 于Duin等的工作[19-28]。Ni/C/H/O的力場參數基于Ni/C/H和Pt/O的力場參數[23,25]。

2 結果與討論

2.1 驗證反應力場參數的可靠性

我們首先對力場參數的可信度進行了驗證。在訓練集中加入了Ni塊體(bulk)的狀態方程，金屬Ni的內聚能(Ecohesive,見式1，EMbulk表示金屬M塊體的能量，NMbulk表示塊體中包含有金屬M原子的數目，EMatom代表氣相 中 金 屬 M 原子的能量)，CH3，CH2，CH，H2和O2等吸附物種在Ni(111)表面上的結合能(Ebinding,見式 2，Esp/Ni(111)和 ENi(111)分別代表吸附物種吸附后的體系和吸附前Ni(111)slab的能量，Espatom表示吸附物種在氣相中的能量)，Ni(111)面的表面能(Esurf，見式 3，Eslab，EMbulk和 NM分別為 Ni(111)slab 模型的能量，塊體中金屬Ni原子的能量和slab模型中Ni原子數目，A是切面的面積)，以及H2和O2在載體上分解的反應熱。采用單參數搜索法[29]，優化了Ni和H，Ni和O的鍵參數，以及Ni和H，O的交叉項。

表1 DFT計算結果和ReaxFF擬合的H 2和O2分別在Ni(111)和Ni1/Ni(111)上分解的反應熱的比較Tabel 1 ReaxFF fits for the reaction energy of the dissociation of H 2 and O2 on Ni(111)and Ni/Ni(111)respectively kJ·mol-1

圖1 不同吸附物種吸附在Ni(111)不同位點以及H2和O2吸附在Ni1/Ni(111)top位點上DFT和反應力場擬合的吸附能結果Fig.1 ReaxFFfits for the binding energy of different adsorbed species(CH3,CH2,CH,O and H)on the different site of Ni(111)surface and H2,O2 molecular on the Ni1/Ni(111)

圖1和表1結果表明力場參數對所研究體系吸附能和反應熱進行了擬合，說明本文所用力場參數是可靠的。

2.2 單原子金屬催化劑(M 1/M(111),M=Cu,Ag,Au,Pd,Pt,Ni,Ru,Fe和Fe1/Fe(100))的穩定性

圖2 從上至下分別為Cu,Ag,Au,Ni,Pd,Pt,Ru,Fe單金屬M1/M(111)模型和Fe1/Fe(100)模型在300,400,500和600 K時動力學構型(100 ps)的側視圖Fig.2 Side view of dynamics configurations(100 ps)Cu,Ag,Au,Ni,Pd,Pt,Ru,Fe single metal model from top to bottom respectively and Fe1/Fe(100)model at 300,400,500 and 600 K

對于單原子金屬催化劑模型，我們使用分子動力學方法研究了具有分散的單原子表面在不同溫度下的熱穩定性。通過DFT對四層擴胞3×3的Cu1/Cu(111)進行結構優化，然后將優化后的構型再次擴胞3×3最終形成一個四層擴胞9×9的結構，在nvt系綜100 K的低溫下模擬25 ps達到平衡。模擬溫度從100 K逐步升溫至300 K，并在300 K溫度下模擬100 ps，對模擬三次后的結果進行統計。發現單金屬原子Cu催化劑會發生聚集(圖2)，并且隨著溫度增加聚集的原子數增加。例如在300 K時有4～5個原子聚集形成五聚體，而400～500 K時9個Cu原子發生聚集，甚至在600 K時表層原子排列變得更加混亂，Cu原子顯著馳豫。這表明單原子Cu在高于300 K時不能穩定存在。

對于大部分過渡金屬 (如 Ag，Au，Pd，Pt和 Ni)，隨著溫度升高，單金屬原子模型變得不穩定并聚集在一起形成金屬團簇。有趣的是對于Fe1/Fe(100)模型，它可以穩定地存在不產生團簇(圖2)，但Fe1/Fe(111)模型的結構不能保持穩定。這說明金屬表面結構和其原子的排列也是影響單原子金屬催化劑穩定性的重要因素。如圖3所示，通過DFT計算可知大多數過渡金屬(除了Fe)均具有內聚能(Ecohesive)高于結合能(Ebinding)的性質(計算方法見公式(1)和公式(2))，說明單原子的聚集是熱力學有利的過程，因此傾向于聚集。然而對于鐵，它的內聚能小于結合能，單金屬鐵原子可以被有效地束縛住，意味著聚集過程不是熱力學有利過程，因此Fe單原子可以穩定地存在。

圖3 金屬單原子和載體的吸附強度以及對應的金屬原子的聚合能兩者之間的關系Fig.3 Relationship between binding strength of single atom and the corresponding cohesive energy

高溫時大部分ad-atom催化劑不穩定，因此研究低溫時M1/M(111)表面單原子團聚成簇現象與溫度的關系(如圖4所示)，同樣溫度下團簇中原子數越少越穩定。可以看出低于200 K穩定性Au(Ru)＞Pt＞Cu＞Pd＞Ag＞Ni，300 K 時 Au＞Ru(Cu)＞Pt(Pd)＞Ag(Ni)。而Pt1/Pt(111)在實驗中160 K時形成二聚體，210 K時形成三聚體[10]，與本工作中模擬結果(165 K時形成二聚體，200 K形成三聚體)比較接近。

圖4 M1/M(111)模型低溫時表面單原子發生團聚現象,溫度與團簇原子數的關系Fig.4 Relationship between temperature and atom number of cluster when the aggregation of M1/M(111)model accurs at low temperature

2.3 Ni1/Ni(111)單原子金屬模型催化劑在H 2(O2)氣氛下的穩定性研究

在實驗上，催化劑一般是在各種氣氛下參與反應，如Ni可以高效地使甲烷部分氧化。因此我們也研究了H2和O2氣氛是如何影響模型催化劑(Ni1/Ni(111))的穩定性。100 K 時(圖 5)，Ni1/Ni(111)催化劑在H2氣氛和O2氣氛下保持穩定，在真空環境中Ni單原子則會聚集。在H2氣氛下，對Ni1/Ni(111)進行動力學研究。發現100 K的模擬溫度下，時間進行至9 ps時出現第1個H2分子分解(圖S1)，至24 ps時出現第2個H2分子分解，之后則無明顯變化，并且Ni1/Ni(111)模型催化劑保持穩定。300 K條件時，有更多的H2分解，同時催化劑無法保持穩定，表面孤立的Ni單原子團聚成簇。在O2氣氛下，負載的Ni單原子移動受到O2分子阻礙，在100 K時保持穩定。當模擬溫度設定為300 K時，升溫過程中出現O2分子向鄰近Ni原子移動的情況(圖S1d)，持續升溫會出現O2分解現象，隨后Ni單原子發生聚集。體系升至300 K，全部Ni單原子團聚成簇。而隨著溫度繼續升高至400 K以上，表面Ni原子會脫出表面，分解的少部分O原子則會嵌入表面空穴中。與真空狀態模擬結果相比較(圖S2)，H2(尤其是O2)的存在，一定程度上提高了Ni1/Ni(111)的熱穩定性。其原因可能是Ni與O2(或O)之間強的相互作用阻礙了聚合現象的發生。

圖5 在H2(a)和O2氣氛(b)下，Ni1/Ni(111)催化模型在100、300、400、500和600 K時的動力學模擬結果Fig.5 Dynamics results of Ni1/Ni(111)catalyst model under H2(a)and O2(b)atomosphere at 100,300,400,500 and 600 K

3 結 論

在本工作中，通過反應性力場分子動力學模擬研究了金屬載體的單原子催化劑(M1/M(111))在300到600 K的熱穩定性。模擬結果表明，對于SAC模型，只有Fe1/Fe(100)催化劑體系可以在高的反應溫度(約600 K)中穩定存在，而別的單原子金屬催化劑會在載體表面形成金屬團簇。原因總結如下：Fe1/Fe(100)的吸附有較強的結合能，高于它的內聚能，而別的過渡金屬相反，與載體的結合能力要弱于內聚能，因此Fe1/Fe(100)可以穩定存在。對于Ni1/Ni(111)催化劑，H2和O2氣氛在一定程度上提高了它的穩定性，氣氛會阻礙Ni單原子在表面發生團聚的過程。

Supporting information is availableat http://www.wjhxxb.cn