尹玉明 趙伶玲 高 騰

(東南大學(xué)能源熱轉(zhuǎn)換及其過(guò)程測(cè)控教育部重點(diǎn)實(shí)驗(yàn)室, 南京 210096)

水和scCO2組成的混合流體(water-bearing supercritical CO2fluids，WBSF)環(huán)境對(duì)Mg2SiO4表面的吸附層特性具有重要影響.已有研究表明，WBSF中水量的增加導(dǎo)致了Mg2SiO4表面薄水膜厚度的增大，薄水膜厚度的增加提高了Mg2SiO4表面與WBSF的反應(yīng)速率和反應(yīng)限度[5].Felmy等[6]通過(guò)實(shí)驗(yàn)觀測(cè)到水的飽和scCO2溶液與Mg2SiO4納米顆粒反應(yīng)產(chǎn)生了MgCO3沉淀，而在scCO2飽和水溶液中僅有水合MgCO3出現(xiàn).Kerisit等[7]采用密度泛函理論(DFT)研究了Mg2SiO4表面的薄水膜特性，結(jié)果顯示影響巖石表面薄水膜性質(zhì)的關(guān)鍵參數(shù)包括H2O與CO2分子在巖石表面的競(jìng)爭(zhēng)吸附能差(即H2O取代巖石表面吸附CO2的驅(qū)動(dòng)力)及H2O分子表面覆蓋率的函數(shù).然而，H2O和CO2分子在Mg2SiO4表面的吸附性質(zhì)差異的微觀機(jī)制尚未明晰.

此外，關(guān)于巖石的溶解規(guī)律國(guó)內(nèi)外學(xué)者們通過(guò)實(shí)驗(yàn)和模擬的方法開(kāi)展了廣泛研究.實(shí)驗(yàn)方面，Giammar等[8]應(yīng)用X射線衍射研究了碳酸溶液中Mg2SiO4的溶解過(guò)程，發(fā)現(xiàn)在一定的溫度(303.15～368.15 K)和壓強(qiáng)(0.1～10 MPa)范圍內(nèi)，Mg2SiO4中的Mg2+和Si4+以化學(xué)計(jì)量比同時(shí)溶解，溫度和壓強(qiáng)的升高均增加了Mg2SiO4的溶解率；但隨著溶液中的碳酸被Mg2SiO4溶解產(chǎn)物中和，溶液pH值不斷升高，Mg2SiO4溶解速率逐漸減緩.Pokrovsky等[9]研究了混合反應(yīng)器中Mg2SiO4在不同濃度的NaCl溶液中的溶解反應(yīng)速率，結(jié)果表明，Mg2SiO4表面的Mg2+和Si4+溶解速率滿足一級(jí)反應(yīng)速率方程，鹽離子濃度的增加可以提高M(jìn)g2SiO4在初期的溶解速率.模擬方面，分子模擬已成為研究固-液界面相間作用的有效手段.Wang等[10]采用分子動(dòng)力學(xué)的方法研究了硅酸鈣在Na2SiO4溶液中的溶解特性，發(fā)現(xiàn)硅酸鈣表面具有親水性，其溶解是由于Ca—O鍵的破壞.Manzano等[11]采用分子模擬的方法研究了常溫常壓下硅酸鈣表面的水合及溶解過(guò)程，通過(guò)對(duì)礦物質(zhì)表面能和水吸附能的計(jì)算，發(fā)現(xiàn)在硅酸鈣晶胞(100)表面，H+離子以游離態(tài)隨H2O滲透到晶體中，在硅酸鈣水合物層中發(fā)生硅酸鈣的溶解反應(yīng).綜上可知，鹽離子對(duì)Mg2SiO4溶解過(guò)程的影響尚需進(jìn)一步研究.

本文采用分子動(dòng)力學(xué)模擬的方法，分別開(kāi)展了H2O、CO2分子在Mg2SiO4表面的吸附以及鹽離子對(duì)Mg2SiO4溶解過(guò)程的影響研究.本文通過(guò)計(jì)算分子吸附能、徑向分布函數(shù)(RDF)、平均力勢(shì)、溶解活化能等，首先分析了H2O、CO2分子在Mg2SiO4表面吸附性能存在差異的微觀機(jī)制，然后從原子尺度探討了溶液中的鹽離子對(duì)Mg2SiO4表面溶解的影響，對(duì)CO2地質(zhì)封存選址和提高其安全性具有理論指導(dǎo)意義.

1 模擬方法及驗(yàn)證

1.1 計(jì)算模型與方法

圖1 Mg2SiO4 (010)表面H2O、CO2吸附和Mg2SiO4溶解的模型示意圖

(1)

式中，kb為化學(xué)鍵的鍵伸縮彈力常數(shù);b為原子間的化學(xué)鍵長(zhǎng)度,b0為其平衡鍵長(zhǎng)；ka為鍵角彎曲的彈力常數(shù);θ為原子間形成的鍵角大小,θ0為其平衡時(shí)的角度；N為系統(tǒng)模型中原子和離子的總數(shù);rij為原子i和j間的距離;εij和σij為原子間混合勢(shì)參數(shù),采用Lorentz-berthelot混合法則[17]計(jì)算;qi、qj分別為原子i和j所帶電荷;ε0為真空介電常數(shù).模擬系統(tǒng)采用周期性邊界條件，利用PME(particle mesh Eward)方法[18]計(jì)算靜電相互作用.范德華作用勢(shì)和靜電作用勢(shì)的位能截?cái)喟霃骄鶠?.2 nm，速度積分采用Velocity-Verlet算法，考慮到模擬的計(jì)算誤差和模型的弛豫時(shí)間，模擬的時(shí)間步長(zhǎng)設(shè)為1 fs.

系統(tǒng)模型首先應(yīng)用定溫定壓系綜(NPT)對(duì)初始模型進(jìn)行了10 ns的平衡模擬，然后在正則系綜(NVT)下運(yùn)行30 ns.模擬采用Parrinello-Rahman恒壓方法控制系統(tǒng)壓力，采用Nose-Hoover熱浴方法控制系統(tǒng)溫度.在模擬運(yùn)行中每隔50 ps保存一次模型原子坐標(biāo)和系統(tǒng)能量，用于數(shù)據(jù)處理.本文采用軟件Lammps開(kāi)展模擬計(jì)算.

1.2 模型驗(yàn)證

本文根據(jù)H2O、CO2密度和Mg2SiO4晶體彈性模量的實(shí)驗(yàn)結(jié)果，對(duì)所建模型和計(jì)算方法進(jìn)行了驗(yàn)證.NPT系綜(壓強(qiáng)20 MPa, 溫度325～400 K)下H2O、CO2密度的模擬計(jì)算結(jié)果和美國(guó)標(biāo)準(zhǔn)與技術(shù)研究所(NIST)的實(shí)驗(yàn)數(shù)據(jù)[19]示于圖2(a).由圖可知，本文模擬得到的H2O、CO2密度隨溫度的升高而降低，與已有實(shí)驗(yàn)值誤差較小；在本文吸附模擬的溫度(350 K)下，H2O、CO2密度的模擬計(jì)算結(jié)果與實(shí)驗(yàn)值間的誤差分別為1.42%和5.97%，與實(shí)驗(yàn)值吻合較好，驗(yàn)證了H2O和CO2模型的準(zhǔn)確性.此外，NPT系綜(壓強(qiáng)20 MPa, 溫度350 K)下Mg2SiO4晶體拉伸應(yīng)力-應(yīng)變(σ-ε)曲線的模擬計(jì)算結(jié)果示于圖2(b).由圖可看出，Mg2SiO4晶體在拉伸過(guò)程中經(jīng)歷了彈性形變(ε≤ 0.07)、塑性形變(0.07<ε≤0.125)和斷裂(ε>0.125)3個(gè)階段；同時(shí)對(duì)彈性形變階段曲線進(jìn)行擬合,可得出Mg2SiO4晶體彈性模量E約為293 GPa,該值與Liu等[20]的實(shí)驗(yàn)值(287 GPa)的誤差僅為2.10%，驗(yàn)證了Mg2SiO4晶體模型的準(zhǔn)確性.綜上可知，本文模型力場(chǎng)和計(jì)算方法的可靠性得到驗(yàn)證.

(a) H2O和CO2的密度

(b) Mg2SiO4晶體的應(yīng)力-應(yīng)變曲線

2 結(jié)果與討論

本文首先分別研究了在CO2地質(zhì)封存環(huán)境(壓強(qiáng)20 MPa、溫度350 K)下H2O、CO2在Mg2SiO4(010)表面的吸附過(guò)程，并基于上述研究結(jié)果，建立了Mg2SiO4表面溶解模型，探討了鹽離子對(duì)Mg2SiO4表面溶解活化能的影響.

2.1 H2O、CO2在Mg2SiO4表面的吸附

(2)

圖3 不同覆蓋度x下H2O、CO2分子在Mg2SiO4(010)表面的吸附能

統(tǒng)更加穩(wěn)定.據(jù)此可知，當(dāng)Mg2SiO4表面附近存在大量的H2O和CO2時(shí)，Mg2SiO4(010)表面更傾向于吸附H2O分子.

為研究H2O、CO2分子在Mg2SiO4表面附近吸附層的微觀結(jié)構(gòu)，本文計(jì)算了Mg2SiO4(010)表面上Mg2+周?chē)鶫2O和CO2分子的RDF，即g(r)Mg-Ow、g(r)Mg-C.同時(shí)，為進(jìn)一步分析H2O、CO2分子在Mg2SiO4(010)表面上解吸附的平均自由能壁壘ΔEH，本文根據(jù)RDF統(tǒng)計(jì)結(jié)果，計(jì)算了Mg2SiO4(010)表面的Mg2+與H2O、CO2分子間的平均力勢(shì)U(r)：

U(r)=-NAkBTln(g(r))

(3)

式中，r為H2O、CO2分子距Mg2SiO4(010)表面上Mg2+質(zhì)心平面的距離，nm；NA為阿伏伽德羅常數(shù)6.022 140 76×1023mol-1；kB為玻爾茲曼常數(shù)1.380 649×10-26kJ/K；T為系統(tǒng)溫度，K；g(r)為Mg2SiO4(010)表面上Mg2+周?chē)鶫2O、CO2的徑向分布函數(shù)在r處的值.

Mg2SiO4(010)表面上Mg2+周?chē)鶫2O或CO2分子的RDF和平均力勢(shì)的計(jì)算結(jié)果見(jiàn)圖4.由圖可見(jiàn)，Mg2SiO4表面上Mg2+周?chē)鶫2O分子的RDF曲線在r=0.222 nm和r=0.365 nm處形成2個(gè)峰，這表明H2O分子在Mg2SiO4表面存在2個(gè)吸附層，具有雙層吸附的特征；而Mg2SiO4表面上Mg2+周?chē)鶦O2分子的RDF曲線僅在r=0.328 nm處形成單峰，且該峰值僅為g(r)Mg-Ow第一峰值的0.23倍，甚至小于g(r)Mg-Ow的第二峰值，表明CO2在Mg2SiO4表面僅存在1個(gè)吸附層，具有單層吸附的特征，在Mg2SiO4(010)表面的吸附能力弱于H2O，與上述吸附能計(jì)算結(jié)果相符.

圖4 Mg2SiO4(010)表面上Mg2+周?chē)鶫2O、CO2分子的徑向分布函數(shù)與平均力勢(shì)

2.2 鹽離子對(duì)Mg2SiO4表面溶解的影響

2.2.1 Mg2SiO4表面溶解的活化能

Mg2+和Si4+的溶解速率采用原子的振動(dòng)半徑R確定，即在固定溫度下，Mg2SiO4晶體內(nèi)原子的振動(dòng)半徑為R，當(dāng)Mg2SiO4表面的原子某時(shí)刻的位置與初始時(shí)刻位置的距離大于2R時(shí)，表明該原子已脫離晶體，發(fā)生溶解.根據(jù)文獻(xiàn)[21]計(jì)算方法，本文將Mg2SiO4表面第1層Mg2+和Si4+的厚度假定為1，則Mg2SiO4表面原子的溶解反應(yīng)速率v為表面單層離子完全溶解所需時(shí)間的倒數(shù).Mg2+、Si4+的溶解活化能EA根據(jù)阿倫尼烏斯公式計(jì)算：

(4)

式中，A為指前因子，ns-1；EA為單層離子的溶解活化能，kJ/mol.對(duì)式(4)中反應(yīng)速率的自然對(duì)數(shù)ln(v)與反應(yīng)溫度的倒數(shù)T-1進(jìn)行線性擬合，即可得到EA.

EA與溫度無(wú)關(guān)，但溫度越高，Mg2+和Si4+的溶解反應(yīng)速率v越大，為提高M(jìn)g2+和Si4+的v，節(jié)約計(jì)算資源，本文分別在溫度1 600、1 700、1 800、1 900、2 000 K下進(jìn)行了溶解模擬.根據(jù)模擬結(jié)果，本文首先統(tǒng)計(jì)了Mg2SiO4表面Mg2+、Si4+在不同溫度下的振動(dòng)半徑R，據(jù)此結(jié)果計(jì)算了相應(yīng)溫度下Mg2+、Si4+的v，然后根據(jù)式(4)擬合得到Mg2+、Si4+在不同溶液中的EA和A.本文以Na2CO3溶液為例，將Mg2+和Si4+在不同溫度下的R及在Na2CO3溶液中的EA擬合結(jié)果示于圖5.由圖可知，隨著溫度的升高，Mg2SiO4晶體中Mg2+、Si4+的R不斷增加，且Mg2+的R大于Si4+的R.此外，隨溫度的升高，ln(v)呈現(xiàn)增大的趨勢(shì)，表明Mg2SiO4表面Mg2+和Si4+的溶解速率隨溫度的升高而加快.Mg2SiO4表面Mg2+和Si4+在Na2CO3溶液中的EA分別為125.8和215.4 kJ/mol.

圖5 Mg2SiO4晶體中Mg2+、Si4+的波動(dòng)半徑及其在Na2CO3溶液中溶解活化能的阿倫尼烏斯公式擬合

2.2.2 Mg2SiO4表面溶解的形態(tài)分析

表1 不同溶液中Mg2SiO4晶體Mg2+、Si4+的溶解活化能、指前因子的自然對(duì)數(shù)及350 K下溶解反應(yīng)速率

為深入分析鹽離子對(duì)Mg2SiO4表面溶解的影響過(guò)程，本文分析了溫度1 700 K、NVT模擬時(shí)長(zhǎng)10 ns時(shí)，Mg2SiO4(010)表面的Mg2+在不同鹽溶液中的排布結(jié)構(gòu)，結(jié)果見(jiàn)圖6.可以看出，在純水溶液中Mg2SiO4(010)表面尚未發(fā)生溶解，而在NaCl、KCl溶液中，Mg2SiO4表面發(fā)生輕微破壞，有少量的Na+、K+離子出現(xiàn)在Mg2SiO4(010)表面，Na+、K+在該表面沒(méi)有明顯的排布規(guī)律(見(jiàn)圖6(b)、(c)).與之相比，在CaCl2溶液中，較多的二

(a) 純水

(b) NaCl溶液

(c) KCl溶液

(d) CaCl2溶液

(e) H3OCl溶液

(a) NaCl溶液

(b) Na2CO3溶液

3 結(jié)論

2) 與在純水中相比，Mg2SiO4晶體的Mg2+、Si4+在強(qiáng)酸弱堿鹽(CaCl2)、弱酸強(qiáng)堿鹽(Na2CO3)及酸性(H3OCl)溶液中的EA均明顯下降，v350均顯著上升，而在強(qiáng)酸強(qiáng)堿鹽(NaCl、KCl)溶液中的EA和v350變化較小.強(qiáng)酸弱堿鹽、弱酸強(qiáng)堿鹽的存在可以加速CO2向碳酸鹽的轉(zhuǎn)化，提高CO2封存的安全性.