一种基于约束正则模式的复杂函数极小值搜索方法技术

本发明专利技术属于计算化学和物理技术领域，具体为一种基于约束正则模式的复杂函数极小值搜索方法。本发明专利技术利用输入的原子坐标，已知的势能面能量函数，能量函数对坐标的一阶导数，求解能量极小值对应的原子坐标。其步骤为，从一个极小值对应的坐标体系开始，随机产生并优化得到一个约束的正则模式，通过在该模式方向上持续添加偏置势函数并重复能量极小值优化，达到越过势能面能量极大值的目的，最终优化得到新的极小值对应的坐标体系。本发明专利技术能快速搜索全局极小值的效果，并适用于复杂函数体系，同时具有搜索最优反应通道的功能。可用于遍历复杂分子、周期性晶体体系的势能面。

【技术实现步骤摘要】

本专利技术属于计算化学和物理
，具体涉及一种快速的复杂函数极小值搜索方法，可用于遍历复杂分子、周期性晶体体系的势能面。
技术介绍
结构预测和反应路径搜索作为当代化学和物理计算模拟研究的核心任务，对于理解和预测材料的热力学性质和动力学性质具有不可替代的作用。虽然分子动力学模拟作为一种搜索势能面和模拟化学反应过程的常规工具已经得到广泛的应用，但是这一方法的预测能力在处理高维度复杂势能面体系、或具有高反应活化能的化学过程体系时会大幅下降。比如分子动力学方法，不能很好处理蛋白质折叠或碳纳米管的生长过程。对于这些具有较高反应活化能的化学过程，一种比较常用的方法是先搜索化学过程所对应的过渡态，再找到经过这一过渡态的反应通道。这一方法最明显的缺点就是它高度依赖于科研人员通过化学直觉所预先猜测的可能的反应通道，而这种预先猜测在高维体系中是几乎不可能做出的。另一类常用的方法是通过在预先定义的反应坐标上增加限制条件(比如添加偏置势函数)来增加高活化能过程的采样概率，比如metadynamics方法(Phys. Rev. Lett. 2003,90, 238302)和伞状抽样方法(umbrella sampling) (J. Chem. Phys. 1998，109，7737)。由于需要根据体系预先定义限制条件，这些方法通常只能用于搜索单步化学反应，不能用于遍历势能面上可能的能量极小点。与前面方法相反,全局最优化方法，比如BH方法(Basin-Hopping方法，J. Phys. Chem. A 1997,101，5111)和遗传算法(Genetic Algorithms in Search, Optimization and Machine Learning, Addison-ffesley Reading, MA, 1989)是以得到函数的最小值为最终目标，也就是势能面上的最稳定结构。这些方法普遍采用了激进的结构形变，来达到跳出势能面上能量极小值区域的目的，所以这些方法完全忽略了势能面上不同极小值间互变的反应通道信息，而且对于复杂势能面的处理效率也大大降低。
技术实现思路
本专利技术的目的是为提供一种简便普适的势能面函数搜索方法，克服以往方法中对于高维度势能面的复杂体系预测能力大大降低的缺点。本专利技术的另一个目的是提供一种方法，可以同时用于能量极小值的结构预测和反应通道的预测。为实现上述目的，本专利技术基于偏置势函数驱动分子动力学思想，通过构造随机的运动模式，并沿着优化得到的约束正则模式，使构型从势能面上一个能量极小值向另一个能量极小值转化，最后利用Metropolis蒙特卡洛方法确定新找到的能量极小值构型是否被接受。其中描述势能面的能量函数的求解方法包括经验力场方法，量子力学从头算方法， 第一性原理密度泛函方法。本专利技术方法示意图如图I所示，具体步骤为利用输入的原子坐标Qi, i=l, 2，…η,已知的势能面能量函数E=Z^qi), i=l, 2，…η,能量函数对坐标的一阶导数F=^(Qi), i=l, 2，…η,求解能量所有极小值对应的原子坐标，以及连接这些极小值的反应通道，具体步骤如下步骤一、输入一个能量极小值坐标Α,并随机产生一个模式NO ；步骤二、在随机产生的模式NO基础上，通过约束双子方法优化得到一个约束的正则模式N;步骤三、在原有的势函数基础上，新添加该N方向的一个偏置高斯势函数，则总势能面L可以表示为真实势能面Vreal和偏置势函数匕的和，即 Vm = K- +I ,沿着该N方向外推当前坐标，到由偏置势函数和原有势函数叠加产生的能量极小值附近；步骤四、能量极小值优化当前坐标；步骤五、判定是否达到最多的偏置高斯势函数个数(NG)的上限，或达到预设的能量阈值4 (满足Emax〈 4);达到则转到步骤六，否则重复步骤二、步骤三、步骤四；步骤六、去掉所有添加的偏置高斯势函数，能量极小值优化得到新的坐标体系B ； 步骤七，采用Metropolis蒙特卡洛方法,判定B是否接受；如果接受，A被B替代；步骤八，重复步骤一到步骤七，寻找下一个能量极小值，或程序结束。本专利技术中，势能面能量函数的求解方法包括经验力场方法，量子力学从头算方法， 或第一性原理密度泛函方法。本专利技术中，步骤一随机产生的模式NO由两部分组成第一部分Ng，为一定温度条件下的Maxwell-Boltzmann速率分布,该部分为全局模式,第二部分N1,为体系内部任意两个不相邻原子相互靠近的运动模式；N0由Ng和N1的构成方法如下N0 = (Nf + m)/|Nf + aNi|，参数λ的取值为O. I到10之间的一个随机数。本专利技术中，步骤二中所述约束双子方法的步骤为定义势能面上的两个坐标点，分别记为R0，Rl,并使得R0，Rl满足关系Rl = HD+ΔΗ_Ν，其中Λ R为双子间距，一般预先设定为O. 005埃。在Rl上作用偏置势函数Vn，之后运用双子方法优化得到N。本专利技术中，步骤三中所述偏置高斯势函数的具体形式如下权利要求1.，其特征在于，利用输入的原子坐标q” i=l, 2，…η,已知的势能面能量函数E=ZXqi)，i=l, 2，…η,能量函数对坐标的一阶导数F=g(qi)，i=l，2，…n，求解能量所有极小值对应的原子坐标，以及连接这些极小值的反应通道，具体步骤如下步骤一、输入一个能量极小值坐标A,并随机产生一个模式NO ；步骤二、在随机产生的模式NO基础上，通过约束双子方法优化得到一个约束的正则模式N ;步骤三、在原有的势函数基础上，新添加该N方向的一个偏置高斯势函数，则总势能面L可以表示为真实势能面Vreal和偏置势函数匕的和，即 ^ = K- +4 ,沿着该N方向外推当前坐标，到由偏置势函数和原有势函数叠加产生的能量极小值附近；步骤四、能量极小值优化当前坐标；步骤五、判定是否达到最多的偏置高斯势函数个数(NG)的上限，或达到预设的能量阈值4 ;达到则转到步骤六，否则重复步骤二、步骤三、步骤四；步骤六、去掉所有添加的偏置高斯势函数，能量极小值优化得到新的坐标体系B ； 步骤七、采用Metropolis蒙特卡洛方法,判定B是否接受；如果接受，A被B替代； 步骤八、重复步骤一到步骤七，寻找下一个能量极小值，或程序结束。2.根据权利要求I所述方法，其特征在于，势能面能量函数的求解方法采用经验力场方法，量子力学从头算方法，或第一性原理密度泛函方法。3.根据权利要求I所述方法，其特征在于，随机产生的模式NO由两部分组成第一部分Ng,为一定温度条件下的Maxwell-Boltzmann速率分布,该部分为全局模式，第二部分N1, 为体系内部任意两个不相邻原子相互靠近的运动模式；N0由Ng和N1的构成方法如下4.根据权利要求3所述方法，其特征在于，参数λ的取值为O.I到10之间的一个随机数。5.根据权利要求3所述方法，其特征在于，所述约束双子方法的步骤为定义势能面上的两个坐标点，分别记为RO，Rl,并使得RO，Rl满足关系IU =HO+ΔΗ·Ν，其中Λ R为双子间距，在Rl上作用偏置势函数Vn，之后运用双子方法优化得到N。6.根据权利要求I所述方法，其特征在于，步骤三中，偏置高斯势函数的具体形式如下其中，i为偏置高斯函数计数。为当前坐标，Η〗为添加当前对应高斯势函数本文档来自技高网...

【技术保护点】
一种基于约束正则模式的复杂函数极小值搜索方法，其特征在于，利用输入的原子坐标qi,？i=1,2,…n，已知的势能面能量函数E=f(qi),？i=1,2,…n，能量函数对坐标的一阶导数F=g(qi),？i=1,2,…n，求解能量所有极小值对应的原子坐标，以及连接这些极小值的反应通道，具体步骤如下：步骤一、输入一个能量极小值坐标A，并随机产生一个模式N0；步骤二、在随机产生的模式N0基础上，通过约束双子方法优化得到一个约束的正则模式N；步骤三、在原有的势函数基础上，新添加该N方向的一个偏置高斯势函数，则总势能面Vtot可以表示为真实势能面Vreal和偏置势函数VG的和，即？？？？？？？？？？？？？？？？？？？？？？？？？？？？？？？？？？？？？？？？？？？？？？？？，沿着该N方向外推当前坐标，到由偏置势函数和原有势函数叠加产生的能量极小值附近；步骤四、能量极小值优化当前坐标；步骤五、判定是否达到最多的偏置高斯势函数个数（NG）的上限，或达到预设的能量阈值Er？；达到则转到步骤六，否则重复步骤二、步骤三、步骤四；步骤六、去掉所有添加的偏置高斯势函数，能量极小值优化得到新的坐标体系B；步骤七、采用Metropolis？蒙特卡洛方法，判定B是否接受；如果接受，A被B替代；步骤八、重复步骤一到步骤七，寻找下一个能量极小值，或程序结束。964079dest_path_image001.jpg...

【技术特征摘要】

【专利技术属性】
技术研发人员：刘智攀，商城，
申请(专利权)人：复旦大学，
类型：发明
国别省市：