梁 簫, 陳 錳

白云母–流體界面鹽離子分布的計(jì)算模擬研究: 對離子地球化學(xué)循環(huán)的啟示

梁 簫1, 2, 陳 錳1*

(1. 中國科學(xué)院 廣州地球化學(xué)研究所礦物學(xué)與成礦學(xué)重點(diǎn)實(shí)驗(yàn)室/廣東省礦物物理與材料研究開發(fā)重點(diǎn)實(shí)驗(yàn)室, 廣東 廣州 510640; 2. 中國科學(xué)院大學(xué), 北京 100049)

黏土礦物是表生土壤、沉積物、碎屑巖結(jié)晶顆粒表面的常見礦物, 常與地質(zhì)流體處于共存狀態(tài)。離子–水–黏土礦物的相互作用影響著表生環(huán)境的流體組成、元素輸運(yùn)、土壤有機(jī)質(zhì)分餾和沉積物的力學(xué)特性等。通過分子動力學(xué)模擬認(rèn)識離子在黏土礦物–流體界面的分布, 是揭示黏土礦物參與離子地球化學(xué)循環(huán)過程的基礎(chǔ)。以白云母為研究對象, 使用分子動力學(xué)模擬和熱力學(xué)計(jì)算方法, 分析了白云母–堿金屬鹽溶液界面的離子吸附結(jié)構(gòu)和靜電作用特性。研究發(fā)現(xiàn)Na+離子在白云母表面形成了兩個(gè)內(nèi)球絡(luò)合吸附層, K+、Cs+離子僅形成一個(gè)內(nèi)球絡(luò)合吸附層。內(nèi)球絡(luò)合吸附層外可以形成多個(gè)外球絡(luò)合吸附層, 表現(xiàn)為離子的層狀吸附結(jié)構(gòu)。增大鹽濃度沒有改變同種離子的吸附結(jié)構(gòu), 僅改變吸附層分布的大小。定量對比了使用傘狀采樣方法計(jì)算的離子吸附自由能曲線和基于離子密度分布得到的自由能曲線, 二者總體一致, 基本確證了模擬采樣的充足性。研究表明界面流體不能視為均一的介電環(huán)境。

白云母; 流體; 堿金屬離子; 自由能; 模擬

0 引 言

離子在礦物–流體界面的分配與自然過程息息相關(guān)。黏土礦物是表生土壤、沉積物、碎屑巖結(jié)晶顆粒表面的常見礦物, 常常與地質(zhì)流體處于共存狀態(tài)(Stoffelsma et al., 2010; Tinnacher et al., 2016; Resasco et al., 2017; Thomas and Peter, 2018; Magnussen and Gross, 2019; Waegele et al., 2019)。黏土礦物由帶負(fù)電荷的片層組成, 并通過層間域或外表面的離子吸附實(shí)現(xiàn)電中和(Rotenberg et al., 2007; Lee et al., 2013)。黏土礦物表面的離子吸附和交換控制著含水層中的流體成分, 并影響有毒和放射性廢物的地質(zhì)儲存(Rotenberg et al., 2007; Sakuma and Kawamura, 2011; Tinnacher et al., 2016; Waegele et al., 2019)。離子–水–黏土礦物的相互作用影響表生環(huán)境的流體組成、元素輸運(yùn)、土壤有機(jī)質(zhì)分餾和沉積物的力學(xué)特性等(Bourg et al., 2017; Brown et al., 2017)。因此, 揭示黏土礦物的層間和外表面離子的競爭吸附、交換機(jī)理十分重要。過去通過結(jié)合分子動力學(xué)模擬和熱力學(xué)計(jì)算, 揭示了黏土礦物層間域的離子交換機(jī)理(Zhou et al., 2020)。針對白云母的實(shí)驗(yàn)研究較多, 許多實(shí)驗(yàn)方法已經(jīng)用于黏土礦物–電解質(zhì)溶液界面離子的研究中, 如X射線反射、原子力顯微鏡等(Lee et al., 2013, 2017, 2021; Zou et al., 2021)。本研究以白云母為對象, 使用分子動力學(xué)模擬和熱力學(xué)計(jì)算方法, 分析白云母–堿金屬鹽溶液界面的離子吸附結(jié)構(gòu)和靜電作用特性, 探討?zhàn)ね恋V物外表面對堿金屬離子的吸附和分配作用。

雙電層結(jié)構(gòu)(electric double layer, EDL)是常用來描述礦物表面離子行為的模型(Shchukin et al., 1958), 帶正電的離子被吸附在黏土礦物表面形成吸附層(Stern層), 吸附層外為擴(kuò)散層(Gouy-Chapman層)。根據(jù)Gouy-Chapman-Stern(GCS)模型, 靜電勢在吸附層內(nèi)線性下降, 而在擴(kuò)散層中呈指數(shù)下降。GCS模型依賴于點(diǎn)電荷假設(shè)和連續(xù)介電環(huán)境假設(shè), 連續(xù)介電環(huán)境假設(shè)擴(kuò)散層中的介電常數(shù)與純?nèi)芤合嗤? 從實(shí)驗(yàn)數(shù)據(jù)和模擬結(jié)果可以看出GCS模型的局限性, 即使是相同電荷的離子, 不同的離子半徑或離子類型會造成吸附層的結(jié)構(gòu)差異(Strmcnik et al., 2009; Ricci et al., 2014; Lee et al., 2017; Resasco et al., 2017)。不同的堿金屬離子, 雖然具有同樣的離子電荷, 但在白云母–流體界面吸附時(shí)與水分子形成的絡(luò)合結(jié)構(gòu)存在差異(Bourg et al., 2017; Kobayashi et al., 2017)。經(jīng)典密度泛函理論(density functional theory, DFT)可用于理論計(jì)算離子的密度分布和礦物–流體界面的靜電勢分布(Voukadinova and Gillespie, 2019; Cats et al., 2021; Wang et al., 2021), 這種方法優(yōu)于傳統(tǒng)的GCS模型, 因?yàn)樗呀?jīng)將離子尺寸的影響考慮在內(nèi)。但是該理論仍然以均一的介電環(huán)境作為假設(shè), 這與實(shí)際情況不符。

分子模擬被廣泛應(yīng)用于表界面的離子吸附性質(zhì)的研究, 包括第一性原理分子動力學(xué)模擬和基于經(jīng)驗(yàn)力場的分子動力學(xué)模擬。分子動力學(xué)模擬不以均一的介電環(huán)境為假設(shè), 較DFT方法更優(yōu)。使用第一性原理分子動力學(xué)方法模擬的固–液界面體系的尺寸比較小(Creazzo et al., 2018), 不能完全反映固液界面的性質(zhì), 因此, 基于經(jīng)驗(yàn)力場的模擬方法被更廣泛地應(yīng)用于這一領(lǐng)域(Rotenberg et al., 2007; Sakuma and Kawamura, 2011; Dewan et al., 2014; Bourg et al., 2017; Kobayashi et al., 2017)。過去使用經(jīng)驗(yàn)力場研究與黏土礦物表面結(jié)合的陽離子(包括Li+、Na+、K+、Cs+、Ca2+、Mg2+、Sr2+等), 這些離子在表面的吸附形態(tài)可以分為外球絡(luò)合和內(nèi)球絡(luò)合。外球絡(luò)合是指離子被水分子組成的水化殼包圍, 然后間接吸附于礦物表面; 內(nèi)球絡(luò)合是指離子直接吸附于礦物表面, 離子與礦物表面間沒有被水化殼隔斷。在堿金屬離子里, Li+、Na+離子傾向于形成外球絡(luò)合, K+、Cs+離子則傾向于形成內(nèi)球絡(luò)合, 堿土金屬族的Ca2+、Mg2+、Sr2+離子均傾向于形成外球絡(luò)合。但是, 當(dāng)黏土礦物的四面體結(jié)構(gòu)層存在較多離子替代配位的時(shí)候, 由于表面電荷增加, 內(nèi)球絡(luò)合的離子會增多(Delville, 1991; Skipper et al., 1991; Chatterjee et al., 2000; Young and Smith, 2000; Hensen et al., 2001; Sutton and Sposito, 2001, 2002; Tambach et al., 2004; Whitley and Smith, 2004; Malikova et al., 2004, 2006; Liu and Lu, 2006; Mignon et al., 2010; Suter et al., 2012; Boek, 2014; Zhang et al., 2014; Dazas et al., 2014, 2015; Vao-soongnern et al., 2015; Yang et al., 2017; Rahromostaqim and Sahimi, 2018)。

本研究通過描述黏土礦物的ClayFF力場(Cygan et al., 2004)、離子的Horinek力場(Horinek et al., 2009)和水分子的SPC/E模型(Berendsen et al., 1987), 研究不同堿金屬離子(Na+、K+、Cs+)在白云母/流體界面的分布, 使用傘狀取樣方法(Torrie and Valleau, 1977)計(jì)算離子從體相流體吸附到礦物表面的自由能, 將該自由能曲線與基于離子分布推導(dǎo)的自由能曲線進(jìn)行定量對比, 并考慮不同離子濃度對結(jié)果的影響。通過模擬計(jì)算離子分布是廣泛使用的方法, 但是, 由于模擬可能存在系綜不遍歷的問題, 模擬統(tǒng)計(jì)的離子分布可能是偏離實(shí)際的。將基于離子分布推導(dǎo)的吸附自由能和直接計(jì)算的吸附自由能做對比, 是確定系綜是否遍歷的重要手段, 也是本研究對今后模擬研究的方法學(xué)貢獻(xiàn)。

1 模擬與計(jì)算方法

白云母的化學(xué)式為K+[Al2(Si3Al)O10(OH)2]–, 本次研究的白云母為2M1多型, 白云母中各個(gè)原子的坐標(biāo)由Richardson使用X射線衍射確定(=0.52 nm,=0.90 nm,2.01 nm,90.0°,95.8°,=90.0°) (Richardson and Richardson, 1982)。根據(jù)Loewenstein規(guī)則, 四面體片的鋁原子對硅原子的替換位不相鄰(Loewenstein, 1954)。本研究構(gòu)建了一個(gè)12×6的白云母超胞, 其在平面內(nèi)尺寸為6.3 nm× 5.4 nm, 被置于一個(gè)長方體的模擬盒子中。3層云母片層堆疊于盒子的方向上, 我們將K+離子插入云母的層間域用于中和電荷。每一層云母的厚度大約為2.7 nm。模擬體系包含3層云母片層, 相當(dāng)于方向上的1.5個(gè)超胞。由于2M1型白云母具有兩個(gè)不同的(001)表面, 兩個(gè)不同表面可能存在差異, 所以本研究構(gòu)建了非整數(shù)個(gè)超胞來制造兩個(gè)不同的表面。在方向上, 厚度約為10 nm的鹽溶液層被放置于白云母片層上。模擬盒子的所有方向均使用周期性邊界條件。白云母與鹽溶液界面系統(tǒng)的示意圖見圖1。本研究模擬NaCl、KCl和CsCl 3種鹽溶液, 用于比較不同濃度、不同類型的鹽溶液的離子在白云母表面的吸附情況, 并使用ClayFF力場(Cygan et al., 2004)描述白云母片層中的原子, Horinek力場(Horinek et al., 2009)描述溶液中的離子, SPC/E模型(Berendsen et al., 1987)描述水分子。這些力場模型均使用勒讓德·瓊斯勢, 在力場形式上是相容的。在過去的研究里, 使用這些力場組合模擬蒙脫石的結(jié)果與實(shí)驗(yàn)認(rèn)識一致(Zhou et al., 2020)。

分子動力學(xué)模擬采用GROMACS 5.1軟件包(Hess et al., 2008)。模擬的時(shí)間步長為1.0 fs, 模擬采用z系綜: 溫度保持在298 K, 使用Nosé-Hoover控溫器(Nosé, 1983; Hoover, 1985)控制溫度; 使用Parrinello-Rahman控壓方法(Parrinello and Rahman, 1982), 控壓過程中調(diào)整盒子在方向的尺寸,平面尺寸不變, 壓力維持為一標(biāo)準(zhǔn)大氣壓。短程非鍵作用(范德華力和靜電力)的截?cái)喟霃綖?.2 nm。使用particle-mesh Ewald(PME)方法(Essmann et al., 1995)描述長程靜電作用。

將NaCl、KCl、CsCl 3種鹽溶液設(shè)置0.5 mol/L、1.0 mol/L和1.5 mol/L 3個(gè)濃度梯度。通過置換溶液中的水分子, 將鹽離子插入到溶液片層中, 并將一定數(shù)量的K+離子以相同的方式插入溶液來平衡白云母片層的負(fù)電荷。對每個(gè)系統(tǒng)進(jìn)行10 ns的NPT系綜模擬, 最后2 ns的模擬軌跡用來收集數(shù)據(jù)。對每個(gè)體系的平衡模擬重復(fù)5次。由于白云母/溶液體系具有2個(gè)界面, 結(jié)構(gòu)分析基于10組數(shù)據(jù)的平均值。

Si、Al、Na、K、Cl、O、H原子分別由黃色、粉色、藍(lán)色、紫色、青色、紅色和白色的小球或點(diǎn)表示。

本次研究計(jì)算了陰陽離子從云母表面到溶液內(nèi)部的平均力勢(potential of mean force, PMF)。隨機(jī)選擇一個(gè)離子并將其沿法線方向從表面拉出(使用彈簧勢, 勁度系數(shù)=500 kJ/mol), 速度為0.01 nm/ps。每20 ps記錄一個(gè)結(jié)構(gòu), 共記錄8個(gè)結(jié)構(gòu)并用它們進(jìn)行傘狀采樣模擬(Torrie and Valleau, 1977)。對每個(gè)傘狀采樣的窗口進(jìn)行0.2 ns的隨機(jī)分子動力學(xué)模擬, 使用彈簧勢(=500 kJ/mol)使離子的位置圍繞目標(biāo)位置波動。上述操作在每個(gè)系統(tǒng)中重復(fù)了10次, 然后使用加權(quán)直方分析方法(weighted histogram analysis method, WHAM)分析數(shù)據(jù)并得到PMF曲線及其誤差(Hub et al., 2010)。

2 結(jié)果與分析

2.1 離子密度分布

通過計(jì)算各種離子在模擬體系方向上的數(shù)量密度(圖2、3), 可以清晰地反映不同鹽濃度條件下, 離子在白云母表面流體中的分布情況。密度分布圖的橫坐標(biāo)為離子的坐標(biāo)位置與白云母表面位置之差, 其中白云母表面位置定義為白云母表面O原子的平均位置。陽離子的密度分布圖(圖2a～c)展示了與其他研究發(fā)現(xiàn)的位置相似的密度峰(Bourg et al., 2017; Kobayashi et al., 2017)。Bourg et al. (2017)和 Kobayashi et al. (2017)研究指出, 距離白云母表面0.3 nm之內(nèi)的密度峰代表了陽離子和水分子、表面基團(tuán)形成了內(nèi)球絡(luò)合物。陽離子密度分布反映了類似的內(nèi)球絡(luò)合結(jié)構(gòu)。不同濃度體系中同種陽離子的分布形式在定性上是一致的(圖3a～c), K+和Cs+離子在表面0.3 nm內(nèi)只形成了一個(gè)顯著的單密度峰; Na+離子則形成兩個(gè)密度峰。Na+的第1和第2個(gè)峰分別代表了位于硅氧骨架表面復(fù)三方孔洞中的Na+和結(jié)構(gòu)表面帶電氧原子正上方的Na+(Bourg et al., 2017)。K+和Cs+離子的單密度峰均代表位于表面復(fù)三方孔洞中的陽離子。在濃度較高的情況下(如1.5 mol/L), Na+離子在距離表面約0.5 nm處形成了一個(gè)寬而小的峰(圖3a), 這代表外球絡(luò)合的離子水化結(jié)構(gòu)。對相同濃度體系中不同陽離子的密度分布(圖2a～c)進(jìn)行比較時(shí)發(fā)現(xiàn), 在0.5 mol/L、1.0 mol/L、1.5 mol/L濃度條件下, 峰值大小順序均為: K+> Cs+> Na+; 濃度增大, Cs+與K+的峰值逐漸接近。隨著鹽濃度上升, Na+和K+的密度峰逐漸遠(yuǎn)離云母表面(圖3a、b), 但在Cs+的密度分布圖中觀察不到這種現(xiàn)象(圖3c)。

圖2 離子在白云母–流體界面的密度分布

圖3 離子在白云母–流體界面的密度分布

通過計(jì)算陰離子(Cl?)的密度分布情況可見(圖2d～f、圖3d～f), 在不同陽離子種類和鹽濃度情況下, Cl?都會呈現(xiàn)出兩個(gè)峰, 出現(xiàn)的位置大致相同。在不同濃度的同種鹽體系中, 隨著鹽濃度的增加, 峰值逐漸增大(圖3d～f)。相同濃度的不同陽離子鹽體系中Cl?的分布情況見圖2d～f, 兩個(gè)Cl?峰出現(xiàn)的位置基本一樣, 峰值的大小有如下關(guān)系: 第1個(gè)峰CsCl > KCl > NaCl, 第2個(gè)峰NaCl > KCl > CsCl。

2.2 離子脫附的平均力勢

陽離子和陰離子從云母表面到流體內(nèi)部的自由能變化, 記為+()和?(), 它們和密度分布的關(guān)系為:

式中:ρ()和?()分別為陽離子和陰離子在位置的密度;b +和b ?為它們在體相的密度;為玻爾茲曼常數(shù);為溫度。利用2.1節(jié)得到的密度分布, 我們可以據(jù)公式(1)、(2)推算+()和?()。但是, 由于密度分布可能受模擬的隨機(jī)性和采樣不足的影響, 可能和實(shí)際情況存在偏差。為此, 同時(shí)使用了傘狀采樣的方法直接計(jì)算+()和?(), 并與以密度分布推算的結(jié)果進(jìn)行對比。這樣, 一方面驗(yàn)證密度分布結(jié)果的可靠性; 另一方面, 確定離子從表面到內(nèi)部的自由能變化。

通過傘狀取樣得到的PMF曲線基本與基于密度分布推導(dǎo)得到的PMF曲線大致吻合(圖4), 說明方法的可靠性。PMF曲線中的最小值點(diǎn)即表示離子移動到此位置時(shí), 體系的能量最低、最穩(wěn)定。圖中的PMF最小值出現(xiàn)的位置與濃度基本無關(guān), 而是取決于陽離子的種類, 這一發(fā)現(xiàn)與前文密度分布的結(jié)果一致, 可以說明不同的陽離子在云母表面有不同的吸附狀態(tài)和偏好。PMF曲線的局部極小值點(diǎn)代表亞穩(wěn)定的狀態(tài)。從同種陽離子不同鹽濃度的PMF曲線中可以看出: Na+的PMF曲線有4個(gè)極小值點(diǎn), 而K+和Cs+則分別只有兩個(gè)極小值點(diǎn)。通過對極小值點(diǎn)出現(xiàn)位置與密度分布的結(jié)果進(jìn)行對比, 可以看出Na+的前兩個(gè)極小值點(diǎn)分別與密度分布中的前兩個(gè)峰相對應(yīng); K+與Cs+的第1個(gè)極小值點(diǎn)對應(yīng)了密度分布中最明顯的峰, 它們代表形成內(nèi)球絡(luò)合結(jié)構(gòu)的陽離子。Na+的PMF曲線上的第3個(gè)極小值點(diǎn)、K+與Cs+的第2個(gè)極小值點(diǎn)對應(yīng)于形成外球絡(luò)合結(jié)構(gòu)的陽離子。對于0.5 mol/L的NaCl的體系, Na+的PMF曲線還明顯出現(xiàn)了第4個(gè)極小值點(diǎn)(圖4a), 該極小值點(diǎn)隨溶液的鹽濃度增大而變得不明顯(圖4b、4c)。這說明了陽離子在白云母表面為多層吸附結(jié)構(gòu), 而PMF曲線較之密度分布圖更明顯地展示了這種多層吸附結(jié)構(gòu)。Martin-Jimenez et al. (2016)研究中, 高分辨原子力顯微鏡也揭示過這種多層吸附結(jié)構(gòu)。

圖4 離子從白云母表面到溶液內(nèi)部脫附的平均力勢

此外, 通過計(jì)算還得到了對應(yīng)體系中Cl?的PMF曲線, 發(fā)現(xiàn)Cl?在不同離子種類和濃度體系中都存在兩個(gè)最小值點(diǎn), 出現(xiàn)在=0.5 nm兩側(cè), 這與密度分布的結(jié)果相符。Cl?離子的PMF曲線也展現(xiàn)了離子多層吸附的特點(diǎn)。

通過傘狀取樣得到的PMF曲線與基于密度分布推導(dǎo)得到的PMF曲線在定量上存在一定的偏差, 尤其是距離表面較近處的陽離子的內(nèi)球絡(luò)合吸附層(圖4)。這種偏差可能是采樣不足所導(dǎo)致的。如果增加傘狀采樣的重復(fù)次數(shù)和通過平衡模擬計(jì)算密度分布的次數(shù), 這種偏差將可以被消除。

2.3 白云母表面流體的靜電勢

雙電層理論和經(jīng)典密度泛函理論對礦物–流體界面離子分布的描述皆以均一的介電常數(shù)()為假設(shè)(Shchukin et al., 1958; Voukadinova and Gillespie, 2019; Cats et al., 2021; Wang et al., 2021)。我們基于分子動力學(xué)模擬分析的結(jié)果來驗(yàn)證這個(gè)假設(shè)與實(shí)際情況的相符程度。

式中:0為真空中的介電常數(shù);為元電荷量;z為原子或離子的帶電荷量;ρ()為原子或離子的密度分布。

圖5為靜電勢的分布圖。圖中可見,()隨著表面距離增大, 先是到達(dá)峰值, 然后以近似指數(shù)方式下降。()則呈現(xiàn)波動式的衰減, 反映了水分子在白云母表面的分層結(jié)構(gòu), 是由離子水合形成的多層吸附所造成的結(jié)果。如果白云母表面的流體為均一的介電環(huán)境,()和()隨距離的變化方式應(yīng)當(dāng)相似。二者在定性上的差別說明介電常數(shù)在界面上的分布為非均一的。本質(zhì)上, 這種非均一的介電環(huán)境, 是由于離子與水分子結(jié)合, 在白云母外表面形成多層吸附。這種多層吸附的特性造成了介電環(huán)境的波動。

圖5 靜電勢(f和ψ)在白云母–流體界面的分布

3 討 論

在方法上, 本研究通過結(jié)合基于平衡模擬計(jì)算得到的密度分布和使用傘狀采樣方法得到的PMF曲線共同揭示了離子在白云母表面的分布狀態(tài)。過去的模擬研究既有通過平衡模擬計(jì)算得到密度分布(Bourg et al., 2017), 也有通過采樣方法計(jì)算自由能(Kobayashi et al., 2017), 但沒有將二者的分析進(jìn)行定量對比。本研究對比利用密度分布推導(dǎo)的PMF曲線與使用傘狀采樣方法得到的PMF曲線, 確證了在模擬采樣充足的情況下, 二者皆可準(zhǔn)確揭示離子在礦物表面的分布狀態(tài)。PMF曲線較之密度分布圖, 更好地展示了離子的外球絡(luò)合吸附層的存在。這種分析方法的結(jié)合, 有望今后繼續(xù)用于礦物/流體界面物質(zhì)分布的分析。

本研究說明礦物–流體界面的介電常數(shù)并非是均一的, 這對雙電層理論和經(jīng)典密度泛函理論所應(yīng)用的均一介電常數(shù)假設(shè)(Shchukin et al., 1958; Voukadinovaand Gillespie, 2019; Cats et al., 2021; Wang et al., 2021)提出了質(zhì)疑。由于離子在礦物表面形成多層吸附結(jié)構(gòu), 每個(gè)吸附層的離子都與水分子在一定程度上形成了水化殼結(jié)構(gòu)。這樣, 將水環(huán)境視為均一的介電環(huán)境與實(shí)際情況不符。今后應(yīng)當(dāng)注意將介電性質(zhì)隨表面位置的變化納入考慮中。當(dāng)然, 仍然需要在理論上探索準(zhǔn)確計(jì)算介電常數(shù)的方法。

通過實(shí)驗(yàn)方法可以認(rèn)識礦物–流體界面離子的分布、離子與表面的結(jié)合狀態(tài)等。分子動力學(xué)模擬的認(rèn)識和實(shí)驗(yàn)認(rèn)識形成對照, 這可以在現(xiàn)象與原理兩方面確立礦物–流體界面的精細(xì)結(jié)構(gòu)認(rèn)識。本研究在分子動力學(xué)模擬基礎(chǔ)上采用熱力學(xué)計(jì)算, 可以用于計(jì)算不同離子在白云母表面的分配關(guān)系, 具有重要前景。

4 結(jié) 論

(1) 對于不同濃度鹽體系Na+、K+、Cs+離子在白云母表面的吸附, Na+離子形成了兩個(gè)內(nèi)球絡(luò)合層, 而K+、Cs+離子形成一個(gè)內(nèi)球絡(luò)合層, 內(nèi)球絡(luò)合層外可以形成多個(gè)外球絡(luò)合層, 表現(xiàn)為離子的層狀吸附結(jié)構(gòu)。

(2) 增大鹽濃度在定性上沒有改變同種離子的吸附結(jié)構(gòu), 僅在定量上改變吸附層分布的大小。

(3) 白云母–流體界面的介電常數(shù)不能以均一的介電常數(shù)來表示, 介電常數(shù)應(yīng)當(dāng)是白云母表面距離的函數(shù)。

致謝:感謝西南科技大學(xué)張紅平副教授和西北農(nóng)林科技大學(xué)張弛副教授給出的專業(yè)評審意見, 對本文質(zhì)量的提升有很大幫助。

Berendsen H, Grigera J R, Straatsma T P. 1987. The missing term in effective pair potentials., 91(24): 6269–6271.

Boek E S. 2014. Molecular dynamics simulations of interlayer structure and mobility in hydrated Li-, Na- and K-montmorillonite clays., 112(9–10): 1472–1483.

Bourg I C, Lee S S, Fenter P, Tournassat C. 2017. Stern layer structure and energetics at mica-water interfaces., 121(17): 9402–9412.

Brown K M, Poeppe D, Josh M, Sample J, Even E, Saffer D, Tobin H, Hirose T, Kulongoski J T, Toczko S. 2017. Theaction of water films at ?-scales in the Earth: Implications for the Nankai subduction system., 463: 266–276.

Cats P, Evans R, Hrtel A, Roij R V. 2021. Primitive model electrolytes in the near and far field: Decay lengths from DFT and simulations., 154(12), 124504.

Chatterjee A, Iwasaki T, Ebina T. 2000. A novel method to correlate layer charge and the catalytic activity of 2∶1 dioctahedral smectite clays in terms of binding the interlayer cation surrounded by monohydrate., 104(35): 8216–8223.

Creazzo F, Galimberti D R, Pezzotti S, Gaigeot M P. 2018. DFT-MD of the (110)-Co3O4cobalt oxide semiconductor in contact with liquid water, preliminary chemical and physical insights into the electrochemical environment., 150(4), 041721.

Cygan R T, Liang J J, Kalinichev A G. 2004. Molecular models of hydroxide, oxyhydroxide, and clay phases and the development of a general force field., 108(4): 1255–1266.

Dazas B, Ferrage E, Delville A, Lanson B. 2014. Interlayer structure model of tri-hydrated low-charge smectite by X-ray diffraction and Monte Carlo modeling in the Grand Canonical ensemble., 99(8–9): 1724–1735.

Dazas B, Lanson B, Delville A, Robert J L, Komarneni S, Michot L J, Ferrage E. 2015. Influence of tetrahedral layer charge on the organization of inter layer water and ions in synthetic Na-saturated smectites., 119(8): 4158–4172.

Delville A. 1991. Modeling the clay-water interface., 7(3): 547–555.

Dewan S, Carnevale V, Bankura A, Eftekhari-Bafrooei A, Fiorin G, Klein M L, Borguet E. 2014. Structure of water at charged interfaces: A molecular dynamics study., 30(27): 8056–8065.

Essmann U, Perera L, Berkowitz M L, Darden T, Lee H, Pedersen L G. 1995. A smooth particle mesh Ewald method., 103(19): 8577– 8593.

Hensen E J M, Tambach T J, Bliek A, Smit B. 2001. Adsorption isotherms of water in Li-, Na-, and K-montmorillonite by molecular simulation., 115(7): 3322–3329.

Hess B, Kutzner C, van der Spoel D, Lindahl E. 2008. GROMACS 4: Algorithms for highly efficient, load-balanced, and scalable molecular simulation., 4(3): 435–447.

Hoover W G. 1985. Canonical dynamics: Equilibrium phase- space distributions.,31(3): 1695–1697.

Horinek D, Herz A, Vrbka L, Sedlmeier F, Mamatkulov S I, Netz R R. 2009. Specific ion adsorption at the air/water interface: The role of hydrophobic solvation., 479(4–6): 173–183.

Hub J S, de Groot B L, van der Spoel D. 2010. g_wham—A free weighted histogram analysis implementation including robust error and autocorrelation estimates., 6(12): 3713–3720.

Kobayashi K, Liang Y F, Murata S, Matsuoka T, Takahashi S, Nishi N, Sakka T. 2017. Ion distribution and hydration structure in the stern layer on muscovite Surface., 33(15): 3892– 3899.

Lee S S, Fenter P, Nagy K L, Sturchio N C. 2013. Changes in adsorption free energy and speciation during competitiveadsorption between monovalent cations at the muscovite (001)-water interface., 123: 416–426.

Lee S S, Fenter P, Nagy K L, Sturchio N C. 2017. Real-time observation of cation exchange kinetics and dynamics at the muscovite-water interface., 8(1), 15826.

Lee S S, Koishi A, Bourg I C, Fenter P. 2021. Ion correlations drive charge overscreening and heterogeneous nucleation at solid-aqueous electrolyte interfaces., 118(32), e2105154118.

Liu X D, Lu X C. 2006. A thermodynamic understanding of clay-swelling inhibition by potassium ions., 118(38): 6448–6451.

Loewenstein W. 1954. The distribution of aluminum in the tetrahedra of silicates and aluminates., 39(1–2): 92–96.

Magnussen O M, Gross A. 2019. Toward an atomic-scale understanding of electrochemical interface structure and dynamics., 141(12): 4777–4790.

Malikova N, Cadene A, Marry V, Dubois E, Turq P. 2006. Diffusion of water in clays on the microscopic scale: Modeling and experiment., 110(7): 3206–3214.

Malikova N, Marry V, Dufreche J F, Simon C, Turq P, Giffaut E. 2004. Temperature effect in a montmorillonite clay at low hydration-microscopic simulation., 102(18): 1965–1977.

Martin-Jimenez D, Chacon E, Tarazona P, Garcia R. 2016. Atomically resolved three-dimensional structures of electrolyte aqueous solutions near a solid surface., 7(1), 12164.

Mignon P, Ugliengo P, Sodupe M, Hernandez E R. 2010.molecular dynamics study of the hydration of Li+, Na+and K+in a montmorillonite model. Influence of isomorphic substitution., 12(3): 688–697.

Nosé S C. 1983. A molecular dynamics method for simulations in the canonical ensemble., 52(2): 255–268.

Parrinello M, Rahman A. 1982. Strain fluctuations and elasticconstants., 76(5): 2662–2666.

Rahromostaqim M, Sahimi M. 2018. Molecular dynamics simulation of hydration and swelling of mixed-layer clays., 122(26): 14631– 14639.

Resasco J, Chen L D, Clark E, Tsai C, Hahn C, Jaramillo T F, Chan K, Bell A T. 2017. Promoter effects of alkali metal cations on the electrochemical reduction of carbon dioxide., 139(32): 11277–11287.

Ricci M, Spijker P, Vo?tchovsky K. 2014. Water-induced correlation between single ions imaged at the solid-liquid interface., 5(1): 1–8.

Richardson S M, Richardson J W. 1982. Crystal structure of a pink muscovite from Archer’s Post, Kenya: Implications for reverse pleochroism in dioctahedral micas., 67(1–2): 69–75.

Rotenberg B, Marry V, Vuilleumier R, Malikova N, Simon C, Turq P. 2007. Water and ions in clays: Unraveling the interlayer/micropore exchange using molecular dynamics., 71(21): 5089–5101.

Sakuma H, Kawamura K. 2011. Structure and dynamics of water on Li+-, Na+-, K+-, Cs+-, H3O+-exchanged muscovite surfaces: A molecular dynamics study., 75(1): 63–81.

Shchukin E D, Pertsov A V, Amelina E A, Zelenev A S. Colloid and Surface Chemistry. 2001. Amsterdam: Elsevier: 165–220.

Skipper N T, Refson K, McConnell J D C. 1991. Computer simulation of interlayer water in 2∶1 clays., 94(11): 7434–7445.

Stoffelsma C, Rodriguez P, Garcia G, Garcia-Araez N, Strmcnik D, Markovic N M, Koper M T C. 2010. Promotion of the oxidation of carbon monoxide at stepped platinum single-crystal electrodes in alkaline media by lithium and beryllium cations., 132(45): 16127–16133.

Strmcnik D S, Kodama K, van der Vliet D, Greeley J P, Stamenkovic V R, Markovi N M. 2009. The role of non-covalent interactions in electrocatalytic fuel-cell reactions on platinum., 1(6): 466–472.

Suter J L, Sprik M, Boek E S. 2012. Free energies of absorptionof alkali ions onto beidellite and montmorillonite surfaces from constrained molecular dynamics simulations., 91: 109–119.

Sutton R, Sposito G. 2001. Molecular simulation of interlayer structure and dynamics in 12.4 angstrom Cs-smectite hydrates., 237(2): 174–184.

Sutton R, Sposito G. 2002. Animated molecular dynamics simulations of hydrated caesium-smectite interlayers., 3: 73–80.

Tambach T J, Hensen E J M, Smit B. 2004. Molecular simulations of swelling clay minerals., 108(23): 7586–7596.

Thomas G, Peter A E. 2018. Simulating donnan equilibria based on the Nernst-Planck equation., 232: 1–13.

Tinnacher R M, Holmboe M, Tournassat C, Bourg I C, Davis J A. 2016. Ion adsorption and diffusion in smectite: Molecular, pore, and continuum scale views., 177: 130–149.

Torrie G M, Valleau J P. 1977. Nonphysical sampling distributionsin Monte Carlo free-energy estimation: Umbrella sampling., 23(2): 187–199.

Vao-soongnern V, Pipatpanukul C, Horpibulsuk S. 2015. A combined X-ray absorption spectroscopy and molecular dynamic simulation to study the local structure potassium ion in hydrated montmorillonite., 50(21): 7126–7136.

Voukadinova A, Gillespie D. 2019. Energetics of counterion adsorption in the electrical double layer., 150(15), 154706.

Waegele M M, Gunathunge C M, Li J Y, Li X. 2019. How cations affect the electric double layer and the rates and selectivity of electrocatalytic processes., 151(16), 160902.

Wang X P, Liu K, Wu J Z. 2021. Demystifying the Stern layer at a metal-electrolyte interface: Local dielectric constant, specific ion adsorption, and partial charge transfer., 154(12), 124701.

Whitley H D, Smith D E. 2004. Free energy, energy, and entropy of swelling in Cs-, Na-, and Sr-montmorillonite clays., 120(11): 5387–5395.

Yang G M, Neretnieks I, Holmboe M. 2017. Atomistic simulations of cation hydration in sodium and calcium montmorillonite nanopores., 147(8), 084705.

Young D A, Smith D E. 2000. Simulations of clay mineral swelling and hydration: Dependence upon interlayer ion size and charge., 104(39): 9163–9170.

Zhang L H, Lu X C, Liu X D, Zhou J H, Zhou H Q. 2014. Hydration and mobility of interlayer ions of (Na, Ca)-montmorillonite: A molecular dynamics study., 118(51): 29811–29821.

Zhou H J, Chen M, Zhu R L, Lu X C, Zhu J X, He H P. 2020. Coupling between clay swelling/collapse and cationic partition., 285: 78–99.

Zou Y C, Mogg L, Clark N, Bacaksiz C, Milovanovic S, Sreepal V, Hao G P, Wang Y C, Hopkinson D G, Gorbachev R, Shaw S, Novoselov K S, Raveendran- Nair R, Peeters F M, Lozada-Hidalgo M, Haigh S J. 2021. Ion exchange in atomically thin clays and micas., 20(12): 1677–1682.

Molecular simulation study of distribution of salt ions at the muscovite-fluid interface: Implications to the ionic geochemical cycle

LIANG Xiao1, 2, CHEN Meng1*

(1. CAS Key Laboratory of Mineralogy and Metallogeny / Guangdong Provincial Key Laboratory of Mineral Physics and Materials, Guangzhou Institute of Geochemistry, Chinese Academy of Sciences, Guangzhou 510640, Guangdong, China; 2. University of Chinese Academy of Sciences, Beijing 100049, China)

Clay minerals are commonly found in soil and sediments. They generally coexist with geological fluid. The interaction between ions, water, and clay minerals affects fluid composition in the environment, transportation of elements, fractionation of organic matter, and mechanical properties of sediments. We investigated the distribution of ions at the clay mineral-fluid interface using molecular dynamics simulations to reveal the geochemical cycle of ions. We studied the ionic adsorption structure and electrostatic properties at the muscovite-alkali fluid interface through molecular dynamics simulations and thermodynamic calculations. We found that Na+ions formed two inner-sphere-complex adsorption layers on the surface of muscovite, whereas K+and Cs+ions only formed one inner-sphere-complex adsorption layer. The inner-sphere-complex adsorption layer was followed by several outer-sphere-complex adsorption layers, resulting in a layered adsorption structure of ions. Increasing the salt concentration did not change the adsorption structure when the salt type was maintained, even though the adsorption amount was altered. We compared desorption free energies derived using density distributions and the umbrella sampling method and found that they were consistent, which confirmed the sampling sufficiency. We found that the fluid on the surface cannot be regarded as a uniform dielectric environment.

muscovite; fluid; alkali ion; partition; simulation

P574; P579

A

0379-1726(2022)03-0261-10

10.19700/J.0379-1726.2022.03.002

2021-12-24;

2022-02-21

中國科學(xué)院青年創(chuàng)新促進(jìn)會資助項(xiàng)目(2019345)資助。

梁簫(1996–), 男, 碩士研究生, 地質(zhì)工程專業(yè)。E-mail: 740631451@qq.com

陳錳(1987–), 男, 副研究員, 主要從事計(jì)算礦物學(xué)研究。E-mail: chenmeng@gig.ac.cn