液体界面现象的分子动力学模拟.doc

液体界面现象的分子动力学模拟 在研究生期间，我的课题跟计算模拟有关，故本文利用分子动力学计算机模拟的方法,研究了液体的界面现象。主要分三部分来说明，这三部分分别是分子动力学模拟方法，分子动力学模拟过程和具体的结果与讨论。 1.分子动力学模拟方法 要研究气液界面性质，在模拟体系中必须要有气相、液相以及气液界面的存在，于是我们采用如图1所示的模拟盒，让液相和气相处于模拟盒的中部，固相处于模拟盒的两侧，这样在整个模拟体系中将有多种界面。 对气液界面现象进行分子动力学模拟研究时，以氩和甲烷为例，其分子间的相互作用采用 L-J12-6 模型，分子 i 和 j 相互作用的势函数表示为： 采用直角坐标系0xyz，分子在初始时刻为饱和液相，分子以边长为L的正方体FCC结构分布于模拟盒中部，其质心在原点（0，0，0）处。模拟过程中，系统采用正则系综，用Verlet算法对运动微分方程求解，并不断调整体系质心使之处于原点位置。其中的运动方程和运动速度分别表示为： 式中：Uij为分子i和j之间的势能，rij为它们之间的距离，fi为第i个分子所受到的作用力，mi为分子i的质量，vi为其速度，t为时间，ri为分子i的位移。由于我们采用的是正则（NVT）系综，需要进行温度修正。模拟盒在x、y方向的长度Lx=Ly=L，为了保证足够的气相空间，取在z方向的长度Lz=3L。本模拟在三个方向上都采用周期性边界条件。在计算机分子动力学模拟过程中，将所用到的各物理量进行无量纲化处理，以确保模拟体系各物理量之间的数量级较为接近，减小模拟过程中出现的误差。 2．分子动力学模拟过程 在进行分子动力学模拟之初，每个分子的具体位置及其速度的大小和方向需要赋予初值。在该研究中，采用直角坐标系0xyz，分子在初始时刻以边长为L的正方体FCC结构分布于模拟盒中部，其质心在原点（0，0，0）处；接着随机地给每一个分子速度赋予初值（本文采用Boltzmann分布对速度赋值，以节约机时），这样就可以跟踪每一个分子的运动状态。由于对气液界面现象的研究需要在气液共存的平衡系统中进行，计算机分子动力学模拟需要让分子经过较长时间的运动，以使整个系统内部达到平衡状态，在达到平衡态之后才对相关物理量进行统计计算。从其他人的分子动力学模拟研究结果可以看出，在不同的温度和分子数条件下，系统内部达到平衡所需要的时间不同，系统的分子数越多，它达到平衡所需的时间越少；而系统的温度越高，系统达到平衡所需的时间越长。 在判断气液系统是否达到平衡态时，需要不断地对系统的能量、气液界面的密度分布、气相和液相的分子数以及压力分布等进行统计。当分子动力学模拟系统达到平衡时，气相和液相的分子数应分别近似保持不变，当然，如果采用的是分子数N不变的系综，在对气液相分子数进行统计时，只需对其中一相的分子数进行统计就行了；而且，密度分布和压力分布也基本上不随时间的变化而变化。另外还需要将模拟中得到的气相和液相各自的热力学量与相关的文献进行对比，如果二者吻合较好，那么就可能认为模拟体系达到了气液平衡态。 3. 结果与讨论 步长取0.01个时间单位,截断距离是2.5,为了研究界面现象和表面张力,在给定了初始位形之后,运行了20000步以达到动态平衡.图2给出了平衡后运行到238000步时全部原子瞬时俯视图,固体原子仍然重迭为简单立方图象,液、气界面已不是初始位形时的明确位置,它具有约三个rm厚度的过渡层,固、液界面也形成一个约有三个rm厚度的过渡层,因此,真正具有液体结构,并具有相应密度的只有初始位形的中间约五、六个rm厚度部分。有趣的是,在气、固交界处,固体“吸附”了大约二个rm厚度的原子层,这与实验事实相符。此外,我们还看到,在固体与液体交界处,液体虽不象固体那样有晶格结构,但它也有类晶格分布。关于这一点,可以从图3中看得清楚,它是在27900步的瞬时位形,在z方向取一个厚度的原子俯视图,图中有规律的三角形结构,可以认为是类晶体结构。 图4给出了2,0000步的平均密度分布.把液、气部分分成40个小区,真正的液沐只有10一12个区,其余3一4个固气“吸附区”,5一6个固液过渡区和5一6个液气过渡区。 如果不考虑原子分布的各向异性,即不考虑角度因素,原子零级径向分布函数Gn(R)曲线如图5(a)一(b),它们是21000步的平均值,n表示图4中的分区位置,Gn( R)/G1(R)是各区的粒子数分布与第1区粒子数分布的比值,可以看到。n=1一3各区均为典型的液体区,n=4一6是过渡区.曲线下面覆盖的面积表示该区的粒子密度,有趣的是,不仅在典型液体区,即使是过渡区,粒子的统计分布也不是杂乱无章的,相反地,它们有规律的分布,告诉我们,液体中的粒子分布仍是有序的. 4. 结论 氩气的固一液、液一气和气