基于Metropolis方法奇偶数碳团簇势能之比较.docVIP

基于Metropolis方法奇偶数碳团簇势能之比较 　　本文基于Metropolis模拟方法一般性地得出了原子数N20的碳团簇的结构,并用Brenner势计算了其势能。发现平均势能曲线当N32时出现锯齿状起伏,与实验观测到的N40时只有偶数碳团簇的现象相吻合。本文结果对于研究奇数碳团簇的结构具有指导意义。 　　 　　一.引言 　　自从1985年提出C60结构的“布基球”猜想以来[1],人们对碳团簇的性质做了大量的实验及理论研究[2,3]。其中引人注目的现象之一是在较大的碳团簇(40个以上的碳原子)中,只发现由偶数个碳原子构成的团簇,而缺乏由奇数个碳原子构成的团簇。对此可以根据碳团簇的几何性质定性地加以解释。设这些团簇都是中空的多面体结构并且碳原子的配位数为3,即每1个顶点有3条棱与之相连。因此顶点数(即碳原子数)N与棱数L(即sp2键数)之间满足关系式3N=2L,由L为整数知N必为偶数。若进一步假设碳团簇由P个五边形和H个六边形构成,则顶点数N=(5P+6H)/3,棱数L=(5P+6H)/2,外表面数为F=P+H,代入欧拉定理N+F-L=2,则有N=20+2H,P=12。满足该条件的最小碳团簇即由12个五边形构成的正十二面体,原子数为20。在此基础上每增加1个六边形,碳原子数必然增加2个。这种结构的碳团簇即通常所称富勒烯[4],显然富勒烯必然含有偶数个碳原子。 　　然而,依照以上的思路,人们必然想知道是否存在稳定的奇数碳团簇?其结构如何?奇数碳团簇的平均势能与偶数碳团簇究竟差多少?这样的定量分析目前尚未见报道。就此问题本文基于Metropolis模拟方法计算出可能的奇数碳团簇结构及其平均势能,并与同样方法得出的偶数碳团簇的结构及平均势能进行系统比较。 　　二.模拟方法 　　为了普遍地得出碳团簇的结构,我们采取Metropolis模拟方法。即将N个碳原子随机地置于空间,并用Brenner势[5,6]计算其势能E0,随机选择一个碳原子将其空间位置游走一个微量δ,并计算其势能E1。若E1E0则接受该碳原子的位置改变量δ,否则只以概率exp[-(E1-E0)/kBT]接受δ。其中kB为Bolzman常数,T为绝对温度。大量重复以上操作,即可得到在环境温度T下“最可几”(Most Probable)的碳团簇结构,该结构的势能接近于某一个局部最小值。 　　为减少弛豫时间和计算量,我们在模拟中使用了球面近似。也就是将N个碳原子的初始位置赋为一个球面,并使游走量δ的切向分量远远大于(100倍以上)其径向分量。这一限制是基于碳团簇的结构都是接近球面的多面体假定之上的,也是基本合理的。初始球面的半径根据C60最稳定结构(即Buckyball)的半径按比例估算得出。已知Buckyball结构的半径约为R0=3.4~3.5A,故N个碳原子组成的笼结构半径约为R=sqrt(N/60)*R0。模拟中环境温度T置为6000K。对于δ我们选择使其切向方向的分量不超过石墨键长的十分之一,即在0.15A左右。 　　用Metropolis方法得到的碳团簇结构还不处于一个势能极小点。为得到在此结构附近的一个势能极小的结构,我们采取分子动力学方法将该结构优化到其附近的一个局部势能凹点。这个势能凹点所对应的结构就是我们最后得到的结构。 　　我们希望看到用以上方法得到的结构就是N个碳原子构成的团簇的势能最小的结构,然而通常得到的都还不是势能的全局最小点,而是一个局部极小点。为了尽可能地逼近势能的全局最小值,我们对于每一个原子数N,用不同的初始位置和随机数,构造了多个样本,每个样本计算20,000,000步,然后从这些样本结果中挑出势能最小的结构。由于碳团簇的异构体很多,例如C60笼结构仅满足纯五边形和六边形约定的,就有1812个之多[7],如果算上非经典的C60异构体就更多了。目前还不能保证我们的模拟一定就能给出势能最低的结构。 　　三.结果与讨论 　　表1给出了原子数N在20~83之间碳团簇的Brenner势。将团簇的总势能除以其相应的原子数N,可以得到图1所示的团簇平均势能曲线。很明显,随着原子数的增加,碳团簇平均势能总体上是呈下降趋势的。当N32以后,这种下降趋势明显地附加了锯齿形的起伏:锯齿的凹点对应于偶数的团簇, 凸点对应于奇数的团簇。表明随着碳团簇逐渐长大,其稳定性并非单调递增的。偶数团簇比相邻的奇数团簇要稳定得多。图2给出了随原子数增长的团簇势能的差分曲线,表示每增加一个碳原子,团簇总势能相应变化几个电子伏特。显然,势能的差分也是一条锯齿形的曲线,奇数与偶数分别对应于曲线的凸点与凹点。粗略的估计,这两者的差异大约有1.5个电子伏特。 　　图3、图4分别给出了奇数及偶数碳团簇的结构图。可以看出,对于绝大部分奇数团簇,都有一个原子配位