（机械工程专业论文）单晶al弯曲变形的分子动力学模拟.pdf

江苏大学工程硕士学位论文 摘要 纳米结构和纳米材料代表着本世纪材料的发展方向，对其性能的研究具有重 大的理论和现实意义，本文利用分子动力学模拟技术对若干纳米结构进行了数值 模拟。 本文综述了纳米尺度下的材料结构和力学行为以及分子动力学数值模拟领 域己经取得的研究成果：然后较为全面的阐述了分子动力学方法的基本内容和实 施细节。重点研究了原子间势函数的选取、固体中常用的势函数以及运动方程算 法的确定等方面的内容。 在此基础上，文章选择f c c 金属中的a 1 作为研究对象，建立了单晶a l 弯 曲变形时的分子动力学模型，目标温度设为2 9 3 k ，仿真原子数3 9 7 个，采用二 维仿真模拟。模型采用m o r s e 势函数描述原子间的相互作用，并以牛顿方程建立 i 力学运动方程，使用v e r l e t 算法求解原子运动轨迹。在求解过程中，为了减少由 于仿真原子数远小于实际情况中原子数所引起的“有限尺寸效应”，对模型底层 边缘施加“刚性边界条件”。为减少计算时间，采用了“截断势”法。文章设计 了二种弯曲模型，包括理想无缺陷情况和表面有裂纹二种情况，通过比较来分析 表面缺陷对材料力学性能的影响。 最后，分析解释了两种模型在平衡阶段、裂纹萌生和扩展阶段的力学现象， 并对单晶a l 这种材料的裂纹萌生和扩展过程进行了描述。通过对单晶a l 的分子 动力学模拟，发现单晶a l 材料的裂纹萌生于空位，由空位的合并形成纳米级裂 纹，随后的微观裂纹的扩展类似于宏观裂纹。而表面裂纹缺陷对裂纹萌生和扩展 过程起到促进作用。 关键词：单晶a l ；弯曲；分子动力学；模拟 江苏大学工程硕士学位论文 a b s t r a c t t h en a n o s t r u c t u r e dm e t e r i a li so n eo ft h em o s tp r o s p e c t i v er e s e a r c hf i e l d s i nt h e 2 1 c n t h e r e f o r et os t u d yt h e i rm e c h a n i c a lb e h a v i o r si sv e r ym e a n i n g f u lb o t hf o r mt h e a c a d e m i ca n da p p l i c a t i o np o i n to fv i e w i nt h i st h e s i st h en u m e r i c a ls i m u l a t i o n so n t h em e c h a n i c a lb e h a v i o r so fs e v e r a ln a n o e l e m e n t sa r ep e r f o r m e db yt h em o l e c u l a r d y n a m i c sm e t h o d f i r s tt h ep r e s e n ts i t u a t i o na n dd e v e l o p m e n t so nn a n o t e c h n o l o g y , n a n o m e c h a n i c s a n da t o m i s t i cs i m u l a t i o na r e b r i e f l yv i e w e d ，t h e n t h eb a s i c t h e o r y a n d a p p l i c a t i o n d e t a i l so fm o l e c u l a rd y n a m i c sm e t h o di si n t r o d u c e d i np a r t i c u l a r t h e i m p l e m e n t a t i o no ft h em o l e c u l a rd y n a m i c ss i m u l a t i o no nt h ec a n o n i c a le n s e m b l ei s d e f i n e d ，b a s e do nt h ee x t e n d e ds y s t e mm e t h o do fn o s r s h o o v e rc o n s t a n t t e m p e r a t u r ec o n t r 0 1 t h eb a s i cc o n c e p ta n dt h e o r ya b o u tm o l e c u l a rd y n a m i c sa r e i n t r o d u c e d t h ee m p h a s e sa r ep u to ns o m ea s p e c ta sf o l l o w s ：t h ee f f e c to fs e l e c t i o no f p o t e n t i a lf u n c t i o na n dt y p i c a lp o t e n t i a lf u n c t i o ni ns o l i d ，t h es e l e c t i o no fa l g o r i t h m b a s e do nt h e s eb a s i ct h e o r i e so fm o l e c u l a rd y n a m i c ss i m u l a t i o n ，t h ep a p e rs e l e c t a l u m i n u mo ff c cm e t a la st h es t u d ys u b j e c t ，t h eb e n d i n gm o d e lo fs i n g l ec r y s t a l a l u m i n u mi sb u i l t t h es i m u l a t i o na d o p tt w o d i m e n s i o n a lm o d e l ，a i m t e m p e r a t u r ei s s e t u p o n19 3 k a n dt h ea m o u n to fa t o m i ci s3 9 7 i nb o t h m o d e l s ，t h en e w t o n s e q u a t i o n so fm o t i o na r es o l v e du t i l i z i n gt h em o r s ep o t e n t i a lf u n c t i o n t h ev e r l e t a l g o r i t h mi se m p l o y e dt o s o l v ea t o mt r a j e c t o r i e s i nt h e s i m u l a t i o n ，t h e ”r i g i d b o u n d a r yc o n d i t i o n ”o ft h em o d e lf i r s tf l o o re d g ei su s e dt or e d u c e ”l i m i td i m e n s i o n e f f e c t ”，”c u t t i n g - o f fp o t e n t i a l ”i su s e dt or e d u c es i m u l a t i o nt i m e t h ep a p e rd e s i g nt w o b e n d i n gm o d e li n c l u d i n gi d e a lc r y s t a lm o d e lw i t h o u tg r o o v ea n di d e a lc r y s t a lm o d e l w i t hg r o o v e ，i no r d e rt of i n dt h ei n f l u e n c eo fd e f e c tt ot h em e c h a n i c a lc a p a b i l i t yb y c o m p a r e d f i n a l l y , i na n a l y z i n gs i m u l a t i o nr e s u l t s ，s o m ep h e n o m e n o ni nt h es t a g eo f e q u i l i b r i u ma n ds o m ep h e n o m e n o ni nt h es t a g eo fc r a c ki n i t i a t i o na n dg r o w t ha r e i l l u s t r a t e d ，a n dt h i sk i n do fm a t e r i a l sa r es u m m e r r i z e d b yt h es i m u l a t i o no fs i n g l e i i 江苏大学工程硕士学位论文 c r y s t a la l u m i n u m ，t h ep a p e rf o u n dc r a c ki n i t i a t i o nw a sf o r mw i t hv a c a n tl a t t i c es i t e ， a n dc r a c kg r o w t hr e s e m b l et ot h em i c r o s c o p i c a lc r a c k t h ec r a c ko nt h es u r f a c e a d v a n c e dt h es t a g eo fc r a c ki n i t i a t i o na n dg r o w t h k e yw o r d s ：s i n g l ec r y s t a la l u m i n u m b e n dm o l e c u l a rd y n a m i c ss i m u l a t i o n 1 1 1 江苏大学工程硕士学位论文 1 1 纳米科技和纳米材料 第一章绪论 纳米科技( n a n os c i e n c e & t e c h n o l o g y ) 是指在纳米尺度( 0 1 n m 到l o o n m 之间， l n m 等于1 0 g m ，大约相当于3 - 4 个原子直径的长度) 上研究物质的特性和相互作 用，以及利用这些特性的多学科交叉的科学和技术。纳米科技的最终目标是基于 原子、分子层次研究，进行产品设计和制造。纳米科技是8 0 年代末、9 0 年代初 才逐步发展起来的前沿、交叉性新型学科领域，它将在2 1 世纪促使几乎所有工 业领域发生革命性的变化。目前，所有发达国家都在对纳米科技的研发进行大量 投入。 纳米器件的最初构想产生于1 9 5 9 年，在美国物理协会召开的年会上，n o b e l 奖获得者r i c h a r dpf e y n m a n 阐述了基于原子、分子层次制造材料和器件的设想 ” ，将科学家的研究目光从宏观范围引向了微观世界，从此微尺度下材料和结构 的研究成为了新的热点。 从1 9 8 1 年d r e x l e r 发表第一篇关于分子纳米技术的论文开始 2 1 ，在纳米结构、 分子微机械等方面的研究取得了很大的进展。尤其是自组装体系【3 、单电子晶 体管7 8 1 等方面的进步让人们看到了微纳米机电系统研制和应用的广阔前景。 1 9 8 2 年b i n n i g 等【9 】首先成功研制出扫描隧道显微镜( s c a n n i n gt u n n e l i n g m i c r o s c o p y ，s t m ) ，为科学家们在原子尺度上研究物质结构开辟了新途径。1 9 8 4 年德国萨尔兰大学( u n i v e r s i t yo f s a a r l a n d e s ) g l e i t e r 掣1 0 1 首次采用惰性气体凝聚 法制备了具有清洁表面的纳米粒子，然后在真空中原位加压成纳米固体，并提出 纳米材料界面结构模型。1 9 9 1 年美国i b m 公司的e i g l e r 等【l i 】利用扫描隧道显微 镜直接操作原子，成功地在镍板上按自己的 意志安排原子，组合成“i b m ”字样。1 9 9 4 年中国科学院真空物理实验室庞世谨等在 s i ( 1 1 1 ) 一7 * 7 表面用s t m 通过移走s i 原子 写成的中国，为此新华社曾发表了“搬动原 子写中国”的报导。 图1 1 用4 8 个f e 原子组成的量子围栏 江苏大学工程硕士学位论文 1 9 9 3 年i b m 公司用原子搬迁法在单晶铜表面用4 8 个原子围成直径1 4 3 纳 米的圆圈，相邻f e 原子间距为1 纳米，由f e 原子组成的围栏，可以看见美丽的电子 波浪，如图l 一1 所示。日本科学家己成功地将硅原子堆成一个“金字塔”，首次实 现了原子三维空间的立体搬迁。1 9 9 1 年i b m 的科学家还制造了速度为2 0 0 亿分 之一秒的x e 原子开关，专家们预计这一突破性的纳米科技研究成果将可能使美 国国会图书馆的全部藏书存储在一个直径仅为0 3 厘米的硅片上。 近年来，具有空间三维结构的纳米机电系统叫a n oe l e c t r om e c h a n i c a l s y s t e m ，n e m s ) 1 2 1 出现，纳米机械制造技术迅速发展。1 9 9 9 年k i m 等1 3 1 研制出 用碳纳米管制成的纳米镊子，成功地抓住了直径约5 0 0 纳米的聚苯乙烯原子团， 并能够抓取和操纵亚微米团簇和纳米线。2 0 0 0 年s o o n g 等【1 4 】制造了一台纳米机 器，通过把生物微型马达与一个金属螺旋桨组装在一起，组合成纳米机器，螺旋 桨能够以每秒8 转的速度旋转。图1 2 ，1 3 分别为美国加州理工( c i t ) 的科学家 c a g i n 用分子动力学模拟的纳米齿轮和纳米转子。 “纳米”是一个尺度的度量，日本在1 9 7 4 年最早把这个术语用到技术上， 但是以“纳米”来命名的材料出现在2 0 世纪8 0 年代，它作为一种材料的定义把 纳米颗粒限制到1 - 1 0 0 n m 范围。在纳米材料发展初期，纳米材料是指纳米颗粒 和由它们构成的纳米薄膜和固体。现在，广义的纳米材料是指在三维空间中至少 有一维处于纳米尺度范围或由它们作为基本单元构成的材料【l ”。按照这个定义， 图1 - 2 纳米齿轮图1 - 3 纳米转子 纳米固体材料、纳米薄膜、纳米丝、纳米杆、纳米管等都是纳米材料。 纳观材料的发展至今大体可以分为3 个阶段1 6 1 。第一阶段为1 9 9 0 年以前， 这一阶段主要是在实验室研究用不同手段制备各种材料的纳米粉体，探索其不同 江苏大学工程硕士学位论文 于常规材料的特殊性能及其原因。由于这一阶段的研究对象局限在单一材料和单 相材料，故一般把这类材料称为纳米晶( n a n o c r y s t a l l i n e ) 或纳米相( n a n o p h a s e ) 材 料。第二阶段为1 9 9 0 1 9 9 4 年，这一阶段的热点是如何利用已发现的纳米材料的 特殊的性能来设计纳米复合材料，即掀起了对纳米材料的复合及其机制研究的热 潮。通常的复合方式有：纳米微粒之间的复合、纳米微粒与常规块体复合以及纳 米复合薄膜。第三阶段从1 9 9 4 年直至现在，这一阶段人们越来越多的关注纳米 尺度的图案材料( p a a e m i n g m a t e r i a l so n t h e m a n o m e t e r s c a l e ) 体系，即展开了对纳 米自组装体系( n a n o s t m c t u r e da s s e m b l i n gs y s t e m ) 的研究。它的主要目的是以纳 米颗粒以及纳米丝、纳米管为基本单元，在一维、二维和三维空间组装排列成具 有纳米结构的体系。 早期对纳米材料的实验研究主要采用x 射线衍射、穆斯堡尔谱、正电子湮灭 等手段 1 7 1 8 、19 1 ，实验发现纳米材料主要由晶粒和晶界组成，晶粒粒径为纳米量 级，晶粒内部由完美的点阵结构组成，粒径尺寸一般呈对数正态分布，晶界结构 在材料中占相当大比例，原子混乱无序，能量较高，密度以及配位数比理想晶体 低，具有类气态特征。 近年来应用高分辨率电子显微镜( h r e m ) 、扫描隧道显微镜( s t m ) 、原子 力显微镜( a f m ) 等手段可以直接观察纳米材料的微观结构，在证实一些早期结论 的同时，也发现了一些新的实验结果。现在知道纳米晶粒结构与理想完整晶体有 很大差别，即存在点阵应变、晶粒畸变、晶格常数变化等等。t h o m a s i 2 1 等对 纳米p d 的界面结构进行了高分辨率电镜的观察，没有发现界面存在扩展的无序 结构，原子排列的有序程度很高，和常规粗晶材料的界面并没有明显的差别。 w u n d e r l i c h t 2 2 1 等对纳米p d 样品的高分辨率电镜观察，发现晶粒中存在着孪晶， 在晶粒与晶界的交接处存在着位错亚结构。2 0 0 0 年2 月s c i e n c e 杂志发表了中国 科学院沈阳金属所卢柯( kl u ) 教授课题组在纳米材料超塑性方向的研究成果，他 们成功地制备出能够压延伸长5 0 倍的纳米多晶铜块体【2 3 1 。 纳米材料的微观结构主要包括晶粒、晶界、孔隙、杂质和点阵缺陷，纳米级 晶粒和大比例晶界是纳米材料的基本结构特点【1 8 19 1 。纳米金属材料因其特殊的 微观结构，而表现出许多奇异的力学行为和性能【2 4 2 副。其中最典型的是强度和粒 径的关系曲线，如图1 - 4 所示。当粒径大于临界尺寸时，材料强度和晶粒尺寸表 江苏大学工程硕士学位论文 现出典型的h a l l - p e t c h 关系，而当粒径小于临界尺寸并降到纳米量级时，材料强 度与粒径的关系表现出反 h a l l p e t c h 特征。 我们知道，几乎材料所有的特 性都具有相关特征长度。当材料尺 寸小于特征长度时，该特性就开始 发生变化。材料在纳米尺度范围内 时，由于量子效应、小尺寸效应、 表面效应、界面效应，物质的许多 性质将发生改变，呈现既不同于宏 观物体、也不同于单个孤立原子的 d cd , g r a i ns i z e 图1 - 4 材料强度与晶粒粒径关系曲线 奇异特性。纳米尺度下大量的新现象势必将物理、化学、力学等众多学科推向一 个新的研究领域层次，并加深、拓展对现有物理规律的理解。 1 2 纳米力学和分子动力学模拟 1 2 1 纳米力学 人们很早就知道材料的力学行为随尺寸发生变化f 2 8 l 。在微机电系统( m e m s ) 中，尺度效应和表面效应对器件的力学性能有显著的影响【2 9 】。毋庸置疑，当尺 度接近纳米量级时，经典连续介质力学失效而无法解释出现的新现象，必须从基 本的原子、分子运动过程来了解体系的力学行为【3 0 】。 纳米力学主要研究原子尺度下材料的结构力学特征以及其变形与运动过程， 是传统力学在纳米尺度领域的延伸。纳米科技的发展为纳米力学开辟了广阔的研 究领域，使纳米力学研究具有现实的理论意义和应用价值。纳米机械、纳米多晶 材料、纳米薄膜、纳米加工技术等都是纳米力学的研究对象。 扫描隧道显微镜( s t m ) 、原子力显微镜( a f m ) 、高分辨率电子显微镜( h r e m ) 等近来发展起来的实验手段使直接观察纳米尺度下原子的运动过程成为可能，这 极大地增加了人们对纳观力学现象的认识，大大推动了纳米力学的发展。如 k i z u k a 3 1 】用动态高分辨率显微镜原位直接观察了纳米金线的拉伸、压缩、剪切的 动态过程。 4 盖采5io露_h焉砖定 江苏大学工程硕士学位论文 随着计算机的发展，计算机模拟已独立于理论分析和实验研究的第三种研究 手段，成为沟通理论与实验的桥梁，计算机模拟技术可实现在实验上很难或根本 无法完成的研究。计算机模拟无需过于简化的假设，能够接近实际的复杂情况， 得到实验无法测量的结果，并深入揭示它们的内在行为机制；同时，计算机模拟 可以比较各种简化理论，评估模型的正确与否。随着计算机硬件、计算技术、模 拟模型等不断进步和发展，原子尺度模拟方法必将在科学研究中发挥愈来愈大的 作用。 1 2 2 分子动力学模拟 借助于计算机科学在原子尺度模拟纳米材料是研究纳米力学的一个有力工 具。原子层次模拟材料力学行为的基本步骤是：根据物理学原理建立原子间相互 作用的模型，针对具体的目标建立简化的原子模型，并根据实际问题的需要，选取 模拟的算法，计算并对模拟结果进行分析，找出可能的趋势或者线索【3 2 1 。目前发 展比较成熟的原子模拟方法主要有分子动力学方法【3 3 3 4 】和蒙特卡罗方法。 分子动力学( m o l e c u l a rd y n a m i c s ，m d ) 是一个重要的原子尺度计算机模拟手 段。分子动力学方法的基本思想是通过原子间的相互作用势，求出每一个原子所 受的力，在选定的时间步长、边界条件、初始位置和初始速度下，对有限数目的 原子( 分子) 建立其牛顿动力学方程组，用数值方法求解，得到这些粒子经典运动 轨迹和速度分布，然后对足够长时间的结果求统计平均，从而得到所需要的宏观 物理性质和力学性质。 1 9 5 7 年a l d e r 等1 3 5 】首次发表了分子动力学模拟文章，采用分子动力学成功 解决了硬球模型系统的固液相变问题。2 0 世纪八十年代，分子动力学方法得到 了飞速发展，方法本身日趋成熟。近年来，分子动力学广泛应用于凝聚态物理、 材料科学、核技术、化学反应动力学、生物化学、药物设计等领域。 目前，人们对单晶、位错缺陷、裂纹扩展、材料的脆韧转化、材料的表面界 面性能等许多现象做了大量工作。c h e u n g 等【3 6 】研究了a - f e 中裂尖的脆韧转化过 程。l y n d e r ”】模拟展示了单晶拉伸破坏的全过程。z h o u 等【3 8 】模拟了零温下晶体 铜( 2 0 0 0 0 0 ，到3 5 0 0 0 0 0 个原子) 中位错相互交截的过程。h o l i a n 等口9 】模拟了单晶体 铜( 1 0 8 个原子) 中冲击塑性的演化过程。张永伟等【4 0 】对裂纹尖端位错发射过程进 江苏大学工程硕士学位论文 行了模拟。 纳米多晶材料是纳米科技的一个重要研究方向。实验过程中，影响纳米多晶 材料实验结果的因素很多，并且许多关键性细节无法通过实验获得。因此，纳米 多晶材料的分子动力学模拟十分必要【4 ”。 采用分子动力学模拟纳米多晶材料的力学行为一般包括两个步骤：纳米多晶 材料初始构型的生成；纳米多晶材料的加载。纳米多晶材料的初始构型可以通过 模拟纳米颗粒粉末的烧结或压结过程获得，也可以通过纯几何方法直接生成。 s c h i o z e 等 4 2 1 模拟了粒径3 3 6 6 n m 纳米多晶铜的低温拉伸性能。文章中首 先采用v o r o n o i 方法直接构造纳米多晶铜的初始几何构型，然后对几何构型进行 弛豫，以获得具有实际物理特征的纳米多晶铜数值模型。此后，在对纳米多晶铜 拉伸模拟过程中发现一种新的变形机制：变形过程中，晶界存在大量相互无关的 只有几个到几十个原子同时运动的小滑移事件，这导致了晶界滑移。模拟结果表 明纳米多晶铜的塑性主要来源于晶界滑移，晶粒内的位错活动对塑性贡献极其微 小。新机制使纳米多晶铜在理想致密的情况下，出现强度随晶粒尺寸减小而发生 软化的现象( 即反h a l l - p e t c h 效应) 。长期以来，由于各种复杂因素的影响，实验 上无法确认是否存在反h a l l p e t c h 效应。s c h i o z e 的模拟结果证实了反h a l l p e t c h 效应是纳米多晶材料的本质特征。 分子动力学存在模拟时间短和模拟系统尺寸小的缺陷。目前并行计算机能够 模拟的粒子数目最多为千万量级，但模型尺寸远未达到微米量级：模拟时间最多 只有纳秒量级，但许多实际物理过程需要达到秒或更长的时间周期。另外，金属 材料本质上是多尺度的，包括空间和时间的多尺度，因此当前的分子动力学技术 无法真实地再现金属材料的全部力学行为。 计算机硬件的发展将允许模拟更大的尺寸、更长的时间和采用更精确的原子 势。但即使如此，在可以预见的将来，能够模拟的粒子数目和模拟的时间仍然是 十分有限的。仅仅依靠计算机的发展难以满足实际需求，因此模拟方法有待进一 步突破，当前出现的多尺度方法是一种有益的尝试。 多尺度方法的思路是将模拟对象按纳观、细观、宏观等不同尺度划分多个区 域，不同的区域采用不同效率、不同精度的算法，这样可以大大消除敷余计算自 由度，克服计算机速度的限制，有效地模拟更大体积材料和模拟更长的时间。目 6 江苏大学工程硕士学位论文 前已经发展了多种有限元与分子动力学嵌套的多尺度方法。一种方式是将模拟对 象简单地分为相互藕合的原子区域和有限元区域【4 3 “4 1 ，原子区域按分子动力学 计算，有限元区域则按常规的有限元方法计算。该方法需要解决原子区域和有限 元区域的边界匹配问题。 t a d m o r 等【4 5 “6 】提出了另一种先进的准连续( q u s c o n t i n i o u s ) 模型。其基本思路 为：将研究对象划分为不同区域，变形剧烈区按原子模型计算，变形均匀区域则 采用有限元分析，这样可以大大降低体系的自由度。计算过程中根据需要自动进 行区域划分，有限元区域的本构方程不再采用唯象的经验公式，而直接通过原子 间作用势计算得到。m i l l e r 掣4 7 1 基于这种准连续模型成功地分析了金属镍中裂纹 的扩展问题。相信在不久的将来，基于这种准连续模型思路的多尺度方法会得到 快速发展。 1 3 本文研究目的、内容和意义 纳米科学是新兴的研究领域，其发展速度日新月异。纳米结构和纳米材料的 尺寸非常小，所以必须从原子尺度来了解纳米结构和纳米材料的力学性能。由于 实验手段的限制，现阶段仅仅通过实验方法并不能彻底了解纳米材料的力学性能 和内在机制。因此，原子尺度的计算机模拟成为了实验研究的有效补充，分子动 力学作为原子尺度模拟的有力工具能够获得相关原子运动细节，有效的澄清实验 现象。综合实验观察和原子尺度的数值模拟的结果，才能真正的了解纳米材料的 微观机制，为最终实验纳米材料和纳米机械的设计、制造提供参考。 本文第二章详细介绍了分子动力学方法的基本思想。从系统粒子运动的基本 方程入手，介绍了势函数、周期性边界条件、动力学差分格式以及宏观物理量的 获取等若干基本概念和实施细节。然后介绍了微正则系统和正则系统分子动力学 方法的实现，并介绍了等温控制方法的相关细节问题。 第三章建立了二种单晶a l 弯曲变形的分子动力学模拟的模型，确定了势函 数，建立了系统运动方程，并且给出了模拟实现的初始化参数。模型包括理想无 缺陷情况和表面有裂纹二种情况，采用m o r s e 势函数描述原子间的相互作用，以 牛顿方程建立力学运动方程后，使用v e r l e t 算法求解原子运动轨迹。在求解过程 中，为了减少由于仿真原子数远小于实际情况中原子数所引起的“有限尺寸效 江苏大学工程硕士学位论文 应”，对模型底层边缘施加“刚性边界条件”。为减少计算时间，采用了“截断势” 法。 第四章进行了单晶a l 弯曲变形的分子动力学模拟实验。给出了实验结果， 并且通过绘制系统总能量随分子动力学步变化曲线，得到弯曲变形时系统总能量 随分子动力学时间步的演化规律。通过对单晶a l 进行了弯曲变形时的分子动力 学模拟，研究这种材料的弯曲变形时裂纹萌生和扩展的本质规律。结果表明材料 的裂纹萌生于空位，空位的合并形成纳米级裂纹，随后微观裂纹的扩展类似于宏 观裂纹。同时通过对无裂纹模型和有表面裂纹模型的裂纹萌生和扩展时间步长进 行对比，表明裂纹缺陷对裂纹萌生和扩展过程起到促进作用。 第五章对全文工作进行了总结。对采用理想无缺陷情况和表面有裂纹二种情 况模拟的结果做了综合说明，同时也分析了利用分子动力学模拟研究材料微观力 学行为的优点，指出了在模拟过程中存在的主要问题。文章最后对下一步要开展 的工作方向进行了展望。 江苏大学工程硕士学位论文 第二章分子动力学模拟方法概述 2 1 引言 从物质的内在微观结构及其原子尺度结构出发来预测材料的宏观性能，并按 预期目标设计材料是材料科学的终极目标。计算材料学( c o m p u t a t i o n a lm a t e r i a l s s c i e n c e ) 是实现这一目标的有效工具。计算材料学研究范围可按尺度划分为4 个 层次，即电子层次、原子层次、细观层次和宏观层次。针对材料不同的层次有不 同的模拟方法：电子层次采用第一原理( f i r s tp r i n c i p l e s ) 的量子力学方法；原子 和细观层次采用分子动力学、蒙特卡罗方法和格子气方法；宏观层次通常采用有 限元法、边界元法和有限差分法等。 根本上而言，材料问题都属于量子力学问题，似乎从第一原理出发就可以计 算得到任何需要的结果。但由于理论上的困难以及计算机条件的限制，目前并不 能够通过求解薛定谔方程获知材料的力学特性。量子计算的理论困难在于：局域 密度泛函近似计算并不总能满足实际问题的需要，理论上有待于新的突破。此外， 材料本质上是多层次、多尺度的，各结构层次具有自身特定的长度尺度、时间尺 度和能量尺度，层次间除了藕合外，还存在脱藕，这种层展现象( e m e r g e n t p h e n o m e n a ) 限制了量子计算的适用范围4 9 棚5 1 3 大量模拟结果表明，对于位错运动、表面界面现象等原子层次的力学现象而 言，采用经验势函数特别是采用镶嵌原子法等多体势的分子动力学可以获得较理 想的结果。分子动力学使我们能够深入地了解复杂的机制，发现本质上崭新的现 象，而且可定量地模拟真实固体中所发生的过程。分子动力学是对诸如表面结构 和扩散中的动力学和稳定性的许多见解的唯一来源，是研究固体、液体结构、表 面界面性能、弛豫过程等凝聚态现象的有力工具。 1 9 5 7 年a d l e r 等口”首次采用分子动力学成功解决了硬球模型系统的固液相 变问题，初步展现了分子动力学处理多体问题的强大能力。1 9 6 0 年g i b s o n 等 5 2 根据连续势函数研究了晶体的辐射损伤，这是计算机模拟第一次应用于材料科 学。7 0 年代，分子动力学开始逐步扩展到不同领域。8 0 年代由于计算机技术的 9 江苏大学工程硕士学位论文 快速发展，分子动力学研究非常活跃，出现了不同系综的模拟方法，镶嵌原子法、 有效介质理论等多体势也同时被提出。近年来，分子动力学业已广泛应用于凝聚 态物理、纳米力学、材料科学、核技术、化学反应动力学、生物化学、药物设计 等领域。 由于分子动力学存在模拟时间短和模拟系统尺寸小的缺陷。目前模拟的粒子 数目最多为千万量级，系统尺寸远未达到微米量级：模拟时间只有纳秒量级，无 法再现许多长时间周期的物理现象。但同时应注意到，纳米器件包含的原子数目 较少，基于现阶段的模拟理论和计算机速度，我们可以较好地模拟纳米器件的力 学性能。 本章从分子动力学基本原理出发，详细介绍分子动力学的理论基础和实施细 节，并给出微正则系统和正则系统的分子动力学方法。 2 2 分子动力学方法的基本概念 分子动力学方法是物理系统的确定性描述，它属于确定性方法。分子动力学 是通过原子间的相互作用势，求出每一个原子所受到的力，在具体的时间步长、 边界条件、初始位置和初始速度下，对由有限个粒子组成的粒子系统建立其牛顿 运动方程，用数值方法求解，得到这些粒子在相空间的运动轨迹，然后运用统计 力学方法求其微观量的系统平均，得到所需要的宏观量。 分子动力学有两个基本假设： 1 ) 所有粒子的运动都遵循经典牛顿运动定律： 2 ) 粒子间的相互作用满足叠加原理。 因此，分子动力学虽然是从原子层次研究问题，但是其本身忽略了量子效应 的影响，这里的原子点实际上属于宏观粒子点的范畴，只是粒子间的作用特性通 过势函数体现。因此本质上，可简单的将分子动力学看作是广义牛顿运动方程的 数值积分，仍然属于近似计算。 1 0 江苏大学工程硕士学位论文 2 2 1i g , 子运动基本方程 原子运动由经典牛顿运动方程描述： m o 等- f 口 ( 2 1 ) m 。为原子质量，r 。为原子空间位置，r 为时间，e 为原子间作用力， 通过确定原子间作用势e ( r ”) 确定e p：塑竺(2-2) p = 一 。 觎。 r ”表示所有原子核空间坐标，尺。表示第a 个原子核空间坐标。考虑到电子结 构，在确定原子间作用势时，需解多体( m u l t i b o d y ) 薛定谔( s c h r o d i n g e r ) 方程。 h q ( r , ，r 。) = e t o t a l g ( r i ，r 。) ( 2 3 ) e 。眦表示系流总能量，因核动能与核坐标无关，故核静止时系统基态能量e ( r n ) 可替代e 。再根据( 2 2 ) 式求得f a 。r 、r 分别表示核、电子的空间坐 标，a ，表示核序号，f ，表示电子序号，( ，r 。) 为系统波函数。日为系统 哈密顿算子( h a m i t o n l i a n ) 。 h = l4 - t4 - u 。4 - u 。 ( 2 4 ) 丁为动能算子，u 为势能算子，n ，e 分别表示原子核和电子，展开如下： 阱卫2 m a3 - 莩瓤否南一j 1 丢南+ 三否尚协s ， p ，p 分别为核和电子的动量算子， m ，m 为核和电子质量，z 为电荷 数。 原子间作用势的确定是分子动力学的关键问题，由于精确求解薛定谔方程的 困难性，近年来，通常要用实验拟合或半经验解法得到原子间作用势，根据对原 子间作用势的不同简化处理，分子动力学可划分为对势、多体势( m b m d ) 、可 变电荷分子动力学模拟( v c m d ) 紧束分子动力学模拟( t b m d ) ，密度泛函分 子动力学模拟( d f m d ) 等。 江苏大学工程硕士学位论文 2 2 2 原子间作用势 1 对势 早期的分子动力学主要用于模拟惰性气体和液体，发展了一些对势模型h 4 、 5 7 。其中应用非常普遍的是l e n n a r d j o n e s 势，最初是为了模拟液态氢的性质而 提出的，其势能曲线如图2 1 所示，其形式如下： 叫”制 p s ， 其中，e 、o 分别为能量参数和长度参数。 相应的参数为= 1 2 0 k 。k b ，o = o 3 4 nm 。玻尔 兹曼常数kb = 1 3 8x1 0 2 3 j 瓜。上式( 2 6 ) 中第 一项代表短程排斥力项，第二项代表远程吸引 力项。 虽然对势能够较好地模拟惰性气体的凝聚 态性能，但在模拟金属性能时出现了问题，如 不可避免地导致c l l = c 4 4 的不正确的c a u c h y 关 系。由于金属晶体中电子云发生相互作用，金 属键具有多体性质，只简单地考虑原子对的作 用必然无法正确计算金属系统的能量。为了克 服对势的困难，人们基于均匀电子气假设引入 对于液氩分子( 单原子分子) ， w ，i 11 盟口， 7厂 图2 - 1l j 势能 了体积项加以修正，但同时却导致图2 - 1 所示l j 势能曲线另外的一些困难5 5 1 。 2 嵌入原子势( e a m ) 上世纪八+ 年代，镶嵌原子法例1 ( e m b e d d e da t o mm e t h o d ，e a m ) ， f i n n i s s i n e l a l ( f s ) 方法【5 7 】、有效介质理论嘲3 等多体势相继提出。这些多体原子 势考虑了原子在局域背景电子云密度环境下的结合能，它们均将系统总能量分解 为统一的形式，特别对于金属原子来讲，嵌入原子势( e a m ) ( 考虑了原子核外 电子及其它原子产生的背景电子的静电作用) 能较好描述金属原子之间的相互作 用，是最常用的一种势函数。 1月nn e 。，= 去妒卢( 月筇) + f o ( p 。) ( 2 7 ) 江苏大学工程硕士学位论文 蝣为原子。和卢相互作用势能，r 妒= i 玛一r 口l 表示原子a 和原子卢之 间距离。 岛= f ( r o o ) ( 2 8 ) j l 岛表示第a 个原子处电子云密度，即周围原子的电子云密度在a 原子处的线 性叠加，f ( r 帮) 表示原子在距离r 邓处产生的电子云密度，镶嵌能函数f 表示将原予a 嵌入到密度为p 。的背景电子云中的嵌入能。函数中采用中心对势 形式，r 郇为原子0 、间距离，式中去表示结合能为两原子共有。镶嵌能函 数f 的形式与背景电子云密度来源无关。 目前基于e a m 等多体势框架： 1nn 巨。= 去妒f l ( r o f l ) + s a p 。) ( 2 - 9 ) 二0 = 1 口= 1 d0 = 1 已经发展了许多适合不同材料的原子势函数巧9 1 需要注意，每种势函数都存在一 定的适用范围。e a m 等原子势虽然比传统对势复杂，但计算量并没有明显增加。 目前采用e a m 等原子势的分子动力学在个人计算机上可模拟数十万原子，在并 行计算机上则己达到千万原子。 如对于面心立方( f c c ) 金属： 妒( r 帮) = a 1 ( r 。一r 筇) 2e x p ( 一c l r 即) ( 2 1 0 ) f ( r o 卢) = 4 2 ( r 旬一r 帮) 2e x p ( 一c 2 r 卵) ( 2 1 1 ) f ( p o ) = d p a ，。岛 ( 2 1 2 ) 岛= f ( r o a ) ( 2 1 3 ) f 蜂计 i 咚j 斟 见，、r c , 分别为对势项和电子云密度项的截断 图2 - 29 c c 结构 半径。即认为当两原子间距离大于相应截断半径时，对势部分和电子密度贡献部 分可以忽略不计。采用截断半径能大大减小计算量而只造成可以忽略的很小的截 断误差，这是一个加速算法的技术。对于面心金属晶体基本长度单位取为原子最 近邻距离d ，( d = 2 5 5 2 1 0 “。r r l ) 面心立方晶格原子排列如图2 2 所示，晶 格常数a = o 3 6 1nm ，最近邻原子间距离d = a j ，次近邻原子间距离为j d ， 江苏大学工程硕士学位论文 第三近邻为i d ，第四近邻为2 d 。对于体心立方金属，f ( ry ) 一般在1 4 d 处截断。 通过确定参数a l 、a 2 、c l 、c 2 和d 可完全确定原子势。 2 2 3 动力学积分算法 牛顿运动方程求解主要通过多种有限差分法，由于分子动力学计算原子数较 多，对时间和空间的要求较高，目前常用时间积分法，通常用v e r l e t 算法和g e a r 预报法。 1 v e f l e t 算法 v e r l e t 算法由v e l e t 由1 9 6 7 提出，是目前运用最为广泛的算法，该算法通过 原子在r 时刻的位置r ( t ) ，加速度a ( t ) 及t a t 时刻的位置计算出t + a t 时刻 的位置r o + f ) 。 对f + 出时刻，t a t 时刻的位置进行t a l a y 展开： r ( h f ) = 川+ h i l 口o ) a t 2 + i b ( r ) a t 3 + o ( 出4 ) ( 2 一1 4 ) 尺( t - a t ) = r o ) 一y ( f ) r + i 1 口p ) a t 2 - l b ( f ) a t 3 + o ( a t 4 ) ( 2 - 1 5 ) 合并( 2 1 4 ) 、( 2 1 5 ) 式可得v e r l e t 算法基本形式： r p + a t ) = 2 r ( t ) 一r o 一f ) + a ( t ) a t 2 + o ( a t 4 ) ( 2 - 1 6 ) = r ( t + a t 瓦) - 广r ( t - a t ) 考虑到( 2 1 6 ) 式无显式速度项 v ( t + 1 f ) = 矿o 一互1 r ) + 口( r ) f r ( t + a t ) = r ( f ) + y ( f + j 1 f ) f y p ) = 三 矿。一j 1 ，) + 矿( + 互1 r ) ( 2 1 7 ) h o c k n e y 提出l e a p - - f r o g 方法即蛙跳法 ( 2 - 1 8 ) ( 2 1 9 ) ( 2 2 0 ) 2g e a r 预报校正算法1 3 3 1 g e a r 预报校正算法的基本思路是 1 4 江苏大学工程硕士学位论文 ( 1 ) 预报：己知当前时刻t 的坐标、速度、加速度、以及加速度的高阶导数， 通过t a y l a r 级数计算这些量在时刻t + d t 的值： ( 2 ) 计算力：根据坐标预测值计算原子间作用力，并计算实际的加速度； ( 3 ) 校正：加速度实际值与预测值存在差别，根据此加速度差值校正得到t + d t 时刻的坐标、速度、加速度、以及加速度的高阶导数。 有关g e a r 预报校正算法的详细操作可以参考文献3 ，文献吲1 有关于v e r l e t 算法和g e a r 算法比较。 分子动力学算法还有许多变种，s a t o h ( 1 9 9 7 ) 6 2 1 就各种分子动力学格式的 稳定性所作的对比研究表明，最优越的算法格式应是v e r l e t 速度格式，之后才依 次是蛙跳格式，b e e m a n 格式等。g e a r 格式比上述这些格式都要差些，因其更复 杂，且稳定性和能量波动性较差。 2 2 4 周期性边界条件 就目前的计算能力而言，分子动力学方法能模拟的原子数目非常有限，大