《计算化学》-计算化学-第8章

第8章分子模拟分子模拟是随着计算机科学的发展逐步形成的一门新兴学科。它以量子化学、统计力学为理论基础.以原子水平的分子模型来模拟分子的结构与行为.进而模拟分子体系的各种物理与化学性质。分子模拟不但可以模拟分子的静态结构(图8.1).也可以模拟分子的动态行为.如蛋自的折叠与去折叠行为(图8.2)。因而该方法能使一般的实验化学家、物理学家和生物学家在电脑屏慕上直观地分析分子的结构及其变化。分子模拟方法也被称作“计算机实验”(计算机+分子模拟软件).即在硅片上完成的实验("InSilico".这与传统的体外实验("InVitro")以及体内实验("InVivo")形成鲜明的对比。随着理论方法和计算技术所取得的成就.以及计算机的飞速发展.分子模拟取得了引入注目的进展.已经成为化学、物理、生物以及材料研究中的有力工具.是人们除实验与理论研究之外.了解和认识微观世界的“第三种方法”(图8.3)。返回8.1分子模拟基本方法8.1.1量子力学量子力学(QuantumMechanics,QM)的基本模型是将原子核和电子看做在势箱中运动的粒子.其状态采用波函数Ψ进行描述.通过一系列简化求解薛定谔(Schrödinger)方程:或其中下一页返回8.1分子模拟基本方法得到能量E和波函数Ψ.从而得到分子的电子结构(包括静电势、电子百密度分布等).而电子结构决定了分子的化学性质。尽管量子力学的出发点是求解薛定愕方程的近似解.但目前有数种不同的量子力学理论和计算方法用于分子体系的计算.如分子轨道理论(MolecularOrbital,MO)、休克尔分子轨道理论(HükelMolecularOrbital,HMO)、从头计算法(AbInitio)以及各种半经验CSemi-empirical)分子轨道的计算方法等.后两者是最为普遍应用的两种方法。上一页下一页返回8.1分子模拟基本方法1.从头算法(AbInitio)量子力学从头计算方法是只利用普朗克常量、电子质量、电量三个基本的物理常数以及元素的原子序数.而不借助于任何经验参数计算体系全部电子的分子积分.求解薛定愕方程。其计算是建立在:①非相对论近似;②玻恩一奥本海默(Born-Oppenheimer)近似(或称绝热近似);③轨道近似(或称单电子近似)三个基本近似基础上的。除此之外.不借助任何参数.因此取名“从头”(AbInitio.拉丁文)。薛定愕方程在引入以上三个近似后称为Hartree-Fock-Roothaan方程(自洽场Self-ConsistentFi1c1d,SCF理论).求解这一方程是从头计算法的核心问题。进行从头计算首先要选定一组适合的原子轨道基函数.最常见的是STO-GTO系基函数:STO即Slater型函数.具有以下形式。式中，ζ为轨道指数;Yim(Θ,φ)护为类氢函数的角度分布;n,l,m对应主量子数、角量子数和磁量子数;r,Θ,φ为球极坐标。上一页下一页返回8.1分子模拟基本方法GTO即Gauss型函数.用直角坐标表示.具有以下形式:其中要注意的是.这里的l,m,n不同于Slater函数中的量子数.为非负整数;a为轨道指数。当l=m=n=0时.即为is型Gauss函数:上一页下一页返回8.1分子模拟基本方法而当l,m,n分别等于1.而其他均为0时.即分别为2px、2py、2pz型Gauss函数:STO-GTo基组是以STO作为自洽场的基函数.而每个STO用若干个GTO来逼近.解出的分子轨道仍用STO的线性组合来表示.而GTO不作为原子轨道.仅是中间过程的数学工具。常用的三种STO-GTO系基组为①STO-NG型极小基组:每个占据的原子轨道只对应一个STO，而每一个STO又用N个GTO来逼近.N一般取1~6.如STO-3G.表示采用三个高斯函数来描述Slater轨道。上一页下一页返回8.1分子模拟基本方法②N-31G型双ξ混合基组:每一个原子内层轨道用N个GTO组合成单ξ的STO基函数表示.而每一个原子的价层轨道则分成指数ξ'和ξ''的两个STO基函数.分别用N'个GTO和N''个GTO组合来表示。例如常用的6-31G基组.第一个数字6表示用6个Gauss型函数逼近内层的一个STO轨道.数字31的意思是用3个Gauss型函数逼近STO内轨.另外用1个Gauss型函数拟合STO外轨。③双ξ扩展基组:即在双ξ基组的基础上.再增加极化函数来描述原子轨道。如用p型函数描述s轨道的极化.用d型函数描述p轨道的极化.这就是6-31G(d)(或6-31G*)和6-31G(d,p)(或6-31G**)等基组的由来。如6-31G(d,p)表示在6-31G基组上.对于重原子增加了d轨道的成分.同时在氢原子轨道中加入了p的成分。上一页下一页返回8.1分子模拟基本方法三种近似的引入给自洽场从头计算带来一定的误差.这些误差主要集中在相对论近似误差和轨道近似误差(即电子相关能的误差)上.引起体系能量的升高。相对论近似误差的校正可采用微扰法(MP)或相对论Hartree-Fock理论.而电子相关能的校正主要通过组态相互作用(CI)或自旋非限制性Roothaan方程。上一页下一页返回8.1分子模拟基本方法2.半经验法(semi-empirical)从头计算法是一种精确的计算方法.但其缺点是计算量大.以致一般的计算机无法完成.其应用范围受到了限制。因此.量子化学工作者不得不寻求进一步的近似计算.即分子轨道的半经验计算方法.借助经验参数来简化从头计算法。通常所说的半经验量子化学计算方法主要指的是建立在零微分重叠(zerodifferentialoverlap,ZDO)近似基础上的计算方法.如MINDO,MNDO以及AM1和PM3等.均是目前很常用的半经验方法。它们比较适合于有机化合物的分子计算.可以得到定性的、有意义的计算结果.但对于无机化合物.尤其主族兀索.即使最好的半经验方法也难以得到好的结果。上一页下一页返回8.1分子模拟基本方法①MINDO(modifiedintermediateneglectofdifferentoverlapmethod):由Dewar在1969年提出.这种方法保留了全部单中心排斥积分.而对双电子积分进行参数化。依据参数化的不同.可以分为MINDO/1,MINDO/2,MINDO/3等.其中MINDO/3对于分子基态计算较为成功。②MNDO(modifiedneglectofdiatomicoverlapmethod):由Dewar在1977年提出.是第一个建立在基于忽略双原子微分重叠(neglectingdiatomicdifferentoverlapmethod,NDDO)基础上的方法。研究表明.多数基态性质的平均绝对误差MNDO比MINDO/3大约可以减少一半.所以前者比后者有更广泛的应用范围。上一页下一页返回8.1分子模拟基本方法③AMl(Austinmodel1):由Dewar在1985年提出.是建立在MNDO方法上的另一种方法。AM1采用高斯函数对MNDO中核一核排斥函数进行了修正.消除了当原子分开到接近VanDeWaals半径时原子间斥力估计过高的缺点。同时.AM1采用大量的实验数据来进行参数量化.因此比MNDO计算结果有显著的改善.主要表现在:a)在氢键处理上.AM1明显优于MNDO;b)AMl对于反应活化能的计算明显好于MNDO;c)对于高价磷化合物.AM1的计算结果与MNDO相比有一定的改进;d)一般AM1计算的生成热比MNDO的计算结果的误差低40%左右。上一页下一页返回8.1分子模拟基本方法④PM3(Parametricmethod3):由Stewart在1989年提出.是基于MNDO(其命名由第三个参数化MNDO而来.AM1被认为是第二个)一种新参数量化的方法.其哈密顿算符与AM1相同.但参数使用了JJPStewart发明的自动参数化程序.对于每一个非氢原子有由18个参数组成的一个参数集.而对于氢原子只有11个参数。PM3方法通常用于分子结构的全参数优化。从测评的结果看.PM3比AM1有一定的改进.表现在:a)PM3计算的生成热误差要小于AM1;b)PM3在处理高价态化合物时优于AM1。上一页下一页返回8.1分子模拟基本方法总之.以上半经验算法通过对大量标样分子的计算.选择合适的半经验公式和优选其中的参数.可以给出相当可靠的计算结果。AM1和PM3是目前应用最广泛的两种半经验量化计算方法。这些方法可以用于下列体系:①非常大的体系.如对于蛋自质等生物分子;②对于一些半经验方法可以得到很好结果的体系;③进行一些分子结构的定性研究。但是.半经验方法也有它的局限性.它只能处理已经有了良好参数的体系。同时.计算产生的误差随意性大.使得结构差异很大的体系依据半经验方法计算的结果来进行性质比较时.往往可靠性不高。上一页下一页返回8.1分子模拟基本方法3.密度泛函理论(DensityFunctionalTheory,DFT)1964年，Hohnberg和Kohn证明了分子基态的电子能量与其电子密度有关.提出了体系基态的电子密度分布完全决定体系性质的定理.奠定了现代密度泛函理论的基础。DFT属于第三类电子结构理论方法.它不用波函数而是通过电子密度来表征多电子体系.其能量可以表述为电子密度的泛函;用变分法对Schrödinger方程进行求解.求得能量的极小值即可确定电子密度。由于DFT融入了统计的思想.不必考虑每个电子的行为.只需计算总的电子密度.因而计算量大大减小。1998年.Kohn因其在DFT的开创性工作与另一位著名量子化学家Pople一起获得该年度的诺贝尔化学奖。上一页下一页返回8.1分子模拟基本方法根据Hartree-Fock理论.系统能量形式为其中.v为核排斥能;P为密度矩阵;<hP>为单电子动能加势能;1/2<PJ(P)>为电子的经典库仑排斥能一1/2<PK(P)>为来自电子的量子(费密子)特性的交换能量。在密度泛函理论中.单行列式的精确交换能量(HF)被一更广泛的公式代替.即交换一相关泛函.该泛函可包含被Hartree-Fock理论省略的电子交换和电子相关能量:其中.Ex[P]为交换泛函;E[P]为相关泛函。上一页下一页返回8.1分子模拟基本方法Hartree-Fock理论是密度泛函理论的特殊情况.Ex[P]由交换积分-1/2<PK(P)>定义.且Ec=0.密度泛函中通常使用的泛函为电子密度和可能的电子密度梯度的某种函数的积分:其中不同的密度泛函的Ex[P]和Ec[P]使用不同的函数f。由于所有的相关泛函必须和交换泛函一起使用.所以一般看到的密度泛函都包含这两部分·如BLYP表示Becke88交换泛函与Lee-Yang-Parr相关泛函进行组合而构成的密度泛函。同时.还有多种由Hartree-Fock交换和DFT交换一相关混合形式的混合泛函。其中.最著名的是Becke于1993年提出的三参数泛函.其能量表达式为上一页下一页返回8.1分子模拟基本方法其中.参数A，B，C是Becke通过拟合G1计算结果确定的常数。类似的比较常用的混合泛函还有B3LYP,B3PW91和G961LYP等。DFT方法由于其相对较高的计算精度和比其他从头计算方法快得多的计算速度.而备受研究者的青睐.其应用的范围日益广泛。但是DFT发展较晚.远不及基于HF和相关能校正的分子轨道法完善.本身尚存在很多限制。例如.从原理上来说.DFT属于仅适用于基态的理论方法·无法应用于激发态研究·近年来发展的含时DFT(Time-DependentDFT,TD-DFT)在单重激发态的研究中取得了实质性的进展。在DFT的发展完善过程中.核心的问题是构建具有普适性、非经验性的交换一相关密度泛函.或通过大量模型化合物优化确立能量泛函中的参数.以对中等和较大体系得到满足化学和物理精度要求的计算结果。上一页下一页返回8.1分子模拟基本方法8.1.2分子力学分子力学(MolecularMechanics,MM)方法.是基于经典牛顿(Newton)力学方程的一种计算分子的平衡结构和能量的方法。与量子力学不同.它求解的是Newton方程.而不是Schrodinger方程。其计算的基本依据是分子中原子之间的作用可以分解为化学键和非键相互作用。计算的核心是构筑力场(Forcefield).把分子处理为一套由经典力学势能函数支配的原子排列.即该势能是原子坐标的函数.因而可以通过分子力学计算和能量最低化操作.对化合物的空间结构或构象加以确定。分子力学忽略电子运动.而把体系能量看作仅仅是原子核坐标的函数.因此可以用于包含有大量原子的分子计算。上一页下一页返回8.1分子模拟基本方法1.分子力场分子力场是分子力学的核心。分子力学的基本理论就是一个分子力场由分子内相互作用和分子间相互作用两大部分构成.即力场的势能包括成键和非键相互作用.所有势能的总和即为分子的构象能。由于分子力学是经验的计算方法.不同的分子力学方法会采用不同的势能函数(PotentialEnergyFunction,PEF)表达式·而且其力场参数值也会不同。一般将分子的PEF分解成键伸缩能Ec(stretchingenergy)、键弯曲能Eb(bendingenergy)、二面角扭转能Et(torsionalenergy)、范德华作用能Ev(vanderWaalsenergy)和静电作用能Ee(electrostaticenergy)等项(如图8.4所示).总能量E总可以表示为其中.前三项为成键作用能.后两项为非键作用能。对于更精确的分子力场.可在以上PEF表达式中增加一些交叉项.如键伸缩一键弯曲相互作用项.也可加入氢键函数项等其他修正项。上一页下一页返回8.1分子模拟基本方法不同的力场包含有不尽相同的能量项.例如:CHARMm22力场的势能函数表达为式中.第一~第三项分别为键的伸缩能、键角的弯曲能和二面角扭转能;第四项为非正常二面角(impropertorsion.平面外弯曲)能量项;第五项Urey-Bradley为1,3非键间的相互作用;最后两项分别为范德华相互作用和静电相互作用。上一页下一页返回8.1分子模拟基本方法又如AMBER力场的函数表达式为式中.前三项表示键、角、二面角的变化引起的势能改变值;第四项表示范德华相互作用;第五项表是静电相互作用;最后一项表示氢键相互作用。上一页下一页返回8.1分子模拟基本方法2.常见的分子力场简介较早出现的分子力场是Allinger等于1973年提出的MM1和MMP.前者用于计算非共扼体系的有机分子.后者用于计算共扼体系分子。后来发展的MM2力场为MMP的改进.扩大到含有离域二电子体系分子.可以计算线形、平面和三维力场;MM3除了可计算以上力场外.还能计算部分化学影响因索.可以用于处理胺、多肤和蛋自质;MM4作了进一步的改进.能较好的处理经典力场和电负性或超共扼等化学影响.其重点集中在烃、碳碳键长的超共扼效应和共扼烃.其中特别重点在振动频率上。随着计算机科学的迅速发展.分子力场也随之得到很大的发展。目前.常见的力场(适用体系)与网址如下。上一页下一页返回8.1分子模拟基本方法①MM2/MM3/MM4(有机小分子):http://europa.chem.uga.edu/allinger/mm2mm3.html②UFF,UniversalForceField(主族元素化合物、有机分子、金属配合物):http://franklin.chm.colostate.edu/mmac/uff.html③AMBER,AssistedModelBuildingwithEnergyRefinement(蛋自质核酸):http://amber.Scripps.edu④OPLS,OptimizedPotentialsforLiquidSimulation(蛋自质核酸):http://zarbi.chem.yale.edu⑤CHARMm,ChemistryatHARvardMacromolecularmechanincs(生物大分子):http://yuri.harvard.edu上一页下一页返回8.1分子模拟基本方法分子力学中有三个重要的概念:函数形式、原子类型和参数组。函数形式即为力场.用数学函数及其参数来描述分子中各原子之间的相互作用。每一种力场包括一套原子类型.划分原子类型是为了更好地使用参数来描述体系的性质。每个力场可以包含多套参数。其中参数可以通过拟和量子力学计算结果以及相关实验数据得到。对于分子力场参数文件一般为.INP格式文件.可以用写字板将其打开。例如:CHARMm22力场的具体参数文件如下。上一页下一页返回8.1分子模拟基本方法上一页下一页返回8.1分子模拟基本方法上一页下一页返回8.1分子模拟基本方法上一页下一页返回8.1分子模拟基本方法上一页下一页返回8.1分子模拟基本方法上一页下一页返回8.1分子模拟基本方法上一页下一页返回8.1分子模拟基本方法上一页下一页返回8.1分子模拟基本方法分子力学的优势是计算量小.比量子力学的计算可快几个数量级.因而在复杂体系中的应用尤为重要.可主要用来研究一些较大分子和柔性较大分子的构型和性质.同时它也是构象搜索、分子动力学模拟、MontCarlo模拟的基础。另外由于它实际上是在原子层次而不是在电子一原子核层次上考虑问题.所以也是连接量子力学和统计力学计算的重要桥梁。然而分子力学也有局限性.主要是不能够提供分子中基于电子分布的相关性质.无法描述化学键的断裂和生成。分子力学的最大优点是计算速度快.其缺点是在构象优化时一般只能达到势能面上最近的能量最低点(局部最低点)。但是.从现在分子力学的广泛应用范围和处理实际问题的快速方面看.仍不失为一种实用的构象分析方法。上一页下一页返回8.1分子模拟基本方法8.1.3QM/MM法近年来人们发展了一种量子力学与分子力学结合的方法(QM/MM法)来研究一些较大体系的性质。例如.先利用分子力学优化结构得到稳定构型.再用量子力学方法计算单电能研究电子的结构性质;或者用现有的QM/MM方法.对核心部位(局域)采用量子化学从头计算方法.对其他部位采用分子力学计算。该方法既包括了量子力学的精确性.又利用了分子力学的高效性。可以看出.QM/MM法的基本思想是:用量子力学处理感兴趣的中心.如酶与底物结合的活性位点.其余部分用经典分子力学来处理。量子力学方法可以是从头算法、半经验法或DFT法;分子力学可以采用各种不同的力场(如CHARMm,Amber,MM系列等)。关键是量子力学与分子相连接部分的处理.目前已发展了多种方法来解决这一问题.如所谓的“连接原子”方法、伪势能与伪原子方法、固定边界的杂化轨道方法等。通过分子动力学或MonteCarlo方法.QM/MM法在研究溶液中的化学反应、生物大分子的结构、性质与功能方面已经取得了很大的成功.有望在分子模拟方面有更广泛的应用。上一页下一页返回8.1分子模拟基本方法以肌红蛋自(Myoglobin,Mb)为例(图8.5),Rovira等应用此方法研究了一氧化碳(CO)结合Mb中血红索空穴(homepocket)构象对home-CO键的结构与构象的影响。结果显示，home-CO结构(QM处理)相当牢固.并不受其远端空穴构象(MM处理)的影响.因而排除了FeCO的扭转与MbCO分子中所观测到CO不同红外吸收之间的相关性。上一页下一页返回8.1分子模拟基本方法8.1.4分子动力学分子动力学(Moleculardynamics,MD)是一种重要的统计物理方法.其基本思想是将分子中的每一个原子均看做牛顿运动方程中的基本质点或粒子.以某起始状态为起点.根据其受力情况依据牛顿第二定律:或上一页下一页返回8.1分子模拟基本方法通过求解牛顿运动方程得到体系中所有原子的轨迹.并从原子轨迹中计算得到分子的各种性质。因此.分子动力学的实质是求解所有原子的运动方程组:式中.ri为第i个原子的移动矢量;mi为第i个原子的质量;Fi为第i个原子所受的作用力;t为时间。上一页下一页返回8.1分子模拟基本方法求解Newton方程通常使用有限差分方法.其基本思想是把整个积分距离分成一个个小的时间段(δt).由于δt通常在飞秒(10-15s)级.所以在每个时间段中粒子的受力被认为是不变的。如果知道t时刻粒子的位置r(t).速度是v(t).那么(t+δt)时刻粒子的位置和速度就可以用泰勒展开式表示:上一页下一页返回8.1分子模拟基本方法其中加速度a(t)可以根据力求算。采用有限差分方法的分子动力学计算可以按照以下步骤进行:①设定粒子的初始位置和速度;②根据粒子的位置计算每个粒子的受力;③根据粒子的位置、速度和受力计算粒子的新位置和新速度;④更新粒子的位置和速度.然后回到第②步。以上各步计算的重要性大小排列依次为积分法(第③步)、初始化(第①步)和受力的求算(第②步)。所谓积分法.就是根据粒子的位置、速度和受力计算δt时刻后粒子的新位置和新速度。目前.MD广泛使用的一种积分法是1967年Verlet提出的Verlet算法:上一页下一页返回8.1分子模拟基本方法此外.还有1970年Hockney等人提出的Leap-Frog算法:Leap-Frog算法与Verlet算法相比有两个优越之处:一是计算时没有δt2项.保证了计算的精度;二是计算速度快。上一页下一页返回8.1分子模拟基本方法初始化工作包括两部分:一是初始化粒子的坐标(即生成初始构象);二是初始化粒子的速度。对于初始构象的生成过程一般对实验数据(如蛋自质机构数据库.ProteinDataBank,PDB)进行预处理.通常包括重新加氢、加抗衡离子、加入水环境等。如图8.6所示为脱辅基(血红索.home)细胞色索b5(apocytochromeb5,apocytb5,PDB编码1IET)置于周期性水分子模型(TIP3P)中的蛋自构象示意图。速度的初始化主要依据Maxell-Boltzmann分布随机生成.在一定温度T下的速度分布如下:式中.P(v)表T一个质量为m的粒子.在温度T的情况下在x方向上的速度为v的概率。上一页下一页返回8.1分子模拟基本方法对于粒子受力的求算一般生物大分子的分子动力学计算中.95%以上的时间都是在计算分子的受力.可见这一过程的重要性。按照经典力学理论一个粒子的受力等于体系势能函数对粒子坐标的一阶导数.而各种势能函数在分子力学中已有描述.其中势能是由力场的形式和参数所决定。MD中所使用的模拟系统分为:①NVE系统(ConstantNumber,VolumeandEnergy);②NVT系统(ConstantNumber,VolumeandTemperature);③NPT系统(ConstantNumber,PressureandTemperature)。上一页下一页返回8.1分子模拟基本方法在NVE系统中.模拟过程中体系的总能量不变.而温度可变.因而可以模拟分子中原子在不同温度下的运动状态.从而可以研究给定温度下分子构象的变化.这在研究蛋自质分子的折叠与去折叠方面有着广泛的应用。MD方法的威力在于它可以跨越较大的能皇.因此.可以通过升温搜索很大的构象空间.尽可能找到真正能量最低的构象。缺点是难以产生大的拓扑重排.因此不适用于分子过大的体系.因为对分子过大的体系.即使很小的能皇.跨越的时间也要超出模拟时间.得不到真实的结果。而最常用的MD模拟系统还是NVT和NPT系统.前者温度和体积保持不变.后者温度和压力保持不变。一般化学反应都是在恒温条件下完成的.如蛋自质的折叠与去折叠.以及高分子的玻璃化相转变等.这些过程用恒温动力学模拟更加合理。上一页下一页返回8.1分子模拟基本方法通过MD模拟.得到的是一组相空间的轨迹.即在某一时刻/体系中每个粒子的三维坐标(x,y,z)和速度(v)。轨迹文件所能提供的重要信息就是体系结构随时间的变化。从轨迹文件中可以得到分子构象的变化过程(蛋自质二级结构α-helix,β-sheet的变化等)、分子结构与某一构象的均方根偏差(RootMeanSquaredDeviation,RMSD)的变化、在某一时间段的平均结构等。例如.蛋自Apocytb5在298K,1atm下的骨架Cα原子与其初始结构(1IET)间的RMSD随分子模拟时间的变化关系如图8.7所示。其中.RMSD的定义如下:上一页下一页返回8.1分子模拟基本方法其中.Nα是用于和其位置相比较的一系列原子;Nt是用来比较的一系列时间步骤内的原子位置;ra(tj)是α原子在时间tj时的位置;而<rα>表示α原子与rα(tj)位置间的平均值。上一页下一页返回8.1分子模拟基本方法第一个报道的MD模拟工作是Alder和Wainwright在半个世纪前(1957年)完成的.他们用硬球势来模拟粒子之间的相互作用(即粒子之间的碰撞是完全弹性碰撞).两次碰撞之间的粒子做匀速直线运动。由于Alder用的模型过于简单.并没有提供什么新的有用信息。1964年.Rahman用MD模拟了氢气体系;1971年Rahman又和Stilliger合作模拟了水分子体系.得到了一些有意义的结论。而第一个将MD应用于模拟生物大分子.则是由哈佛大学(HarvardUniversity)的Karplus教授在30多年前完成的.这一工作于1977年发表于Nature。尽管当时只是在真空体系中模拟了8.2ps.却为我们展现了一个运动的蛋自分子体系.从而改变了人们传统认为蛋自分子是相对刚性结构的思维模式。随着计算机水平的飞速发展.目前MD可以模拟蛋自质的水溶液体系长达1μs，进而使得从原子水平上研究蛋自质的结构、性质与功能成为现实。上一页下一页返回8.1分子模拟基本方法8.1.5分子蒙特卡洛法蒙特卡洛(MonteCarlo,MC)法因利用“随机数”对模型系统进行模拟以产生数值形式的概率分布而得名.美国数学家冯·诺伊曼(vonNeumann)和乌拉姆(Ulam)在1945年引进这个用语.将其定义为“使用随机数处理确定性数学问题的方法”。MC方法也称为随机模拟Crandomsimulation)方法.或称为随机抽样(randomsampling)技术或统计试验Cstatisticaltesting)方法。它的基本思想是.为了求解数学、物理、化学、生物及工程技术等方面的问题.首先建立一个概率模型或随机过程.使它的参数等于问题的解.然后通过对模型或过程的观察或抽样试验来计算所求参数的统计特征.最后给出所求解的近似值.其中解的精度可用估计值的标准误差来表示。上一页下一页返回8.1分子模拟基本方法假设所要求的量二是随机变量的数学期望E(ξ).那么近似确定二的方法是对ξ进行N次重复抽样.产生相对独立的ξ值的序列ξ1,ξ2,...,ξN,并计算其数学平均值根据柯尔莫戈罗夫(Kolmogorov)加强大数定理有因此.当N充分大时.下式成立的概率等于1.亦即可以用作为所求量x的估计值。上一页下一页返回8.1分子模拟基本方法MC方法与一般计算方法的主要区别在于:它能比较简单地解决多维或因索复杂的问题。其突出特点是.它的解是由试验得到的.而不是计算出来的.不像常规数理统计方法那样通过真实的实验来解决问题.而是抓住问题的某些特征.利用数学方法建立概率模型.然后按照这个模型所描述的过程通过计算机进行数值模拟实验.以所得的结果作为问题的近似解。因此.

MC法是数理统计与计算机相结合的产物。如果所要求解的问题是某种事件出现的概率.或者是某个随机变量的期望值.就可用MC法得到这种事件出现的概率.或者以这个随机变量的平均值作为问题的解。这就是MC法的基本原理。其解题过程一般分为三个步骤:建立概率模型、产生随机数.并根据概率模型进行计算、统计计算结果。上一页下一页返回8.1分子模拟基本方法在化学中.MC方法常用来计算一个分子或分子体系的平均热力学性质.并可扩展到研究分子结构以及液体、溶液的平衡性质.也可用来在非平衡态的条件下进行构象搜索。同时.由于生物大分子链由大量的重复单元构成.聚合反应存在着随机性。高分子化合物中分子量的大小分布、共聚物中的序列分布、蛋自质与DNA分子的构象变化、折叠与去折叠.都存在着随机性问题.这些都无疑成为MC法研究的最佳对象。上一页下一页返回8.1分子模拟基本方法以蛋自质的折叠为例，Okamot。等应用退火-淬火(annealing-quenching)的MC方法模拟了特殊的蛋自分子片断BPTI(16-36)的结构与折叠的过程。这种蛋自分子的氨基酸序列为Ala-Arg-His-His-Arg+-Iyr-Phe-Tyr-Asn-Ala-Lys+-Ala-Gly-Leu-Cys-Gln-Thr-Phe-Val-Tyr-Gly(其中残基Arg与Lys带有正电荷)。模拟过程中将系统的温度升至1000K.再逐渐降至250K。模拟由随机选取的结构开始.重复这样的退火一淬火过程20次。图8.8为模拟的结构与X射线衍射实验的比较。图中显不在水溶液中模拟的最优化结构形成R折叠状·与观测的结构一致。但真空系统的模拟则仅形成弯曲转折(turn)，而不形成β折叠的结构。图8.9显示由随意选取的结构开始.逐步形成S结构的过程。图中虚线表示形成氢键的残基对。综合实验结果，Okamot。等认为BPTI蛋自质形成R折叠结构的原因主要是受溶液中水分子的影响。这些水分子迫使BPTI分子的带电侧链相互吸引.形成氢键.而导致分子折叠。上一页返回8.2分子模拟软件简介从20世纪旱期开始.化学家就试图通过建立分子模型来表达分子的结构信息.通过分子模型可视化技术.使一些复杂的信息变成容易并直接为人们所理解的图形结构(二维或三维图形模型)。对于分子模型的可视化.可以分为分子结构可视化和分子属性可视化。对于小分子或晶体分子的结构模型一般有金属线结构模型(WireFrame)、棍状模型(Sticks)、球棍模型(Ball&Sticks)和空间填充模型(SpaceFilling).如图8.10所示;对于生物大分子的结构模型一般有球棍模型(Ball&Sticks)、管状模型(Tube.表征分子骨架)、带状模型(Ribbon.表征蛋自分子的二级结构)以及CPK模型(Corey,Paining,Koltum模型.表征空间填充)等.如图8.11所示。下一页返回8.2分子模拟软件简介仅仅知道一个分子的空间结构对于详细了解一个复杂分子的内部作用来说并不充分。同时.还需要考虑一些分子的性质.如静电势、分子轨道、疏水性与亲水性等.这便要求分子属性可视化。所有这些分子模拟的实现均需要借助一定的分子模拟软件。目前.常见的分子模拟软件主要有ChemOffice,HyperChem,Gaussian,WinMOPAC,Gamess以及MaterialsStudi。等。本节着重介绍前三者.其他软件可以查阅相关Internet站点。上一页下一页返回8.2分子模拟软件简介8.2.1ChemOffice简介ChemOffice是由CarbridgeSoft公司开发的功能强大的化学绘图及分子模拟计算程序包.一般有Std,Pro,Ultra等多种版本.其功能依次增强.网址为http://www.cambridgesoft.comp目前常用的ChemOfficeUltra2004包含ChemDrawUltra8.0(化学结构绘图),Chem3DUltra8.0分子模型及仿真)、ChemFinderPro8.0化学信息搜寻整合系统)。此外还加入了E-NotebookUltra8.0;BioAssayPro8.0;量化软件MOPAC,Uaussian和GAMESS的界面;ChemSARServerExcel,CLogP,CombiChem/Excel等。ChemOfficePro还包含了全套ChemInfo数据库.有ChemACX和ChemACX-SC,Merck索引和ChemMSDX。当前采用的版本是ChemOfficeUltra2006(10.0)。我们常用的是其中2D绘图程序ChemDrawUltra10.0和分子模拟程序Chem3DUltra10.0。上一页下一页返回8.2分子模拟软件简介1.ChemDrawChemDraw是世界上最受欢迎的化学结构绘图软件.是各论文期刊指定的格式。可以绘制各种分子的化学结构与立体结构、分子轨道与电子云模型、纸层析不意图.以及实验室常见玻璃仪器等。其程序界面如图8.12所示。这里着重介绍几种常用的功能。①Templates:与以前的版本相比.ChemDrawUltra10.0又增加了很多的分子模板.使绘制分子图形更为方便。如增加了碳纳米管(图8.13(a))与环糊精(图8.13(b))等模型分子。上一页下一页返回8.2分子模拟软件简介②AutoNom:这是Bei1steiny最强的软件.现已内含在ChemDrawUltra内.它可自动依照IUPAC的标准命名化学结构。选定所绘的分子后选择Structure菜单中的ConvertStructuretoName命令.可给出分子的系统命名。相反.使用ConvertNametoStructure命令.可以将所输入的IUPAC化学名称自动产生ChemDraw结构。③ChemNMR:在ChemDraw内预测13C,和1H的NMR光谱.节省实验的花费。选定所绘的分子后选择Structure菜单中的Predict1H-NMRShifts或Predict13C-NMRShifts.程序会相应给出分子中各个H或C的化学位移。例如.咪III的1H-NMR光谱预测结果如图8.14所示。上一页下一页返回8.2分子模拟软件简介2.Chem3DChem3D同ChemDraw一样.是ChemOffice的组成部分.它能很好地同ChemDraw一起协同工作.ChemDraw上imp的二维结构式可以自动地在Chem3D中转化为正确的三维结构。它提供工作站级的3D分子轮廓图及分子轨道特性分析.并和数种量子化学软件结合在一起。由于Chem3D提供完整的界面及功能.已成为分子仿真分析最佳的前端开发环境。Chem3DUltra10.0版本不但包括分子力学MM2程序.还包含GAMESS计算模块.无需另外中请.并能与著名的量子化学程序Gaussian98/03和Mopac连接.同时还增加了量子化学计算程序Jaguar的图形接日.作为它们的输入、输出界面.能够以三维的方式显T量子化学计算结果.如分子轨道、电荷密度等。其程序界面如图8.15所示。上一页下一页返回8.2分子模拟软件简介这里着重介绍几种常用的功能①ChemPropPro:预测BP、MP、临界温度、临界气压、吉布斯自由能、lgP、折射率、热结构等性质。绘好分子结构后.如图8.15中的脯氨酸.选择菜单Calculatio二中的ComputePropoties...，在其PropertyPicker对话窗口中选择ChemPropPro.执行计算.结果如下所示。上一页下一页返回8.2分子模拟软件简介上一页下一页返回8.2分子模拟软件简介②MM2分子力学计算:能够进行常见的有机分子构型优化(如脯氨酸的优化构型如图8.10所示).计算结果以Output形式给出。上一页下一页返回8.2分子模拟软件简介另外，还可以进行MM2分子动力学模拟(Calculations/MM2/MolecularDynamics...),其对话窗口如图8.16所示.可以进行相关参数的设定。上一页下一页返回8.2分子模拟软件简介③分子结构数据的测定:使用Structure/Measurements/GenerateAllBondLengths(BondAngles,DihedralAngles)命令.可以分另of显不分子的几何参数键长、键角和二面角。例如.CH3CN在用MM2分子力学优化后的分子结构几何参数如图8.17所示。④分子轨道的计算与显示:以咪唑(Imdazole,Im)为例.先在ChemDraw中画出其化学分子式.粘贴到Chem3D中.使用Surfaces/ChooseSurface/MolecularOrbital(HuckelCalcutaion)命令进行表面计算.然后使用Surfaces/SelcetMolecularObital命令分别显示其最高占据分子轨道(HighestOccupiedMolecularOrbital,HOMO)(图8.18(a))与最低末占据分子轨道(LowestUnoccupiedMolecularOrbital,LUM())(图8.18(b)),同时可以看出其计算结果分别为-10.559eV和7.191eV。此外.也可以使用同样的方法按编号选择和显示其他的分子轨道图。上一页下一页返回8.2分子模拟软件简介3.ChemFinderChemFinder化学信息搜寻整合系统.可以建立化学数据库、储存及搜索.或搭配ChemDraw,Chem3D使用.也可以使用现成的化学数据库。ChemFinde:是一个智能型的快速化学搜寻引擎.所提供的ChemInfo是目前世界上最丰富的数据库之一包含ChemACX,ChemINDEX,ChemRXN,ChemMSDX.并不断有新的数据库加入。ChemFinder可以从本机或网上搜寻Word,Excel,Powerpoint,ChemDraw,ISIS格式的分子结构文件。还可以与微软的Excel结合.可连结的关联式数据库包括()racle及Access,输入的格式包括ChemDraw,MDLISISSD及RD文件。其程序界面如图8.19所示。上一页下一页返回8.2分子模拟软件简介8.2.2HyperChem简介HyperChem是HyperCubeInc.公司开发的一款以质量高、灵活、易操作而闻名的分子模拟软件。网址为。它可以构建复杂的分子模型.具有强大的综合计算和分析功能.是优秀的分子图形和分子设计的工具软件之一。HyperChem自带了很多计算方法.可直接进行量子化学(半经验、从头算和密度泛函)、分子力学、分子动力学、蒙特卡罗模拟等计算。其程序界面如图8.20所示:这里着重介绍几种常用的功能。①Display(分子的显示):以Cytb5为例.使用File/Open导入Cytb5的蛋自数据pdb文件后.再将HyperChem程序主界面的背景颜色改为自色:File/Preferences/WindowColor/white。然后选择Display/Rendering/RendringMethod.分别将AtomRendering选为，Tubes以及将SecondaryStructureRendering选为ThickRibbons。即可将Cytb5显示为如图8.21所示的结构模式。上一页下一页返回8.2分子模拟软件简介②Databases(数据库):使用此命令.可将残基库中的氨基酸和核酸以特定连接规则构成多肤和聚核茸等生物大分子。例如构建L-型a-Helix多肤(10-Ala)分子.先在标准氨基酸库的窗口(图8.22)中选择Conformation/Alphahelix及Isomer/L.连续点击10次Ala.程序会自动生成所要求的多肤.然后再使用DatabasesiMakeZwitterion命令.将其转化为两性离子.即将多肤的末端基团-NH和-Co转化为-NH3+和-COO-分别带上正负电荷。最后选择Display的模式.如AtomRendering选为Tubes;SecondaryStructureRendering选为ThinRibbons(图8.23)。同时.我们可以对所构建的多肤分子(或其他蛋自分子)进行氨基酸残基的突变。例如.要将10-Alat突变成N-端和C-端.即1,10位氨基酸残基Ala均突变成Phe.先使用命Sclcct/Residues.选择所要突变的氨基酸残基.再使用命令Databases/Mutate.在弹出的窗口中选择Phe.点击OK即可如图8.24所示。上一页下一页返回8.2分子模拟软件简介③Compute(计算):HyperChem能进行SinglePoint单点计算(即计算分子结构或被选定原子的总能量)、GeometryOptimization(几何构型优化)以及MD,MC分子模拟等。以优化dA(Adenine)分子为例.使用Databases/NucleicAcids绘制目标分子式.再使用Setup设置优化方法.我们选择MolecularMechanics,在弹出的窗口中选择AMBER.点击OK。然后使用Computer/GeometryOptimization命令.程序自动弹出的分子力学优化的参数设置窗口.我们使用默认的设定参数.即不作任何改动.直接点击OK。优化后的分子构型如图8.25所示.在界面的左下方会显示分子的总能量(-20.485183kcal/mot).右下方显示所采用的分子力学类型Amber99。上一页下一页返回8.2分子模拟软件简介此外.进行分子轨道(Orbital)、振动光谱(VibrationalSpectrum)以及电子光谱(Electronic印ectrum)的计算也非常重要。以NH3分子为例.要计算分子轨道.先要使用量子化学(如选择Setup/Semi-empirical/AM1)进行构型优化.然后再使用命令Computer/Oribitals计算其分子轨道.结果如图8.26所示;要计算其振动光谱.在完成优化后.再进行振动频率计算，使用命令Compute/Vibration,RotationAnalysis.最后再使用命令Compute/VibrationalSpectrum来显示其振动光谱(如图8.27(a));而要计算其电子光谱.则在完成优化后.进行单激发的组态相互作用(CI)计算.使用命令Compute/SinglePointCI.在弹出的窗口中选择SinglyExcited作为CI方法.以及OribitalCriterion指定被占有轨道和空轨道的数目。最后再使用命令Compute/ElectronicSpectrum来显示其电子光谱(如图8.27(b))。上一页下一页返回8.2分子模拟软件简介④MoelucularDynamisc(分子动力学):这一功能包含在Compute菜单内.下面以二肽Cys-His两性离子的MD模拟为例.熟悉这一方面的相关操作。a.首先按照在Databases中介绍的方法构建Cys-His.将其转化为两性离子.使用Build/SetCharge命令.将N-端Cys的N原子设定1个单位正电荷.为C-端His的两个O原子分别设定0.5个单位负电荷。b.使用Setup/MolecularMechanics.选择AMlBER方法.接着使用Build/AddH&ModelBuild,再使用Compute/GeometryOptimization.将RMSgradient设为0.1.进行真空环境的几何构型优化。c.使用Setup/Periodic13ox命令.在输入框中键入X,Y,Z的值分别为10,10,20,点击OK.即为二肤两性离子加上周期性水盒。上一页下一页返回8.2分子模拟软件简介d.使用Setup/MolecularMechanics.选择AMBER方法.将其Option选项中的Dielectric(介电常数)设定为Constant(常数)。使用Compute/GeometryOptimization(默认值)命令.进行周期边界条件的几何构型优化·优化收敛后的能量值约为-1282.6kcal/mol。上一页下一页返回8.2分子模拟软件简介8.2.3Gaussian简介Gaussian98/03就是目前国际上公认的功能最为强大、使用最为广泛的量子化学标准程序库。Gaussian软件是以1998年诺贝尔化学奖得主JohnPople为代表的研究集体经过30多年不断创新计算方法.在Gaussian70,76,80,82,86,88,90,92及94的基础上于1998年推出Uaussian98版本.2003年又推出了最新的版本Gaussian03,目前这两个版本的软件广泛地用于物理学、化学、生物学和材料学的研究中。其内容涵盖了分子力学、半经验方法、从头计算方法及密度泛函理论等精度要求各不相同的计算程序.并配备有GaussView图形输入/输出软件.可以在Unix,Windows及NT等操作平台上运行。其站点为http://www.gaussian.com。程序界面如图8.29所示。上一页下一页返回8.2分子模拟软件简介Window系统下的Gaussian03W程序的安装和运行相对简单.关键是准备输入文件和分析运行结果。以水分子的计算为例.来了解如何运用Gaussian程序进行计算并获得所需要的分子信息。对于初始输入文件的生成.利用Gaussian的配套软件GaussView(也可以利用Chem3D软件).先在GaussView的主窗口中问出水分子模型.然后选择菜单Calculate/Gaussian.程序弹出如图8.30所示的参数设定窗口.选择其中DFT的B3LYP方法.在6-31G基组水平上对水分子进行计算。同时在JobTypc中选择Opt+Freq.即进行水分子结构优化和频率计算。将标题，Title取名为H2O.同时在Link0中将检查点(Checkpoint)文件也取名为H2O.即%chk=D:\GaussianScratch\H2O.chk.点击Retain.使用命令File/save.保存为H2O.gjf文件。上一页下一页返回8.2分子模拟软件简介准备好输入文件后.便可在Gaussian03W中读入此输入文件(如图8.31所示).点击界面右栏的RUN按钮.即可开始进行计算。Gaussian在编写中采用了层状结构将程序分成几大块.每层内又分为很多小块.每个小块成为一个Link.对应一个程序.执行一种功能。为了节省内存.Gauss