分子动力学模拟橡胶

1、 摘要随着计算机的性能提高和计算技术的长足发展，利用计算机对物质的微观结构和运动进行数值模拟计算得到了迅速发展，并在此基础上发展了运用数值运算统计求和的方法:分子动力学法。本文通过分子动力学方法，应用Accelrys公司的软件MaterialsStudio来模拟二甲醚在不同橡胶中的扩散行为，用Einstein法和Green-kubo法求得二甲醚在6中橡胶不同温度下的扩散系数，并对两种算法进行比较。基于COMPASS力场进行建模。为了使模拟体系接近真实体系，模型的构建采用三维周期性边界条件，选取能量最低、体系最稳定的周期边界作为初始构型。然后对模拟体系进行能量优化和预平衡，进行NVT系综的分子动

2、力学模拟以便得到二甲醚在橡胶内的扩散过程和运动轨迹，最后采用两种算法计算二甲醚在橡胶内的扩散系数。通过本研究的计算数据，发现分别通过两种算法得到的二甲醚在橡胶中的扩散系数在数量级上保持一致，除去个别点外，两者的相对误差在一定稳定范围内。Green-Kubo法计算需要积分，较为繁琐，误差较大，Einstein法计算简便，计算结果与实验值比较接近。用Green-Kubo法求得二甲醚在343k温度下的氟橡胶中的模拟扩散系数最大，可推测此温度下在所有研究的橡胶中氟橡胶的耐二甲醚的能力最差。基于本文研究得到的数据，进行综合分析后认为本文所研究的橡胶中丙烯腈含量为20%的氢化丁腈橡胶的耐二甲醚性能最好。关

3、键词：二甲醚，分子动力学，橡胶，扩散系数ABSTRACTWiththerapiddevelopmentofcomputerperformanceandcomputingtechnology,theuseofcomputersandthemovementofmaterialmicrostructurenumericalsimulationhavebeendevelopingrapidly,andthedevelopmentofstatisticalsummationofnumericalcalculationmethodbasedonthis:MolecularDynamics.Intheart

4、icle,useMaterialsStudioofAccelrysCompanytosimulationDMEdiffusionbehaviorinthedifferentrubber,andobtainthediffusioncoefficientofDMEin6typesofrubberthathasdifferenttemperaturebyusingEinsteinandGreen-kuboalgorithms,andcomparethetwoalgorithms.CreatemodelbasedontheCOMPASSforcefield.Inordertosimulatethesy

5、stemclosetotherealsystem,buildamodelofthree-dimensionalperiodicboundaryconditions,andselectthelowestenergy,themoststablesystemwithperiodicboundaryastheinitialconfiguration.Thenthesimulationsystemforenergyoptimizationandpre-equilibrium,moleculardynamicssimulationswithNVTensembleinordertogetDMEdiffusi

6、onprocessesandtrajectoriesintherubber,finally,calculateDMEdiffusioncoefficientintherubberwithtwoalgorithms.Afteranalyzingthedataofthisstudy,foundthattheDMEdiffusioncoefficientsintherubberwereobtainedbythetwoalgorithmsconsistentinmagnitude,exceptfromafewpoints,therelativeerrorofthetwoalgorithmswithin

7、astablerange.Green-Kubomethodrequiresintegration,morecumbersome,errorislarge,Einsteinmethodissimple,andtheresultsareclosetotheexperimentalvalue.Underthe343k,simulationofDMEdiffusioncoefficientinfluorinerubberobtainedwithgreen-kubomethodisthebiggestoftheotherrubberinthestudy,itcanbespeculatedthatthea

8、bilityresistantingDMEoffluorinerubberistheworstoftheotherrubberintherubberunderthetemperature.Inthispaper,basedonthedataobtained,afteracomprehensiveanalysis，ithinkthattheabilityresistantingDMEofHNBR20isthebestoftheotherrubberstudiedinthispaper.KEYWORDS:Dimethylether,Moleculardynamics,Rubber，Diffusio

9、ncoefficient目录TOC o 1-5 h z HYPERLINK l bookmark2 第一章绪论1 HYPERLINK l bookmark4 .1.1分子模拟方法的分类1:1.1.1量子力学法1；1.1.2分子力学法1 HYPERLINK l bookmark6 i1.1.3蒙特卡洛(MC)法2I1.1.4分子动力学(MD)法2 HYPERLINK l bookmark8 i1.2分子动力学的发展3 HYPERLINK l bookmark10 I1.3本文的研究意义5 HYPERLINK l bookmark12 :1.4本文研究内容6 HYPERLINK l bookmar

10、k14 I第二章分子动力学方法介绍8 HYPERLINK l bookmark16 I2.1模拟软件MaterialsStudio8 HYPERLINK l bookmark18 I2.1.1MaterialsStudio软件简介8I2.1.2使用模块简介8 HYPERLINK l bookmark26 装2.2分子动力学基本原理9I2.2.1牛顿运动方程和数值积分算法9 HYPERLINK l bookmark64 I2.2.2分子动力学模拟的基本步骤11 HYPERLINK l bookmark66 I2.3分子动力学模拟的关键技术12 HYPERLINK l bookmark68 I2.

11、3.1系综12订2.3.2初始条件12I2.3.3势函数12 HYPERLINK l bookmark108 I2.3.4常见力场14I2.3.5周期性边界条件与最近镜像理论15I2.3.6截断距离与电荷基团15 HYPERLINK l bookmark110 线2.3.7长程力的计算16 HYPERLINK l bookmark116 调温控制17 HYPERLINK l bookmark141 能量最小化19 HYPERLINK l bookmark143 2.4扩散系数的计算20 HYPERLINK l bookmark161 第三章分子动力学模拟过程22 HYPERLINK l boo

12、kmark163 3.1建立二甲醚与橡胶分子模型22 HYPERLINK l bookmark173 3.2构建二甲醚与橡胶分子的晶格体系24 HYPERLINK l bookmark175 3.3进行分子动力学模拟27 HYPERLINK l bookmark177 3.4模拟分析计算均方位移32 HYPERLINK l bookmark179 3.5模拟分析计算速度自相关函数33 HYPERLINK l bookmark181 第四章模拟结果计算与分析35 HYPERLINK l bookmark183 I4.1Einstein方法模拟结果计算与分析35 HYPERLINK l bookm

13、ark187 4.2Green-Kubo方法模拟结果计算与分析37 HYPERLINK l bookmark189 4.3对两种算法进行比较39 HYPERLINK l bookmark191 结论及展望41致谢错误！未定义书签。 HYPERLINK l bookmark195 参考文献42二甲醚在橡胶中扩散特性的分子动力学模拟共44页第 页第一章绪论分子模拟方法的分类分子模拟方法主要包括量子力学(QM)法、分子力学(MM)法、蒙特卡洛(MC)法和分子动力学(MD)法。1.1.1量子力学法量子力学(QuantumMechanics)是描写微观物质的一个物理学理论，和相对论一起被认为是现代物理学

14、的两大基本支柱，是研究微观粒子的运动规律的物理学分支学科，它主要研究原子、分子、凝聚态物质，以及原子核和基本粒子的结构、性质的基础理论，它和相对论一起构成了现代物理学的理论基础。量子力学不仅是近代物理学的基础理论之一，而且在化学等有关学科和许多近代技术中也得到了广泛的应用。量子力学的基木原理包括量子态的概念，运动方程、理论概念和观测物理量之间的对应规则和物理原理。在量子力学中，一个物理体系的状态由波函数表示，波函数的任意线性叠加仍然代表体系的一种可能状态。状态随时间的变化遵循一个线性微分方程，这个方程预言体系的行为，物理量由满足一定条件的、代表某种运算的算符表示;测量处于某一状态的物理体系的某

15、一物理量的操作，对应于代表该量的算符对其波函数的作用;测量的可能取值由该算符的木征方程决定，测量的期待值由一个包含该算符的积分方程计算。波函数的平方代表作为其变数的物理量出现的几率。根据这些基本原理并附以其他必要的假设，量子力学能够解释原子和亚原子的各种现象。随着计算机硬件的发展，量子力学法已可用于大分子包括聚合物的低聚物在内的体系，尽管如此，将量子力学法应用在多达数千的聚合物体系仍然是不现实的。1.1.2分子力学法分子力学(MM)法，也叫做ForceFieldMethed,就是忽略了电子运动而只计算与原子位置相关的体系能量，是一种在原子分子水平上用力场解决问题的非量子力学方法，目前广泛地用于

16、计算分子的构象和能量。分子力学从本质上说上是能量最小值方法，就是在原子间相互作用势的作用下，通过改变粒子分布的几何位型，以能量最小为判据，从而获得体系的最佳结构。分子力学主要应用分子力场方法计算分子的势能，分子力学使用解析经验势能函数来描述分子的势能面。每个分子都有固定的力场。在分子力学实际计算时，将力场分解成不同的组分，使用理论计算和实验拟合的方法建立力场参数，力场参数要有可移植性，它要适应于同类分子的计算，即同一类分子也要有很高的计算精度。分子力学认为分子体系的势能函数是分子体系中原子位置的函数。分子力学将分子体系作为在势能面上运动的力学体系来处理，求解的是经典力学方程,而不是量子力学的薛

17、定愕方程。所以分子力学方法可以求得分子的平衡结构和热力学性质等,但不能得到子体系与电子结构有关的其他性质。I现在，分子力学是模拟蛋白质、核酸等生物大分子结构和性质以及配体和受体相i互作用的常用方法,在生命科学领域得到了广泛的应用。随着分子图形学的不断发展，I分子力学已经广泛应用于分子模型设计。当今优秀的分子设计程序都将分子力学作为II初始模型优化的主要方法，分子力学和分子图形学己经充分地揉合在分子设计中，分I子模型的构建也是分子力学为主，分子力学方法是九三级辅助分子设计中常用的方法,I尤其是在药物设计中，已经离不开分子力学计算和模拟方法。应用分子力学方法，可1以非常迅速得到分子的低能构象，通过

18、构象分析可以获得合理的药效构象和药效集II团。如果己知受体的三维结构，可以用分子力学模拟药物与受体的相互作用。在分子!的定量结构活性关系研究中，也需要用分子力学方法进行计算。I装1.1.3蒙特卡洛(MC)法IIIi蒙特卡洛(MonteCarlo)方法，或称计算机随机模拟方法或随机抽样方法或统计试I验方法，属于计算数学的一个分支，是一种基于“随机数”的计算方法oMonteCarloII方法也称为统计模拟方法，是二十世纪四十年代中期由于科学技术的发展和电子计算订机的发明，而被提出的一种以概率统计理论为指导的一类非常重要的数值计算方法。IIMonteCarlo方法可以很好地用来对付维数的灾难，因为这

19、种方法的计算复杂性II不再依赖于维数，以前那些本来是无法计算的问题现在也能够计算。科学家们提出了I许多“方差缩减”技巧来提高方法的效率。MonteCarlo模拟适用于研究复杂体系，线研究具有多得数不清的结构、状态的体系，对此我们可以采用蒙特卡洛模拟，以统计I的方法寻找出现几率最高的结构、状态，或相应的有关数据。IMonteCarlo方法作为一种概率性统计方法在相空间中形成Markov链，每一尝试步移动的结果只依赖于上一步，也就是相空间中的随机行走。它局限于平衡态热力学I量的计算，一般不能预测体系的动力学特性，平衡态物理量通过系综平均得到。II1.1.4分子动力学(MD)法III自五十年代中期开

20、始，分子动力学方法得到了广泛的应用。它和蒙特卡洛方法一起己经成为计算机模拟的重要方法。应用分子动力学方法取得了很多重要成果，例如气体或液体的状态方程、相变问题、吸附问题等，以及非平衡过程的研究。它的应用I已经从化学反应、生物学的蛋白质到重离子碰撞等广泛的学科研究领域。分子动力学是在原子、分子水平上求解多体问题的重要的计算机模拟方法。通过求解所有粒子的运动方程，分子动力学方法可以用于模拟与粒子运动路径相关的基本过程。在分子动力学中，粒子的运动行为是通过经典的牛顿运动方程所描述。分子动力学的发展1956年，Alder和Wainwrightl通过研究硬球的集聚首次报道了正确的分子动力学模拟，采用分子

21、动力学研究气体和液体的状态方程，结果显示这种方法在解决实际问题中能够起到重要的作用，从而开创了利用分子动力学模拟方法研究物质宏观性质的先例。后来，人们对这一方法作了许多改进，并运用它对固体及其缺陷以及液体作了大量的研究。真实材料模型的首次分子动力学模拟在1959年被报道了（发表于1960年），该研究由位于Brookhaven的Vilneyard2领导的小组开展，模拟了钻晶体的辐射危害O1963年Rahman3采用连续势模型对真实流体（氢）进行了首次分子动力学模拟。1967年Verlet给出了著名的Verlet算法，就是在分子动力学模拟中对粒子运动的距离、速度和加速度的逐步计算法，这种算法后来被

22、广泛采用，为分子动力学模拟计算做出了很大贡献。1973年Rahman和Stillinger对水（液态多原子分子）做了第一原理分子动力学模拟计算。1980年，Andersen做了恒压状态下的分子动力学研究，提出了等压分子动力学模型；在同一年，Hoover7进行了非平衡态分子动力学的研究。1981年，Parrinello和Rahman给出了恒定压强的分子动力学模型，他们将等压分子动力学推广到允许元胞的形状可以随其中粒子的运动而发生改变的范围，为之后的发展做出了里程碑的贡献。1983-1984年间，Daw和Baskes9，10根据准原子概念和密度泛函理论，提出了一个能较好描述过渡金属各种性质的新方法

23、一镶嵌原子法（EmbeddedAtomMethod，简称EAM）。后来经过Johnson11的改进，EAM法在固态材料的诸多领域己经取得了很大成功12-15。1984年，Nose提出恒温分子动力学方法14。1985年,Car和Parrinello16又提出了共价键系统的分子动力学从头计算法,此算法需在高性能计算机上进行。1986年，Frisch,Hasslacher和Pomeau17提出晶格气体自动学。1989年，Frenkel和Emst提出晶格气体自动学和长时尾的处理方法。因为计算机速度运行速度慢以及内存不够大，早期模拟的空间尺度和时间尺度都受到很大限制。多体势函数的提出和发展，为分子动力学

24、模拟技术注入了新的活力。分子动为学模拟不但能得到原子的活动细节，还能像做实验一样进行各种观察。对于平衡系统，可以用分子动力学模拟作适当的时间平均来计算一个物理量的统计平均值。对于非平衡系统，发生在一个分子动力学观察时间内的物理现象也可以用分子动力学计算进行直接模拟。特别是许多在实际实验中无法获得的微观细节，而在分子动力学模拟中都可以方便地观察得到。随着计算机技术的飞速发展和分子模拟技术的不断进步，利用计算机模拟粒子在分子、原子微观水平上的运动，并对粒子运动轨迹进行分析，可以计算出粒子的扩散系数。用分子模拟方法求取扩散系数，一方面节省了大量的时间和成本，另一方面能够模拟在极端环境，比如在高温高压

25、以及临界条件下的扩散系数。对扩散过程的分子模拟研究，国内外己经有很多报道。周健等18,19运用分子动力学方法模拟计算了甲烷、氧气等气体在常温水中的扩散系数，其模拟计算结果与实验值吻合较好，验证了模型的准确性。他们用同样的方法计算了在不同温度下水的自扩散系数，除了在低温条件下模拟值与实验值差异较大外，其它温度下模拟的扩散系数值与实验值吻合较好。赵长伟等20测定了298.15K328.15K不同温度下葡萄糖水溶液的积分扩散系数，考察了温度对扩散系数的影响，提出一个有温度关联项的半经验模型。朱春英等21以Fick定律和不可逆过程热力学为基础，提出了氨基酸在水溶液中扩散系数模型，并利用该模型对多种氨基

26、酸及尿素在水溶液中的扩散系数进行了计算。孙炜等22采用分子动力学模拟方法，计算了L-J流体氩及氩/氪溶液的自扩散系数和互扩散系数，并研究了扩散系数与温度的关系。在扩散系数的数据处理方面，陈敏伯等23针对模拟计算的均方位移曲线，在原来Einstein方法的基础上有所改进，提出了微分一曲线变分法以及线性回归-区限变分法，能够有效的剔除无效数据，使计算结果更加准确。Pavel等24通过分子动力学模拟计算氧气和二氧化碳在无定形聚对苯二甲酸乙二酯、芳环聚酯及其共混物中的扩散过程，考察了温度、高分子动态性质、密度和自由体积分布对扩散系数的影响。Karayiannis等25模拟小分子在PEI和PET中的扩散

27、过程，研究发现透过物分子在高分子中的扩散主要取决于能容纳渗透物分子的自由体积及其分布。小分子在聚合物或其它多孔介质内的扩散较之本体扩散要复杂得多，这是当前理论研究的热点之一。有关小分子在聚合物中扩散的模拟研究方面，黄宇等对聚二甲基硅氧烷(PDMS)和聚三甲基硅烷基苯乙烯(PSI)体系进行模拟，发现COMPASS力场对于计算氧气和氮气在PDMS和PSI中的扩散非常接近实验值,氧气运动范围大于氮气，小分子的运动轨迹基本与聚合物的自由体积分布相对应。Hofmann等先后采用PCFF和COMPASS等力场对小分子气体在PDMS、聚酰亚胺(PI)等聚合物中的扩散进行了模拟，并分析了聚合物的结构对扩散的影

28、响，他们发现，由于玻璃态聚合物链段结合紧密而难于形成分子扩散的通道，小分子气体在其内部的扩散系数要小于其在橡胶态聚合物中的扩散系数。Boyd等研究了甲烷在聚对苯二甲酸乙二酯(PET)、聚2,6-萘二甲酸乙二酯(PEN)、聚乙烯(PE)和聚异丁烯(PIB)等聚合物内的扩散，分析了体系温度对扩散的影响；研究发现随着温度的升高，甲烷跳跃频繁，扩散速率也增大；由于PIB的自由体积小于PE,甲烷在PIB内的扩散系数要小于其在PE中的扩散系数。在这些气体扩散的模拟研究中，模拟体系中气体分子的数目是任意给定的，这势必会影响气体扩散性质计算的准确性，因为多孔介质内吸附质的扩散和其在吸附剂内的密度紧密相关，如果

29、任意给定小分子在聚合物微孔中的数目，经常不能反映实际条件下的吸附量而导致得到的模拟结果不可靠。2010年，Farkhondeh等26在很大温度范围内对氩、甲烷、二氧化碳、氮气等多种气体在聚苯乙烯中的渗透和扩散进行了分子动力学模拟，计算得到的扩散系数与实验数据以及以前的模拟结果吻合很好。2010年Garruchet等27应用分子动力学方法对在950K-1600K时氧原子在镍中的扩散进行了研究，发现氧气扩散符合空位扩散机理。这些研究都对以后的科学研究提供了宝贵的经验。然而，进一步扩大分子动力学模拟的应用面仍然是一个热门的发展方向。本文的研究意义能源作为重要的战略资源，关乎国家的军事安全、经济发展、

30、乃至国民生活。在当今能源消费结构中，石油因其广泛的使用与不可再生性,成为最重要的能源之一。随着我国汽车保有量的增加和经济持续的高增长,我国能源不足的问题日益严重。尤其近年来，我国石油的对外依存度逐年攀升，这严重威胁到国家的能源安全。随着人民生活水平的提高和车价的不断降低，人们积累的汽车需求不断得到释放。预计我国汽车保有量在今后一段时间内将持续增加，因此国内车用燃料油需求会进一步增大。为了减少对石油的消耗，寻找车用替代燃料势在必行。二甲醚是一种理想的车用替代燃料。常温常压下，二甲醚是一种淡淡香味的无色无毒气体。同时二甲醚不会腐蚀金属。二甲醚在光化学反应中不会产生甲醛，因此不会破坏臭氧层，是一种环

31、保燃料。天然气、甲醇、液化石油气的十六烷值都不超过10只能用于点燃式发动机。而二甲醚的诸多性质决定了它适用于压燃式发动机，可以做柴油发动机的替代燃料。首先，二甲醚的十六烷值大于55，自燃温度低于柴油；其次，二甲醚的蒸发压力低，容易被液化。二甲醚做为柴油发动机代用燃料的优势。二甲醚的化学结构式为CH3-0-CH3，它是最简单的一种醚结构，其中不含C-C键，只有C-O键和C-H键，氧气质量分数为34.8%，因此燃烧后生成的一氧化碳、碳氢化合物和微粒比较少，发动机可以承受较高的废气再循环率以降低氮氧化合物的生成和排放。二甲醚的低热值比柴油小，因此为了达到和柴油相当的动力水平，必须增大每循环供油量。二

32、甲醚理论混合气热值是3066.7KJ/Kg，而二甲醚的理论混合气热值是2911KJ/Kg，因此柴油机燃用二甲醚的升功率会升高。二甲醚的沸点低，这使得它在喷入气缸后就可以气化，它的油束雾化特性将明显优于柴油，有可能在低喷射压力下就能满足燃烧要求。二甲醚作为柴油机代用燃料还可以降低发动机噪音。现有的汽车发动机只要略加改造就可以使用二甲醚燃料。二甲醚是重要的甲醇衍生物，传统上采用气相甲醇催化脱水法生产，这种方法的优点是简单，缺点是产量低，成本高。HaldorTopsoeA/S公司开发了一种由合成气生产二甲醚的低成本技术，省去了中间制甲醇步骤，而这种合成气几乎可以从任何含碳的原料转化而得到，如煤、天然

33、气、木材等，从而大大降低了生产成本。所以，二甲醚是一种可以不依赖于石油的可再生燃料。利用我国丰富的煤炭资源，再结合二甲醚生产技术与工艺的改进和批量的增加，二甲醚的价格将会降低，可以加快二甲醚用作柴油机的代用燃料的推广。城市公交汽车由于行驶里程多，运行又集中于市区，重型大客车和汽油车因汽缸容量大，排气量大，排放的氮氧化物较多，因而城市公交汽车已成为城市机动车尾气污染治理的重中之重。由于二甲醚具有清洁燃料性能，所以二甲醚燃料汽车将有望成为新一代城市客车的替代燃料。二甲醚做为车用燃料有助于缓解我国柴油供应紧张的局面，有利于环境保护，市场前景十分广阔。但DME是一种溶剂，能够溶解多种橡胶材料。传统柴油

34、机如果不经过改造直接使用DME燃料，燃油供给系统中的耐油橡胶密封件会因长期接触DME发生溶胀和老化等现象，力学性能会逐渐变差从而失去密封作用，影响发动机的可靠运行。国内外关于DME作为柴油机代用燃料的研究时间还很短，对橡胶材料耐DME的性能研究还不足。因此研究开发耐DME的橡胶密封材料具有重要的实际意义，有利于DME汽车的推广。通过研究二甲醚在橡胶中的扩散系数可以了解橡胶的耐二甲醚性能。预测扩散系数的方法可以分为理论方法、经验关联和实验方法三类.理论方法主要有Stokes-Einstein方法、Darken方法、Hartley-Crank方法、Eyring方法和各种摩擦系数方法。它们一般不能精

35、确地预测液相的扩散系数。经验关联方法因为实验数据有限，应用范围也很窄。实验方法包括物理实验和计算机实验.由于物理实验代价昂贵，所以计算机实验方法受到广泛的重视。本文运用分子动力学方法对二甲醚在橡胶材料的扩散系数进行研究。随着计算机模拟技术的发展，运用分子动力学方法获得的数据和实际试验数据比较接近。分子动力学模拟通过对研究物系中粒子的运动方程的求解从而得到粒子的运动速度和运动轨迹，通过统计计算即可求得扩散系数，由于是计算机实验，它的成本比较低，并且可以缩短试验时间，因此它是研究扩散系数等传递性质的有效手段。并且运用分子动力学模拟可以看到很多在实际实验中无法获得的微观细节，这可以在橡胶材料的研究和

36、开发中发挥重大作用。本文研究内容本文分别用Einstein方法和Green-Kubo方法计算二甲醚在氯丁橡胶、氟橡胶、丙烯腈含量为20%的丁腈橡胶和丙烯腈含量为15%、20%、30%的氢化丁腈橡胶在三种不同温度下的扩散系数，并对两种算法进行比较。第二章分子动力学方法介绍模拟软件MaterialsStudioMaterialsStudio软件简介美国ACCELRYS公司在2000年初推出了新一代的模拟软件MaterialsStudio,将高质量的材料模拟带入了个人电脑（PC）的时代。本文计算采用的软件正是Accelrys公司推出的基于PC平台的材料模拟软件MaterialsStudiooMate

37、rialsStudio是Accelrys专门为材料科学领域开发的可以在个人电脑上运行的材料计算软件,能够帮助研究人员解决现在化学及材料工业中的许多重要问题。支持Windows,Unix以及Linux等多种操作平台的MaterialsStudio使化学及材料科学的研究者们能更方便的建立三维分子模型,深入的分析有机、无机晶体、无定形材料以及聚合物。该软件可以生成高质量的图片,同时可以处理各种不同来源的图形、表格和文本。因为MaterialsStudio可以综合运用多种先进算法,所以它的模拟能力非常强大。不管论是性质预测、聚合物建模还是X射线衍射模拟,我们都可以通过一些简单易学的操作来得到切实可靠的

38、数据。该软件模拟的内容包含了催化剂、聚合物、固体化学、结晶学、晶粉衍射以及材料特性等材料科学研究领域的主要课题。2.1.2使用模块简介MaterialsStudio包含了很多模块,每个模块都有不用的功能。本文研究二甲醚在橡胶中的扩散系数,主要用到其中的AmorphousCell,Discover和Forcite三个模块。（1）AmorphousCell模块AmorphousCell是一个可以建立复杂无定型系统的代表性模型,并且可以预测它的主要性质的模块。AmorphousCell创建结构采用的基于建立的产生包含链分子的疏松无序系统。它的特点是：任意混合体系的建模方法、特殊的产生有序向列型中间相

39、以及层状无定型材料的能力。它通过观察体系结构和性质的关系,对分子的一些重要性质有更深入的了解,从而设计出更好的新物质的配方。它可以预测并研究的性质包括状态方程行为、内聚能密度、链堆砌和局部链运动。（2）Discover模块Discover是MaterialsStudio中非常重要的模块,它提供了功能强大的原子模拟方法,适用于大部分的分子和材料。Discover集成的很大范围已经被证明适用于分子设计的动力学方法,它以PCFF、COMPASS、CVFF等多个力场为基础,能够准确地计算出分子的最低能量构像，同时能够给出不同系综下体系结构的动力学轨迹oDiscover为AmorphousCell等模块

40、提供了计算的基础。周期性边界条件的引入使得它可以对固态体系进行研究，如晶体、非晶和溶剂化体系。Discover还提供强大的分析工具，能对模拟结果进行分析，从而得到各类结构参数、力学性质、动力学量和热力学性质等。(3)Forcite模块Forcite是分子力学模块,能够使用经典力学对任意分子和周期系统进行势能和几何优化计算。Forcite支持UFF,、COMPASS和Dreiding力场。因为支持力场的广泛性，Forcite原则上可以处理任何材料。几何优化算法提供了最速下降法、共扼梯度法和牛顿法，并且提供了连续使用这些方法的灵巧方法。这使得研究人员可以准确地进行能量最小化计算。2.2分子动力学基

41、本原理2.2.1牛顿运动方程和数值积分算法分子动力学是对物理系统确定的微观描述，该系统既可以是多体系统，也可以是少体系统。该方法的基本原理是:建立一个粒子系统，对所研究的微观现象进行仿真，各粒子间的相互作用根据量子力学来确定。对于符合经典牛顿力学规律的大量粒子系统，通过粒子运动学方程组的数值求解，得到粒子在相空间的运动规律和轨迹，然后按照统计物理原理得出该系统相应的宏观物理特性。以下是分子动力学模拟的基本原理:在包括N个粒子的分子动力学模拟系统中，系统的总能量是系统中粒子动能与势能的总和。其中，势能是系统中粒子位置的函数：U(r,r,.r)。12N通常认为势能为各粒子对间势能的加成，即:2-1

42、)Uu+u+.+u+u+.U+.U/、因石(r)121323242N(N-1JNijiji1ji+1公式中，uij按照经典力学的作用力为：表示第i、j个粒子的势能；r是粒子i、j间的距离。ij2-2)系统中任何一个粒子所受的力为势能的梯度，则粒子j对粒子idU()FFj_j.jjiorij公式中F是粒子j对粒子i的作用力；F是粒子i对粒子j的作用力。ijji根据叠加原理和牛顿运动定律，可得到第i个粒子所受的作用力F和加速i度a的表达式为：iF=Fiijj=lj幻Fa二limi2-3)2-4)根据加速度a与粒子的位移和速度的导数关系，可以对加速度a积分从而得到第iii个粒子的速度和位移:d2rd

43、vi=i=adt2dt2-5)公式中，rv分别是第i个粒子的位置和速度。ii然后，一般采用各种有限差分法来求解运动方程，从而得到粒子的运动，实现对粒子的运动规律和运动轨迹的模拟，最终实现对系统的宏观特性模拟计算的目的。有很多可以求解牛顿运动方程的有限差分算，但是分子动力学计算的原子数目比较多，因此对空间和时间的要求都比较高，需要给出适合用于分子动力学的时间积分算法。时间积分算法是基于有限差分法，对时间进行离散，At是时间步长，根据t时刻的变量和它的时间导数，可以推演得到t+At时刻的物理量，逐步计算就可以得到整个历史过程的物理量。这个过程是近似的，因此必然会有误差，包括截断误差和舍入误差。截断

44、误差主要是因为有限差分法用Taylor展开截断某些高阶项引入的，舍入误差是由于计算机的本身数值精度引入的。两种误差都可以通过减小时间步长At而减小。对于比较大的At，主要是截断误差起作用，但是随着At的减小截断误差会迅速降低。舍入误差随着At下降相对来说比较慢，在较小的At时起主导作用，采用双精度数据类型(64位)有助于把截断误差降到最低。在分子动力学计算领域主要采用Gear算法和Verlet算法。(l)Verlet算法在1967年被提出，之后有对它进行了一些变化。该算法基于Taylor展开：r(/+At)=r(/)+v(tkt+a(t)At2+b(t)At3+OCt4)26r(t-At)=r

45、(t)-v(t)At+a(t)At2-b(t)At3+oCt4)262-6)两项相加，可以得到/、rC+At)二2rC)-r(/-At)+a(tkt2-OAt4丿2-7)这就是Verlet算法的基本形式。从此式能得到粒子在相空间的轨迹，而且该式计算的截断误差是oCt4)。公式(2-7)中没有v(t)，可以通过下面的公式计算：v(t)=厂C+山)厂C一山)(2-8)2At得出粒子的速度就可以计算系统的动能，从公式(2-8)可以看出速度计算的截断误差是OCt2)。在下一步的位置没有得到之前，难以得到速度项。而且，它不是一个自启动算法。在t=0时，只有一组位置，因此需要运用别的方法得到t-At的位置

46、。为了解决这个问题，对Verlet算法进行了一些改变。Verlet算法的速度形式具有的形式如下：在一定距离(spline宽度)ij2-3图所示内能够让能量的变化变的缓慢，同时计算量又不会过大，其作用原理如U=工wSC丄NBijijiji、jR丫2ijrIj丿R丫ijrIj丿丄yX)q.q.+wSijijij4k88ri、j0ij2-15)公式内的w是加权因子(0w1),S(r)的函数形式如下方程所示ijijS(r)=ij1ra1+y(rbby(r)-3】arb2-16)电荷基团通常采用基于基团的截断距离来减小直接采用基于原子截断距离而引起的人为的单极一单极相互作用和偶极分裂。电荷基团是由临近原

47、子组成中性基团。一般将电荷基团中临近几何中心的原子定义为中心原子。如果中心原子在截断距离之内，电荷基团中所有的原子都会包括在相互作用的计算中。电荷基团通常是化学官能团。计算时，引入电荷基团能够减小非键相互作用的计算量，可以提高计算的速度同时计算的精度不会受到影响。长程力的计算处理短程力一般用周期性边界条件和截断半径，但是对于静电引力等长程有效作用力，简单的位能截断已经不再适用。需要使用别的处理方法，经常使用的处理方法如下：(1)Tree-Based方法基于树结构的方法是一种被广泛应用的计算算法，这种方法是基于树的方法计算效率非常高，在处理长程作用方面的优势非常明显，但是不经常被用在计算模拟上。

48、(2)EwaldSum加和法EwaldSum加和是比较常用于处理长程力作用的方法，它的作用原理是在空间里处理截断半径内的作用，周期性边界条件下处理截断半径以外的长程力。反应场方法反应场方法也是一种常见的计算长程力的方法，它基本思想是作用偶极子的场分成两个部分。一部分是由位于截断半径内的短程作用，另一部分由截断半径外的作用贡献，这部分堪称是介电常数的反应场。N(4)PPMM法PPMM法由Eastwood等人在1980年提出，短程力处理PP(particle-particle)法，与长程力计算处理方法比较相似；长程力处理采用PM(particle-meth)法。PP和PM两种方法结合得到的结果比较

49、精准。调温控制因为外部或者摩擦力引起的过热，积分误差和力的截断，这些情况都会引起能量的飘移，所以控制系统的温度就显得非常必要。经常使用的控温方法有:速度标度法、Gaussian热浴、Nose-Hoover方法和Berendsen弱耦合方法。(1)速度标度法改变温度最简单和最直观的方法是直接对速度进行标度，如果t时刻的温度是T(t)，速度乘以因子九后温度的变化是八盔(2-控制温度最简单的方法是在每一步乘上标度因子九。上面公式中，T是期望的req参考温度。(2)Gaussian热浴这种方法的思路是在运动方程中加上一项“摩擦力”项，这个力又和系统的速度联系在一起。运动方程如下(2-19)md2rid

50、t2mv.v，用来保持总动能守恒，为了把温度调到期望的温度下，需要把初始的温度设定为期望值。(3)Nose-Hoover方法1984年Nose提出了“Nose动力学”，这种方法最有力的方法是扩展体系Hamilton量的做法，不但能够用于恒温体系，而且可以启发延用到恒压体系或恒温一恒压体系。Nose动力学的直接结果是实现了正则分布。之后，Hoover对这个方法进行了修正，在运动方程中引入了一个摩擦项和热储项来延展系统的哈密顿量，摩擦力和每个粒子速度和摩擦参数g(热浴参数)的乘积成比例，这个摩擦参数是一个动态量。Hoover的运动方程能够用“实时间”的等间距抽样得到正则分布，这对非平衡的模拟显得更

51、加重要。Hoover方程中，粒子的运动方程被替代为:d2rdridt2dt(2-20)其中，热浴参数g的运动方程是:dgdt(T-T)(2-21)T为系统的当前瞬时温度，T为参考温度。藕合的力度由常数Q(储热的“质量0参数”)和参考温度共同决定。用质量参数描述耦合强度有些困难，为了保持藕合强度，可以使Q和T相0联系起来:Q二to(2-22)4兀2T是系统和热储之间的动能振荡周期，它独立于系统的大小和参考温度。T(4)Berendsen弱耦合方法设体系放置在一个温度为T的大热浴内，体系温度为T(t)。体系和热浴之间的导o热过程在唯象上服从Fourier定律:dT)二-TTC)(2-23)d&丿九

52、0比例系数是(1/九,)的量纲是时间，叫做弛豫时间。导热弛豫时间越小表明导热速度越快。把Tt+绎2)在t-卜2处Taylor展开，再根据Fourier定律就可以得到：1+T(tAt2)T(tAt2)1At=n2九i(2-24)或者(2-25)(2-26)化简后可得到in是标度因子，上式在推导过程中认为温度是体系粒子动能之和。即:3NKTG)仝2B2iiii=1当体系粒子速度的绝对值V厂偏小时，贝I相对应的温度T于把Tt+，反之亦然。所以通过这种方法可以调节温度。从公式(2-25)可以得到n1，把n1代到公式(2-26)得到的V+A%就会变大，等同Berendsen弱耦合算法把一个系统弛豫到目标

53、温度的效果很好。但是当系统达到平衡后，弱耦合算法不能探索到一个正确的正贝系综。能量最小化为了使模拟计算结果更准确需要对分子结构进行优化。当一个无规贝晶胞生成时，分子有可能不是等价地分布在晶胞中的，这样就会造成真空区。为了矫正这个，也要对晶胞优化。对分子的结构和晶胞进行优化的过程实质上就是对两者的能量最小化的过程。对于能量的优化有很多种方法，如最陡下降法、共轭梯度法、牛顿一拉斐孙方法等。最陡下降法是最简单的一种结构优化算法。在最低能量搜索的过程中，最陡下降法对能量函数进行微分，计算梯度，每次沿能量下降最多的方向前进。当搜索位置离能量极小点比较远时，用这种方法能够迅速向极小点靠近。当接近极小点时，

54、会产生振荡，收敛速度变慢。共轭梯度法在最低能量搜索过程中和最陡下降法一样对能量函数进行微分并计算梯度，但是在选择搜索方向时，不但考虑当前的梯度，而且还考虑原来的搜索方向，经过综合决定下一步的搜索方向。共轭梯度法的收敛速度快，但是容易陷入能最局部极小点，而且构象能远离极小点时共轭梯度法稳定性较差。牛顿一拉斐孙方法利用一阶微分确定搜索方向，同时用二阶微分确定梯度在什么地方改变方向。这种方法计算精度高，可以迅速收敛，但是计算量非常大。最陡下降法和共轭梯度法可以用于优化较大的系统，但是计算精度较差。而牛顿一拉斐孙方法不适合用于处理较大的系统。本研究使用的能量最小化方法MaterialsStudio中的

55、SmartMinimizer方法，它集成了最陡下降法、共轭梯度法和牛顿法。在进行能量优化时，如果系统初始构型离能量极小点较远时，首先采用最陡下降法进行优化，使体系迅速向极小点靠近。当能量迭代达到1OOOkcal/mol时，采用共轭梯度法。最后当能量迭代小于l0kcal/mol时采用牛顿一拉斐孙方法，可以提高计算的精度。SmartMinimizer方法结合以上三种方法的优势，能够大大加快构型优化的速度和精度。扩散系数的计算分子扩散系数是表征物质分子扩散能力的物理量。根据斐克定律，组分A在组分B中的分子扩散系数，等于该物质在单位时间内，单位浓度梯度下、经单位面积沿扩散方向传递的物质量。分子动力学模

56、拟（MD）近年来成为模拟传递过程的一种重要手段。分子动力学计算扩散系数有两种方法:一种是通过计算体系的各种速度相关函数并对其积分，即Green-Kubo法;一种是通过计算体系平均平方位移并对其微分，即Einsitein法。（1）Green-Kubo法从微观角度，自扩散系数可以通过粒子的速度相关函数获得工v(0)工vC)ii2-27)如果扩散物质的浓度比较低，粒子之间的相互作用力为短程力，方程（2-27）可以改写为方程（2-28）2-28)D二1辿（0）臥）”t0公式（2-27）,（2-28）即为平衡系统时间相关函数求得的扩散系数Green-Kubo式，其中vC）表示原子i在t时刻的速度。i（2

57、）Einstein法动力学将所模拟分子进行标记，计算标记粒子的浓度分布，假定在时间t=o,标记粒子被集中在坐标原点。为了计算浓度分布的时间演绎，须联立Fick定律和物质守恒方程得到:-DV2c（t）=0（2-29）dt2在边界条件c（r,0）=5（r）下，（2-29）得到在式(2-30)中，令jc(r,t丄2dr=1(2-30)(2-31)将式(2-31)乘以r2并在全空间上积分，得到方程：jr2c(r,tr=Djr2V2c(r,t)drdt认r2(t一(2-32)左边等于。右边进行分布积分，得到(2-33)，它将扩散系数D与dt浓度分布联系起来，D是宏观传递系数，r2G)是微观解释。(2-3

58、2)=DJr2V2c(r,t)dr=2dD(2-33)dt对于每个粒子i,测量在时间t内经过一定的距离，得到均方距离对时间的函数关系如式(2-34)/Ar()2=Ar()2t=1基于Einstein方程得到无限稀释扩散系数为:Dg=lim/|r(t)-r6tzz自扩散系数为：(2-34)(2-35)D=lim丄:迄|r(t)-r(0)|selfiitTg2iii=1(2-36)6Nt式中极限的部分近似用均方位移对时间的变化率来代替，就是说用均方位移随时间变化的曲线斜率a来代替。因此上式简化为D=(2-37)6理论上，通过Green-Kubo法和Einstein方法得到的扩散系数是相等的，但对于

59、不同的体系两种方法的使用范围有所不同。本文使用这两种算法分别计算二甲醚在橡胶内的的扩散系数。第三章分子动力学模拟过程3.1建立二甲醚与橡胶分子模型首先构建二甲醚的分子模型，二甲醚由一个氧原子，两个甲基构成，其中与氧原子相连的两个甲基之间的键角为109.403。如图3.1：图3.1二甲醚分子模型图中红色球体代表氧原子，灰色球体代表碳原子，白色球体代表氢原子。然后构建橡胶分子的模型。本次模拟对象是极性橡胶，包括丁腈橡胶、氯丁橡胶、氢化丁腈橡胶和氟橡胶。丁睛橡胶分子模型的构建丁腈橡胶，简称NBR,是由丁二烯（CH2=CH-CH=CH2）和丙烯腈（CH2=CH-CN）经过乳液共聚而制得的一种高分子弹性

60、体。耐油性极好，耐磨性较高，耐热性较好，粘接力强。其缺点是耐低温性差、耐臭氧性差，电性能低劣，弹性稍低。丁腈橡胶主要用于制造耐油橡胶制品。其结构式为：CH=CHCN图3.2.丁腈橡胶的结构式丁腈橡胶的中的丙烯腈含量一般在15%50%,因此其按照丙烯腈的含量可以分成极高丙烯腈丁腈橡胶（丙烯腈含量43%以上）、高丙烯腈丁腈橡胶（丙烯腈含量3642%）、中高丙烯腈丁腈橡胶（丙烯腈含量3135%）、中丙烯腈丁腈橡胶（丙烯腈含量2530%）、低丙烯腈丁腈橡胶（丙烯腈含量24%以下）五种类型。本文采用丙烯腈含量约为20%的丁腈橡胶来进行模拟。通过对丙烯腈的含量的计算，能够得到丁二烯与丙烯腈的分子比例，从而