微纳尺度的受限水体系中, 受限水大面积接触受限腔内壁, 受到腔壁势场影响发生强弱程度不同的电极化[1], 与此同时水分子的排列受到空间约束和腔壁微观结构势场调制, 致使受限水体系呈现丰富多彩的有序度和物相结构[2].界面水分子局域极化能力也影响着局域水分子之间、水分子和其它带电粒子之间的相互作用, 进而影响诸多重要的界面物理与化学过程[3~5], 如表面水化、离子溶剂化、分子输运、颗粒自组装等.受限水系统兼具界面水和受限物相的物理特性, 毫无疑问其介电性质的探索具有重要的科学意义和应用价值, 近年来也已经获得越来越多的科研关注.例如, 受限水分子在其平动与转动自由度部分缺失、周边电场环境的状态下如何影响静态与动态介电性质, 至今尚未有定量的认识.掌握受限水介电性质的变化规律, 将帮助微纳尺度药物输运技术和高效海水淡化盐离子技术取得进步[6, 7], 也将促进新型纳米储能与传感器件设计[8, 9]领域的发展.
界面水分子分层有序伴随着电极化[10], 施加电场难以进一步驱使已极化界面水分子转动作出线性介电响应, 因此实验测量介电常数极具挑战.近二十年来受限水介电性质的研究集中在理论与计算模拟领域[11~14].例如, Ballenegger和Hansen[11]基于线性响应理论推导出受限腔体极性分子液体介电常数空间分布函数的理论表达式; Galli[13]和Netz[14]小组基于其开展的分子模拟研究, 分别报道了受限水介电常数张量和介电弛豫性质的强烈各向异性; Bonthuis等[12]和Schlaich等[14]计算了受限腔壁亲疏水性和柔性腔壁对受限水介电常数的影响等. 2018年, Geim团队[15]在受限水介电性质的测量方法上取得突破, 他们利用先进的微纳加工技术结合扫描介电显微镜与原子力显微镜, 首次精密测量了受限尺寸小至0.8 nm的受限液相水的介电常数. Geim等[15]提出, 对受限水的精密实验测量结果, 将对计算模拟研究提出更高的要求.这是因为已有受限水介电性质的分子模拟研究大多应用恒定偶极矩的简单水分子模型, 并不考虑水分子极性的变化.而实际受限水的分子极性、介电性质以及限域条件之间的物理关联都亟待探究.
本工作的另一个研究动机涉及受限水的固-液结构相变. Algara-Siller等[16]运用透射电镜直接观测到了范德瓦尔斯作用下两片石墨烯中二维正方晶格的受限方冰, 并监测到电子束移动可以驱动受限方冰的无序化.这个实验发现意味着应用外场控制限域水不同物相间结构相变的可能性, 及利用受限冰相变定量调控介电储能性质的潜在应用价值.虽然上文介绍了受限水介电性质的实验测量技术上的突破, 但目前仍无法做到同时定量控制结构相变和定量测量介电性质.本文作者已陆续开发了低维固-液界线的平衡态分子动力学模拟技术, 模拟了单分子层受限冰-水(固-液)两相平衡[17], 并预言了二维方冰晶界的预熔化相变液膜宽度的定量物理规律[18].在这些模拟研究基础上, 近一步精确计算单分子层受限冰和受限水的介电性质, 对于评估受限冰相变介电器件应用的可行性至关重要.
因此本工作的目标是计算单分子层受限水和受限冰体系的介电性质.此前, 受限水的模拟研究大多使用经典恒定偶极矩刚性水分子模型, 以牺牲分子极性描述精度换取计算速度优势.但当关注的问题涉及介电性质时, 恒定偶极矩刚性水分子模型不再适用, 需要使用极化水分子模型以动态适应外界电场环境变化[19], 从而更加准确地预测介电常数.到目前为止, 尚未有研究系统比较过简单的恒定偶极矩水分子模型和极化水分子模型在描述受限水体系介电性质的具体优劣表现.
基于上述研究背景和研究目标, 本论文将利用平衡态分子动力学模拟和受限体系介电性质计算方法, 详细比较两种水分子模型(恒定偶极矩模型和极化模型)在描述受限单分子层液相水和固相冰的介电性质上的优劣势.我们分别用两种水分子模型模拟了在0.65 nm限域尺寸、5×108 Pa限域压强条件下的四个不同温度的平衡态单分子层水和冰.计算了每个模拟体系的静态结构因子、静态偶极空间关联函数、静态介电常数、体系偶极时间关联函数和德拜弛豫时间.另外, 我们统计了极化水模型描述单分子层受限水和受限冰的水分子偶极矩量值概率密度分布.得到了关于极化水模型在描述受限水/冰体系介电性质具体优势较为全面系统的认识, 并首次发现了单分子层受限水、冰的分子极性区别于块体系统的奇异特性.本工作获得的计算结果能够为受限水领域研究提供重要数据, 对进一步开展基于受限水系统储能、传感、输运的理论设计具有重要指导意义.
采用两种水分子模型开展本论文中的模拟与计算工作.一种是最常用的SPC/E水模型[20](具体作用势参数参考引文[17]).另一种是2006年Lamoureux等开发的Drude极化水模型[21](全称SWM4-NDP, 4-site simple water model with new drude polarizability).如图 1所示, Drude极化水模型中引入一个假想的极小质量带电荷粒子.这个粒子通过一个假想的弹簧连接着氧原子粒子, 带电粒子的电荷、弹性振荡的瞬时位置、其所在空间点位的电场强度与自身极化率决定了水分子的动态诱导偶极矩数值.
除了Drude极化水分子模型之外还有两套著名的极化水分子模型, 分别是诱导点偶极水分子模型[22]和基于电负性平衡原理的涨落电荷水分子模型[23].这两套模型, 前者相较于Drude模型缺点在于矢量化长程库伦能力场的修改过于繁杂, 后者相较于Drude模型的缺点在于极化率张量描述能力的欠缺(如线性分子只能够获得轴向的极性). Drude模型很好弥补了这些欠缺, 在计算效率上有了显著提高.除此之外, Drude极化水分子模型能够成功再现常温常压条件下液相块体水的静态介电常数(见表 1)、蒸发焓、密度、自扩散系数、剪切粘度和水化自由能等一系列关键热力学性质.需要指出的是, 本工作采用的经典Drude模型(基于经典力学)[24]仍具有一定的缺陷和不足.例如, 它所给出的原子极化率是各向同性的, 无法精准描述高阶极化率, 无法完美再现强场梯度和高阶极矩作用下的原子间色散力.针对这些不足, 量子Drude模型近年来浮出水面[28], 但尚未被广泛应用, 计算效率也有待提升.
下载:
导出CSV
| SPC/E | Drude | Expt.c | |
| μGas./(Debye) | 2.35a | 1.85b | 1.85 |
| μLiq./(Debye) | 2.35a | 2.46b | 2.95 |
| ε0 | 71 (3) | 78 (2) | 78.4 |
| a SPC/E model, gas-phase dipole moment | |||
单分子层受限冰和受限水的限域尺寸和受限压强采用参考文献[17]和[18]中的参数设定, 限域尺寸为0.65 nm, 限域压强为5×108 Pa.这个受限条件借鉴了Algara-Siller等[16]实验观测结论, 高限域压强由双层石墨烯间范德华力致使, 0.65 nm是当前实验上可获得单分子层受限水的最小尺寸.另外, Algara-Siller等还证实过单分子层方冰结构不依赖于受限腔壁的结构(简单的刚性墙或是柔性石墨烯)、墙-水的相互作用强度以及疏水性.
与参考文献[17]和[18]中做法一致, 我们将水分子置于两平行刚性壁内, 刚性腔壁无原子结构[29, 30], 限域腔壁与受限水分子相互作用采用10-4-3势能函数[31], 下标中w代表受限墙, o代表氧原子. z表示水分子中氧原子到腔壁的垂直距离.对于SPC/E水分子模型,
|
$ \begin{array}{c} {\varPsi }_{\mathrm{w}\mathrm{o}}\left(z\right)＝2\mathrm{\mathsf{π} }{ϵ}_{\mathrm{w}\mathrm{o}}\left(\frac{4}{5}\right)\left[\frac{2}{5}{\left(\frac{{\sigma }_{\mathrm{w}\mathrm{o}}}{z}\right)}^{10}-{\left(\frac{{\sigma }_{\mathrm{w}\mathrm{o}}}{z}\right)}^{4} \\ -\frac{\sqrt{2}}{3{\left(\frac{z}{{\sigma }_{\mathrm{w}\mathrm{o}}}+0.61\sqrt{2}\frac{z}{{\sigma }_{\mathrm{w}\mathrm{o}}}\right)}^{3}}\right] \end{array} $ |
(1) |
分子动力学模拟盒子的构建基于直角坐标系.盒子中包含单层水分子以及上下各一个刚性腔壁提供限域.定义z方向垂直于单层水分子和刚性腔壁表面, 水分子水平分布在xy平面方向.模拟盒子在xy方向上应用周期性边界条件, z方向上没有应用周期性边界条件. z方向盒子的尺寸为0.65 nm, 等同于两个刚性腔壁表面间距. x和y方向盒子尺寸均为4.4 nm.盒子中有256个水分子.为了精确计算长程静电相互作用, 我们使用Particle-Particle Particle-Mesh算法[33], 并应用了Yeh等[34]提升Particle-Particle Particle-Mesh算法精度的技巧, 在z方向上引入一个尺度约为5 nm的真空区域用以去除模拟盒子与其各个镜像间的偶极-偶极相互作用.
本工作中所有分子动力学模拟均通过美国Sandia国家实验室开发的LAMMPS软件实现[35].模拟过程中, 应用SHAKE算法[36]来实现SPC/E和Drude水分子(Drude弹簧振荡粒子之外部分)运动的刚性约束条件.氢氧键键长和键角维持不变. Drude振荡粒子由单独的弹性约束运动方程来进行时间积分.两种体系模拟积分时间步长取1 fs (10－15 s). SPC/E水分子模型体系利用Nosé-Hoover热浴恒温法控温, 利用Anderson恒压法控压[17]. Drude水分子模型体系则相对复杂一些, 需要采用双朗之万控温和控压方法[37]才能够避免虚拟的Drude振荡粒子与系统真实粒子达到热平衡状态, 而进行非物理的能量交换导致错误的动力学演化行为.在具体控温中, 非Drude粒子温度控制在目标温度, 弛豫时间为0.1 ps; Drude粒子温度控制在远小于目标温度的低温数值, 弛豫时间为0.02 ps.
针对两套水分子模型, 本工作分别模拟了4组不同的温度的固相(100 K、200 K、300 K、400 K)和液相(400 K、450 K、500 K、550 K)受限水体系.调控限域水x和y方向上压强维持在5×108 Pa.对每组体系先使用NPxyT系综模拟10 ns使之达到结构、热力学和静电力学平衡.之后转为NVT系综模拟2 ns, 每0.1 ps记录一次所有粒子坐标信息用以后续统计分析与计算.需要说明的是, 在Drude水分子模型体系会观察到众所周知的(由扩展拉格朗日方法造成)的飞冰(flying-ice-cube)问题, 为了解决这个问题, 需要每5 ps去除一次系统整体平动与转动动量.
为了和受限水体系比较, 我们利用模拟了768个SPC/E或Drude模型水分子组成的液体水和六角相冰结构.控温、控压、动力学积分步长和平衡态模拟时长均与受限水体系保持一致.因为模拟盒子是三个方向应用周期性边界条件, 仅应用Particle-Particle Particle-Mesh算法来处理长程静电相互作用. x, y, z三个方向盒子尺寸大约在2.7 nm、3.1 nm、2.9 nm附近(随温度变化).
此外, 本工作额外模拟了不同温度下的单层石墨烯作为受限腔壁的(Drude模型)单分子层受限水和受限冰, 用以检查受限条件对Drude模型水受限水分子极性变化的影响.采用两组720个碳原子的单层石墨烯(上下两层相距0.65 nm)代替无原子结构的刚性墙.水分子数目在模拟中可改变以维持压强恒定.碳原子在模拟中始终固定, 仅与水分子中氧原子通过Lennard-Jones势函数相互作用[32].长程静电相互作用的处理、模拟温压数值、平衡态模拟时长均与刚性墙受限水体系一致.
电介质的介电常数是定义极化矢量和外加电场之间关系的一个宏观物理量.当体系具有非均质性或各向异性时, 介电常数也不再是一个简单量值.例如本论文中所研究的准二维受限系统, 介质在z方向对称性缺失, 介电常数需要写成张量形式[11],
|
$ \varepsilon =\left(\begin{array}{ccc}{\varepsilon }_{x}& 0& 0\\ 0& {\varepsilon }_{y}& 0\\ 0& 0& {\varepsilon }_{z}\end{array}\right) $ |
(2) |
其中xy和z分别表示平行于和垂直于限域腔壁的分量方向, 极化密度矢量和受限腔体中介质层内电场强度的响应关系也分为平行和垂直两种分量,
|
$ {\boldsymbol{P}}_{x}=\frac{{\varepsilon }_{x}-1}{4\pi }{\boldsymbol{E}}_{x} $ |
(3a) |
|
$ {\boldsymbol{P}}_{y}=\frac{{\varepsilon }_{y}-1}{4\pi }{\boldsymbol{E}}_{y} $ |
(3b) |
|
$ {\boldsymbol{P}}_{z}=\frac{1}{4\pi }\frac{{\varepsilon }_{z}-1}{{\varepsilon }_{z}}{\boldsymbol{E}}_{z}$ |
(4) |
Ballenegger和Hansen[11]利用统计物理方法(配分函数求平均值)获得了极化密度矢量分量的表达式,
|
$ {\boldsymbol{P}}_{\alpha =x, y, z}=\frac{1}{{k}_{\mathrm{B}}T}\sum \limits_{\gamma =x, y, z}\left[\langle {\boldsymbol{m}}_{\alpha }{\boldsymbol{M}}_{\gamma }\rangle -\langle {\boldsymbol{m}}_{\alpha }\rangle \langle {\boldsymbol{M}}_{\gamma }\rangle \right]{\boldsymbol{E}}_{\gamma } $ |
(5) |
其中m和M分别代表电解质的瞬时局域极化密度和偶极矩总和, 尖括号代表时间平均.将式(5)带入(3)和(4), 即可获得受限偶极分子电介质的介电张量分量的计算式:
|
$ {\varepsilon }_{x}=1+\frac{4\pi }{{k}_{\mathrm{B}}T}\left[\langle {\boldsymbol{m}}_{x}{\boldsymbol{M}}_{x}\rangle -\langle {\boldsymbol{m}}_{x}\rangle \langle {\boldsymbol{M}}_{x}\rangle \right]$ |
(6a) |
|
$ {\varepsilon }_{y}=1+\frac{4\pi }{{k}_{\mathrm{B}}T}\left[\langle {\boldsymbol{m}}_{y}{\boldsymbol{M}}_{y}\rangle -\langle {\boldsymbol{m}}_{y}\rangle \langle {\boldsymbol{M}}_{y}\rangle \right] $ |
(6b) |
|
$ {\varepsilon }_{z}^{-1}=1-\frac{4\pi }{{k}_{\mathrm{B}}T}\left[\langle {\boldsymbol{m}}_{z}{\boldsymbol{M}}_{z}\rangle -\langle {\boldsymbol{m}}_{z}\rangle \langle {\boldsymbol{M}}_{z}\rangle \right] $ |
(7) |
从分子动力学模拟中每0.1 ps记录的所有粒子的坐标信息, 可以获得式(6)和(7)中极化密度分量
介电常数之外, 我们还计算了xy平面内径向偶极-偶极关联函数[13]用以分析受限水分子的偶极方向结构性质,
|
$ F\left(r\right)=\frac{1}{N}\langle \sum \limits_{i=1}^{N}\frac{1}{{n}_{i}\left(r\right)}\sum \limits_{j=1}^{{n}_{i}\left(r\right)}{\widehat{\boldsymbol{m}}}_{i} \cdot {\widehat{\boldsymbol{m}}}_{j}\rangle $ |
(8) |
表达式中
本工作研究的最后一个介电性质是受限水的德拜(偶极)弛豫时间
|
$ {\varPhi }_{\alpha }\left(t\right)=\frac{\langle {M}_{\alpha }\left(0\right) \cdot {M}_{\alpha }\left(t\right)\rangle }{\langle {M}_{\alpha }^{2}\rangle } $ |
(9) |
使用指数函数
图 2中记录了NVT系综模拟中某一个瞬时时刻单分子层受限冰(a), (b)和受限水(c), (d)的xy平面坐标快照. T＝400 K. 图 2(a)和(c)为SPC/E水分子模型结果, 图 2(b)和(d)为Drude极化水分子模型结果. 图 2中用红色代表受限冰, 蓝色代表受限液相水, 下文图中同样使用红蓝颜色区分受限冰和水.恒定偶极矩与极化水分子模型模拟的受限冰均呈现与Algara-Siller等工作中(分子模拟部分结果)非常一致的结构, 即长程有序的平面晶格结构, 氢氧键的方向集中在xy平面内.液相限域水中, 两种水模型体系均呈现长程无序性结构, 平动与转动型扩散远大于受限冰.在SPC/E水模型描述的受限水中, 观察到皮秒量级短时间尺度, 纳米空间尺度局域有序的团簇结构(10到20个水分子)集体形式地参与密度涨落与质量输运.在Drude极化水模型体系中, 局域瞬态有序的团簇结构仍可被观察到, 但似乎尺寸要小一些、生存时间更短一些.
SPC/E水分子模型的偶极矩为恒定数值2.35 Debye, 它不随物相、温度、受限条件发生任何变化.而Drude极化水模型的显著优势即能够展现物相对分子极性的强烈影响, 前言中介绍的其它两类极化水模型也都能够很好地描述液、汽两相中水分子不同的偶极矩数值.汽相水分子间距远大于液相水中分子间距, 因此分子极性的明显变化是容易理解的. 图 3(a)中, 我们计算了Drude水分子模型模拟常压下块体水和六角晶格块体冰相的分子极性统计结果.其中插图是分子偶极矩平均值随温度的变化结果. 图 3(a)中的结果展现出水分子偶极矩的分布呈正态分布, 温度越高分布中值逐渐左移变小, 分布范围变宽.可以理解为水分子瞬态极性数值的涨落和其在凝聚态物相中平动和转动自由度的扩散能力关联.从50 K到400 K共350 K的温度变化范围内, 偶极矩的平均值从2.62 Debye减小到了2.34 Debye, 减小幅度近0.3 Debye.固相和液相区域的温度梯度变化斜率不等, 液相温度变化梯度更大.在固-液相变温度处, 固-液两相水分子的偶极矩平均值相等, 并未像密度和势能等热力学量呈现跳变(块体水-冰的固-液结构相变是一级相变), 而是展现出连续变化形式.
然而, 在单分子层受限体系中的计算结果让我们非常吃惊. 图 3(b)中两处结果与块体体系明显的差异: (1)相同温度变化范围(100 K到450 K, 350 K温区)的水分子偶极矩平均值从2.84 Debye降到2.34 Debye, 减幅0.5 Debye明显大于块体水体系, 说明温度对受限水体系的极性变化影响更大. (2)更意想不到的是, 在400 K受限冰和受限水两相共存温度下, 水分子偶极矩平均值出现跳变, 与块体体系的过渡方式完全不同.
此区别意味着单分子层受限水-冰结构相变过程中, 固-液界面处可能伴随着偶极矩分布的不均匀性, 偶极矩数值变化与结构相变可能是协作发生的.为了验证这一猜测, 利用本文作者开发的低维受限水-冰固-液平衡态分子动力学模拟技术[17], 参考文献[17]中的分子动力学模拟细节参数设定并结合本工作所使用的Drude极化水分子模型参数和腔壁-水分子间相互作用形式, 我们实现了Drude极化水分子模型描述的单分子层受限冰-水(固-液)两相平衡态模拟.如图 4(a)所示, 色谱颜色变化标注的水分子/氧原子的结构序参量(计算方法详见参考文献[17]).红色代表高序参量数值的固相, 蓝色代表低序参量数值的液相.结构序参量中色彩变化可直观表征出固-液界面的位置.在5 ns的平衡态模拟中, 固-液界面始终稳定存在, 说明当前单分子层限域冰-水体系经历了一级(固-液)相变.通过统计分析Drude水分子偶极矩沿着界面法向方向上的分布函数(参照参考文献[17]中列热力学参量沿x轴分布函数的计算方法), 见图 4(b), 我们发现受限冰和受限水偶极矩数值处于不同的平台值, 偶极矩数值变化与结构相变是协作发生的.由此可以推断单分子层受限冰和水的分子偶极矩的跳变, 相较SPC/E等恒定偶极矩水分子模型体系, 将会导致受限冰、水两相的化学势差别增大, 两相形成界面自由能要更高[38], 从而影响到受限冰形核形貌和晶体生长速度.
需要说明的是, 国内外不乏关于受限水(或冰)体系水分子重要物理化学性质异于块体水(或冰)的报道, 例如, Luo和Zeng等[39]发现的铁电极化受限冰, Wang等[40]预言的受限水分子O—H键振动频率在受限腔壁离域π键作用下发生红移效应等.本工作中关于水分子极性随温度和固-液结构变化方式在受限条件下显著区别于块体体系的理论预言, 和Luo[39]和Wang[40]等的研究工作一样, 一定程度上丰富了对受限水体系奇特结构物性的认识.显而易见, Drude极化水分子模型描述的受限维度对水分子极性的强烈影响效应, 是恒定偶极矩的SPC/E模型无法企及的.为了进一步证实这个发现不依赖于受限条件的模拟细节, 我们额外模拟了单层石墨烯作为受限腔壁的(Drude模型)单分子层受限水和受限冰, 观察到450 K温度下两相水分子偶极矩平均值更显著的跳变, 见图 3(b)插图中三角形数据点.最后, 我们需要指出的是, 经典的Drude极化水分子是针对液体结构和动力学进行优化设计的, 并不能保证其具备精准描述单分子偶极矩的能力.偶极矩的变化本质上来源于电子结构的变化, 要严格研究这一单分子层受限冰-水相变带来的偶极矩跳变的问题, 需要考虑量子效应的影响.
图 5比较计算了单分子层受限冰氧原子结构因子(b~e, g~i), 左右两列(SPC/E和Drude水分子模型)结构因子结果具有非常相似的强度分布结构, 高波数峰随温度升高而退化的行为(温度接近熔点, 受限冰长程有序逐渐消失)也非常接近.值得注意的是, 结构因子第三个峰位(
图 6比较了SPC/E和Drude水分子模型模拟单分子层受限水xy平面内径向偶极-偶极关联函数
通过图 5和图 6的比较分析, 我们可以获得一个结论, 即极化水分子模型在描述受限冰和受限水的物相的结构性质时, 和SPC/E等恒定偶极矩水分子模型效果相当, 和水分子极性的涨落变化无明显关系.因此我们建议如果开展仅限于受限水物相结构性质的模拟研究的话, 选用传统恒定偶极矩水分子模型进行模拟足以胜任, 且更节省计算资源.
我们用公式(6)和(7)计算了单分子层受限冰和受限水的介电张量分量, 5×108 Pa受限压强、不同温度下SPC/E水分子模型和Drude极化水分子模型体系的比较结果都在表 2中详细列出. 表 2前4行数据比较了受限冰的介电张量分量.由于处于固相的水分子几乎没有平动和转动的能力, 100 K到400 K介电张量分量数值都非常接近1, 但最大数值只有1.7, 低于质子有序块体冰的介电常数数值(3~4)[44].低温下
下载:
导出CSV
| T/K | Phase | SPC/E | Drude | |||||
| εx | εy | εz | εx | εy | εz | |||
| 100 | s | 1.140(1) | 1.143(2) | 1.47(2) | 1.156(2) | 1.156(3) | 1.53(1) | |
| 200 | s | 1.177(2) | 1.172(3) | 1.50(1) | 1.167(2) | 1.166(4) | 1.49(2) | |
| 300 | s | 1.232(5) | 1.234(3) | 1.60(2) | 1.255(7) | 1.257(4) | 1.54(1) | |
| 400 | s | 1.386(5) | 1.399(12) | 1.70(1) | 1.52(2) | 1.54(3) | 1.62(1) | |
| 400 | l | 66(4) | 63(4) | 2.24(4) | 54(1) | 54(2) | 1.95(3) | |
| 450 | l | 59(3) | 58(6) | 2.17(6) | 47(1) | 47(1) | 1.90(2) | |
| 500 | l | 50(1) | 49(2) | 2.15(3) | 38(1) | 38(1) | 1.86(2) | |
| 550 | l | 40(1) | 41(1) | 2.15(5) | 29(1) | 29(1) | 1.85(2) | |
表 2中后4行数据是单分子层受限水介电张量分量的计算结果.受限液相水的
两种水分子模型预言的单分子层受限液相水
就可极化水分子模型对于恒定偶极矩水分子模拟的优势而言, 通过上述介电性质的计算结果, 我们获得一个重要认识, 模型改进(偶极矩变化能力的引入)对于介电常数描述的优势在分子偶极方向转动未被限制的情况下可以展现, 当水分子的偶极转动涨落在固相晶格和受限方向被大幅限制时, 极化模型描述静态介电性质的优势也被大幅削减.此外, 400 K单分子层受限水固-液两相介电张量
最后, 我们计算了偶极矩时间关联函数
下载:
导出CSV
| T/K | Phase | SPC/E | Drude | |||||
| 400 | l | 10.7(9) | 9.7(9) | 0.154(3) | 6.4(3) | 6.7(3) | 0.055(2) | |
| 450 | l | 5.6(2) | 5.5(5) | 0.126(2) | 3.4(3) | 3.4(2) | 0.042(1) | |
| 500 | l | 3.6(2) | 3.5(2) | 0.112(2) | 2.1(2) | 2.1(3) | 0.031(4) | |
| 550 | l | 2.4(1) | 2.5(2) | 0.096(4) | 1.3(1) | 1.4(2) | 0.029(3) | |
比较两种水分子模型在描述单分子层受限水动态介电响应性质上的表现, 我们发现SPC/E水分子模型预言的
在本工作中, 我们利用平衡态分子动力学的模拟和受限体系介电性质计算方法, 详细比较两种水分子模型——经典恒定偶极矩的SPC/E水分子模型和可极化Drude水分子模型, 在描述(0.65 nm限域尺寸、5×108 Pa限域压强、4个不同温度)单分子层受限水和受限冰结构与介电性质上的优劣势.统计了Drude水模型描述单分子层受限水和受限冰的水分子偶极矩量值概率密度分布, 并且计算了每个模拟体系的静态结构因子、静态偶极空间关联函数、静态介电常数、体系偶极时间关联函数和德拜弛豫时间.我们发现了: (1) Drude极化水分子模型能够描述出受限水分子极性的奇异特性, 受限冰-水结构两相共存的同一温度下, 水分子偶极矩平均值不相等; (2)两种模型描述受限冰结构和受限水空间偶极分布结构没有明显区别; (3)受限冰的三个分量和受限水z分量介电张量数值因水分子转动自由度的晶格和受限腔约束而被大幅限制; (4)受限水极化弛豫动力学行为也呈现强烈各向异性.两种水分子模型体系z分量偶极矩时间关联函数呈周期性振荡衰减结构.
通过本工作系统性比较研究, 我们获得了极化水模型在受限水体系结构与介电性质描述具体优势较为全面的认识.我们首次发现了极化水模型在描述低维度受限水分子极性特殊变化的绝对性优势, 并提出了极化水分子模型在受限水介电极化弛豫动力学性质描述上的潜在优势.我们也认识到极化水分子模型对于静态介电常数描述的优势在受限水系统中会因受限条件的增强而被大幅削减.因为发现了恒定偶极矩水分子模型和极化水分子模型在描述受限水和受限冰的静态结构性质的效果相当, 我们建议受限水物相结构的模拟研究无需选用极化水分子模型.综上, 我们认为Drude极化水分子模型在探索受限冰形核形貌、受限水结构相变动力学、介电性质的调控以及受限体系离子输运溶剂化动力学过程的模拟研究中是比常用的恒定偶极矩水分子模型更好的选择.
关于受限水z方向介电张量分量物理贡献组成, 本工作应用的分子模拟技术尚不予以澄清, 在未来的工作中, 我们将尝试采用量子Drude水分子模型结合能够描述恒定电势且伴随电荷涨落的受限腔壁[47]的分子模拟技术进一步深入探索.另外, 我们发现的受限条件下水分子在冰水两相中极性的跳变现象的物理机理, 也需要考虑量子效应的影响, 开展独立的研究提供更加充足的科学证据予以澄清.
Eijkel, J. C. T.; Berg, A. Microfluid Nanofluid 2005, 1, 249. doi: 10.1007/s10404-004-0012-9
Granick, S. Science 1991, 253, 1374. doi: 10.1126/science.253.5026.1374
Israelachvili, J. N. Intermolecular and Surface Forces, 3rd ed., Academic Press, Burlington, 2011.
Leikin, S.; Parsegian, V. A.; Rau, D. C.; Rand, R. P. Annu. Rev. Phys. Chem. 1993, 44, 369. doi: 10.1146/annurev.pc.44.100193.002101
Honig, B.; Nicholls, A. Science 1995, 268, 1144. doi: 10.1126/science.7761829
Cohen-Tanugi, D.; Grossman, J. Nano Lett. 2012, 12, 3602. doi: 10.1021/nl3012853
Szymczyk, A.; Fatin-Rouge, N.; Fievet, P. J. Colloid Interface Sci. 2007, 309, 245. doi: 10.1016/j.jcis.2007.02.005
Lin, Y.; Shiomi, J.; Maruyama, S.; Amberg, G. Phys. Rev. B 2009, 80, 045419. doi: 10.1103/PhysRevB.80.045419
Mikami, F.; Matsuda, K.; Kataura, H.; Maniwa, Y. ACS Nano 2009, 3, 1279. doi: 10.1021/nn900221t
Toney, M. F.; Howard, J. N.; Richer, J.; Borges, G. L.; Gordon, J. G.; Melroy, O. R.; Wiesler, D. G.; Yee, D.; Sorensen, L. B. Nature 1994, 368, 444. doi: 10.1038/368444a0
Ballenegger, V.; Hansen, J. P. J. Chem. Phys. 2005, 122, 114711. doi: 10.1063/1.1845431
Bonthuis, D. J.; Gekle, S.; Netz, R. R. Phys. Rev. Lett. 2011, 107, 166102. doi: 10.1103/PhysRevLett.107.166102
Zhang, C.; Gygi, F.; Galli, G. J. Phys. Chem. Lett. 2013, 4, 2477. doi: 10.1021/jz401108n
Schlaich, A.; Knapp, E. W.; Netz, R. R. Phys. Rev. Lett. 2016, 117, 048001. doi: 10.1103/PhysRevLett.117.048001
Fumagalli, L.; Esfandiar, A.; Fabregas, R.; Hu, S.; Ares, P.; Janardanan, Q. Y.; Radha, B.; Taniguchi, T.; Watanabe, K.; Gomila, G.; Novoselov, K. S.; Geim, A. K. Science 2018, 360, 1339. doi: 10.1126/science.aat4191
Algara-Siller, G.; Lehtinen, O.; Wang, F. C.; Nair, R.; Kaiser, U.; Wu, H.; Geim, A.; Grigorieva, I. Nature 2015, 519, 443. doi: 10.1038/nature14295
杜涵, 梁洪涛, 杨洋, 化学学报, 2018, 76, 483. doi: 10.6023/A18040128Du, H.; Liang, H. T.; Yang, Y. Acta Chim. Sinica 2018, 76, 483 (in Chinese). doi: 10.6023/A18040128
Liang, Z.; Du, H.; Liang, H. T.; Yang, Y. Mol. Phys. 2019, 117, 2881. doi: 10.1080/00268976.2019.1593532
Vega, C.; Abascal, J. L. F. Phys. Chem. Chem. Phys. 2011, 13, 19663. doi: 10.1039/c1cp22168j
Berendsen, H. J. C.; Grigera, J. R.; Straatsma, T. P. J. Phys. Chem. 1987, 91, 6269. doi: 10.1021/j100308a038
Lamoureux, G.; Harder, E.; Vorobyov, I. V.; Roux, B.; MacKerell, A. D. Chem. Phys. Lett. 2006, 418, 245. doi: 10.1016/j.cplett.2005.10.135
Lybrand, T. P.; Kollman, P. A. J. Chem. Phys. 1985, 83, 2923. doi: 10.1063/1.449246
Stuart, S. J.; Berne, B. J. J. Chem. Phys. 1996, 100, 11934. doi: 10.1021/jp961076d
Dequidt, A.; Devèmy, J.; Pádua, A. A. H. J. Chem. Inf. Model. 2015, 56, 260.
Shepard, A. C.; Beers, Y.; Klein, G. P.; Rothman, L. S. J. Chem. Phys. 1973, 59, 2254. doi: 10.1063/1.1680328
Gubskaya, A. V.; Kusalik, P. G. J. Chem. Phys. 2002, 117, 5290. doi: 10.1063/1.1501122
Fernandez, D. P.; Mulev, Y.; Goodwin, A. R. H.; Sengers, J. M. H. L. J. Phys. Chem. Ref. Data 1995, 24, 33. doi: 10.1063/1.555977
Jones, A. P.; Crain, J.; Sokhan, V. P.; Whitfield, T. W.; Martyna, G. J. Phys. Rev. B 2013, 87, 144103. doi: 10.1103/PhysRevB.87.144103
Kimmel, G. A.; Matthiesen, J.; Baer, M.; Mundy, C. J.; Petrik, N. G.; Smith, R. S.; Dohnálek, Z.; Kay, B. D. J. Am. Chem. Soc. 2009, 131, 12838. doi: 10.1021/ja904708f
Giovambattista, N.; Rossky, P. J.; Debenedetti, P. G. Phys. Rev. Lett. 2009, 102, 050603. doi: 10.1103/PhysRevLett.102.050603
Magda, J. J.; Tirell, M.; Davis, H. T. J. Chem. Phys. 1986, 84, 2901.
Werder, T.; Walther, J. H.; Jaffe, R. L.; Halicioglu, T.; Koumoutsakos, P. J. Phys. Chem. B 2003, 107, 1345. doi: 10.1021/jp0268112
Hockney, R. W.; Eastwood, J. W. Computer Simulation Using Particles, CRC Press, 1988, p. 55.
Yeh, I. C.; Berkowitz, M. J. Chem. Phys. 1999, 111, 3155. https://www.researchgate.net/publication/234978223_Ewald_Summation_for_Systems_with_Slab_Geometry
Plimpton, S. J. Comput. Phys. 1995, 117, 1. doi: 10.1006/jcph.1995.1039
Ryckaert, J. P.; Ciccotti, G.; Berendsen, H. J. J. Comput. Phys. 1997, 23, 327.
Jiang, W.; Hardy, D. J.; Phillips, J. C.; MacKerell, A. D.; Schulten, K.; Roux, B. J. Phys. Chem. Lett. 2011, 2, 87. doi: 10.1021/jz101461d
Wu, A. K.; Lin, S. C.; Karma, A. Phys. Rev. B 2016, 93, 054114. doi: 10.1103/PhysRevB.93.054114
Luo, C. F.; Fa, W.; Zhou, J.; Dong, J. M.; Zeng, X. C. Nano Lett. 2008, 8, 2607. doi: 10.1021/nl072642r
王新华, 冯莉, 曹泽星, 化学学报, 2014, 72, 487. doi: 10.6023/A13121223Wang, X. H.; Feng, L.; Cao, Z. X. Acta Chim. Sinica 2014, 72, 487 (in Chinese). doi: 10.6023/A13121223
Zhou, W.; Yin, K. B.; Wang, C. H.; Zhang, Y. Y.; Xu, T.; Borisevich, A.; Sun, L. T.; Idrobo, J. C.; Chisholm, M. F.; Pantelides, S. T.; Klie, R. F.; Lupini, A. R. Nature 2015, 528, E1. doi: 10.1038/nature16145
Wang, F. C.; Wu, H. A.; Geim, A. K. Nature 2015, 528, E3.
Chen, J.; Schusteritsch, G.; Pickard, C. J.; Salzmann, C. G.; Michaelides, A. Phys. Rev. Lett. 2016, 116, 025501. doi: 10.1103/PhysRevLett.116.025501
Petrenko, V. F.; Whitworth, R. W. Physics of Ice, Oxford University Press, Oxford, 1999.
Hill, N. E. Trans. Faraday Soc. 1963, 59, 344. doi: 10.1039/tf9635900344
Kindt, J. T.; Schmuttenmaer, C. A. J. Phys. Chem. 1996, 100, 10373. doi: 10.1021/jp960141g
杨鹏里, 王振兴, 梁尊, 梁洪涛, 杨洋, 化学学报, 2019, 77, 1045. doi: 10.3866/PKU.WHXB201905058Yang, P. L.; Wang, Z. X.; Liang, Z.; Liang, H. T.; Yang, Y. Acta Chim. Sinica 2019, 77, 1045 (in Chinese). doi: 10.3866/PKU.WHXB201905058
图 1 Drude水模型结构及参数示意图. D代表drude振荡粒子, M代表模型中刚性虚拟无质量电荷点.键长为0.09572 nm, HOH键角为104.52°
Figure 1 The schematic diagram of the Drude water model. "D" is drude oscillator particle, "M" is rigid virtual massless charge point. The bond length is 0.09572 nm and the HOH angle is 104.52°
图 2 单分子层受限冰和受限水的平衡态分子动力学模拟坐标快照. (a), (c)分别是SPC/E水分子模型模拟的受限冰和受限水. (b), (d)分别是Drude极化水分子模型模拟的受限冰和受限水.温度400 K, 压强5×108 Pa.白色球为氢原子, 红色球为受限固相冰中水分子内氧原子, 蓝色球为受限液相水分子中氧原子.绿色球为Drude水分子模型中M点, 黄色球为Drude水分子模型中的Drude振荡粒子.
Figure 2 Snapshot of monolayer square ice and water confined within two from the equilibrium molecular dynamics simulation. Panel (a) and (c): simulation results using SPC/E model potential, panel (b) and (d): simulation results using SWM4-NDP polarizable model. Confined pressure is 5×108 Pa and simulation temperature is set as 400 K. Hydrogen atoms are represented with white spheres, oxygen atoms are depicted as red (crystalline phase) and blue (liquid phase), respectively. Masslees site "M" and Drude oscillator particle in SWM4-NDP model are represented with the green and yellow spheres.
图 3 Drude极化模型水分子偶极矩随温度和物相的变化统计结果. (a)不同温度下, 常压下块体水(蓝色)和六角晶格块体冰(红色)的水分子瞬时偶极矩概率密度分布函数. (b) 5×108 Pa受限压强, 不同温度下, 单分子层受限水(蓝色)和受限冰(红色)中水分子瞬时偶极矩概率密度分布函数. (a)和(b)中插图为水分子平均偶极矩和温度的变化关系, (b)插图中三角形为石墨烯受限腔壁模拟结果
Figure 3 The calculated polarity variation of the water molecules, described by SWM4-NDP polarizable model, with changing the temperature and (confined) structural phases. Panel (a): measured probability distributions of the instantaneous water molecular dipole moments, blue lines represent bulk water phase at ambient pressure while red lines represent bulk hexagonal ices. Panel (b): under 5×108 Pa confinement pressure, blue lines represent probability results of confined monolayer liquid water phases, red lines represent results of confined monolayer ice phases. Inset figures in (a) and (b) show average dipole moments as the function of the temperature. Triangle symbols in the inset figure in (b) represent simulation results with graphene-confinement.
图 4 (a) Drude极化水分子模型单分子层限域冰-水相平衡分子动力学模拟分子坐标快照.压强5×108 Pa、温度400 K.氧原子颜色变化遵循结构序参量数值, 红色代表固相, 蓝色代表液相. (b)单分子层限域冰-水相平衡平面内细粒氧原子密度分布函数(实线)和粗粒偶极矩分布函数(虚线). x轴0点对应固液界面位置.结构序参量和粗细粒度分布函数计算方法参见文献[17]
Figure 4 (a) SWM4-NDP polarizable model, snapshot of the monolayer confined ice-water equilibrium interface obtained from a molecular dynamic simulation. Lateral pressure 5×108 Pa, under a coexistence temperature of 400 K. Oxygen atoms in (a) is colored coded based on the value of the structural order parameter, red for crystal and blue for melt. (b) In-plane fine-grained area-density of oxygen atoms profile (solid line) and coarse-scaled profiles for dipole moment of SWM4-NDP water (dashed line) for the monolayer confined ice-water coexistence. The zero point of x corresponding to the time averaged position of the crystal-melt coexistence interface. The calculation methods for the structural order parameter as well as the fine/coarse-scaled profiles can be found in the Ref. [17]
图 5 单分子层受限水(a, f)和受限冰(b~e, g~j)氧-氧原子结构因子.左侧SPC/E水分子模型, 右侧Drude极化水分子模型
Figure 5 The calculated oxygen-oxygen atomic structural factors for the monolayer confined water (a, f) and ices (b~e, g~j). Left column: results of SPC/E model systems, right column: results of polarizable SWM4-NDP polarizable model systems
图 6 单分子层受限水的径向偶极-偶极关联函数.左侧SPC/E水分子模型, 右侧Drude极化水分子模型.虚线、浅灰实线和蓝色实线分别代表xy方向、z方向分量和所有方向平均结果
Figure 6 Radial dipole-dipole correlation function (blue solid lines) and its xy (dashed lines) and z (grey lines) components of monolayer confined water. Left column: results of SPC/E model systems, right column: results of SWM4-NDP polarizable model systems
图 7
单分子层受限水偶极矩(分量)时间关联函数. (a), (b):
Figure 7
Components of the auto-correlation function. Panel (a), (b):
表 1 比较SPC/E和Drude两种水分子模型所预言的块体物相的分子偶极矩与介电常数.压强为105 Pa, 温度为300 K
Table 1. Comparison of the predictions (using SPC/E and Drude water models) for the bulk phases properties (dipole moments and dielectric constants) of water. At ambient pressure and temperature (105 Pa and 300 K)
| SPC/E | Drude | Expt.c | |
| μGas./(Debye) | 2.35a | 1.85b | 1.85 |
| μLiq./(Debye) | 2.35a | 2.46b | 2.95 |
| ε0 | 71 (3) | 78 (2) | 78.4 |
| a SPC/E model, gas-phase dipole moment | |||
下载: 导出CSV
表 2 4个温度下单分子层受限冰(用s表示)和受限水(用l表示)xyz三个方向的介电张量分量的计算结果.第3~5列是SPC/E水分子模型模拟体系的结果, 第6~8列是Drude极化水分子模型模拟体系的结果.括号内数字代表 95%置信度的统计误差
Table 2. Calculated dielectric tensor data of the monolayer confined ice (labeled with "s" for the solid phase) and water (labeled with "l" for the liquid phase), for the four temperatures and xyz components. The 3rd to 5th column present results of the SPC/E model systems, the 6th to 8th column present results of the SWM4-NDP polarizable model systems. Numbers in parentheses are 95% confidence errors for the last digits shown
| T/K | Phase | SPC/E | Drude | |||||
| εx | εy | εz | εx | εy | εz | |||
| 100 | s | 1.140(1) | 1.143(2) | 1.47(2) | 1.156(2) | 1.156(3) | 1.53(1) | |
| 200 | s | 1.177(2) | 1.172(3) | 1.50(1) | 1.167(2) | 1.166(4) | 1.49(2) | |
| 300 | s | 1.232(5) | 1.234(3) | 1.60(2) | 1.255(7) | 1.257(4) | 1.54(1) | |
| 400 | s | 1.386(5) | 1.399(12) | 1.70(1) | 1.52(2) | 1.54(3) | 1.62(1) | |
| 400 | l | 66(4) | 63(4) | 2.24(4) | 54(1) | 54(2) | 1.95(3) | |
| 450 | l | 59(3) | 58(6) | 2.17(6) | 47(1) | 47(1) | 1.90(2) | |
| 500 | l | 50(1) | 49(2) | 2.15(3) | 38(1) | 38(1) | 1.86(2) | |
| 550 | l | 40(1) | 41(1) | 2.15(5) | 29(1) | 29(1) | 1.85(2) | |
下载: 导出CSV
表 3 4个温度下单分子层受限水(用l表示)xyz三个方向的德拜弛豫时间计算结果.第3~5列是SPC/E水分子模型模拟体系的结果, 第6~8列是Drude极化水分子模型模拟体系的结果.括号内数字代表 95%置信度的统计误差
Table 3. Calculated Debye relaxation time of the monolayer confined water (labeled with "l" for the liquid phase) for the four temperatures and xyz components. The 3rd to 5th column present results of the SPC/E model systems, the 6th to 8th column present results of the SWM4-NDP polarizable model systems. Numbers in parentheses are 95% confidence errors for the last digits shown
| T/K | Phase | SPC/E | Drude | |||||
| 400 | l | 10.7(9) | 9.7(9) | 0.154(3) | 6.4(3) | 6.7(3) | 0.055(2) | |
| 450 | l | 5.6(2) | 5.5(5) | 0.126(2) | 3.4(3) | 3.4(2) | 0.042(1) | |
| 500 | l | 3.6(2) | 3.5(2) | 0.112(2) | 2.1(2) | 2.1(3) | 0.031(4) | |
| 550 | l | 2.4(1) | 2.5(2) | 0.096(4) | 1.3(1) | 1.4(2) | 0.029(3) | |
下载: 导出CSV
扫一扫看文章
扫一扫关注我们
