（粒子物理与原子核物理专业论文）纳米碳团簇力学性能的模拟研究.pdf

摘要 自从c 6 0 等小富勒烯和碳纳米管的发现以及实验制备技术的发展，人们对于 这些新型纳米结构材料的物理、化学性质进行了大量的探索研究。在固体表面散 射的研究中发现，小富勒烯具有相当高的结构稳定和可逆形变的特点；小富勒烯 沉积在固体表面以产生具有某种单一物理特性的薄膜。近年来，随着碳纳米管内 充c 6 0 等纳米新材料的涌现，纳米金属、单壁碳纳米管以及内充c 6 0 的单壁碳 管等纳米结构材料的力学性能的研究引起了众多学者的兴趣。这些研究具有工业 应用的前景。虽然小富勒烯的性质与合成薄膜的性质密切相关，但对于其力学性 质我们却知之甚少。本文采用分子动力学模拟，用b r e n n e r 多体势以及长程对势 相结合来描述原子间的相互作用，通过建立小富勒烯和内充小富勒烯的碳纳米 管的压缩模型，以研究其纳米力学性能。主要内容和结论如下： 1 我们建立了石墨层对小富勒烯( c 2 0 ，c 3 6 ，c 6 0 和c 7 0 ) 压缩和恢复的微观 模型，采用分子动力学方法研究小富勒烯表现出的力学性能。我们发现，在富勒 烯被压缩至只有原来1 2 到1 ，3 半径( 压缩方向的半径) 范围并释放的过程中， 所有的富勒烯都表现出相似的弹性行为。体系的总势能与石墨层所受的压强均随 着应变可逆变化，而这正是富勒烯结构变化可逆的表现。当富勒烯被压缩至一个 碟状结构时，压缩几乎饱和，而在释放外界压强后，富勒烯重新返回初始的笼状 结构。为了比较不同富勒烯的弹性性能，我们采用应变能对应变的二次导数作为 体系的近似杨氏模量，发现随着富勒烯尺度的增大，杨氏模量也会增大。 2 ，我们采用分子动力学模拟对填充小富勒烯( c 2 0 ，c 3 6 和c 6 0 ) 的单壁碳纳米 管( a r m c h a i r 型) 进行轴向压缩，研究其力学性能。研究表明，对同一种填充物 的碳纳米管压缩时，在较小的应变范围内，纳米管的应力与应变呈线性关系，改 变填充富勒烯的种类和密度，其杨氏模量的值几乎相同。在纳米管压缩至某一临 界应变值时，出现了塑性扭曲，其应力也达到了一个临界值，称为扭曲力。我们 发现，对同一种填充物，随着填充密度增大，纳米管的临界应变增加，而且扭曲 力会增随之增大，表示其弹性范围扩大。另外，我们还发现了一些反常现象，当 纳米管直径相对填充小富勒烯的直径差别较大时，应力的临界值几乎不受填充密 度的影响，甚至会随填充密度的增大而减小。 a b s t r a c t s i n c et h ed i s c o v e r yo ft h ef u l l e r e n e sa n dc a r b o nn a n o t u b e sa n dt h et e c m o l o g i c a l d e v e l o p m e n ti nt h es y n t h e s i so fm a c r o s c o p i ca m o u n t s ，av a r i e t yo f s t u d i e so nt h e i r c h e m i c a la n dp h y s i c a lp r o p e r t i e sh a sb e e nm a d e t h eh i g hs t a b i l i t ya n dt h er e v e r s i b l e d e f o r m a t i o no fs m a l lf u l l e r e n e sw e r eo b s e r v e di nt h es c a t t e r i n go ff u l l e r e n e so nt h e s u r f a c e s m a l lf u l l e r e n e sw e r ed e p o s i t e do nt h es u r f a c et op r o d u c ef i l m so fu n i q u e p r o p e r t i e se x p e r i m e n t a l l y f u r t h e r m o r e ，e n c a p s u l a t e dc 6 0 一f u l l e r e n e si nn a n o t u b e s ，a n e wk i n do fn a n o s t r u c t u r e dc a r b o nm a t e r i a l ，h a v eb e e ne x p e r i m e n t a l l yo b s e r v e d i n r e c e n t y e a r s ，t h e m e c h a n i c a l p r o p e r t i e s o fn a n o s t r u c t u r e d m a t e r i a l s ，s u c h a s n a n o m e t a l s ，s i n g l e w a l l - c a r b o n n a n o t u b e s ( s w c n t s ) ，a n ds w c n t s f i l l e dw i t hc 6 0 ， h a v eb e e ns t u d i e db e c a u s et h e i r p o t e n t i a li n d u s t r i a la p p l i c a t i o ni se x p e c t e d h o w e v e r , m u c hl e s si sk n o w na b o u tt h ee l a s t i cp r o p e r t yo fs m a l lf u l l e r e n e s ，w h i c hi sr e l a t e dt o t h em e c h a n i c a lp r o p e r t yo ff u l l e r e n e - f i l m s i nt h i sp a p e r , t h em o l e c u l a rd y n a m i c s ( m d ) s i m u l a t i o nw a su s e dt oi n v e s t i g a t et h e n a n o m e c h a n i c a lp r o p e r t i e so fs m a l l f u l l e r e n e sa n ds w c n t sf i l l e dw i t hs i n a i lf u l l e r e n e s ( c 6 0 ，c 3 6a n dc 2 0 ) w i t ha c o m p r e s s i o nm o d e l t h ei n t e r a c t i o n b e t w e e nc a r b o na t o m sw a sd e s c r i b e d b ya c o m b i n a t i o no ft h em a n y - b o d yb r e n n e rp o t e n t i a lw i t hat w o - b o d yp a i rp o t e n t i a l o u r m a l nr e s u l t sc a nb es u m m a r i z e da st h ef o l l o w i n g ： 1 t h en a n o m e c h a n i c a lp r o p e r t i e so fs m a l lf u l l e r e n e s ( c 2 0 ，c 3 6 ，c 6 0a n dc 7 0 ) c o m p r e s s e db e t w e e ng r a p h i t ep l a r i e sw e r ei n v e s t i g a t e db ym d s i m u l a t i o n w ef o u n d a l lf u l l e r e n e se x h i b i ts i m i l a re l a s t i cb e h a v i o re v e ni ft h em i n i m u mr a d i io fc l u s t e r si n t h ed i r e c t i o no f c o m p r e s s i o n a r er e d u c e dt o1 3o r1 2t h a to ft h ef r e ec l u s t e r sb o t h p o t e n t i a le n e r g y o ft h e s y s t e m a n dr e s t o r i n gp r e s s u r eo nt h e g r a p h i t ep l a n e s a s f u n c t i o n so ft h ea x i a ls t r a i nw e r ef o u n dt ob er e v e r s i b l e t h ef u l l e r e n ew a sd e f o r r n e d s e v e r e l yt ob ea d i s c - l i k es t r u c t u r ed u et oc o m p r e s s i o n t h ed e f o r m a t i o nw a sa l m o s t s a t u r a t e d t h e n ，i nt h ep r o c e s so fd e c o m p r e s s i o nt h ef r e e f u l l e r e n es t r u c t u r ew a s r e t u r n e d t oc o m p a r et h ee l a s t i c i t yo fd i f i e r e n tf u l l e r e n e s ，t h es e c o n do r d e rd e r i v a t i v e o fs t r a i ne n e r g yt os t r a i nw a sc a l c u l a t e da say o u n g sm o d u l u sa p p r o x i m a t e l y t h e m o d u l u sb e c o m e s l a r g e rw i t hi n c r e a s i n gf u l l e r e n es i z e 2 t h em e c h a n i c a lp r o p e r t i e so fs w c n t sf i l l e dw i t hs m a l lf u l l e r e n e s ( c 2 0 ，c 3 6a n d c 6 0 ) w e r ei n v e s t i g a t e du s i n gm o l e c u l a rd y n a m i c s ( m d ) s i m u l a t i o n ，w ef o u n dt h a t t h es t r e s so fs w c n ti n c r e a s e sl i n e a r l yw i t ht h es t r a i na n dt h ey o u n g sm o d u l u so f s w c n tw i t hd i f f e r e n tf i l l i n gd e n s i t i e si sa l m o s tt h es a m ef o rs m a l ls t r a i n t h es t r a i n i n c r e a s e su n t i li tr e a c h e sac r i t i c a lv a l u e ，w h e r es w c n t b e g i n sb u c k l i n g i tw a s o b s e r v e dt h a ta tr o o mt e m p e r a t u r et h eb u c k l i n gf o r c e ，w h i c hi sd e f i n e da st h es t r e s s c o r r e s p o n d i n gt ot h ec r i t i c a ls t r a i n ，b e c o m e sh i 曲e ra st h ef i l l i n gd e n s i t yo fs w c n t i si n c r e a s e di ng e n e r a l ，w h i c hi m p l i e sh i g h e re l a s t i c i t yo fs w c n t h o w e v e r , i nt h e c a s eo fs w c n tw i t hl a r g e rr a d i u sf i l l e dw i t hs m a l l e rf u l l e r e n e s ，t h ed e p e n d e n c eo f t h eb u c k l i n gf o r c eo nt h ef i l l i n gd e n s i t ym i g h tb ed i f f e r e n t h i 第一章引言 第一章引言 碳在元素周期表中位于第一短周期中的i v 族位置，有利于形成共价键类型 物质。碳的同素异构体的多样性也是该元素的一个重要特点。出于碳原子可以形 成s p 、s p 2 和s p 3 杂化轨道而与其它原子成键，所以自然界中存在着众多碳的化 合物。碳元素捌有大量的同素异构体，如最常见的石墨和金刚石晶体。石墨晶体 是自然界大量存在的碳单质，是由s p 2 杂化的碳原子形成的具有层状六边形结构 的一种碳晶体。它具有柔软，耐高温的特性，是一种良好的固体润滑剂，同时它 也是一种良导体。而自然界当中存在的另一种碳单质金刚石确具有与石墨截然不 同的性质。它是由s p 3 杂化的碳原子之间形成的具有面心立方结构的晶体。这种 结构导致金刚石独特的力学、电学和光学性质。它是自然界中最坚硬的物质，是 一种良好的耐磨材料，被广泛的应用于钻探和切割行业。随着以c 6 0 为代表的一 大批富勒烯的发现和其独特的物理化学性质【1 1 ，再一次引起人们对碳元素研究 的高潮。近来，对碳纳米管的力学、电子学等特性的研究成为现在纳米研究领域 的热点。 1 1c 6 0 和小富勒烯的实验发现及其特性 在1 9 8 5 年9 月，美国r i c e 大学和英国s u s s e x 大学的研究人员在s m a l l e y 和 k r o t o 两位教授的领导下，在r i c e 大学用镭射蒸发石墨的气相生长法发现了c 6 0 分子，并命名为布基球或富勒烯( f u l l e r e n e ) ，图1 1 a 。布基球有六十个碳原子对 称分布于截角二十面体的顶角处，共有1 2 个五边形面和2 0 个六边形面组成，具 有i h 对称性的足球结构，它具有出乎寻常的稳定性( 结合能为7 0 1 e v a t o m ) 。富 勒烯的发现是科学发展史上一个激动人心的事件。 1 9 9 8 年，加里福利亚大学的z e t t l 研究小组使用弧光放电的方法制各了固态 形式的c 3 6 2 】。发现它与c 6 0 完全不同。如c 6 0 溶于甲苯，但c 3 6 不溶；c 6 0 膜较 软，c 3 6 膜较硬等。他们采用净化的方法把c 3 6 同非晶碳和其他富勒烯成功地分 离开来。测得其直径为0 5 n m ，结构不是球状。固态核磁共振( n m r ) 测量显示 c 3 6 分子具有d 6 h 对称结构，如图1 1 b 所示。另外，在c 3 6 固体中内部分子问间 距为6 6 8a 当它掺杂碱金属后，发现其电导率会显著的增力1 1 3 4 1 。 2 0 0 0 年9 月，德国的h o r s tp r i n z b a c h 小组在d o d e n c a h e d r a n e 形式的c 2 0 h 2 0 中用b r 原子替代h 原子，然后在气相中将b r 原子解吸，从而在实验室产生了 最小的富勒烯。同时他们也用类似的方法产生了碗状的c 2 0 分子异构体【5 【6 】。发 第一章引言 现c 2 0 三种异构体的电子亲和性和振动特征彼此之间差别很大。图1 1 c 为c z o 富 勒烯的结构特征。 图1 r a ) 自然祭志的封面( a n n o u n c i n g t h ed i s c o v e r yo f f u l l e r e n e ) ，c 6 0 富勒烯的结构。 图1 1 ( b ) d 6 h 对称的c m 富勒烯结构图( 1 2 个五边形面和8 个六边形面，c 1 、c 2 和c 3 代表三 不同的原子) 图】1 ( c ) 具有c 2 ，对稼的c 2 0 富勒烯结构( 只有1 2 个五边形面) 最小的富勒烯由1 2 个五边形面组成。小富勒烯结构中五边形相连在起， 使其具有较高的曲率( c u r v a t u r e ) ，从而增加了应变能，因而呈现出较高的化学 反应活性，这预示小富勒烯具有独特的化学和物理性质【7 【8 】。 1 2碳纳米管的实验发现及其特性 1 9 9 1 年在日本n e c 实验室由s u m i ol i j i m a 通过电子显微镜发现了多壁碳纳 米管 9 ，1 9 9 3 年n e c 实验室i i j i m a 1 0 和i b m 研究中心【1 1 分别同时发现了单 壁碳纳米管。它们是由很多微小的石墨层卷成的圆柱体组成。从几个特殊的纳米 管结构看，可将其看成是把一种富勒烯切成两个半球然后拉伸得至l j 1 2 ，见图1 2 。 缀 c 5 舟| l r m c h n l r 礴燃 馔 ( 1 n 却内h i 图1 2 ：单壁磺纳米管a r m c h a i r 型、z i g z a g 型、c h i r a l 型结构示意图。 。礁固 第一章引言 因为碳纳米管和石墨在原子排列方式上的相似，一般可以将单壁碳纳米管的 结构设想为将石墨层沿某一个方向卷起来而得到，卷的方向可以用，= r i g + m b 来描述，其中口、b 为石墨平面格点基矢，n 、m 为表征纳米管结构的一对整数， 称为t ( n ，m ) 型的纳米管。当m 和n 相等时，纳米管的结构被称为a r m c h a i r 型； 当m 等于0 时，结构被称为z i g z a g 型( 见图l _ 3 ) ；而当m 与n 不相同时，则一 般性称为c h i r a l 型。由纳米管的c c 键长和矢量，的长度可以得到，纳米管的直 径为o 7 8 3 x n 2 + m + m 2a 1 3 。如果是多个单壁纳米管嵌套在一起，具有同心轴， 则称为多壁碳纳米管。 ( a ) g r a p h i t e ( b ) a r m c h a i rt ( 5 ，5 ) ( c ) z i g z a gt ( 1 0 ，0 ) 图l _ 3 ：单壁碳纳米管构成示意图。( a ) 为石墨结构；( b ) 为a r m c h a i r 型( 5 ，5 ) 单壁碳管 结构；( c ) 为z i g z a g 型( 1 0 ，0 ) 单壁碳管结构。 在实验上最先发现的碳纳米管是多壁碳纳米管，见图1 4 ，其中左上图是单 个多壁碳纳米管的高清晰透射电子显微镜( h r t e m ) 图像；左下图是单个单壁 碳纳米管束在无定型碳上的h r t e m 图像；右上图是一个多壁碳纳米管管端结构 的h r t e m 图像：右下图是一个多壁碳纳米管由于出现5 边形和7 边形碳环而导 致的集合变形【1 4 】。 由于碳纳米管的碳原子与石墨的碳原子s p 2 杂化轨道相似，又具有特殊的结 构，使得碳纳米管表现出异常的力学特性和导电、导热特性，并引起了众多科研 第一章引言 工作者的兴趣，例如用单壁碳纳米管做纳米器件和纳米电路，做人造肌肉纤维， n a s a 也利用纳米管做太空应用的材料，甚至还有公司用它做场发射显示器的纳 米器件等等 1 4 】。 图1 4 ：多壁碳纳米管的高清晰度t e m 图像 1 3纳米粒子的力学性能研究 1 3 1小富勒烯的力学性能研究 c 6 0 的n m r 实验以及光电子谱分析的结果预示了c 6 0 薄膜将保持c 6 0 分子 本身具有的特性，在团簇组装纳米薄膜的研究中也发现富勒烯小团簇在与固体表 面的碰撞过程中能够在剧烈形变之后恢复初始构型，这说明c 6 0 等小团簇有充 当很有效的工业用润滑剂的应用前景。 本文将呈现对小富勒烯力学性能的模拟研究，并对不同的富勒烯作了简单比 较，请见第三章。 第一章引言 1 3 2 碳纳米管的力学性能研究 在碳纳米管被发现之初，已有理论研究指出纳米管比一般纳米材料要坚固很 多，比如美国m i c h i g a n 大学的y o m a n e k 小组通过理论计算预测碳纳米管的杨氏 模量在1 5 0 0 5 0 0 0 0 p a 之间【1 5 】，而后来的理论研究得出的结论有9 7 0 g p a 1 6 】 和1 2 4 0 g p a 1 7 1 等等，n o r t hc a r o l i n as t a t e 大学的研究小组还研究了在剧烈形变 下的碳纳米管的力学性能 1 8 1 1 9 。实验上运用t e m 研究碳纳米管的杨氏模量在 4 1 0 4 1 5 0 g p a 的范围内，其平均值为1 8 0 0 g p a 1 5 。在对多壁碳纳米管压缩强度 的t e m 研究中发现其模量在1 0 0 1 5 0 g p a 之间 2 0 1 ，这已经要比任何已知的纤维 管的压缩强度要高1 0 0 倍。之后很多研究使用原予力显微镜( a f m ) 来测量碳 纳米管的伸缩特性，发现纳米管在外力作用下出现扭曲，并且在一定应变范围内 这种扭曲可以恢复至初始管状结构。另外，随着填充了c 6 0 的纳米管在实验上 的发现【2 l 】，填充小团簇和气体分子原子的纳米管力学性能也在最近得到了很多 关注。 本文将集中讨论单壁纳米空管和填充小富勒烯的纳米管的力学性能，请见第 四童。 第二章计算机模拟方法 第二章计算机模拟方法 本世纪以来，人们逐渐认识到，物质的各种物理性质及其性能是与它们的 微观结构密切相关的。由于真实体系中粒子间相互作用十分复杂，使传统的分析 方法对大多数问题都不能直接获得解析解。计算机的诞生为人们对复杂问题的研 究提供了新的方法和工具。人们针对具体问题设计一个合理的模型，由计算机完 成大部分的计算工作，从而获得复杂体系的微观信息，并从中分析得到体系的宏 观物理性质与行为，这就是计算机模拟方法的基本思路。 新近发展的计算机模拟方法介于传统的实验和理论分析这两种方法之间， 起着双重角色的作用。图2 1 1 给出了这三者之间的关系。 一方面，计算机模拟和解析理论都依赖于某个基本的理论模型。解析理论 通常要在这个基本模型的基础上再做出更多的假设和近似，然后才能进行理论推 导和预言。而计算机模拟则是对这个基本的理论模型直接进行模拟与数值计算， 从而获得这个基本模型的准确解。所以，只要把计算机模拟的结果与解析理论预 言的结果进行比较，就可以判断解析近似的正确与否。从这个意义上说，计算机 模拟可以起到检验理论的作用。而且在计算机模拟中，可以包含一些较复杂却又 比较符合实际情况的因素，从这一意义上说，计算机模拟具有比解析理论方法更 加优越的特点。 另一方面，计算机模拟的结果可以和实验相比较，从而有助于解释新的实验 现象和为实验中所遇到的疑难问题提供有价值的启示。另外，当物理实验尚不能 或难以进行时，计算机实验可以替代部分实验，或优化实验参数，从而得到实验 上所不能得到的结果。当然由计算机模拟得到的结果最终应该接受实验的验证。 最早的计算机模拟实验【2 2 】奠定了现在称之为蒙特卡罗( 简称m c ) 模拟的基 础，它是一种随机方法。在m c 模拟方法产生后不久，又诞生了另一种确定性的模 拟方法，这就是分子动力学( m o l e c u l a r d y n a m i c s ，简称m d ) 模拟方法。 1 9 5 7 年和1 9 5 9 年，a l d e r 和w a i n w r i g h t 2 3 ，【2 4 】通过对原子经典运动方程( 也 就是牛顿运动方程) 求数值解，首次实现分子动力学方法对硬球相互作用体系的 模拟。1 9 6 4 年，r a h m a n 2 5 把m d 方法应用至u l e r m a r d j o n e s ( i 互作用体系中，并 在模拟中采用在每一微小时间步长内计算每个粒子所受到的作用力的方法，称为 s t e pb ys t e p 6 第二章计算机模拟方法 r e a l s y s t e m p e r f o r m e x p e r i m e n t s e x p e r i m e n t a l r e s u l t s m a k e m o d e l s t e s t s o f m 0 0 e l s c a r r yo u t c o m p u t e r s | m u l a t i o n s e x a c t r e s u l t sf o r 。m o d e l m o d e l s y s t e m t e s t s o f t h e o r le s 图2 1 1 计算机模拟与试验和解析理论方法之间的关系 c o n s t r u c t a p p r o x i m a t e t h e o r l e s t h e o r e t i c a l p r e d l c t i o n s 的程序步骤，以研究体系的时间演化。从而使得分子动力学方法对实际体系的应 用成为可能。现在的分子动力学模拟基本上都采用r a h m a n 的s t e p 。b y s t e p 模拟步 骤。 早期的计算机模拟所依据的模型都是非常简单的，模拟的内容主要集中在 一些很简单的单原子流体或固体整体性质研究。近三十年来，随着微电子技术的 迅速发展，计算机的功能越来越强，速度也越来越快，这就为计算机模拟的发展 提供了非常有利的条件。现在计算机模拟所针对的体系越来越复杂，模拟应用的 范围也越来越广。固体表面结构、液体结构、离子注入、辐照损伤、异质催化、 酶的行为、薄膜生长等，都可以通过计算机模拟进行研究。 第二章计算机模拟方法 2 1计算机模拟方法简介 近年来，计算方法的发展以及计算速度快而价格低廉的计算机的相继出现， 使得第一性原理的计算( a bi n i t i o ) 成为研究表面现象的可行并且可靠的工具【2 6 】 这些计算不仅能给出有关结构方面的信息，同时还能给出有关势能面的信息。但 是在一定程度上要受到体系大小的限制。 分子动力学方法尤其适合于研究薄膜生长的动力学过程。分子动力学方法 既可以提供扩散进行中所涉及的过程的定性知识，也可以给出不同机制的扩散系 数、扩散势垒的定量信息。早在1 9 8 5 年，s c h n e i d e r 2 7 等人就首次运用分子动力 学方法研究了l e n n a r d j o n e s 体系的f c c 点阵的生长。随着各种各样的多体势的出 现 2 8 ，这类模拟很快被运用于研究s i ，c u 以及其它金属的扩散 2 9 1 1 3 0 、分子束 外延和一些薄膜生长过程【3 l 】但分子动力学的主要缺点是模拟的时间标度和长 度标度不可能延伸到介观或宏观范围。这主要是因为m d 模拟的时间步长不能比 声子频率的倒数长，一般只能取原子振动周期的几十分之一甚至几百分之一。近 来，v o t e r 和g o n g 等人发展了一种超分子动力学的方法，可以大大地延长m d 模 拟的时间标度。 本文中研究富勒烯与填充富勒烯的纳米管的力学性能，采用的是分子动力 学与经验势相结合的方法，下面是对分子动力学的介绍。 2 2 分子动力学m d ( m o l e c u l a rd y n a m i c s ) 方法模拟 m d 模拟的开创性工作是1 9 6 0 年j b g i b s o n 等人【3 2 关于低能辐照损伤的研 究。虽然四十多年来，m d 模拟随着一代又一代新计算机的出现而有了长足的发 展，但其基本原理并未改变。m d 模拟已经用于固体中原子碰撞中的众多问题， 如缺陷的形成和迁移、离位阈能的计算、移位级联的空间构型、级联发展过程、 溅射机理、高能量密度级联以及薄膜生长、表面重构、界面混合等计算机模拟的 新领域。随着原子间相互作用势的发展，m d 模拟已被广泛应用到对c ，s i ，g a ，a s 等重要的半导体材料的模拟上。它已成为材料科学研究领域里对微观结构和微观 过程进行研究的一个重要方法。 固体中原子碰撞现象的m d 模拟程序一般都是在下述近似基础上进行设计 的： 第二章计算机模拟方法 体系完全可以由经典力学处理：即粒子的运动服从经典运动方程，不考 虑相对论效应和量子效应。在低能情况下，经典近似是成立的。 粒子间的相互作用可用经验势函数近似描述。 所有的势函数都和粒子的动能无关。 势能函数的形式不因粒子在碰撞过程中可能发生激发或电离而改变。 2 2 1 分子动力学方法原理 我们知道体系的微观状态可以由组成体系的n 个粒子的坐标和动量完全确 定。在原子层次的模拟中，原子为组成体系的基本粒子，体系的h a m i l t o n i a n 量为 原子坐标和动量的函数。如果认为体系中粒子的运动服从经典运动规律，那么体 系的h a m i l t o n i a n 量就可以进一步表示为所有原子的动能k 和原子间相互作用的 总势能之和，电子对原子运动的平均效应在势能中反映。如果用 i ”= ( i ，露) ，f ”= ( i ，砭，i ) ，p ”= ( p ，p ，p ) 分别表示体系中n 个粒子的坐标、速度和动量，那么，体系的h a m i l t o n i a n 量为： 妒，p ) - 麒p ) + o ( r ) i 蛳协善等f 2 2 1 】坝r ”，”) = 麒”) +一) ， 智2 珥 f 2 1 1 体系的运动方程可以有多种表达形式。最基本的一种表达方式是拉格朗f | 运动方程： 要晏一鲁：0 讲西- 谚 f 2 2 2 1 其中l ( i “，i “) 是拉格朗日函数，它和体系的h 锄i l t o n i a l l 量的关系是： ( i ”，i ”) = e 曩一h ( r ”，卢“ 一 ( 2 2 ，3 ) 把( 2 4 3 ) 式代入( 2 4 2 ) 式，就可以得到体系的h a m i l t o n i a n 正则方程： 巫：垡：盈 d t 泵码 鲁2 署_ _ v 删 口：川 v ，。善+ 善歹+ 署e 。由于盟d t 一，所以弋。( f 一) 是体系中所有其他粒 础砂， 噶。由于“，所以一v t o ( ，”) 是体系中所有其他粒 子对i 粒子的总作用力。 m d 方法的核心是数值解正则运动方程( 2 2 4 ) ，求出体系的相空间轨道 i “( ，) ，多”( f ) 或i ”( f ) ，i ”( n 。显然方程( 2 2 4 ) 能否求解，求得的解是否正确，关 键在于能否给出正确的粒子间相互作用势函数如( f ”) 。 9 第二章计算机模拟方法 2 2 2 数值积分方法 m d 模拟采用的势函数都是连续势。在连续势作用下，粒子各个时刻都受到 力的作用，运动轨迹一般不再是直线。系统运动方程组( 2 4 4 ) 的数值积分可以有 多种方法。算法选择的关键在于差分方程能否稳定收敛，计算是否准确可靠，在 此基础上再看算法是否简单。分子动力学常用的算法有中心差分法( c e n t r a l d i f f e r e n c e s c h e m e ) 3 2 】、平均力法( a v e r a g ef o r c em e t h o d ) 3 3 】、预测修正法 ( p r e d i c t o r - c o r r e c t o rs c h e m e ) j 3 4 1 、v e r l e t 算法( t h e v e r l e ts c h e m e ) 3 5 等。我们在计 算中使用的是n o r d s i e c k g e a r 预测校正算法d 6 ，因为这一算法在每一步的时间步 长很小的情况下会提高计算精度，下面对5 值n o r d s i e c k g e a r 算法作简要介绍。 假设r 代表了粒子坐标集合，对r 做泰勒展开 ，( t + 8 0 = 心) + 国| v ( 岫2 1 - 6 f a ( f ) + 吉非琊) + 去矿c + - ( 2 2 5 ) 定义= e c ( & o a t ) ，乞= 三2 ( d 2 疵2 ) ，吩= 吉矿( d 3 r o l d t3 ) ，_ = 丢盘4 ( d 4 ，衍4 ) n z , 就得到： 学( f + & ) r ( t 十西) 哆( f 十毋) 口( f + b e ) 学( h 一舀) 11 o1 oo 0 o o o ll 23 l3 01 oo ( f ) ( f ) ( f ) ( f ) ( f ) 这就是第一步预测。在预测了下一步的各项值以后，采用校正算法 4 ( t + s t ) 彳( f + 国) 巧o + 函) 学( f + 国) 4 ( t 十西、 + a r = 弓一哆 在校正步骤中，首先根 然后代入( 2 2 6 ) 。 据曙得到，再计算a r = 巧一譬 1 0 ( 2 2 6 ) ( 2 2 7 ) 坛哇吒 ，f_t_-_j_，l-_f_f-_1_，ll n叫川到 ， lll险if 鄢白。如乱 留甲譬擎 第二章计算机模拟方法 根据g e a r 的计算，得出了校正算法中的各项因子，见表2 2 1 。 l c 。c 1c 2c 3 c 4 l1 9 1 2 03 4l1 21 9 9 0 2 2 3 边界条件 材料的宏观性质是大量粒子的统计行为。模拟体系的粒子必须非常庞大，才 能准确再现真实体系的行为。另一方面，受计算机内存及容量的限制，m d 模拟 能够处理的粒子数是有限的。m d 模拟通过在无限大的体系中取出一块小晶体作 为模拟体系，模拟体系的边界原子要受到外界原子的作用力，在平衡构型中，此 力正好抵消晶体内部对边界原子的作用。而且边界原子在受到扰动时，它要向外 界传递能量，甚至离位。这就要求对模拟体系的边界提出合适的计算模型。一个 最简单的方法是增加体系中的原子的个数，随着计算机的飞速发展，这是个趋 势。目前，根据具体的研究问题不同，模拟体系的边界通常有以下几种处理方法： 1 周期边界条件 通常情况下，为了解决真实体系粒子数的庞大和计算机处理能力有限之间的 矛盾，模拟中一般采用周期性边界条件。 一般取一立方体为模拟体系，体系中含有几百一几万个粒子。如果不模拟表 面，就使中心原胞在三维空间上重复排列，如要模拟表面，就使模拟立方体在二 维水平面上重复排列。于是体系中的像粒子将在三维和二维空间周期性的重复出 现，称为赝粒子。则i 粒子的像坐标为： ( x i + l l x ，y i + m l y ，z i + n l z ) l , m ，n = 0 + 1 ，2 ， 其中k ，l ，l z 为立方中心原胞三边边长。体系的势能包括粒子与像粒子之间 的相互作用在内。 应用周期边界条件后，整个体系的赝粒子数无限多。当中心原胞中的一个粒 子移动时，其像粒子也会发生相应地移动。这样，若一个粒子离开了中心原胞， 必定有它的一个像粒子从另一面进入中心原胞。因而，中心原胞中的粒子数总是 保持不变的。 应用周期性的边界条件隐含着一个要求，即相距为l ( l = m i n ( l 。，l 。l z ) ) 的 两粒子的空间相关性是可以忽略的，也就是说，粒子间的相关长度应远小于模拟 计算的立方体的边长。这就要求模拟体系的粒子数要足够大。 第二章计算机模拟方法 2 ，固定边界条件 对于正则系综( n v t ) 的分子动力学模拟，有时采用固定边界条件，比如对 于粒子与固体表面相互作用体系，在模拟中将靶底层或者下面几层原子固定，以 模拟衬底对模型靶的相互作用。 3 自由边界条件 如果模拟原胞包含的粒子对于所研究的问题足够大，或者为了减小周期性边 界条件或固定边界条件对所研究问题的影响( 冲击波的反射等) ，m d 模拟中也 常采用自由边界条件，即对边界原子不加任何限制。 2 2 4时间步长( a t ) 及其它 我们从前面介绍的体系运动方程组的数值积分可以看出，从计算精度考虑， 时间步长t 的越小越好。但计算效率决定t 不能太小。t 的选择要协调考虑这两 方面的要求。我们在m d 模拟中采用变步长的方法。该方法先设置一个比较小的 时间步长，在计算过程中每隔一定的时间步自动调整步长。调整的原则是该时刻 体系中运动最快的粒子在一个时间步内所走的路程不超过晶格长度a 0 的2 、即 a t = 0 0 2a d v m ，v m 为该最快粒子的速度。上面这个百分比还可调整，尽可能达到 在满足计算精度的要求( 如体系总能量误差不超过0 0 2 ) 的条件下得到最大的 计算效率。典型的积分步长从0 1 f s 1 0 f s ，一个完整的模拟计算通常需要1 0 3 10 4 时间步。 模拟计算中最费时间的是计算每个时刻每个粒子受到的所有其它粒子对它的 作用力。为了提高计算效率，除了调整时间步长，我们在计算中还采用t v e r l e t 3 5 1 提出的在计算过程中每隔n 个时间步长更新一次近邻原子表的方法。即对体系中 的任一个粒子i ，首先对以i 为中心，r m 为半径的圆球内所有其它粒子的编号建立 一个近邻原子表。距离r 。取值要比粒子间的相互作用势程r 口( r d 为m d 模拟中所选 定的相互作用势函数的截止距离) 大，保证没有包含在表内的粒子在下面的1 1 个时 间步内不会进入i 粒子的相互作用势程内。通常r 。的取值可由下式给出 r m r 。mxojcm一；cmo cl一叮一oj- 第三章受压缩小富勒烯的纳米力学性能的模拟研究 b o n dl e n f l t h ( 0 1n m ) 图3 4 ：c 6 0 团簇在压缩过程中的键长分布。d 1 ，d 2 ，d 3 ，i ) 4 分别是在l l ， l 2 l 3 ，h 的刚性距离时得到的f 见图3 1 、。 势能也是变化很小。比如l ，从1 3 4 8 a 到1 1 a ( 压缩终点) ，c 6 0 体系总势能变化 大约为2 0 0 e v ，而仅有1 0 “是来f l c 6 0 内部的势能变化上的。这显示了在这一阶 段中的团簇已达到变形饱和点，很难在被压缩。 图3 4 是c 6 0 团簇键长的分布图，这是固定住刚性层原子使其驰豫1 0 0 0 个时间 步长，统计其键长在不同长度的分布得到的。分布曲线d l 和d 2 分别是在刚性层 距离为l 1 和l 2 的得到的，图中可以看键长峰值出现在1 4 1a 1 4 2 a 和1 4 4a 1 4 5 a 处。d 2 和d 3 分别是在c 6 0 剧烈形变时得到，在这一时刻的分布中，会发现 在1 6 a 和1 7 a 之间出现一个新的峰值，这显示出c 6 0 的一些键很明显变长。这些 co；ods|p口niiee古z 第三章受压缩小富勒烯的纳米力学性能的模拟研究 妄 一 吾 墨 商i 耋 焉 芷 1 01 21 41 61 82 02 2 2 4 d i s t a n c eo fr i g i dl a y e r s ( a ) 图3 5 ：团簇形变因子随压缩距离的变化 键就是图1 中加深的黑线所标出的，可以看出，它们环绕富勒烯形成了一个封闭 的圈。 为了描述团簇形变，我们采用了形变因子描述，见图3 5 。形变因子为z 方 向的最大半径与x y 平面方向的最大半径之比。从图中可以看出，在压缩开始之 初，各团簇的形变因子基本维持在一个值附近，这个阶段团簇的变形很小。随着 压缩的进行，形变因子开始减小，并且在中间某一时刻出现大幅下降，这预示着 团簇开始出现较大的形变，是一个扁平状的结构。然后在压缩后期几乎维持在 个值附近，说明这个阶段的团簇形变已经达到饱和。从图3 5 还能看出，越大的 团簇最后达到的形变因子相对于初始时的越小，说明越大的团簇形变程度越严 重。 为了能够比较不同富勒烯的特性，图3 6 显示了应变能与轴向应变的关系， 其中每一条曲线在小应变范围中都可以近似拟和为一个抛物线。所以，我们利用 下面的公式计算出在抛物线的二次导数作为体系的近似杨氏模量( 见表3 2 ) 。其 中e 是杨氏模量，w 表示应变能，表示应变，确表示系统初始体积。 。 la 2 i , v 胙百万 第三章受压缩小富勒烯的纳米力学性能的模拟研究 岁 3 x p m 山 g 巴 历 s t r a i n 图3 6 ：系统应变能与小应变( 0 0 1 到0 1 ) 的关系。 通过对c 6 0 薄表面声波的激光分析得到的。 表3 ，2 ：从抛物线拟和计算得到的杨氏模量 富勒烯 c 2 0c 3 6c 6 0c 7 0 模量( g p a ) a t3 0 0 k 7 61 0 9