1.2 模擬細節

本研究基於LAMMPS 軟件開展分子動力學模擬研究。模擬采用的混合勢函數包括Buckingham 勢(截斷半徑為10 Å)、L-J 12-6 勢(截斷半徑為 10 Å)以及長程庫侖相互作用(截斷半徑為 15 Å)。分子內鍵長與鍵角通過SHAKE 算法進行約束。長程庫侖力采用PPPM(particle-particle particle-mesh)方法計算，精度設置為 10-5。為避免液膜整體平移對應力計算產生影響，模擬過程中每25 000 步對體係動量進行重整化以消除整體漂移。初始構型經控壓(NPT)弛豫後獲得平衡結構；生產計算采用等容過程(NVT)，時間步長均為1.0 fs。

相較於均勻體係，混合物界麵性質的計算需排除界麵間的相互幹涉。文獻明確了界麵體係的建構準則，即模擬盒子沿液膜法向設為長方體構型，液相體係模擬應包含2 000~5 000 個粒子，液膜厚度不低於30 Å。本文針對氫氧化鈉液膜體係的模擬，選取5 000 個氫氧化鈉分子為研究對象，初始模擬盒子為 140 Å ×100 Å ×100 Å 的長方體 (140 Å 為液膜法向尺寸)，在控壓弛豫過程中保持盒體高度 x 不變，通過 y 、 z 平麵的距離調控體係壓力，弛豫後模擬盒子尺寸變為 140 Å ×80 Å ×80 Å，對應液膜厚度為40 Å。為探究水分子含量對熔融氫氧化鈉表麵張力的影響，還分別構建了含250、125 和30 個水分子的液膜體係。水分子擴散係數的計算體係具有各向同性特征，因此采用正方體模擬盒子建模，體係由5 000 個氫氧化鈉分子和250 個水分子組成。初始模擬盒子尺寸為 80 Å ×80 Å ×80 Å，控壓收縮後，尺寸變為 57 Å ×57 Å ×57 Å。圖1 為兩個體係在生產計算前的初始構型。所有模擬體係均采用周期性邊界條件。

圖 1 用於生產計算的初始構型

a--弛豫後的液膜體係(計算表麵張力)；b--弛豫後的體相體係(計算擴散係數)

1.3 計算模型

基於Green-Kubo公式[23]計算黏度（η）：

其中：V為體係體積；KB為玻爾茲曼常數；T為熱力學溫度；t為時間，取t=0為初始時刻；m為原子質量；τ為對應t時刻（式（3））的時間延遲；α、β為方向；riα和viα為第i個分子在α方向上的位移和速度；rijα為粒子i和j在α方向的位置矢量分量；fijβ為粒子i和j在β方向的相互作用力分量；〈·〉表示體係內所有原子的平均。

本研究采用體係應力張量法[24]計算表麵張力係數。該方法利用體係內法向應力和切向應力之間的差異性計算表麵張力。對於以z方向為界麵法向的體係，應力與表麵張力係數的關係如下：

(5)

其中：σxx、σyy和σzz為應力張量的分量；Lx、Ly和Lz為模擬盒子的邊長。

表麵張力係數是表麵張力的定量表征，反映的是單位長度界麵上分子間的作用力。本文采用應力張量法計算得到的即為表麵張力係數，後續對表麵張力變化規律的分析均以此為基礎。

需要指出的是，在高溫條件下，水分子可能從氫氧化鈉液膜中逸出至真空區域，若對整個模擬盒子進行積分統計，將導致計算結果失真，因此僅對液膜所在區域進行積分處理。由於水分子在氫氧化鈉液膜界麵發生了聚集，為合理納入其對液膜性質的貢獻，本文以液膜法向方向上的總分子數密度作為劃分液膜區的判據：將數密度高於液膜中心密度10%的區域定義為液膜，以確保積分範圍涵蓋液膜主體及其附近水分子的影響。

1.4 模型驗證

為驗證勢場精度，對比了不同溫度下熔融氫氧化鈉密度、黏度的計算值與實驗值的誤差。根據文獻[25]的實驗可得知密度ρ（kg/m3）和黏度η（mPa·s）與溫度T（K）的關係式為：

根據式（6）、（7）可知，氫氧化鈉的密度與溫度呈負線性關係，黏度隨溫度呈指數衰減。熔融氫氧化鈉密度和黏度及其與實驗值的對比如表3所列。

表  3  熔融氫氧化鈉密度與黏度計算值與實驗值的對比

由表3可知，熔融氫氧化鈉密度計算值與實驗值的相對誤差小於5%，黏度計算值與實驗值的相對誤差小於9%，該誤差處於可接受區間。可見，Buckingham勢能準確描述熔融氫氧化鈉密度和黏度性質，對該體係具有良好的適用性。

2 結果與討論

2.1 溫度和水分子對熔融氫氧化鈉表麵張力的影響

表麵張力係數是液體界麵分子間相互作用的宏觀表現。不同體係液膜的表麵張力係數如圖2 所示。由圖2 可見，對於純熔融氫氧化體係，600~750 K 溫度區間，表麵張力的計算值與實驗值吻合良好；750 K 以上的高溫區，二者出現偏差，計算值高於實驗值。但總體來看，該混合力場可以預測熔融氫氧化鈉表麵張力係數的變化趨勢。與純氫氧化鈉體係對比，含水體係的表麵張力係數隨水分子含量的升高呈梯度式降低。600 K 時，5% 含水體係的表麵張力係數僅為54 mN/m，相較於純氫氧化鈉體係(143 mN/m)降幅超過 60%，這表明在低溫區間，水分子對液膜張力有不可忽視的影響。所有含水體係的表麵張力係數均隨溫度的升高持續上升，當體係溫度接近900 K 時，不同含水體係的表麵張力係數趨於一致，水分子對表麵張力的影響顯著減弱。

圖 2 不同體係液膜的表麵張力係數

為揭示水分子影響液膜界麵張力的微觀機製，本研究采用徑向分布函數gAB(r)（RDF）評價各粒子間的相互作用強度。RDF是描述係統中i、j兩粒子間距離為r處相對分布的函數。其定義為：

其中：ρA和ρB分別為粒子A和粒子B的平均密度；NAB（r）為在距離r範圍內，粒子A周圍出現粒子B的數目；V為體係體積；nA和nB分別為粒子A和粒子B的總數。

通過RDF曲線特征峰的高度與峰位，可判斷粒子間的相互作用強度與空間分布規律：特征峰越高，表明粒子間相互作用越強；峰位則反映體係在不同分子間距下的結構有序性。對RDF曲線積分可計算以A粒子為中心的B粒子的平均配位數MAB：