MOPAC2012——最强大的半经验计算软件包！

文号t65472

43587

博丽灵梦2014/05/31化学 IP:浙江

本帖最后由 718281828kc 于 2014-6-1 15:12 编辑

一、MOPAC简介 　

最著名的半经验分子轨道程序包。MOPAC从1981年由J. J. P. Stewart开始开发，最早基于Dewar和Thiel的NDDO方法。在发布了MOPAC 7.0之后，MOPAC成为Fujitsu公司的商业软件Fujitsu MOPAC。但是15年后，作者在MOPAC 7.0的基础上重新开始了非商业版本MOPAC的开发，用FORTRAN 90/95语言进行了改写，并加入了很多新的功能。MOPAC 2007是目前唯一由J. J. P. Stewart维护的MOPAC版本。MOPAC 2007（即MOPAC 7.2）功能相当于商业程序Fujitsu MOPAC 2002，增加了PM6，但缺少MINDO/3，PM5，解析导数，Tomasi溶剂模型和系间穿越。

二、功能 　

（7.0标准版）
   1. MNDO，MINDO/3，AM1和PM3哈密顿量
   2. 限制Hartree-Fock (RHF)和非限制Hartree-Fock (UHF)方法
   3. 组态相互作用
      (a) 100个组态
      (b) 单重，双重，三重，四重，五重和六重
      (c) 激发态
      (d) 指定态的几何优化，等
   4. 单自恰场计算
   5. 几何优化
   6. 梯度最小化
   7. 过渡态定位
   8. 反应路径坐标计算
   9. 力常数计算
10. 简正坐标分析
11. 跃迁偶极计算
12. 热力学特性计算
13. 局域化轨道
14. 共价键顺序
15. 键分解为sigma和pi分量
16. 一维空间聚合物计算
17. 动力反应坐标计算
18. 内反应坐标计算

MOPAC 2007在MOPAC 7.0的基础上增加了以下功能：
1. 用FORTRAN 90改写。
2. 增强了周期体系（晶体、表面、聚合物）计算方面的功能。
3. 支持的方法有：MNDO，AM1，PM3，用于镧系元素的Sparkle/AM1，RM1，PM6，实现对所有主族元素及过渡金属的参数化。PM6是对NDDO方法的最近一次参数化，对近似做了三处修改，定义了对芯-芯相互作用的主要影响。除了12种镧系元素和所有的锕系元素之外，PM6参数已对大多数元素进行了优化。有机分子元素(H，C，N，O，F，P，S，Cl，Br，I)的参数已经稳定，其它主族元素、某些过渡元素以及1A和2A元素的参数尚在测试中。
4. 用PM6半经验模型模拟有机固体、无机固体的结构。
5. 用PM6计算氢键。
6. PM6对于仅由H，C，N，O，F，Si，P，S，Cl，Br和I构成的化合物计算极化率，平均误差为2.1%。
7. PM6更精确地计算生成热，熵，焓，热容，进行结构优化。对简单的有机分子精度最高，对复杂有机体系和无机体系精度稍低。PM6尚无法很好地计算频率。
8. 修改了AM1和PM3的一些严重错误。
9. 计算酸碱性(pKa)。
10. 产生动画。
11. 结构用JMol或CHIME显示。

MOPAC 2009的新功能：
1. 模拟酶的反应，优化整个蛋白质、DNA片断以及其它的大体系，可以进行直到15,000个原子的优化。
2. 线性标度算法MOZYME可以把半经验量子化学方法处理的极限从1,000个原子提高到15,000个（几何优化）或18,000个（SCF）原子，对100个原子以上的体系显著增加了计算速度。一个由1,000个原子构成的分子要比普通的MOPAC计算快60倍，并且只需要7%的内存。MOZYME可以与MOPAC的大多数功能结合使用，包括最新的PM6方法。

MOPAC 2012的新功能：
1. 加入新的半经验方法PM7和PM7-TS。

三、PM7方法简介

对于巨大的体系，DFT无论结合什么法子，是无论如何也处理不动的（当然，SCC-DFTB这样的半经验层面上的DFT是另外一回事）。如果用力场，准确度往往很有限，普适性差，而且经常缺乏合适的参数。半经验方法往往是唯一的选择，尽管其重要性已经很大程度上被人们忽视和遗忘。
传统的半经验方法，诸如主流的AM1、PM3，以及后来重新参数化而成的RM1和PM6，对于弱相互作用都比较差劲，就连氢键这样带有主要静电成份的弱相互作用尚不好，就更别提范德华作用了。自从色散校正在从头算、DFT领域流行开来后，人们纷纷把半经验方法也结合上色散以及其它分子间作用的校正项，使得半经验方法对弱相互作用的表现能力有了质的飞跃。包括AM1-D、PM3-D、OM1/2/3-D、PM6-DH1、PM6-DH2、PM6-DH2X、PM6-DH+、PM6-D3H4和PM7。其中2007年由McNamara和Hillier提出的AM1-D和PM3-D用的人比较少。2008年提出的OMx-D系列只有其作者Walter Thiel的MNDO2005程序才能算，这东西不公开给大家用（笔者鄙夷这种行为），所以自然也流行不起来，据说OM2-D的弱相互作用计算准确度能达到DFT-D的水准。2009年提出的PM6-DH1是在PM6基础上加上了色散、氢键校正项，后来很快被2010年提出的改进版PM6-DH2所取代。2011年的PM6-DH2X则是在PM6-DH2基础上又增加了对卤-氧和卤-氮型的卤键校正。2012年又提出了PM6-D3H4，除了氢键、色散校正项外还引入了Pauli互斥校正项，可能这是目前最好的基于PM6对弱相互作用进行校正的方法。前述几种修正PM6的方法都是Hobza等人提出的，2010年另一个人Korth也提出了一种PM6加上氢键和色散校正的形式，名为PM6-DH+，比PM6-DH2好，但未必比PM6-D3H4好，另外此方法还存在势能面在一些过程中表现得不平滑的问题。PM7是PM3、PM6和MOPAC的作者XXXXXXXXewart于2013年提出的最新的半经验方法，改进了其上一代PM6存在的一些缺陷，从统计数据上看热力学性质计算并没什么显著提高，但是这次在训练参数的过程中着重地考虑了弱相互作用体系，并且引入了弱相互作用的校正项，使得弱相互作用计算比起PM6有了长足进步。PM7由于比较新，目前来说和PM6-D3H4、PM6-DH+相比较谁对于弱相互作用更准确还不好说，但应该在伯仲之间。目前来说，PM7是整体上最好的半经验方法，巨大体系弱相互作用计算个人最推荐PM7。不过即便从统计数据上看，PM7以及各种修正PM6的方法在弱相互作用计算误差上都很小，和DFT-D系列相比还是有差距，可靠性差距更大，所以对于DFT-D能算得动的体系最好先用DFT-D，但是半经验方法可以作为初算，比如估计一下大致的作用能、初步优化一下结构，或者做动力学、构象搜索之类的很耗时的任务。

四、MOZYME简介

MOZYME是XXXXXXXXewart提出的一种基于定域化轨道加快半经验计算的技术，具体细节在这里不谈，可以参阅其原文IJQC,58,133(1996)。不像许多基于定域化轨道的方法，也就是仅限于作者自己发paper晒晒数据而无法真正投入广泛应用，MOZYME是有巨大实用意义的，加速效果是极其显著的，这将在后面的测试中看到。而且MOZYME是个黑箱方法，并不需要用户额外去考虑什么。MOZYME的存在使得半经验能方法处理的体系尺度大了两个数量级以上，能达到几万个原子的体系。而且MOZYME不光能用于单点计算，还能用于几何优化、搜索过渡态、动力学等需要导数的场合，不过只能用于闭壳层基态。对于大体系弱相互作用问题，MOZYME技术引入的误差通常是很小的（起码对于体系的尺度而言）。但对于极少数电子结构无法通过Lewis式表现的体系，有可能SCF过程会收敛到激发态上而不正确，不过这个问题对于生物分子这类成键方式典型的体系肯定不会出现。MOPAC的后期的商业版本以及免费版本中从2009版开始都支持MOZYME。MOZYME也仅在MOPAC中实现，是MOPAC的一个重要亮点。

五、安装方法

Windows版MOPAC2012安装方法是：把下载到的程序压缩包解压到C:\Program Files里，并且把目录名改为MOPAC。进入此目录，建立一个文本文件，名字为你向Stewart索取的密码，后缀为.mop，例如XXXXXXXXXXXXXXXXXXXp。然后右键点它选“打开”，并用MOPAC2012.exe来打开。在MOPAC2012的程序窗口里激活Graphic1窗口，敲回车，输入Yes，然后选“是”，激活即完成。

六、实例

计算条件就是普通的家用机的条件，CPU为Intel i7-2630QM(4核，最高频率2.5G)，双通DDR3-1333 16GB内存。系统为Win7-64bit。
1. H2O Box，左图12426原子，右图25095原子，左图PM7 MOZYME 4线程 14.4min，-282616.10434kcal/mol，内存用了1.2GB；
右图PM7 MOZYME 4线程 98.6min，-576587.48001kcal/mol，内存用了4.5GB。

可见即便是一两万个原子的体系，在普通计算条件下，利用MOPAC2012结合它特有的MOZYME技术也照样能应付得了。如果内存足够大，做到几万个原子的体系也是有可能的，完全能计算纳米级别体系的弱相互作用。

2. ETK蛋白，3944个原子

PM7 MOZYME 1/4threads 323.8/325.5s -16765.64836 kcal/mol，用了300MB内存。
可见对于几千个原子的普通大小的蛋白，用MOPAC2012来计算易如反掌。同时也看到在目前的MOPAC2012 13.331W版里，对于MOZYME计算，用单线程和4线程实际上速度基本没有差异。由于并行化目前还没带来什么效果，为了最大程度利用计算资源，如果你有n个核，可以同时算n个任务。MOPAC原本是不支持并行的，在撰写本文时MOPAC2012才刚刚支持并行不久，估计是程序现阶段还没把MOZYME的并行化搞好。
其它半经验方法的耗时：
PM6-DH+ MOZYME 4threads 349.6s -16405.40799 kcal/mol
PM6 MOZYME 4threads 342s -14168.08555 kcal/mol
PM3 MOZYME 4threads 330s -10619.64239 kcal/mol
RM1 MOZYME 4threads 295s -13358.60313 kcal/mol
AM1 MOZYME 4threads 271s -10033.36276 kcal/mol
可见，不同半经验方法耗时大致都在一个水平上。从生成焓来看，PM6-DH+和PM7的结果很接近，然而AM1和PM3则大大低估了生成焓，在很大程度上和它们表现弱相互作用很差有关。而RM1和PM6的结果则介于AM1/PM3和PM7之间，在侧面上说明RM1、PM6对弱相互作用的表现能力可能在AM1/PM3和PM7之间。

3. DNA片段

这是个DNA片段，左图共380个原子，含6个碱基对；右图共754个原子，含12个碱基对。

754个原子：
PM7 4threads 超过7小时计算尚未完成
PM7 MOZYME 1/4threads 58.88/57.8s -1161.85336 kcal/mol
可见MOZYME极其重要，不用MOZYME根本没法算。
389个原子：
PM7 1/4threads 56.2s/34.6s -2463.90001 kcal/mol
PM7 MOZYME 4threads 10.9s -2463.88277 kcal/mol
*BLYP D3 GCP(DFT/SVP) def2-SVP def2-SVP/J 124min
这样几百个原子的体系在MOPAC中不用MOZYME终于也能算了，并且并行化终于对计算速度有提升，但是并行效率还是很低。开了MOZYME后计算时间能省好几倍。MOZYME带来的误差可见非常之小，仅有0.117kcal/mol，对于几百个原子的体系完全可以忽略不计。

4. DNA小片段-小分子复合物，196个原子

PM7 MOZYME 4threads 2.5s
*BLYP D3 GCP(DFT/TZ) def2-TZVP def2-TZVP/J 91.4min
BLYP D3 GCP(DFT/SVP) def2-SVP def2-SVP/J 21.7min
对这样的约200个原子的体系，半经验计算只需要一眨眼的功夫，做几何优化也非常容易。BLYP在def2-TZVP级别也能很容易地计算，仅花了一个半小时，虽然比起用def2-SVP长了几倍，但还是值得的。

七、下载~~~

1. MOPAC2012：

MOPAC2012.zip 4.15MB ZIP 314次下载

2. 免费赠送一个可视化程序：

MOPAC2XYZ 1.0.2.zip 339.37KB ZIP 158次下载

注：此版本为14.083W For 32bit，其它版本(Linux or Win7 64bit)请在下文的官网中下载~
14.083W已经对并行作了改进，使并行效率大大提高。

序列号可以免费获得。。。见官方网站，制作人会亲自给你回复邮件~
XXXXXXXXXXXXXXXXXXXXXXXt/XXXXXXXXXXXXml
MOPAC2012独家支持PM7算法，是半经验算法精度最高的，还可以解决例如-SO₃Na之类使PM3报错的问题。
利用独家开发的MOZYME方法，可以在仅引入微量可忽略误差的情况下，成数量级加速运算。
而且。。。MOPAC2012的大小仅为11.1M(压缩后仅4.2M)！真是很小很强大！！

八、参考

1. 大体系弱相互作用计算的解决之道：XXXXXXXXXXXXXXXXXXX/sobereva/item/3cbf8904defdc4803d42e2e9

2. 将MOPAC2012的几何优化过程转换为xyz轨迹文件的小工具：XXXXXXXXXXXXXXXXXXX/sobereva/item/02416ec7f5676010505058e2

[groupid=337]化学[/groupid]

+800 科创币

bh，王者归来 2014/05/31 高质量发帖

+800 科创币

QQQQQQQ 2014/06/01 好东西当然加分

+1000 科创币

dracula1429 2014/06/01 有了这些个软件，个人条件就可以进行研究所做的HEDM筛选预测工作了！！

+899 科创币

qg2010 2014/06/02 高质量发帖

+1 学术分

QQQQQQQ 2014/09/19 刚下载试用，确实十分强大。

来自：数理化 / 化学

本帖最后由 718281828kc 于 2014-5-31 23:10 编辑

那个Example原来是甲醛啊。。。
话说一运行就直接秒杀了(0.04s, by PM7)

*******************************************************************************
** Site#: 19838        For non-commercial use only     Version 14.083W 32BITS**
*******************************************************************************
** Cite this work as:                                                        **
**                                                                           **
** (a) MOPAC2012, James J. P. Stewart, Stewart Computational Chemistry,      **
**  Version 14.083W web: XXXXXXXXXXXXXXXXXXXt         Days remaining: 297    **
**                                                                           **
**                                 AND                                       **
**                                                                           **
** (b) Maia, J.D.C., et al., J. Chem. Theory Comput. 2012, 8, 3072-3081.     **
**                DOI: 10.1021/ct3004645.                                    **
*******************************************************************************
**                                                                           **
**                                MOPAC2012                                  **
**                                                                           **
*******************************************************************************

                            PM7 CALCULATION RESULTS

*******************************************************************************
*  CALCULATION DONE:                                Sat May 31 21:25:57 2014  *
*  SYMMETRY - SYMMETRY CONDITIONS TO BE IMPOSED
*  T=       - A TIME OF 172800.0 SECONDS REQUESTED
*  DUMP=N   - RESTART FILE WRITTEN EVERY  7200.000 SECONDS
*******************************************************************************

     PARAMETER DEPENDENCE DATA

        REFERENCE ATOM      FUNCTION NO.    DEPENDENT ATOM(S)
            3                  1               4
            3                  2               4

             DESCRIPTIONS OF THE FUNCTIONS USED

   1      BOND LENGTH    IS SET EQUAL TO THE REFERENCE BOND LENGTH
   2      BOND ANGLE     IS SET EQUAL TO THE REFERENCE BOND ANGLE
SYMMETRY
Formaldehyde

  ATOM    CHEMICAL      BOND LENGTH      BOND ANGLE     TWIST ANGLE
NUMBER    SYMBOL       (ANGSTROMS)      (DEGREES)       (DEGREES)
   (I)                     NA:I           NB:NA:I       NC:NB:NA:I       NA    NB    NC
     1       O          0.00000000       0.0000000       0.0000000        0     0     0
     2       C          1.21081526  *    0.0000000       0.0000000        1     0     0
     3       H          1.09749962  *  122.0235013  *    0.0000000        2     1     0
     4       H          1.09749962     122.0235013     180.0000000        2     1     3

          CARTESIAN COORDINATES

    NO.       ATOM         X         Y         Z

     1         O        0.0000    0.0000    0.0000
     2         C        1.2108    0.0000    0.0000
     3         H        1.7928    0.9305    0.0000
     4         H        1.7928   -0.9305    0.0000

     ===================== shared-memory parallel CPU Section ======================
     **                                                                           **
     **  MOPAC calculations will run in shared-memory parallel CPUs               **
     **                                                                           **
     **  Number of Threads =  2                                                   **
     **                                                                           **
     **                                                                           **
     **  For information on this work, see:                                       **
     **                                                                           **
     **  Maia, J.D.C., et al., J. Chem. Theory Comput. 2012, 8, 3072-3081.        **
     **        DOI: 10.1021/ct3004645.                                            **
     **                                                                           **
     ===================== shared-memory parallel CPU Section ======================


General Reference for PM7:
"Optimization of Parameters for Semiempirical Methods VI: More Modifications to the
NDDO Approximations and Re-optimization of Parameters", J. J. P. Stewart, J. Mol. Mod., 1:32, 19 (2013)
XXXXXXXXXXXXXXXXXXXXXXXXXXX/XXXXXXXXXXp?genre=article&id=doi:10.1007/s00894-012-1667-x

           Empirical Formula: C H2 O  =     4 atoms

      MOLECULAR POINT GROUP   :   C2v

      RHF CALCULATION, NO. OF DOUBLY OCCUPIED LEVELS =    6

      Geometry optimization using EF

          DIAGONAL MATRIX USED AS START HESSIAN

CYCLE:     1 TIME:   0.008 TIME LEFT:  2.00D  GRAD.:    15.739 HEAT: -25.44830
CYCLE:     2 TIME:   0.000 TIME LEFT:  2.00D  GRAD.:    16.864 HEAT: -25.47569
CYCLE:     3 TIME:   0.000 TIME LEFT:  2.00D  GRAD.:     0.919 HEAT: -25.54241

     RMS GRADIENT =  0.91931  IS LESS THAN CUTOFF =  1.00000

-------------------------------------------------------------------------------
SYMMETRY
Formaldehyde

     GEOMETRY OPTIMISED USING EIGENVECTOR FOLLOWING (EF).
     SCF FIELD WAS ACHIEVED

                              PM7 CALCULATION
                                                       MOPAC2012 (Version: 14.083W)
                                                       Sat May 31 21:25:57 2014
                                                       No. of days left = 297

          FINAL HEAT OF FORMATION =        -25.54241 KCAL/MOL =    -106.86944 KJ/MOL

          TOTAL ENERGY            =       -444.38126 EV
          ELECTRONIC ENERGY       =       -840.31531 EV  POINT GROUP:     C2v
          CORE-CORE REPULSION     =        395.93406 EV
          COSMO AREA              =         60.70 SQUARE ANGSTROMS
          COSMO VOLUME            =         42.52 CUBIC ANGSTROMS

          GRADIENT NORM           =          0.91931
          IONIZATION POTENTIAL    =         10.352991 EV
          HOMO LUMO ENERGIES (EV) =        -10.353  0.402
          NO. OF FILLED LEVELS    =          6
          MOLECULAR WEIGHT        =         30.026

          MOLECULAR DIMENSIONS (Angstroms)

            Atom       Atom       Distance
            H     3    O     1     2.00299
            H     4    O     1     1.65067
            O     1    C     2     0.00000

          SCF CALCULATIONS        =          4

          WALL-CLOCK TIME         =          0.031 SECONDS
          COMPUTATION TIME        =          0.031 SECONDS

  ATOM    CHEMICAL      BOND LENGTH      BOND ANGLE     TWIST ANGLE
NUMBER    SYMBOL       (ANGSTROMS)      (DEGREES)       (DEGREES)
   (I)                     NA:I           NB:NA:I       NC:NB:NA:I       NA    NB    NC
     1       O          0.00000000       0.0000000       0.0000000        0     0     0
     2       C          1.20614565  *    0.0000000       0.0000000        1     0     0
     3       H          1.09115836  *  121.2760970  *    0.0000000        2     1     0
     4       H          1.09115836     121.2760970     180.0000000        2     1     3
                             CARTESIAN COORDINATES

   1    O        0.000000000     0.000000000     0.000000000
   2    C        1.206145652     0.000000000     0.000000000
   3    H        1.772634262     0.932586309     0.000000000
   4    H        1.772634262    -0.932586309     0.000000000

           Empirical Formula: C H2 O  =     4 atoms

      MOLECULAR POINT GROUP   :   C2v

                  EIGENVALUES
-31.15115 -23.10300 -16.82179 -14.94648 -14.14073 -10.35299   0.40195   4.45980
   4.85299   5.06188

              NET ATOMIC CHARGES AND DIPOLE CONTRIBUTIONS

  ATOM NO.   TYPE          CHARGE      No. of ELECS.   s-Pop       p-Pop
    1          O          -0.402102        6.4021     1.87897     4.52313
    2          C           0.202234        3.7978     1.15118     2.64659
    3          H           0.099934        0.9001     0.90007
    4          H           0.099934        0.9001     0.90007
DIPOLE           X         Y         Z       TOTAL
POINT-CHG.     2.873    -0.000     0.000     2.873
HYBRID        -0.207    -0.000    -0.000     0.207
SUM            2.666    -0.000     0.000     2.666

          ATOMIC ORBITAL ELECTRON POPULATIONS

    Atom      s        px        py        pz
    1  O   1.87897   1.26277   1.94677   1.31360
    2  C   1.15118   0.89165   1.06853   0.68640
    3  H   0.90007
    4  H   0.90007

TOTAL JOB TIME:             0.04 SECONDS

== MOPAC DONE ==

结果还是挺准的。。。甲醛的生成焓是-116kJ/mol，误差对于半经验方法来说是可以接受的，大概10kJ/mol。。。