基于cvff的分子间氢键作用的研究

基于cvff的分子间氢键作用的研究

1氢键的再生体学模拟氢键是与电负性强的原子（f、cl、o、s、n等）组成的氢原子，与相邻电负性强或带孤对电子的原子（f、cl、o、s、n等）之间的强非限制作用，大小介于成限制作用和非限制作用之间。分子间氢键在水、醇、酸等溶液中广泛存在,是产生高沸点、溶剂化等作用的重要因素;同时,生命体中分子内的氢键起着举足轻重的作用,如蛋白质的螺旋结构、部分生化分子的活性等,都与氢键有密切的关系。不仅如此,最近,许多X射线衍射实验证据表明,晶体中碳原子也参与氢键,形成所谓弱氢键作用,在许多金属有机分子的结晶过程中,扮演着重要的角色。到目前为止,氢键的性质究竟如何反映、如何界定以及其能量如何计算等尚无定论,对氢键本质的研究尚在进行,仍为化学和物理领域的研究热点之一。迄今为止,取得的共识是,氢键是一种非键作用,形成氢键的原子间既有范德华作用又有静电作用,没有统一标准界定氢键形成的条件,也无法采用一致的函数形式计算氢键的能量。在学术界,氢键的从头算早已为科学家所关注,已有许多报道,对氢键的本质有了进一步的认识,但利用分子力学对氢键本质的研究却很少见报导。本课题组正在将已掌握的晶体中各种氢键的几何数据汇总,希望利用分子力学模拟技术对氢键的本质进行更深入的研究。随着计算机科学的迅速发展,与研究体系的尺度的增加,分子力学已成为当前化学家经常使用的研究手段,也是分子模拟方法的重要组成部分。分子力学基本用途是寻求分子的平衡构型及能量,它是分子中原子间存在化学键,非键原子之间有范德华及静电相互作用这一经典理论的自然扩展。分子力学方法基于经验力场,忽略分子中电子的运动,仅把系统的能量作为各原子核的函数即经验势函数,利用力场方法模拟分子结构,计算分子的性质。在分子力学中,最常用的方法是能量最小化(Minimization)和分子动力学(MolecularDynamics)模拟,这些在MaterialsStudio(简称MS)软件中应用方便。MS软件为世界最大的分子模拟制造商——美国的Accelry公司,将其工作站上使用多年的材料科学模拟软件Cerius2移植到微机上的版本,具有可操作性强、使用方便等特点,已为全世界科研机构、公司、医院、学校等广泛使用。目前国内的正式用户也超过数百家,尚有快速增长的趋势。多年来,本课题组已成功地以MS软件为工具,开展了多领域分子模拟(主要在分子动力学模拟方面)的研究工作,既取得了一些成果[8,9,10,11,12,13,14,8,9,10,11],也积累了不少经验。由于国内使用MS软件起步较晚,用户对软件中许多问题的认识还不充分,比如遇到的一个很重要的问题就是:在MS的力场中如何处理氢键?鉴于氢键在某些系统中的重要性,因此,在使用MS软件时,对力场是如何界定氢键和计算氢键能量的问题,必须有较清醒的认识。本文以乙醇分子为例,通过参与形成氢键的给体和受体间作用能的计算和讨论,对MS软件中常用力场如何处理氢键的问题作一初步说明,为利用MS从事分子力学,特别是分子动力学模拟的科研工作者提供方便。2氢键形成条件通常把与氢原子成键的原子称为给体(donor),以D表示;将与氢形成氢键的原子称为受体(acceptor),以A表示。rAH和rAD分别表示受体与氢或给体间的距离;θ和θhb为两种不同的键角的表示法,互为补角。一般认为氢键的形成条件是:,θhb<90°或θ>90°。氢键势能的计算公式通常采用LJ12-10的形式:式中,D0代表氢键强度,r0为受体原子与氢间的平衡距离。这一形式在Dreiding力场中采用。此外,YETI力场还考虑受体孤对电子与氢键间的夹角ω,在LJ12-10公式中增加一项cos4ω,而GRID程序中则应用一种依赖于孤对电子方向的6-4函数:式中,m通常设定为4。考虑截断的势能公式为:一般设定,θcut=90°。目前,使用的多数力场并未包括精确的氢键势能项,而是依靠静电和范德华作用来重现氢键,如pcff和COM-PASS力场声称可以计算氢键作用,实际上是通过综合考虑静电和范德华作用来实现的。3计算机模拟3.1氢键作用研究的对象以两个乙醇分子间简单的氢键作用为例,主要考虑两乙醇分子间O-H…O的范德华作用和静电作用。3.2力场模块的划分为考察不同力场是否能反映和计算氢键的问题,选用MS中的两个力场pcff(polymerconsistentforcefield)和cvff(consistent-valenceforcefield)以及独立的力场模块COM-PASS(Condensed-phaseOptimizedMolecularPotentialsforAtomisticSimulationStudies)。力场文件中,不同力场对各原子的电荷、杂化类型及力场参数等的考虑列于表1中。3.3岩石能量变化以H和受体0间的距离rHA及氢键夹角θ为变量,范德华能Evdw、库仑作用Ecoul及总非键能Et=Evdw+Ecoul为应变量,计算能量随距离及角度的变化。计算非键能时固定羟基的O-H键长(3种力场中分别为0.977、0.962和),考虑受体O原子与H及给体O原子间的非键作用之和。4结果与讨论4.1分子的构象优化分别采用pcff、cvff和COMPASS力场,再用Discover模块中的Minimizer工具,将两个乙醇分子形成的分子对进行构象优化,结果见表1。3种力场优化后的rHA,均介于1.7～之间,远小于,说明3种力场均允许氢键的存在。优化后的氢键夹角θ均接近180°,十分合理。从总非键能的大小分析,pcff和COMPASS力场中的氢键,应较cvff中的稳定得多,在后面的结果和讨论中,将进一步论证。4.2pc杰夫培训场4.2.1维能量图的形成采用pcff力场,计算不同的rHA和θ下的总非键能Et,得到的结果以三维形式见图2。从图中可看出,三维能量图出现一明显的波谷,当氢和受体氧间的距离及形成的角度处于波谷内时,则形成的氢键较稳定。为了能观察到范德华作用能和静电作用能(库仑力),对总非键能的贡献及对氢键形成稳定性的影响,下面将非键能进行分解。4.2.2能量平衡距离从不同的角度,以范德华作用能对rHA作图,呈典型的LJ曲线(图3),但由于范德华作用势,采用了LJ9-6的形式(表1),故短距离排斥力随距离变化而平缓。将rHA从2.2至的区域放大,可看出,当角度变小时,势能阱处的“平衡距离”(能量最低点对应的rHA)逐渐增大,从2.62增至,但势能阱深度(能量最小值)则几乎不变(从-0.233到-0.242kcal/mol)。说明范德华作用势可体现出角度因素对氢键稳定性的影响:角度越大,当H原子和受体氧原子接近,且两原子间距离越短时,范德华作用能越不起作用,形成的氢键越稳定。4.2.3rha的能量从不同的角度,以库仑力对rHA作图,基本呈双曲线形状(图4)。结果显示,在任何角度,库仑作用能均随rHA的增加而增大,但由于库仑力为氢原子和受体氧原子间的静电引力,及给体氧和受体氧原子间的静电斥力之和,所以,当角度很小时(110～118°)、发现当rHA介于之间,其能量最大(图4中放大部分)。对较弱的范德华作用,rHA较小时具较强的库仑作用,应该是氢键稳定的主要原因。4.2.4rha有氢键且间距在不同角度下以总非键能对rHA作图,呈类似LJ曲线形状(图5),说明总非键能的变化基本上决定于范德华作用。从图中可得出以下几点结论:(1)各个势能阱的深度,与角度成正比,角度越大,势能阱越深,非键“平衡距离”越小,氢键越稳定。(2)当角度小于约130°时,总非键能均为正。假定氢键形成的条件为:总非键能小于零,则只有当θ>130°时,才有可能“形成”氢键。从图中可以看出,角度越大,可能形成氢键的rHA的范围越宽。当角度为180°时,按假定,原则上在(最大值与由pcff力场中力常数计算出的非键“平衡距离”大体相当)之间均能形成氢键;但当rHA超过氢氧原子的范德华接触(pcff力场中氢、氧原子的范德华半径分别为1.2和)时,此时的氢键无意义。因此,图示的氢键能稳定存在的距离rHA和角度θ范围,分别为和130～180°,比较合理。角度和距离间具体的对应关系,可通过定义一“临界角”加以说明。(3)若将氢原子和受体氧原子间距离处于某一rHA时,形成氢键(按(2)中的假定)所需的最小角度定义为“临界角”θC,则θC为rHA函数(图6)。图中阴影部分的区域表示当rHA和θ处于该区域时,按假定可以形成较稳定的氢键。可以看出,这一距离和角度范围还是相当宽,说明pcff力场确能较好地反映氢键作用。(4)若在不同rHA下以总非键能对角度θ作图,也能得到类似的结论(图略)。(5)以上结果表明,通过将氢键作用视为氢及给受体氧原子间的范德华作用和库仑力的综合体现,pcff力场确实反映好的氢键作用。此作用已经在本课题组的前期工作中得到证实:应用pcff力场和分子动力学模拟,研究环糊精和对甲基苯酚所形成的主客体包结物在纳米水滴中的性质,确实观察到较强的水分子间、水分子或客体与环糊精的端口羟基间等的分子间氢键以及环糊精的端口羟基间的分子内氢键作用。(6)虽然在势能曲线中,并未在观察到突变,事实上也不可能出现突变,但若希望力场能更为精确地反映氢键作用,尚需进一步调节力场中相关的参数(详见COMPASS力场部分)。4.3总非键能曲线(1)采取IJ12-6的经典形式计算cvff力场中的范德华作用能,力常数与pcff有所不同(表1)。计算结果表明,能量曲线形状与pcff基本相同,但相同角度下的势能阱“平衡距离”较pcff明显偏小,为,而势能阱深度几乎相同且完全不变化(均为-0.228kcal/mol)。(2)由于cvff中对各原子特别是氢键相关的氢和氧原子的部分电荷分配不同,导致库仑力在范围内始终为负值,即静电吸引作用始终占据主导地位,这与pcff中库仑力随着rHA的增大逐渐由静电吸引转化为排斥的趋势明显不同。这一不同使总非键能曲线的形状明显改变。(3)总非键能曲线的变化趋势和pcff类似(图7),非键“平衡距离”同样随角度减小而变大。不同的是,rHA经非键“平衡距离”变化至时,总非键能始终为负数,且能量变化较平缓,因此总非键能曲线未能显出氢键的特征,所成的氢键很不稳定或说不能很好地反映氢键,这在我们的前期工作中已得到证实。(4)以上结果表明,若将氢键作用视为氢及给受体氧原子间的范德华作用和库仑力的综合体现,cvff力场并不能很好地反映氢键的作用,这与软件销售商在其力场介绍中所声称的相同。4.3.2comper力场的表现(1)COMPASS力场与pcff力场均是基于cff91力场,非键能函数形式相同,只是力场参数和原子电荷分布有所不同。为了更准确地描述氢键作用,COMPASS力场对涉及形成氢键的氢原子和氧原子的力场参数和原子部分电荷作了一些调整,使计算结果更为合理。(2)由于采用和pcff相同的非键势函数形式,COMPASS范德华作用能的曲线和pcff非常相似,不过为了突出静电对氢键的作用,COMPASS力场选用较小的力常数ε(表1)使范德华作用变弱,势能阱明显变浅:能量最小值由pcff时的约-0.24kcal/mol增加至COMPASS时的-0.095kcal/mol。(3)库仑作用能曲线呈现和pcff相似的趋势:随着rHA增大,静电作用由吸引转化为排斥,临界rHA在130°时为,180°时为;角度较小(110～130°)时,能量曲线也出现极大值。但由于COMPASS力场中计算得到的原子部分电荷分布,与pcff中有所不同(表1):羟基中氢原子的部分电荷数值稍小,而氧原子的部分电荷数值稍大,所以在同一角度时的临界rHA,比pcff(图4)小。下面的讨论显示,这一调整显得非常重要。(4)从图8a中的总非键能曲线和图8b的临界角分布图中可见,由于COMPASS力场适当调整了某些非键参数,使曲线在大角度时,氢键形成区域呈现得更为合理:①只有当角度大于138.5°时才能形成氢键;②角度越大,能形成氢键的rHA的范围越宽,可以形成氢键的区域较宽也较合理;③180°时原则上在间,均能形成氢键,此最大值正好为范德华接触距离,也与公认的形成氢键的临界距离rHA完全一致,显示出:以本文4.2.4(2)中氢键形成的假定为基础,COMPASS力场较pcff和cvff力场能更好地反映氢键作用。5模拟结果及分