到目前为止,关于液态水偶极矩的“正确”数值仍然是一个有争议的问题。在实验上,Soper小组利用X-射线衍射方法第一个给出受到周围环境影响的液态水平均偶极矩为2.9D。另一种归纳方法计算了冰Ih晶格的偶极矩为3.09D,这种方法利用比Coulson和Eisenberg更加准确的输入数据和更高级的多项展开式,因而比Coulson和Eisenberg给出的偶极矩2.6D更加准确。密度泛函理论给出冰Ih晶格的单体偶极矩从2.33D到2.97D不等。利用abiniito分子动力学模拟计算给出298K,1.013×105Pa下液态水的偶极矩为2.66D,2.43D和2.95D,这三种模拟利用不同的方法把电子密度分配到分子上。不同的分配电子密度的方法可以明显地改变偶极矩的数值。
ABEEM-7P方法是把每个水分子的密度划分到原子、键和孤对电子上,计算的平均偶极矩为2.80D,与Soper利用X-射线衍射方法得到的偶极矩(2.90D)非常接近。TIP5P方法把电子密度划分在两个氢原子和氧原子的两个孤对电子上给出的平均偶极矩为2.29D;TIP4P-FQ方法利用TIP4P的几何结构,把电子密度划分在两个氢原子以及HOH键角朝向H原子方向的角平分线上,计算的偶极矩为2.6D;POL5模型电子密度的划分类似于TIP5P方法计算的平均偶极矩为2.71D;MCDHO方法利用的是abinitio的偶极极化模型,给出了较大的平均偶极矩3.0D,与这些经验分子力场相比较,ABEEM-7P模型得到的平均偶极矩与POL5模型计算的结果最为接近,同时两种非刚体模型(ABEEM-7P和MCDHO)虽然都给出了较大的偶极矩,但与最新的实验值非常接近。利用单个水分子电负性均衡和电荷守恒的ABEEM-7P-2模型计算的水分子团簇(H2O)n(n=1~6)的平均偶极矩分别为:1.855D、2086D、2.316D、2.465D、2.528D和2.560D,六聚体水的书状、棱状和笼状低能构象的平均偶极矩分别为2.522D、2.462D和2.483D。与这些数值相比较,ABEEM-7P的液态水平均偶极矩2.80D要大一些,这说明周围环境的不同对于极性分子电荷的重新分布有很大的影响。
讨论液态水偶极矩的另一个重要问题就是单体水偶极矩及其分布。在298K,1.013×105Pa下ABEEM-7P力场模拟的液态水的单体偶极矩数值列在表8-5中,相应的分布图画在图8-5中。ABEEM-7P模型给出较宽的单体水偶极矩的分布范围:
2.30~3.39D,最大值和最小值之间相差了1.09D,与abinitio和极化经验力场模型给出的液态水体系中单体偶极矩分布范围较宽这样一个普遍存在的事实相符合。图8-5更加清晰地说明了298K液态水中单体偶极矩的分布情况。ABEEM-7P模型计算的偶极矩分布的半峰宽为0.422,和TIP4P-FQ模型计算的室温下偶极矩分布的半峰宽0.42非常接近,而SPC-FQ模型的半峰宽为0.49[72],这是由于SPC模型中静电相互作用位置距离较大。因而,SPC-FQ模型比TIP4P-FQ和ABEEM-7P模型具有更大的极化性。另外,正如我们前面所提到的,液态水体系中每个单体水周围电场的改变,使得O-H键长增大,HOH键角减小,从而也导致ABEEM-7P等一些非刚体极化模型计算了较大的偶极矩。
介电常数依赖于偶极矩的大小、单位体积偶极子的个数以及偶极子的方向。介电常数为液态水中单体水偶极矩的正确计算提供了另一个判断依据。在所研究的立方体系中,介电常数和体系总偶极矩浮动之间的关系取决于系统中长程相互作用力的计算。介电常数
0可以利用中心模拟盒子总偶极矩M的波动计算:
其中,是温度T时体系的密度,k是Boltzmann常数,Nmol是体系中总的分子个数。一般来说,如果经验水分子模型计算出298K和1.013×105Pa的偶极矩在2.4D到2.7D之间,那么它得到的介电常数会接近80。ABEEM-7P、TIP5P、TIP4P-FQ和POL5水分子力场模型计算介电常数以及实验值都列在表8-4中。从表中可以看出,与实验值78.3相比较,ABEEM-7P模型计算的介电常数761可以很好地模拟298K液态水的介电性质。与其它水分子力场比较,ABEEM-7P与TIP4P-FQ模型计算的结果798较为接近,与实验值符合得都很好,而POL5和TIP5P模型则较高地估计了介电常数,分别为988和822。因而我们可以得到这样的结论,基于密度泛函理论的浮动电荷模型ABEEM-7P和TIP4P-FQ都能够正确反应出在不同环境下液态水体系电子密度的重新分布。
利用液态水中单个分子间的相互作用能Uliquid可以计算摩尔汽化热Hvap:
H(T)U(T)PVU(T)RTvapliquidliquid(8-2)其中,P是压力,V是气态和液态之间摩尔体积改变量,R是气体常数,T是相应的温度。表8-4给出了298K,1.013×105Pa下ABEEM-7P模型计算的分子间相互作用能和汽化热为:-10.26kcal.mol和10.85kcal/mol,实验值分别为:-9.92kcal/mol和10.51kcal/mol,ABEEM-7P的结果与实验值符合得很好,相互作用能和汽化热比实验值大0.34kcal/mol,较大的相互作用能和汽化热主要是由于ABEEM-7P模型的参数是利用气态水分子团簇的各种性质(如优化的几何结构、偶极矩、能量和振动频率等)拟合的,并且没有经过任何修改直接转移应用到模拟液态水体系的性质。从表8-5中我们也可以看出,POL5模型计算的汽化热(或者是水分子间相互作用能量)与实验值几乎完全符合,这是由于POL5模型的参数主要来自于对液态水分子间相互作用能量的调节。TIP4P-FQ和TIP5P计算的结果也和实验值比较接近,TIP4P-FQ的分子间相互作用能和汽化热分别为:-9.89kcal/mol和10.48kcal/mol,相应的,TIP5P结果为:
-9.87kcal/mol和10.46kcal/mol。MCDHO则给出比ABEEN-7P结果更大的分子间相互作用能和汽化热:-10.40kcal/mol和10.99kcal/mol,比实验值大0.48kcal/mol。
径向分布函数是描述体系结构的一种非常有用的方法,尤其是对液体结构的描述。考虑与原子距离为r,厚度r的圆形壳层,该壳层的体积为:
如果单位体积内的粒子数为,则在这个壳层中总的粒子数为4πr2r。
对径向分布函数g(r)给出了一个原子(如果模拟溶液体系则为一个分子)距离另一个原子为r时找到它的可能性,这也是和理想气体分布比较而言的。g(r)是无量纲的,高级径向分布函数(也就是三级径向分布函数)也是可以被定义的,但是很少计算高级径向分布函数,所以通常所指的“径向分布函数”就意味着对径向分布函数。在晶体中,径向分布函数中有一个无限大的尖峰,它的位置和高度就是晶格结构的特征。
液态的径向分布函数是位于固态和气态之间的中间体,液态的径向分布函数有一些较小位置、较小数量的峰,在较大的位置逐渐衰退为一个常数。为了从模拟中获得对径向分布函数,每一个原子或分子的原子都被归类为“distancebins”或直方图(histogram)。在每一个bins中,邻近的分子或原子数目在整个模拟过程中被平均化。
比如说在模拟过程中,对每个原子或分子在其周围2.52.75或2.753.0范围内的邻近原子或分子分别计数。计数或者是在模拟过程中直接得到,或者是通过分析结果组态得到。
从X-ray衍射实验中可以得到径向分布函数,晶体中原子分布的特征X-ray衍射图是明、尖的一些点,分析X-ray衍射分布可以计算实验上的分布函数,而这些分布函数可以用来和模拟结果相对比。利用径向分布函数能够得到许多热力学和动力学性质,如果假设力是成对叠加性,这些性质包含了理想气体和真实气体部分,如真实气体的能量和压力等。对于分子,如果要研究径向分布的真正本质,就必须考虑分子的取向。分子的径向分布函数通常是利用两个固定点之间的距离来计算的,比如质量中心之间的距离,如果考虑分子的方向就要用到方向分布函数来补充。对于线性分子,这些方向分布函数是利用分子轴方向之间的夹角来计算的,角度的变化从-180到180之间。对于复杂的分子,通常计算的是一些点点之间的分布函数,比如三点水分子模型,就定义了三个gOO,gOH和gHH方程。点点模型的优势在于可以直接和实验的X-ray散射相关。O-O,O-H和H-H的径向分布函数对于验证液态水的各种势能模型的正确与否有很大的作用。
液态水分子体系的径向分布函数无论在实验还是理论上都有大量的报导。Soper等人指出,液态水分子体系径向分布函数的峰高在实验上有很大的不确定性,主要是由于实验上利用了不同的方法除去自身原子散射的贡献部分,但是峰的位置却应该是非常确定的,因而峰的位置可以为建立模型势能提供可靠的依据。图8-6给出了ABEEM-7P势能函数的O-O,O-H和H-H的径向分布函数gOO,gOH和gHH,以及和最新实验的比较。与固定的点电荷模型相比较,浮动电荷模型计算的gOO第一峰的位置都较大,这是由于浮动电荷模型中增大的电荷给出了分子间较大的长程序列(long-rangedordering)。O-O之间的径向分布函数gOO是和液态水体系中形成氢键的两个水分子密切相关的,从图8-6a中我们可以看到,ABEEM-7P模型gOO的第一峰位置出现在大概3.0左右,虽然比实验位置稍微大一点,但还是符合得很好,这说明,ABEEM-7P模型能够较好地模拟液态水体系的氢键,毫无疑问,这应该归功于我们引进的描述氢键相互作用区的参数klp,H(Rlp,H)。与实验相比,ABEEM-7P模型gOO的第一峰略显高一点、宽一点,正如前面所述这并不能影响利用ABEEM-7P势能函数对液态水体系结构的描述。径向分布函数gOO的第二峰是和体系中水分子的四面体结构相关连的,ABEEM-7P模型的第二峰位置出现在4.75,与Soper的第二峰位置4.50非常接近,说明ABEEM-7P模型能够充分体现单体水分子的四面体结构。利用实验中第一波谷的位置3.36作为积分上限计算体系中每个水分子周围平均氢键的个数,实验值为4.5,ABEEM-7P模型计算值为4.75,TIP4P-FQ为4.4,SPC-FQ为4.2,可见,ABEEM-7P模型的计算值与实验值非常接近。
转移(transport)是指分子从一个区域流向另一个区域的性质,比如溶液处于一种溶质分布不平衡的状态,这样溶质就要发生扩散直至浓度达到均一的状态。如果体系中有一定的热力学梯度,则能量就要流动直到温度处处相等,同样动量梯度会产生黏度(viscosity)。转移的存在意味体系并没有达到平衡。根据Fick的第一扩散定律,物质Jz的扩散速度(也就是单位时间经过单位体积的扩散量)等于扩散系数(D)乘以浓度梯度:
JD(dN/dz)z(8-4)其中,N是粒子密度(单位体积内原子的个数),(8-4)式是沿着z轴方向的扩散,负号表示在负的浓度梯度方向流动增加。当模拟体系是纯液体时,系数D就是自扩散系数(selfdiffusioncoefficient)。扩散系数是和均方位移(meansquaredistance)2r(t)r(0)相关的,Einstein方程指出均方位移等于2Dt。在三维体系中,均方位移表示为:
上式中,只有当t才是严格成立的,()irt是指分子i的质量中心在时刻t的位置向量,平均值是指要计算体系中所有的分子。
作为一种动力学性质,扩散系数是和较短时间相关联的,无论在实验和模拟过程中都可以非常准确地获得,因而扩散系数是评价经验力场模型的一个非常重要和基本参数,另外扩散系数不但能够反应出分子间势能函数的短程相互作用,还能够表现出长程相互作用。ABEEM-7P,POL5,TIP4P-FQ和TIP5P模型计算的扩散系数列在表8-5中,分别为1.8,1.81,1.9和2.62109m2/s。相对于固定点电荷模型,一般来说,极化模型的扩散系数较小,这可能是由于极化模型中较大的静电相互作用导致的。与实验值2.3109m2/s相比较,极化模型都给出较小的结果,而固定电荷模型的扩散系数较大。扩散系数还和模型势能函数的非刚体性相关,但是目前的研究还存在一定的争议,一些计算结果表明加入键长和键角的振动可以增大自扩散系数,反之一些结果认为会降低自扩散系数。利用ABEEM-7P模型的计算结果表明,加入键长和键角的振动降低了扩散系数。
8.3液态水温度依赖性性质
许多经验力场模型都曾经研究过液态水体系的性质和温度的依赖关系。由于不同的温度和环境可以使水分子体系的各种性质发生一定的改变,能够在不同条件下正确地描述水分子体系对于建立水分子力场模型是非常重要的。大多数的非极化和极化模型,比如TIP4P、SPC/E和TIP4P-FQ等模型都能够给出不同温度条件下水分子体系的各种性质,如径向分布函数、单体水平均偶极矩、介电常数和汽化热等。
由Bordhole,Sampoli和Vallauri建立的专门计算冰Ih体系的极化BSV势能函数,虽然BSV模型能够很好得计算冰Ih的各种性质,但是它不能够模拟液态的水分子体系。到目前为止,没有一个势能函数能够同时描述液态水和固态冰体系的各种性质。