聚合物分子模拟中的力场

V01．15 l 9g9年 7月 高分子材料科学与工程 POLYMER MATERIAIS SCIENCE AND 日 GINE匝RING No．4 Jul l99 聚合物分子模拟中的力场 黄世强朱申敏u 时远 (夕 7、 f 湖北大学化学 幕，武汉．430062) 摘要 舟绍丁力场的基本原理度其发展史+综述丁各种力埸的应用范田，说明了在分子模拟中要堪得到可靠而卫有 意jL的结果．一个好的力埽是首要而基本的条件。在高分子模拟中有两个通用的力场，即 Dreiding力场扣Lhiversal 力场 Dreiding力场参敷基于最简单的杂化原理；Lhiver~l力场参敷的设定基于元素、元素的杂化应#化告性的最基 本的原则。并详细描连了lMiversal力场的奉敷设定。 关键诃 抵 高分子 通用撕 场参数 喂 古 1 前言 计算机辅助分子设计的价值很早就为人们所认 识。在结构化学、生物化学[‘ 以及化学教育[ 中取得 了迅速的发展。尽管在3O年前就报导了有关小分子 的非常好的理论计算。分子模拟技术作为一种工具 源于分子力学的发展[ 。而分子力学 ‘ 是在分子水 平上解决问题的非量子力学计算技术。其原理是分 子内部应力 ，在一定程度上反应被计算的分子结构 的相对位能大小，从而可用来确定分子结构的相对 稳定性，因此可广泛地用于计算各类化合物的分子 构象，热力学参数和谱学参数 分子力学从几个主要 的典型结构参数和作用力上来讨 论分子的结 构变 形 ，即用表征键长、键角和二面角变化以及非键相互 作用位能函数来描述分子结构改变所引起的分子内 部应力或能量的变化。在过去的 20年内，在开发适 合于预测分子的结构和动态性质的力场方面取得了 很大的成绩LsJ。基于分子性质的经验力场的理论计 算首先由Westheimer提出并加以运用 』。还有些力 场由Warsher和Boyd提出 ]。～linger及其小组开发 的 MMa、MMP2和 MM3力场口 ]，对各种有机和无 机系统非常有用，这种力场具有与量子化学交互使 用 的 能 力，由 此 成 为分 子 模 拟 技 术 的界 碑uI 。 F,~ltman及其研究小组开 发的 Al衄 t力场对蛋 白 质 和 核 酸 有 用 ] 。Karrplus及 其 小 组 开 发 的 CHARMM 力场也对 蛋 白质和 核酸有用 。1979年 Iifson、Hagler和 Dauber开发的力场适合于形态学模 型 ；1 988年 Yashonath等开发 的力场适合于沸石 收稿日抽：1998—12—16 联系^及第一作者 ：萤世强 ．男．s0岁．萌I教授． r 0 ， ／ ’ 结构 的吸附作用，对 C、H、O和 Na原子的参数进行 了设定[132；1990年 Van Beest等开发的 bks1．01对硅 石和铝磷酸的原子问 的相互作 用进行 了设 定 ； 1 992年由 Erik de Vos Burchart基于 ab initio和实验 数据开发的力场描述了硅石和铝磷酸盐 的能力学 ， 这种力场包括了4种原子Si、()、S、P。。 在每一种力 场中，对相似的原子在不同的环境中，力参数和几何 参数常有显著的区别，所以对怎样产生新原子或新 的键型就不清楚，为了完成对那些没有或很少有实 验数据 的分 子进行结 构的预测，由 Mayo、Olafc．on、 Goddard等于 1990年开发的 DR玎n1N(；力场 ，对 结构的预溺及动态计算(包括有机物+生物分子和主 族非金属分子)有着非常好 的效果。1 992年．由 A． K，Rappe及其小组开发的 UNIVERSAl 力场 ”。”一是 一 种针对周期表 中的所有元素 的力场，所有的力学 参数基于元素，杂化和化合的规则而产生，并为许多 结构形式所证实。 随着高分子科学的发展，人们已不满足于仅仅 用实验的手段来研制新型高分子材料和提高现有材 料的性能，开展了以计算机模拟为主要内容的高分 子设计。80年代以来，利用计算机对有机分子和药 物分子的设计得到了很大的发展．英国剑桥大学的 晶体结构数库，美国国家布鲁海文实验室的蛋白质 数据库以及其它众多的化学信息和数据．也为高分 子的设计提供了丰富的资源 I 。以晶体数据库为基 础，将分子力学和量子化学，分子图形技术结台在一 起，对高分子进行模拟以期得到结构和性能的关系。 维普资讯 http://www.cqvip.com 第 4期 荫世强等 ：聚音物分子模拟 中的力场 在 高 分 子 模 拟 中 通 用 的 是 DRIELDING 和 1JNⅣERSAI 力场 同时利用随机抽样技巧或称统计 法的“蒙特卡罗”法在高分子模拟方面也取得了广泛 的应用，但这种方法与一般的力场方法有着很大的 区另U 基于力场和动态算法发展起来的分子动态模拟 克服了“蒙特卡罗”法仅能描述不同温度下分子结构 的特征，却不能描述不同温度下从一种热力学状态 向另一种热力学状态演变的过程+并能提供量值，使 我们观察到非零温度下系统 的行 为和性质，并可确 定热力学运动的结果。由此可见，力场的发展将分子 模拟带人 了一个崭新的境界 ，并起着举足轻重的作 用 由于以前开发的MM2、MM3、AMBER、CHARMM 等力场只对周期表中的某几种元素的参数进行 了设 定，对于含有未参数化元素的聚合物分子的设计显 然是 不适用的，DRⅡDING 和 Ⅷ Ⅵ S A【 力场 的 出 现 解 决 了这 一 难 题 ，两 种 力 场 互 为 补 充。 DREIDING力场对分子量太的结构和具有较少实验 数据的物质具有非常好的预测作用，对各种有机系 统导出精确的势垒，但在计算 Ⅵ)W 相互作用方面 受到限制，不能预测电荷分布。而 UNIVERsAI 力场 考虑了静电能。下面对运用到高分子模拟的通用力 场 Universal参数的设定予以介绍。 2 IINTⅥ琅SAI 力场 UNIVERSAl 力场是针对整个周期表的分子力 学和动态模拟的力场，力场参数是依靠元素、元素的 杂化及化合性而计算 出来的，与电荷平衡计算法相 结合 ．对用 UNIVERSAI 力场进行的计算极 力主张 运用静电电荷计算。 参数和函数形式是分子力学和动态力场的关键 分支结构，分子动态学和能量优化的一个重要 的应 用是用来预测新分子的结构+以前的力场 ，基于经典 的工作+仅限于某些原子的化合物 ，如蛋 白质、有机 物、以及核酸 ．不能提供取得参数的系统 的过程，在 用键角弯曲函数描述非线性分子 当 0趋于 1 80。时， 角扭变能得出的结论是不正确的。为了完成对各种 原子的研究 ，开发了一种新的力场，它用最通用的规 则基于简单的相互作用来计算力场参数 ，其基本参 数只与元素 +杂化及其相互作用有关 ，角扭变函数形 式对大的取代体系很有用。 用于产生 UFF的参数依靠原子的成键半径 的 杂化、杂化角、范德华参数，扭变和键角面外弯曲势 垒以及有效核电荷数。 2．1 原子类型 力场元素即原子类型。同种类型的原子有物理 和化学上的相似性 ，通常在力场中被认为是等同的 UFF有 l26种原子类型，每种用 5位来表示：前两位 是元素符号；第三位是杂化类型或几何形状 ；第四位 表示含的氢数 ；第五位表示变更参数 ．比如氧化态 2．2 力场形式 一 个体系的能量可表示为各种成键相互作用能 和非键相互作用能。 E一 + B + 巴 + + ‰ &． (1) 成键相互作用由键伸缩 、角弯 曲、二面角扭转和键角 面外弯曲能等子项组成 非键相互作用由范德华能 和静 电效应部分组成。 2．3 键伸缩 Harmonic形式 ：E 一1／2·K ．‘ 。一 ．) Morse形式：&一岛[exp~--aO‘一 ))一1] 式中 K，为力常数{ ．为键长；日．为离解能 以下 讨论 Harmonic形式的键伸缩 ： 2．3．1 键半径：平衡键长 ‘．定义为所指的原子单键 的半径加上键级修正和电负性修正值： ri J一 ^一 + 一 【2) 说明：H、C、N、O的半径从小 的有机分子得到．氢 的 半径从丙烷中甲基的 C—H距离求得为 0．1112 nm， C一3、C一2、C一1从而烷、丙烯、丙炔 中求得单键、单 键、单键、驭键、三键 相应 地 为 u．1526 DiTI、0．1 5O] fllTl、O．1458 nm、0．1336 nm 和 0．1207 D1TI；(：一R 的 R从苯中求得为 0．1 3390 D1TI；N一3、N—R、N一2、 一 1从二 甲基胺、N一甲基 甲酰胺 、二 甲基二氮烯．￡ 腈中的 C—N单键、C—N单键、N兰N、碳氯 三键 的 键距 求得 分别 为 0．1 426 D1TI、O．1 4 59 1711YI、0．1_917 D1TI、0．11 57 nlTl；氧的半径 O--3()--R()一2从甲基 醚、甲基乙烯醚、丙酮的C一()单键、C一()单键、(’一 O双键而得+分别为 0．1 410 Elm、0．1 428 nm、0．1222 nm；0—1从 C的 C—O键距而得 ．为 0．1128 nm；(] 一 3—2从(a。Sj) 0中的 一()单键而得为 O 1 592 nm；乙硼烷中的桥氢 H—b从 B—H而得为 0．1 326 nln；第一旗元素的原子半径可从相应的气相同饺聚 合物得到 ；Ca、Sr、Ba的半径从相应的羧酸盐的X刳 线结构得来；()一3的半径定为 0．0557 niTl；xc4+4 半径从平面正方形 XeR 得来 ；F一的半径 为 0．0668 nm；Kr4+4半径从 KrR 外推而得；基于 Xe—F的键 长从 Xe 到 Xe 减少 0．007[1ift：其余的稀有气体 半径可 由此外推得到。 对其它元素 ，单键半径可从包含此元素与 形 式的单键的化合物的实验结构数据而得或用类推 法；气体实验结构数据从印刷的书籍可查；X射线结 维普资讯 http://www.cqvip.com 商分子材料科学与工程 构从牛津数据库而得 ；以前化学上的许多重要的结 构的理解是通过标准键距与实验键距相比较而得来 的。结构的电子效应包括 ：电负性、共轭、金属配位 键、金属配位反键和离域作用。一个能够预测全部分 子结构的力场必须重新产生这些结果 。为了说明这 些结果，用共价键半径和键级来探测数据的真实性， 用名为PAULING的式子来修正单键半径： B：I一 一 (^^ 一 t"i)ln (n) (3) 常数 一0．1332(丙烷、丙烯、丙炔 的 C一3、C一2、C 一 1半径)，n为键级。(YKeeffc—Rrese的电负性校正 因子为 ： ‘ =， [(X) 一( )“ ] ／(五 + rj)(4) 2．3．2 力常数：键伸缩力常数基于原子的Badger 规则而产生。结合能曲线： & =E — FR—G(Zi Zj"／R) (5) 从结合能魄线可得出 F；zi’z 是有效荷电荷；G= 332．O6 设结合能 由离子部分的相互吸引与短程 pauli排斥之和组成，所以从平衡条件可得： K． 一 2G·zi z ／ 一 664．12·五 Zi ／r (6) 互 是有效荷 电荷。 2．4 角扭变 傅立叶扩 展式(式 7)被 用于 UFF来 描述角扭 变。 E = K c0s (7) 2．4．1 键 角弯曲： 一 K k三cⅡcos( ) (8) 为满足适 当的边界条件而选择的系数 cⅡ包括 选 1个在 自然键角 处有最小值的函数 ，傅立叶扩 展式 中的 co 因更好地描述大幅运动，所以比 Harmonic更好一些，并且扩展式基 于数学平衡，有 充足的力学基础 对线性系统、三角平面，平面正方形、八面体的 坐标环境 ，式(8)用 n一0时的 ，相应地用 n一1，3， 4，这样式(8)简化为 岛一K ／ [I—cos(n )] (9) 对非线性系统 ，如水，当 一 =104、5。时 岛= 0；乳的二次导数为力常数； 为 180。时有最大值，为 了正 确选择 角度部 分，考虑 6组 ()、 S、 Se、 Te、 Po。O、S、Se、Te、Po用正交不成对 P轨道形 成共价键，氢化物的平衡键角为 9O。，cos20角描述角 扭变，然而 0中的氢相距太近，由于排斥而使键 角 为 104．5。，从键距来考虑 ．这种空 间排斥在零度 有最大值，在 180。有最小值，扩展式的3部分和形式 为 ： 吕=臣 [G +CL c0 十 COS(20)] 【10) c=!=1／(4sin )， Cl：一 4C!~OS 0． c。=C：(2cos~00+ 1) (11) (1)标准键角：主族元素的键角可从参考的母体 氢化物的标准参考结构而得 ：0～3从 Hz()而得、为 104．5。；S--3从硫化氢而得 ，为 92．2。；例外的有 ()-- 3一Z、0一R和 N一2，它们分别从(CI~Si) O、甲基乙 烯醚和二 甲基二 氮烯而得 ，分 别为 1 6c、l18 3 、 112．3。；其余键角可从以上元素外推而得。 (2)力常数 ：用 Badser规则来产生键角弯曲力常 数 。 岛 =蜀 一 F0一 ( ／r, ) (1 2) 矗= 十r 一 2rilr~cos0 Ki =[ 五’ ／ 嗑 · (1一cos200)一r~ocs0ol 2．4．2 扭转U KL乏问的扭转能：扭转能 B 定义 为： — Ki COS“ (1 3) K．“和系数由旋转势垒、势能的周期性和平衡 角决 定，对所给的中心为J—K键的体系，考虑所有的键 的扭转，每一扭转势垒被扭转数所平分。由 一1， 一一cosr~o ，K 一(1／2)· 代^式 (1 3)得： 日 一 (t／2) [1一c∞， c0s婶。] (1 d) 2．4．3 键 角面外弯曲： B —K讪( +Cl COS + COS 2 ，u)(1 5) 是力常数 ， 是 1I 和 [JK面的角。对所给定的 中心原子1，3个独立的轴(IJ⋯II IK)，每～面外弯曲 势被中心原子 1的面外弯曲数目平分 ，当 =90 时 有最小值，当 一O。时有最大值(适 台 PH J；对于乙 烯而言，当 ∞ 0 o时，COSO~有最小值 ；当 一180 时有 最大值 ；对 3取代的C～2、C—R原子类型 c。=1、C． =一1、G一0；如 C与 0—2成键 ，力常数 为 209．2 kJ／tool，其他的力常数为25．1 kJ／tool。 2、5 VDW 非键相互作用 臣 =岛{一2[x．／硼 十[X ／_Yj ) (16) 式中日为离解能；x 为Ⅵ)w 键长．日．、X．可通过运 用化台规则从同棱原子参数得到。在 U、F中分子参 数用Lennard—Jones几何化台规则 ： 日 ： (日日) 。“ ‘17) 式中日为原子的范德华能。目前开发一个力场的最 大原障碍就是如何取得范德华参数，林纳德一琼斯范 德华参数的取得是通过 6形式，函数式基于物理上 的短程相斥 ，包括保里正交和长程的偶极诱导，偶极 分散。 维普资讯 http://www.cqvip.com 弟 4期 黄世强等 ：聚台物抒 子模拟中的力场 w： Axp(一 ／iX)一 G／ (1 8) 相斥指数 B和偶极分散部分 c6对整个周期表都适 用，式 (1 8)中的第 三个 自由度 A决 定于 Lennard- Jones距离文献值和相斥指数 B之间的经验关系 基于波动函数的长程近似联系为： ～exp[一，‘(2IP)“ ] (1 9) 一 对指数的函数的重叠为： S= 1+}·r+ (1／3)(}·r) ]exp(一 ) (20) 一 对原子问的相互斥力 ： S ～ exp(一 25r) (21) 从 ～exp[一，。(2 ]得出： }= (2IP) (22) B ： 2 一 2(2IP)“ (23) c6= 。Kurr~r和Meath报导的稀有气体的半经 验FKS值用于计算S：s—l_166+0．01626 z。等式 (18)的第 三 自由度 可通过将式 (18)转换成势 阱 深岛、距离部分Xii和形态参数 f： 目nw一 {日1[6／(f一6)·exp(f)]exp(一ti(x／Xi，)) 一 [岛弼f／(f一6)]1／ (24) 式(1 8)与式(24)相比较得距离 ： X f／胃 (25) 势阱深 日=(=6 (}一 6)／U霹 (26) 相斥预指数部分： 4一[日6／(}一6)] (f) (27) 如下选择数据，根据经验来计算参数}： (1)用等式(23)计算出 B，以及文献 中报导的有 关指数一6或 6一l 2式的 X值，成键电子数目和周 期表中旗的函数的 f(即／3／ )成一线性关系，第二 周期元素的散射图大部分归因于用于决定 X的 x 射线晶体结构的部分电荷的选择。 (2)给出电荷排布，选择 运用两 中性过程 ，f一 12．73(对石墨中的 C)；f一1 5．44(对 Ne，以得到线性 关系)； f— l0．02+ 0．6775 (对第二周期元素)(28) (3)对第三周期元素，用最小二乘法处理 X 数 据，导出： f一 8．587+ 0．897 n (29) 由对 作圈显示第三和第四周期间的斜率存在 不连续的现象，这归因于固体中的特殊的分子内键 相互作用，这些特殊的键相互作用不应包含在 Ⅵ)w 参数中．对其余的主族元素，用式(29)表示。 (4)f一8十‰设斜率为 1，从等式(28)的斜率外 推而得稀有气体的值为 1 5，相斥指数 B和分散系数 、形态参数 f用于计算 Lennard—Jones距离 X。日用 等式(25)和式(26)来计算 注 ：对于处于正氧化态的金属，以上的计算不适 用 对阳离子金属： B+一 2(2I ．) t30) jP2．是第二电离能 ～ IPc~ 』 ．= X + 3／2，． c 31 J ．+一C ． (， 口 l，ciF,一d ] L 02， 、 是原子 i的第一、第二电离能 +是金属 的中性和阳离子 的 FKS极化 度导出 的 参数． 在文献中可当作 VIEW参数用 明显地与晶体的性质 相符 目前，氢的半径为 0．2886 nm．D-=0．1 84 M／ tool与 I~nnards-Jones的半径和势讲深相似 ．分别为 0．29267 nm和 0．1 40 kJ／rno1．适用于聚乙烯 对碳 ： ，‘一 0．3851 nm，D一 0．439 M／mo[，在 Lennard—Jones 中计算的 r=0．3805 nrn，D=0，289 kJ／tool(适用于 石墨 )。 2．6 静电相互作用 &l一 332．0637tQQ／E R．) t33) 式中Q和Q为电子单元的电荷：R 为距离；E为介电 常数 UFF中电荷可从 ∞ 电荷平衡式求得 2．7 非键排斥 按一般惯例在 IYFF中所计算 的排斥 和静 电相互作用，不包括 1，2成键作用以段 1．3相互作 用 以上的力场参数已被许多有机和无机分子的结 构所证实 总之，L】NIVERSAI 从简单的规则和原子 参数来构造一种力场，这种力场能够就整个周期表 重新产生化台物的结构因子，键距的误差小于 0，01 rim，角伸展误差为 5 ～10。，所以说 UNIVF~SA[ 力 场与电荷平衡的计算方法有关．用此力场最好用静 电电荷计算(包括所有的键级矫正)，要得到正确的 结果，必须在指定原子类型之后。建立能量表达式之 前，计算分级键级(依靠单个原子的化台)。 参 考 文 献 】 Cal&~ll J Kolknan PA PioPolyme) ． 埘B．25 249 2 J0 I Casanova J Clam Educ ．19 1} ．70【11 J 、90 3 DhmrUr1．HaBter&todd r Rev L m L．Chum ．1991．： 99 ～ l6,1 4 DorigoA E．HoukR N J Am Chem So： tgB7．】09 3698 5 Evooks BR J Comput ( 札 ．1’坤z-l3 3】 25 6～ 6 6 W h m r FH，Ma~．er J E J Chem Pl1y^ -j 46．j 4．77 7 Mav。S[ _ af∞n BD-G 矗如rd W 、 I Phys Chem 1 9’⋯ g4 8897 维普资讯 http://www.cqvip.com 高舟 子材料科学 与工程 1神 年 8 ＆0e Jeffrey L-Hoffmann Rdn~rd W 』 Am．Ch,ml So：．． 1991．1】3【l3)：5006 9 K∞ J．AlEnger N L J Am．Chem Soc．1977．99-975 10 SpTa~e J T．Tai J C．Yuh Y．e【a【』 Cotn~tt Chem l鲥 7，8 58． 11 Weiaer S J．KollM n PA et a1．j Chem Soc ．1984，106：765 1 2 ljf∞n S．Hagler A T-Dau~r P 』 An1 Chem．‰ ， 1979． i01：5111 13 Ya出Dnath S．1l帕ma J M．Nomak A K．e【al Nature．1988． 331 60 l4 Van Be猷 BW H Phy5 Re k tc．1990．64 l955 l5 Burchart E de v％ P}lD．n esis ”Smdi~,on Zotites．Mol~ular Mec[maics Framework S Ktity nd Crvstal Growth ． rec hJ J he Universlteit De Lft．1992 l6 RappeA K．Casewi【C』．Colwell K S．e L a【J Am Clam ． 1992．U4：10024～l(1035 l7 Ca．~,mt C L，Alliager A K．et a1 J Arm Chum ‰ ．1 q2． 1l4：l0035 l8 张 福 贵 (Zhang Fugui) 计算 机 与应 用 化幸 Compu and Ap led Cheamist：y) l989-6∞ 】 0RCE Fm IJ)F0R M0I ECIILAR S玎ⅥULAT10NS OF M_ACR0lM0I ECIIU Huang Shiqiang，Zhu Shenmin，Cheng Shiyuan (Department ofChevm'stry．H如f L~versity) ABSTllACr e general principles and the history of molecular mechanics in molecular simulations were summerized．Each force field had its own range．In order to obtain reliable and meaningful results，a good force field was essential "Ilaere are two generic force fields that aI℃ ctsefuI for modeling structures and dynamics of macromolecules．One s the DRⅡD1NG force field whose parameters are calculated based on Nmple hybridizations 1he other is the U VⅡ AI force field where the force field parammers are estimated using general rules based on the elemem ，its hybridlzafion and its conneclivity．Finally，describie in detail the parameters used in Universal force field． Keywords force field．generic force field．parameters of UFF 维普资讯 http://www.cqvip.com