水化蒙脫石拉伸力學特性的分子動力學模擬研究

中圖分類號： TU443 最近更新：2024-11-22 DOI： 10.11835/j.issn.2096-6717.2022.088

摘要

了解蒙脫石在拉伸狀態下的力學行為在地球科學、巖土力學等領域至關重要，但現有理論和方法難以在小間距范圍內預測其水化力學性質及關鍵機理。通過編寫施加應力-計算應變Perl語言腳本，進行不同水化量蒙脫石拉伸應力下的分子動力學計算模擬與應力-應變分析，確定其不同應力階段的力學特性、相互作用機制和微觀結構演化。結果表明：蒙脫石內層水化對極限應力和拉伸模量的弱化效應明顯，且在水化初期弱化幅度會更大；體積水化膨脹主要源于晶格長度c的線性增長。Z方向的拉伸模量遠小于平面內，即應力對表面Z方向的力學行為影響最大，達到極限拉應力后，會出現整層分離的破壞現象；內層是大部分形變的主要原因，并且支配著蒙脫石的拉伸力學性能；Z方向拉應力主要造成晶格長度c和晶格角β的增大，而在X和Y方向拉應力下主要發生β的減小和增大。層電荷密度越高，結合水膜越密實，形成的氫鍵數目越多，體積和晶格長度c越小，抗拉力學性能也越強。

關鍵詞

水化蒙脫石; 分子動力學; 拉伸應力狀態; 應力-應變分析; 微結構力學

正確理解黏土礦物的力學行為，在地球科學[1]、巖土力學[2]等領域至關重要，其不僅為巖石物理建模所必需[3]，也適用于土體性質評價、開挖穩定性分析以及沉降、滑坡等地質災害的監測[4-5]。其中，最為典型的膨脹性黏土蒙脫石，在地層中大量富存，是許多土壤和巖石的主要成分，也是土木工程領域中最為重要的地質材料之一。該黏土礦物具有較強的親水性，并且內層水化后容易軟化，土壤中即使含量不多，也會對其宏觀物理力學性能產生明顯影響。而在巖土基礎設施的建造及運營過程中，這些黏土體會受到拉伸加載路徑的影響。一些強動力荷載，如地震、爆炸和沖擊等，也會帶來強大的拉應力作用[6]。

傳統的黏土力學性質研究，主要通過理論和實驗的手段進行。但是，目前的一些微結構力學理論還處于起步階段[7]，或者不適用于一些特殊的巖土問題分析，如動態拉伸荷載下的黏土層力學特征。此外，黏土礦物的晶層溶脹是短程的（層間距低于20 ?（10-1 nm）），代表了初始的水化機制，但現有的理論和方法，難以描述這一小間距范圍內的膨脹行為[8]。而對于大部分均勻和壓實的黏土來說，水卻主要存在于該區域[9]，并且呈現水化量低、高致密性和膨脹性的特點。盡管實驗的方法，可以幫助獲得不同應力狀態下的力學響應及結構變化，但卻缺乏定量的解釋[10]。納米壓痕和原子力顯微鏡可以被用于獲取黏土的微觀力學性質[11]，但可用的數據卻非常有限，黏土礦物納米級的尺寸和層狀的幾何構型導致樣品制備困難。此外，在這些測量中，水化量并沒有被控制，而在黏土層尺度，力學性能卻對其十分敏感。最為嚴重的是，實驗方法會忽略黏土力學性質改變的關鍵機理，計算機模擬可以解決這些現有理論和實驗方法中的不足。

巖土領域常用的計算機仿真方法，如有限元法和離散元法等，盡管可以較好地預測宏觀屬性，但在遇到小尺度巖土問題的求解，以及其他機理性問題時，會受到諸多限制，例如，這些方法難以預測膨脹黏土的水化力學行為。分子動力學（MD）方法可以較好地解決這些問題[12]。21世紀以來，得益于分子理論、統計力學等學科的快速發展，分子動力學已成為繼實驗和理論手段之后，從分子水平了解和認知客觀世界的第3種手段[13]。由于MD方法具有小間距預測巖土行為等諸多優勢，其已被廣泛應用于黏土的性質研究。Karaborni等[14]應用分子模擬方法較早地開展了蒙脫石水化問題的研究，Cygan等[15]開發了專門針對黏土的CLAYFF分子力場，Katti等[16]應用MD方法研究了蒙脫石的初始水化機制。近年來，中國也展開了相關研究。張亞云等[17]討論了溫壓條件對蒙脫石水化的影響，李勤等[18]采用MD方法模擬了蒙脫石層間的陽離子交換，楊微等[19]基于該技術分析了黏土礦物的微觀行為，Yang等[20]探索了伊-蒙混合層黏土的溶脹機理。楊亞帆等[21]揭示了高溫下鈣蒙脫石的膨脹特性。但是，目前關于水化蒙脫石拉伸力學特性的MD研究仍非常有限，而獲得這些基本的力學性質和內在機理卻十分重要。

鑒于膨脹性黏土水化后，往往伴隨著力學性質的改變，并且其與受力狀態、自身成分以及水化量等因素密切相關。筆者通過編寫施加應力-計算應變MD計算腳本，開展采用傳統手段難以進行的不同水化量蒙脫石拉伸應力下的分子動力學計算模擬與應力-應變分析，結合微結構演化確定其不同應力階段的力學特性和微觀結構機理，并定量分析了水化量、加載方向、層電荷密度等重要因素的影響。

1 模型和計算模擬方法

1.1 水化蒙脫石模型的構建

初始模型的建立和計算腳本的運行均在具有復雜體系建模優勢和后處理優勢的Materials Studio軟件中進行。蒙脫石的基本結構單元，來自美國礦物學家晶體結構數據庫[22]，其晶格取代僅發生在中間八面體片層，化學式為Na[Si8][Al3Mg]O20（OH）4·nH2O。通過Supercell工具將基本單元沿X、Y、Z方向擴展為4a × 2b × 2c的超晶胞，并進行晶格同構取代和內層鈉離子平衡，便可以獲得攜帶負電荷的蒙脫石黏土片層，初始晶胞尺寸為20.92 ?×18.12 ?×19.2 ?，見圖1。其中，圖1（a）展示了蒙脫石模型的構建過程，圖1（b）顯示了蒙脫石的初始水化結構，圖1（c）則為其結構參數和加載方向的示意圖，包括晶格長度a、b、c，晶格角α、β、γ以及加載方向X、Y、Z。

首先，對初始模型進行幾何優化，以更好地開展后續的水分子吸附及分子動力學模擬。其次，執行水分子插入腳本來實現內層中固定數目的水分子的吸附。接下來，使用正則NVT系綜（粒子數、體積和溫度恒定）來使系統保持在常溫環境（298 K）。然后，進行等溫等壓NPT系綜（粒子數、壓力和溫度恒定）下的MD平衡計算來實現充分的溶脹和最為穩定的狀態。此時，體系的體積、密度和能量均能達到收斂平衡。并且，三維周期性邊界條件被用來模擬近似無限堆疊的宏觀黏土層結構。最后，繼續進行更多數目的水分子的吸附與MD溶脹平衡，如此循環，來獲得不同水化量的蒙脫石模型。由于在晶層溶脹階段，黏土內層中大約形成1～3層水化物[23]，通過MD溶脹模擬，水化量Nw被控制在0～240之間，分別對應于內層空間中，從絕對干燥到含有3個水化層的不同水化狀態。

除了水化量的影響，蒙脫石自身的層電荷密度，也會嚴重影響其水化力學行為。層電荷密度是每單位（Si，Al）4O10的層間電荷數，能夠直接反應黏土晶體的物理化學特征。通常情況下，蒙脫石八面體片層中的Al3+會被Mg2+取代，導致晶層表面攜帶負電荷。不同產地因成礦環境的差異，造成片層結構中類質同象取代量的不同，從而具有不同的層電荷密度，溶脹能力和力學性能也存在差異。不同層電荷密度的蒙脫石模型[24]，根據Loewenstein不相鄰取代規則，以八面體中6個Al3+被Mg2+晶格取代作為低層電荷密度（-0.75元電荷/原胞），8個Al3+被取代作為中層電荷密度（-1元電荷/原胞），10個Al3+被取代作為高層電荷密度（-1.25元電荷/原胞），并且內層中分別吸附6個、8個和10個Na+來平衡負電荷。

通過上述方法得到了水化蒙脫石的結構模型。而不同層電荷密度的蒙脫石模型，還需經歷水分子的吸附、幾何優化、NVT系綜溫度平衡以及NPT系綜溶脹模擬，來獲得不同的水化程度。計算過程中，NVT系綜下的MD計算時間為500 ps，NPT系綜下的MD計算則用1 ns來實現充分的溶脹，時間步長為0.5 fs（1 ns = 103 ps = 106 fs），溫度T為常溫的298 K（25 °C），壓力p為1個標準大氣壓（101 325 Pa）。并且，所有原子和晶胞參數的約束被釋放，黏土體系可以自由地膨脹或收縮，內層中的物質也可以自由移動。溶脹模擬之后，輸出結構便被用于后續的拉應力狀態下的MD計算。

1.2 拉應力狀態下的分子動力學計算模擬

編寫的施加應力-計算應變Perl語言腳本，會控制應力沿晶胞X、Y、Z中的某一方向施加拉應力（其他方向應力為零），然后對體系進行MD平衡，并測量體系中的應變。首先，它會導入一個施加應力序列（由用戶輸入）。這是以GPa為單位的應力列表，對應著一系列的MD計算周期，且在每個應力下（包括零應力下的初始平衡相），計算程序都會存在一個MD平衡階段和一個MD生產階段，分別由平衡步數和生產步數來控制（大小為1 000 000，時間步長為0.5 fs），以保證在每次施加應力后，體系均能達到充分的動力學平衡，并獲得穩定、可靠的應變值。其次，在每個應力周期開始時，晶胞的尺寸和形狀都會被重置到前一個周期的平均尺寸，這是為了通過停止晶胞振蕩來加速平衡[25]。最后，將對應的應變值寫入研究表格。

細節的設置中，將最新開發的、專門針對黏土礦物及其水化物的CLAYFF力場嵌入MD計算腳本中[15]來精確地描述蒙脫石原子間的相互作用，并為每個原子分配相應的電荷。該力場將大多數原子間的相互作用視為點電荷和非鍵作用，這將允許其被用于各種各樣的相，以及框架內完全平移的自由度。CLAYFF方法可以正確地預測水化黏土體系的物理力學性質以及黏土與水之間的相互作用[26]，并且其與程序腳本間的協調性較好。此外，選取有限變形下的Cauchy常應力模式，以便依照體系在應力狀態下的構形改變來進行相關應變張量的提取；采用具有較好溫度穩定性的Nosé-Hoover-Langevin（NHL）恒溫器方法控制溫度；使用可控制各向異性加壓的Souza-Martins控壓函數實現體系的快速平衡和晶胞尺寸的改變；應用高效率的particle-particle particle-mesh（PPPM）求解器算法進行體系中靜電相互作用的計算，該方法可以明顯提高大規模帶電體系的計算效率；采用基于原子截斷的atom-based求和算法加速范德華相互作用的計算，截斷半徑則被設置為10.5 ?（近似于晶胞尺寸的一半），以保證較高的計算精度。每個模型的生產計算都重復3次并求取平均值，以減少隨機誤差。通過這些適宜方法的選取和高精度的參數設置，便可以獲得高效準確的計算結果[27]。

輸入參數方面，需要指定一個應力施加序列，來實現拉應力狀態下的MD模擬。通過文獻調研[1]和計算驗證發現，合理的分子尺度應力范圍在0～3.0 GPa，應力增加步為0.1 GPa。完成了參數的設置后，計算程序便會按照預先輸入的應力序列和加載方向，沿體系某一方向循環拉伸應力，并進行長時間的分子動力學平衡計算。在此過程中，晶胞的尺寸和形狀也會發生改變。分子動力學計算體系中，周期結構的平均形狀可以由3個晶胞向量a0、b0和c0來描述，并通常表示為矩陣h0中的列，而h0的行列式則對應于體積。當結構發生形變時，3個晶胞向量會變為a、b和c，而晶胞矩陣則變為h。兩個矩陣h0和h共同定義了應變張量ε[1，28]

ε=12（h0?TGh0?1?I）

（1）

式中：G為度規張量hTh；I為單位矩陣；-T為逆矩陣的轉置。

根據施加應力與測量應變之間的關系，還可以獲得體系在某一應力狀態下的彈性模量[29]

cv=σvεv

（2）

式中：cv為彈性模量；σv為施加的應力；εv為測量的應變。

不僅如此，在拉伸模擬過程中，還會每隔2 000步（1 ps）收集一次原子軌跡和結構數據。這樣便可以獲得蒙脫石在每一個應力狀態下的動態平衡過程，從而為進一步分析其結構參數，包括體積、晶格長度、晶格角、氫鍵數目等，以及觀察其變形破壞過程奠定基礎。

2 拉應力下的力學特性分析

2.1 水化量的影響

編寫的MD計算腳本會通過控制應力的方式，對不同水化程度的蒙脫石沿X、Y、Z中的某一方向施加拉伸應力，以獲取相應的應力-應變曲線。分子尺度施加應力范圍為0～3.0 GPa，主要目的為，研究在一定拉伸應力作用下，土木工程領域最主要的膨脹性黏土蒙脫石各向異性的水化力學性質，并獲得最主要的拉伸破壞方向及破壞機理。圖2顯示了不同熱力學穩態水化量的蒙脫石（層電荷密度為-0.75元電荷/原胞，水化量Nw為0～240）在不同拉應力條件下體系中應變的變化情況。其中，圖2（a）、圖2（b）和圖2（c）分別對應于Z軸、X軸和Y軸方向的加載情況。可以看出，整體上，隨著某一方向拉應力的增加，體系中應變值的大小都在不斷增大。然而，相較于較為接近線性和難以發生變形的干燥狀態（Nw=0），水化黏土的拉伸現象更明顯（Nw=80至Nw=240），也更不穩定，并且水化量越高，拉伸應變值越大，且這一現象在高應力下會更為明顯。其中，當Z方向的加載達到極限應力時，水化蒙脫石體系會發生拉伸破壞，應變值急劇上升，直至整層的剝落分離。圖3顯示了Z方向的極限拉應力隨水化量的變化情況，以及已有研究中蒙脫石Z方向的極限抗拉應力隨含水量增加的改變情況[30]。可以看出，在相同含水量范圍內（0～30%），兩者極限拉應力的變化趨勢較為接近，都隨著水化量的增加而減小，并且，它們的弱化曲線都符合ExpDec1單指數函數模型的變化規律，在水化初期弱化幅度會更大。可見，內層水化程度的加深，減弱了蒙脫石的抗拉力學性能，其進入屈服破壞階段的極限應力值也更小，即從黏土基體中拉出一定數量的黏土層所需的拉應力更小。

然而，在細節方面，它們的應變值隨應力的變化情況卻有著明顯不同。最為顯著的是，蒙脫石在法線Z方向的拉伸應變，遠遠大于平面內其他兩個方向。這說明在拉伸應力狀態下，層間水主要影響礦物表面法線方向的力學行為，即Z方向的拉伸模量和拉伸強度最小，并且水化程度越高，就會越早地出現大幅度的拉伸變形，以及黏土層間的分離破壞。為了探索Z方向的抗拉力學性能與水化量之間的關系，還求取了水化蒙脫石在該方向線彈性階段的平均彈性模量。

圖4顯示了Z方向彈性階段的彈性模量隨水化量的變化情況以及Wei等[30]模擬獲得的鈉基蒙脫石Z方向的彈性模量隨含水量增加的改變情況。可以發現，在相同含水量范圍內（0～30%），兩者的彈性模量值較為接近，都隨著水化量的增加而減小，且它們的弱化曲線都擬合于ExpAssoc函數模型，即在水化初期弱化幅度會更大。這不僅說明了MD計算結果的可靠性，并且提供了分子尺度的細節解釋。圖5則顯示了MD模擬的體積和晶格長度，隨水化量增加的演化情況（圖5（a））以及Mahsa等[31]模擬獲得的基礎層間距的變化情況（圖5（b））。可以發現，體積和晶格長度c都有著不同程度的增長。其中，體積膨脹了74%，晶格長度c增長了72%，兩者的膨脹曲線，都與基礎層間距的變化趨勢較為接近。相比之下，其他兩個方向的晶格長度a和b，則變化幅度較小。晶格角α、β和γ則不發生改變。而黏土內層吸附水后，不僅會造成層間距的增加，并且會抵消黏土層間固有的相互作用，弱化黏土基體的抗拉性能，其也是軟巖吸水軟化誘發大變形的主導因素。

注：" （a） 體積和晶格長度的變化情況；（b） 基礎層間距變化情況。

為了更好地對Z方向的應力-應變曲線進行分析，并獲取其內在機理，不同應力狀態下MD計算平衡后的輸出結構被展示，來捕捉黏土-水體系中的相互作用機制和微觀結構演化。圖6顯示了水化蒙脫石（Nw=240）在Z方向不同拉應力條件下的變形和拉伸破壞過程。其中，在0.1 GPa的低應力狀態下，見圖6（a），一方面，水化結構非常規整，鈉離子集中在內層的中間區域，且周圍環繞著水分子，形成外球表面復合物（OSSCs）。另一方面，大量水分子吸附在了黏土表面，形成引力作用和氫鍵作用的強結合水有序結構。而圖6（b）則顯示了0.3 GPa時的分子結構拉伸，可以看出，內層空間被拉大，結合水膜密實度降低，原子空隙增加，加之層間距增大后的幾何效應以及黏土層間非鍵相互作用的減弱，共同導致了抗拉力學性能的降低。此時，黏土表面與水分子間的距離，也在一定程度上被拉長，一些理論將這一區域定義為硅-水相互作用區。但整體上，蒙脫石的力學行為仍處于彈性變化階段。

圖6（c）顯示了0.5 GPa時的變形和拉伸破壞過程。可以發現，在高應力狀態下，內層空間被進一步拉開，層間微孔開始局部聯通，部分水化層的氫鍵網絡被破壞，一些鈉離子去水化后，又重新吸附到了蒙脫石表面，拉伸模量持續減小。此時，表面仍然牢固地吸附著水分子，但兩者之間的距離卻進一步拉大，底部片層也出現了結構變化。而隨著拉伸進行，鈉離子周圍的水化網絡被持續破壞，并逐漸顯現一些縱向的裂縫。接下來，黏土內層進入屈服階段，中間區域出現了一些貫通的空隙，水分子層間的連接作用減弱。同時，硅-水相互作用區間的水分子數目也開始減少，拉伸進入快速發展階段，模量也達到了最低點。

圖6（d）顯示了0.6 GPa時的破壞失效和剝離過程。此時，一些水分子間的氫鍵完全斷裂，貫穿整個內層的大空隙形成，并開始整層地剝落分離。只有部分強結合水，仍然與黏土表面保持著吸附的相互作用，其他內層水則轉化為自由水，游離于體系之外。而當拉伸破壞進入最終階段，黏土層間實現了整體的分離破壞，剝離片層吸附著強結合水，逐漸遠離于黏土基體。表1統計了該水化蒙脫石在Z方向不同拉應力下的氫鍵數目，可以發現，隨著拉應力增加，氫鍵數目會不斷減少，其中高應力狀態下的0.4～0.6 GPa下降幅度急劇增大，表明其進入破壞失效階段，這與微結構分析得到的結論相一致。

這些微結構的分析工作，從細節剖析了蒙脫石拉伸力學行為的內在機理和變形破壞過程；蒙脫石水化后會發生明顯的內層膨脹，晶格長度c的線性增長是其體積增大的主導原因；Z方向的極限應力和拉伸模量，與水化量呈負相關，即水化弱化了堆垛方向的抗拉性能，達到極限拉應力后，會出現整層分離的破壞現象；黏土片層有著較強的剛度，作為較軟的組分，內層是大部分變形的主要原因，并且支配著蒙脫石的力學性能；內層空間拉長后的幾何效應，結合水膜密實度的降低（層間空隙的形成和擴展），硅-水相互作用區和水化結構網絡的拉伸破壞，以及氫鍵數目的減少是蒙脫石抗拉力學性能減弱的本質機理。

2.2 加載方向的影響

對于水化后的蒙脫石，不同加載方向也會對其力學行為有著重要影響。蒙脫石在Z方向的拉伸應變遠大于平面內，并且其應變變化近似多項式的形式，而平面X和Y方向則接近線性。因此，Z方向在實際工況中也更容易發生拉伸變形，拉伸強度也更低。相較于水化蒙脫石（Nw=240）在Z方向0.5 GPa時便接近屈服破壞，X和Y方向在高達3.0 GPa的拉伸應力下，卻仍然只產生非常有限的應變變化，分別為0.054和0.038，彈性模量分別為56 GPa和78 GPa。蒙脫石在平面內這兩個方向的彈性模量非常大，與Chen等[32]獲得的實驗值較為接近，受水化量的影響也較小，3個方向的力學行為呈現明顯的各向異性，這也與已有的研究報道一致[30，33]。此外，還可以發現，X方向的拉伸應變要略大于Y方向，這是因為在依照Loewenstein規則進行晶格取代工作時，取代位點在平面內的分布不均勻[34]。通過觀察黏土結構也可以發現，沿X方向的晶格取代數要多于Y方向，又由于蒙脫石在平面X、Y方向不同的拓撲結構和離子共價鍵強度[24]，X方向的拉伸模量相對較小，也產生了更大的拉伸應變值。

為了進一步對Z方向和平面X、Y方向拉應力作用下的力學特性進行解釋和說明，一些重要的結構參數被提取，包括體積V，晶格長度a、b和c以及晶格角α、β和γ。圖7顯示了水化蒙脫石（Nw=240）的體積和晶格參數隨Z、X、Y三個方向拉應力增加的變化情況。其中，圖7（a）、（b）對應于Z方向的拉應力，可以發現，體積和晶格長度c都在非線性增大，并且兩者的變化趨勢較為一致，而其他兩個方向的晶格長度a和b則增長幅度較小。此外，晶格角β由95.29°大幅增大到107.21°，晶格角γ則增長幅度很小。這些特征參數的變化趨勢，與上文中的動態拉伸過程相一致，共同描述了黏土層堆垛方向的拉伸破壞機理。

圖7（c）、（d）和圖7（e）、（f）分別對應于X和Y方向拉應力下的結構參數變化。可以看出，在平面內拉應力作用下，晶格參數僅有輕微的增加，而不會發生明顯的內層結構改變。其中，X方向拉應力下，晶格長度a和c輕微增加，晶格長度b輕微減小，晶格角α輕微減小，而晶格角β卻大幅減小。而在Y方向拉應力下，晶格長度b和c輕微增加，晶格長度a輕微減小，晶格角α輕微減小，而晶格角β卻大幅增大。值得注意的是，在平面X和Y方向拉應力作用下，體積局部會出現上下波動。由于此時晶格長度并沒有明顯改變，而晶格角卻變化較大，這一現象可以被解釋為當平面拉應力達到一定值時，蒙脫石黏土片層間會發生錯位，從而造成晶格角和體積的波動變化。通過觀察水化蒙脫石在平面X和Y方向拉應力下的變形過程也可以發現，應力對結構和內層的影響不大，只有輕微的體積增加，而形狀卻發生了較為明顯的改變。這是由于在平行于表面方向施加的應力荷載，主要由晶體結構而不是內層水支撐。并且，平面內成鍵聯結的力學剛度，要遠遠高于黏土層堆疊方向的非鍵作用聯結（包括庫倫力和范德華力等）。因此，Z方向的抗拉性能要遠遠小于平面內。

可見，Z方向拉應力主要造成體積和晶格長度c的非線性增長，以及晶格角β的大幅度增大；而在平面X和Y方向拉應力下，則主要發生晶格角β的大幅減小和增大；蒙脫石在不同水化和拉應力狀態下，會發生不同程度的結構參數改變。

2.3 層電荷密度的影響

除了上述兩個因素，蒙脫石自身的層電荷密度不僅會影響其水化行為，也會改變其在應力加載下的力學性能。雖然實驗和模擬方法中已有一些關于層電荷密度對蒙脫石水化膨脹影響的研究，但目前仍缺少對其水化力學性質作用機理的了解。為此，圖8顯示了通過MD模擬獲得的不同層電荷密度的蒙脫石，在Z方向不同拉應力下，體系中應變的變化情況。其中，圖8（a）、圖8（b）和圖8（c）分別對應于低層電荷密度蒙脫石（L-MMT）、中層電荷密度蒙脫石（M-MMT）和高層電荷密度蒙脫石（H-MMT）的加載情況。

可以看出，整體上，隨著拉應力的增加，體系中應變值的大小都在不斷增大。并且，相較于較為接近線性的干燥狀態（Nw=0），水化黏土的拉伸現象更明顯，也更不穩定。相同應力下，水化量越高，拉伸應變值越大，且這一現象在高應力下會更為明顯。說明內層水化程度的加深，會導致蒙脫石Z方向抵抗拉伸荷載的能力變弱，其進入屈服破壞階段的極限應力值也越小。然而，在細節方面，3種層電荷密度的蒙脫石，其應變隨應力的變化情況，卻有著明顯不同。整體上，拉伸應變值的大小關系為L-MMT gt; M-MMT gt; H-MMT，即層電荷密度越高，拉伸應變值越小，抗拉力學性能越強，應力應變曲線受水化量的影響也越小。

通過觀察3種層電荷密度的水化蒙脫石，在Z方向拉應力作用下的變形過程，也可以發現，L-MMT、M-MMT和H-MMT的內層拉伸現象依次減弱。并且，層電荷密度越高，晶層表面對水分子和鈉離子的吸附作用越強，形成的氫鍵數目越多，內層也更難以發生拉伸變形；而層間陽離子的水化效應則相對變弱，形成的OSSCs更少，因而體積溶脹和基礎層間距更小；此外，拉應力作用下形成的層間微孔會越少，結合水膜密集程度越高，抗拉力學性能也越強。表2統計了上述3種層電荷密度的水化蒙脫石（Nw=240），在Z方向拉應力（0.4 GPa）作用下的物理化學參數。可以發現，層電荷密度越高，體積和晶格長度c越小，水化膨脹能力越弱；然而，彈性模量和形成氫鍵的數目卻越大，蒙脫石也更難以進入屈服破壞階段。這與微結構分析得到的結論相一致，并定量闡釋了層電荷密度的作用機制。

一些學者曾通過實驗方法對不同層電荷密度蒙脫石的水化性能進行研究，并比較了MMT-1（層電荷密度為0.342）和MMT-2（層電荷密度為0.439），兩種層電荷密度的鈉基蒙脫石溶脹能力的差異見表3[35]。可以發現，層電荷密度越高，膨脹容和膨脹指數越低，水化能力越弱。這說明高層電荷密度蒙脫石具有更弱的水化膨脹能力，因而基礎層間距更小，層間引力作用更強，也更難以發生拉伸破壞。這一實驗結論，與MD模擬結果相一致，但后者提供了更多機理性的闡釋。

上述工作從分子尺度揭示了蒙脫石在不同水化量、拉應力狀態、層電荷密度下的基本力學性質及內在結構機理。實際工程中，應控制水化量盡可能小以提升抗拉性能，并預防蒙脫石在動態荷載作用下的強度衰減；也應重點關注內層Z方向的力學行為和相互作用機制，采取必要的物理化學手段（鈍化劑等）來改變黏土的表面結構，增強層間憎水性和增加氫鍵數目；最后，可以結合本工作中結構參數的演化規律對蒙脫石在不同條件下的形變特征進行判斷，來更好地預測土體或者軟巖的結構力學性質，評估其對工程建筑物的適應性，從而“對癥下藥”采取合理的措施，來提升其力學性能和工程行為。

3 結論

通過編寫施加應力-計算應變腳本和嵌入CLAYFF力場開展了采用過去方法難以實現的水化蒙脫石在拉伸應力狀態下的分子動力學模擬研究，確定了其不同應力階段的力學特性、相互作用機制和微觀結構演化，得到以下主要結論：

1）蒙脫石內層水化對力學性能的弱化效應明顯，包括極限拉應力和拉伸模量，且在水化初期弱化幅度會更大；體積水化膨脹主要源于晶格長度c的線性增長；這些力學和結構參數隨水化量的變化情況，與已有的研究結果相一致。

2）蒙脫石的力學性質各向異性明顯，Z方向的拉伸模量遠小于平面內，并且其應變變化近似多項式的形式，而平面X和Y方向則接近線性；應力對表面Z方向的力學行為影響最大，達到極限應力后，會出現整層分離的破壞現象。

3）黏土層有著較強的剛度，作為較軟的組分，內層是大部分形變的主要原因，并支配著蒙脫石的拉伸力學性能；Z方向拉應力主要造成體積和晶格長度c的非線性增長，以及晶格角β的增大，而在X和Y方向拉應力下，則主要發生β的減小和增大。

4）層電荷密度越高，結合水膜越密實，形成的氫鍵數目越多，體積和晶格長度c越小，抗拉力學性能越強，也更難以屈服破壞。實際工程中應控制水化量盡可能小，并采取有效措施增強層間憎水性和Z方向的抗拉性能，實現巖土災害宏（微）觀防治；也應結合結構參數的演化規律，對水化蒙脫石在不同條件下的結構力學進行判斷，為土體性質及工程適應性評價奠定基礎。

參考文獻

1

KATTI D R， SCHMIDT S R， GHOSH P， et al. Molecular modeling of the mechanical behavior and interactions in dry and slightly hydrated sodium montmorillonite interlayer [J]. Canadian Geotechnical Journal， 2007， 44（4）： 425-435. [百度學術]

2

SAYERS C M， DEN BOER L D. The elastic properties of clay in shales [J]. Journal of Geophysical Research： Solid Earth， 2018， 123（7）： 5965-5974. [百度學術]

3

MUKERJI T， DVORKIN J. The rock physics handbook [M]. 2ed edition. Cambridge： Cambridge University Press， 2009. [百度學術]

4

VYZHVA S， PRODAYVODA G， VYZHVA A. Elastic properties of some clay minerals [J]. Nafta-Gaz， 2014， 70（11）： 743-756. [百度學術]

5

楊文鈺， 鄭俊杰， 章榮軍， 等. 考慮黏土土性參數與支護壓力變異性的盾構掌子面穩定性分析[J]. 土木與環境工程學報（中英文）， 2021， 43（6）： 27-37. [百度學術]

YANG W Y， ZHENG J J， ZHANG R J， et al. Face stability analysis of shield tunnel considering variability of soil parameters and support pressure in clay [J]. Journal of Civil and Environmental Engineering， 2021， 43（6）： 27-37. （in Chinese） [百度學術]

6

ZHOU J K， CHEN X D. Stress-strain behavior and statistical continuous damage model of cement mortar under high strain rates [J]. Journal of Materials in Civil Engineering， 2013， 25（1）： 120-130. [百度學術]

7

劉珊. 結構性黏土力學特性與微觀形態試驗研究[D]. 江蘇 徐州： 中國礦業大學， 2014. [百度學術]

LIU S. Experimental study on mechanical properties and micro-structure of structural clay [D]. Xuzhou， Jiangsu： China University of Mining and Technology， 2014. （in Chinese） [百度學術]

8

PRADHAN S M， KATTI K S， KATTI D R. Evolution of molecular interactions in the interlayer of Na-montmorillonite swelling clay with increasing hydration [J]. International Journal of Geomechanics， 2015， 15（5）： 04014073. [百度學術]

9

HONTY M， FREDERICKX L， WANG L， et al. Boom Clay pore water geochemistry at the Mol site： Experimental data as determined by in situ sampling of the piezometers [J]. Applied Geochemistry， 2022， 136： 105156. [百度學術]

10

張明洋， 楊曉松. 蒙脫石的彈性性質實驗[J]. 地震地質， 2020， 42（5）： 1229-1239. [百度學術]

ZHANG M Y， YANG X S. Experimental studies on elastic properties of montmorillonite [J]. Seismology and Geology， 2020， 42（5）： 1229-1239. （in Chinese） [百度學術]

11

BOBKO C P. Assessing the mechanical microstructure of shale by nanoindentation： The link between mineral composition and mechanical properties [D]. Cambridge： Massachusetts Institute of Technology， 2008. [百度學術]

12

周繼凱， 金龍， 陳徐東. 水化硅酸鈣動態力學性能的分子動力學模擬研究[J]. 防災減災工程學報， 2014， 34（3）： 277-282. [百度學術]

ZHOU J K， JIN L， CHEN X D. Molecular dynamic simulation for mechanical properties of C-S-H [J]. Journal of Disaster Prevention and Mitigation Engineering， 2014， 34（3）： 277-282. （in Chinese） [百度學術]

13

LAU D， JIAN W， YU Z C， et al. Nano-engineering of construction materials using molecular dynamics simulations： Prospects and challenges [J]. Composites Part B： Engineering， 2018， 143： 282-291. [百度學術]

14

KARABORNI S， SMIT B， HEIDUG W， et al. The swelling of clays： Molecular simulations of the hydration of montmorillonite [J]. Science， 1996， 271（5252）： 1102-1104. [百度學術]

15

CYGAN R T， LIANG J J， KALINICHEV A G. Molecular models of hydroxide， oxyhydroxide， and clay phases and the development of a general force field [J]. The Journal of Physical Chemistry B， 2004， 108（4）： 1255-1266. [百度學術]

16

KATTI D R， SRINIVASAMURTHY L， KATTI K S. Molecular modeling of initiation of interlayer swelling in Na-montmorillonite expansive clay [J]. Canadian Geotechnical Journal， 2015， 52（9）： 1385-1395. [百度學術]

17

張亞云， 陳勉， 鄧亞， 等. 溫壓條件下蒙脫石水化的分子動力學模擬[J]. 硅酸鹽學報， 2018， 46（10）： 1489-1498. [百度學術]

ZHANG Y Y， CHEN M， DENG Y， et al. Molecular dynamics simulation of temperature and pressure effects on hydration characteristics of montmorillonites [J]. Journal of the Chinese Ceramic Society， 2018， 46（10）： 1489-1498. （in Chinese） [百度學術]

18

李勤， 陸現彩， 張立虎， 等. 蒙脫石層間陽離子交換的分子模擬[J]. 南京大學學報（自然科學）， 2019， 55（6）： 879-887. [百度學術]

LI Q， LU X C， ZHANG L H， et al. Molecular simulation of interlayer cation exchange of montmorillonite [J]. Journal of Nanjing University （Natural Science）， 2019， 55（6）： 879-887. （in Chinese） [百度學術]

19

楊微， 陳仁朋， 康馨. 基于分子動力學模擬技術的黏土礦物微觀行為研究應用[J]. 巖土工程學報， 2019， 41（Sup1）： 181-184. [百度學術]

YANG W， CHEN R P， KANG X. Application of molecular dynamics simulation method in micro-properties of clay minerals [J]. Chinese Journal of Geotechnical Engineering， 2019， 41（Sup1）： 181-184. （in Chinese） [百度學術]

20

YANG Y， ADHIKARI S， XU G Y. Molecular dynamics simulation in the interlayer of mixed-layer clays due to hydration and swelling mechanism [J]. Crystals， 2021， 11（6）： 586. [百度學術]

21

楊亞帆， 王建州， 商翔宇， 等. 高溫下鈣蒙脫石膨脹特性的分子動力學模擬[J]. 物理學報， 2022， 71（4）： 043102. [百度學術]

YANG Y F， WANG J Z， SHANG X Y， et al. Molecular dynamics simulation of swelling properties of Ca-montmorillonite at high temperatures [J]. Acta Physica Sinica， 2022， 71（4）： 043102. （in Chinese） [百度學術]

22

VIANI A， GUALTIERI A F， ARTIOLI G. The nature of disorder in montmorillonite by simulation of X-ray powder patterns [J]. American Mineralogist， 2002， 87（7）： 966-975. [百度學術]

23

CASES J M， BEREND I， BESSON G， et al. Mechanism of adsorption and desorption of water vapor by homoionic montmorillonite. 1. The sodium-exchanged form [J]. Langmuir， 1992， 8（11）： 2730-2739. [百度學術]

24

MORODOME S， KAWAMURA K. Swelling behavior of Na- and Ca-montmorillonite up to 150℃ by in situ X-ray diffraction experiments [J]. Clays and Clay Minerals， 2009， 57（2）： 150-160. [百度學術]

25

WANG Y H， YANG G， WANG W H， et al. Effects of different functional groups in graphene nanofiber on the mechanical property of polyvinyl alcohol composites by the molecular dynamic simulations [J]. Journal of Molecular Liquids， 2019， 277： 261-268. [百度學術]

26

CARRIER B， VANDAMME M， PELLENQ R J M， et al. Elastic properties of swelling clay particles at finite temperature upon hydration [J]. The Journal of Physical Chemistry C， 2014， 118（17）： 8933-8943. [百度學術]

27

宋華苗， 金祖權， 王攀， 等. 水環境下摻鋁硅質基體樹脂黏結界面性能演變的分子動力學模擬與驗證[J]. 土木與環境工程學報（中英文）， 2021， 43（5）： 123-131. [百度學術]

SONG H M， JIN Z Q， WANG P， et al. A molecular dynamics simulation and validation of bonding behavior evolution of aluminum-doped silica matrix-resin bonding interface under water environment [J]. Journal of Civil and Environmental Engineering， 2021， 43（5）： 123-131. （in Chinese） [百度學術]

28

RAPAPORT D C. The art of molecular dynamics simulation [M]. Cambridge： Cambridge University Press， 2004. [百度學術]

29

SCHMIDT S R， KATTI D R， GHOSH P， et al. Evolution of mechanical response of sodium montmorillonite interlayer with increasing hydration by molecular dynamics [J]. Langmuir， 2005， 21（17）： 8069-8076. [百度學術]

30

WEI P C， ZHENG Y Y， XIONG Y， et al. Effect of water content and structural anisotropy on tensile mechanical properties of montmorillonite using molecular dynamics [J]. Applied Clay Science， 2022， 228： 106622. [百度學術]

31

RAHROMOSTAQIM M， SAHIMI M. Molecular dynamics simulation of hydration and swelling of mixed-layer clays [J]. The Journal of Physical Chemistry C， 2018， 122（26）： 14631-14639. [百度學術]

32

CHEN B Q， EVANS J R G. Elastic moduli of clay platelets [J]. Scripta Materialia， 2006， 54（9）： 1581-1585. [百度學術]

33

SAYERS C M， DEN BOER L D. The elastic anisotropy of clay minerals [J]. Geophysics， 2016， 81（5）： C193-C203. [百度學術]

34

ZARTMAN G D， LIU H K， AKDIM B， et al. Nanoscale tensile， shear， and failure properties of layered silicates as a function of cation density and stress [J]. Journal of Physical Chemistry C， 2010， 114： 1763-1772. [百度學術]

35

高仁波， 趙云良， 陳立才， 等. 蒙脫石層電荷密度對其二維納米片剝離的影響[J]. 硅酸鹽學報， 2021， 49（7）： 1420-1428. [百度學術]

GAO R B， ZHAO Y L， CHEN L C， et al. Effect of layer charge density on the exfoliation of montmorillonite to prepare two-dimensional nanosheets [J]. Journal of the Chinese Ceramic Society， 2021， 49（7）： 1420-1428. （in Chinese） [百度學術]