官方introduction文档: http://www.ethz.ch/content/dam/ethz/special-interest/chab/physical-chemistry/csms-dam/doc/CSCBP-res/CSCBP_gro_man_v1_HS17.pdf Should not be confused with softwares, which begin in 1987 年: http://gromos.net/gromos87/GROMOS87_manual.pdf 目前都是基于的gromos96: 最初版本叫 gromos43a1 softwares: https://pubs.acs.org/doi/pdf/10.1021/jp984217f parameters: https://onlinelibrary.wiley.com/doi/abs/10.1002/%28SICI%291096-987X%2819980415%2919%3A5%3C535%3A%3AAID-JCC6%3E3.0.CO%3B2-N trp: CA-C bond type: gb_26 b0 = 0.1530, k_b=7.1500e+06 但这些参数来自 1996 年的 manual , which I cannot find the pdf file. Gromos43a2 is seldomly seen in papers, and 43系列已经不用了。 gromacs package refer it as improved alkane dihedrals to gromos43a1, so, if there's a need, use this one instead of gromos43a1. 2001 年: gromos45a3 : L.D. Schuler, X. Daura, and W.F. van Gunsteren. An Improved GROMOS96 Force Field for Aliphatic Hydrocarbons in the Condensed Phase. J. Comput. Chem., 22:1205–1218, 2001. I. Chandrasekhar, M. Kastenholz, R.D. Lins, C. Oostenbrink, L.D. Schuler, D.P. Tieleman, and W.F. van Gunsteren. A consistent potential energy parameter set for lipids: Dipalmitoylphosphatidylcholine as a benchmark of the GROMOS96 45A3 force field. Eur. Biophys. J., 32:67–77, 2003. ( called gromos45a3_C95 in this paper ) 2003 年: gromos45a4 : 45A3 reparameterised to improve DNA representation. T.A. Soares, P.H. Hu ̈nenberger, M.A. Kastenholz, V. Kr ̈autler, T. Lenz, R.D. Lins, C. Oostenbrink, and W.F. van Gunsteren. An Improved Nucleic-Acid Parameter Set for the GROMOS Force Field. J. Comput. Chem., 26:725–737, 2005. For 芳香环: R.D. Lins and P.H. Hu ̈nenberger. A new GROMOS force field for hexopyranose-based cardohydrates. J. Comput. Chem.,26:1400–1412, 2005. 2004 年: Gromos53a5/6: 一般用 gromos53a6 C. Oostenbrink, A. Villa, A.E. Mark, and W.F. van Gunsteren. A biomolecular force field based on the free enthalpy of hydration and solvation: the GROMOS force-field parameter sets 53A5 and 53A6. J. Comput. Chem., 25:1656, 2004. (影响的 Martini 力场的做法。) Oostenbrink C, Soares TA, van der Vegt NFA, van Gunsteren WF (2005) Validation of the 53A6 GROMOS force field. Eur Biophys J 34:273–284 该力场有问题: dihedral-angle parameters of the backbone transferred from the earlier version of the force field were no longer appropriate. 2011 年: gromos54a7 ,54b7 D. Poger, W.F. van Gunsteren, and A.E. Mark. A new force field for simulating phosphatidylcholine bilayers. J. Comput. Chem., 31:1117–1125, 2010. N. Schmid, A. P. Eichenberger, A. Choutko, S. Riniker, M. Winger, A.E.. Mark, and W.F. van Gunsteren. Definition and testing of the GROMOS force-field versions: 54A7 and 54B7. Eur. Biophys. J., 40:843–856, 2011. https://link.springer.com/article/10.1007%2Fs00249-011-0700-9 【目前境况】 Since the last official release of the GROMOS software and manual in 1996, called GROMOS96, no comprehensive release occurred till 2011. Yet the GROMOS software has seen a steady development since 1996. The programming language has been changed from FORTRAN to C++, the documentation has been put into electronic form, and many new features have been included in the software. 该力场在膜蛋白和大蛋白体系还是有一定市场的: A. Choutko, A. Gl ̈attli, C. Fern ́andez, C. Hilty, K. Wu ̈thrich, and W.F. van Gunsteren. Membrane protein dynamics in different environments: simulation study of the outer membrane protein X in a lipid bilayer and in a micelle. Eur. Biophys. J., 40:39–58, 2010. 使用的还是 gromos45a3. Water 用的是 spc! 可能是因为体系太大了吧。 2014 年 cell: https://www.sciencedirect.com/science/article/pii/S0092867414014287 gromos45a3 和 spc water model. 体系为: Muscle α-Actinin 大蛋白。 这大概就是 gromos45a3 的处境 。 2016 年 gromos53a6: https://www.nature.com/articles/nature17199 用的是 Martini 然后回到全原子。大体系。 2014 science: http://science.sciencemag.org/content/344/6185/1249783.full gromos53a6 和 gromos54a7 都用了,水为 SPC water ,大体系。 2013 jacs ,也是 Martini back mapped 的研究: gromos54a7 https://pubs.acs.org/doi/pdf/10.1021/ja310577u 2015 nature, 力场方面比较细致,先用 gromos54a7 prepare, 然后用 charmm27. https://www.nature.com/articles/nature14158 【后续发展】 gromos54a8: M.M. Reif, P.H. Hu ̈nenberger, and C. Oostenbrink. New interaction parameters for charged amino acid side chains in the GROMOS force field. J. Chem. Theory Comput., 8:3705–3723, 2012
\01988 年: https://pubs.acs.org/doi/pdf/10.1021/ja00214a001 分子间 coulomb 和 vdw 项: 分子内 parameter 取自 amber ,所以可以称作 Amber/OPLS force field. 第一次提出 hydrogen bond 项可以去掉。 Amber94 提出之后: 1996 版本 : OPLS-AA https://pubs.acs.org/doi/abs/10.1021/ja9621760 Gromacs package 提供的是 2001 版本:但电荷参数等基本取自 1996 版本,这里只是小优化。 https://pubs.acs.org/doi/pdf/10.1021/jp003919d force field torsion parameters have been optimized E tot =E bond +E angle +E torsion +E nb 其中 E nb : 此外: 注意这里的 torsion 跟 amber 不同,但似乎后面的并没有用到。 C-CA 键参数:分别对应原子类型 opls_235--opls_224B 怎么和 ffbond.itp 里对应,根据 ffnonbonded.itp 文件。 b0=0.15220 kb=265265.6 果然跟 amber94 力场里是一致的。 CA 的电荷为 0.140 ,而 amber94 的为 -0.02750 ,在氨基酸文件里,差很多呢。 【后续发展】 OPLS-AA/M: Improved peptide and protein torsional energetics with the OPLS-AA force field J Chem Theory Comput, 11 (2015), pp. 3499-3509 2017 年: OPLS3 https://pubs.acs.org/doi/pdf/10.1021/acs.jctc.5b00864 作者在文中坦诚的指出了事实: Although, improvement of the OPLS protein force field has received less attention than for CHARMM and AMBER over the past decade, the major revision of the protein force field in OPLS3 reported here appears to display performance that is competitive with the state of the art. 而且用该力场去做 MD ,发高影响因子的文章很少,更多的是测试。 特点是 used empirical charges , amber 和 charmm 是 QM fitting. 近期应用: Rosenman, D. J.; Connors, C. R.; Chen, W.; Wang, C.; Garcia, A. E. J. Mol. Biol. 2013 , 425, 3338−3359. They say: “Comparison of simulation data published in this study with OPLS-AA/TIP3P and that generated by Sgourakis et al. with AMBER99sb/TIP4P-Ew showed that the former described the experimental data better in longer simulations.” 配合 TIP3P 还是不错的。
gromacs软件包top文件夹下:7种 amber94 amberGS.ffamber96 amber99amber99sbamber99sb-ildn amber03 各个版本: http://ffamber.cnsm.csulb.edu/ffamber.php amber94 就是 jacs1995 文章的参数。 AMBER-96 is a variant of AMBER-94 in which peptide f/y torsional potentials have been adjusted to yeild better agreement between molecular mechanical and quantum mechanical energetics (disfavoring helices, favoring extended conformations). AMBER-GS is identical to AMBER-94 with peptide f/y torsional potentials removed (set to zero). 这里说了, In our ff99SB paper we cite and discuss it, and show that it has some not very good properties. http://archive.ambermd.org/201203/0517.html AMBER-99 is the 3rd generation update to AMBER-94(2 nd generation), including updated parameters for both amino and nucleic acids. The primary changes for peptides/proteins is again in the f/y torsional potentials (again to disfavor helical conformations). 所以可以放弃 amber94, amber96 了。 AMBER-99SB is the AMBER-99 variant developed by Simmerling and co-workers with updated backbone f/y torsions fit to ab initio calculations. 第三方修正的。 AMBER-03 is a variant of the AMBER-99 potential in which charges and main-chain torsion potentials have been rederived based on QM+ continuum solvent calculations and each amino acid is allowed unique main-chain charges. 也就是 amber99 的官方修正版。( 注意 : Instead of relying on the HF/6-31G* approach to provide aqueous-phase charges, a low-dielectric continuum model corresponding to an organic solvent environment was included directly in the QM calculation of the dihedral parameters and electrostatic potential (from which the charges are obtained). Because of these differences, ff03 should be considered a distinct force field model rather than extension of previous Amber force fields ) amber99sb-ildn 针对 I,L,D,N 四种氨基酸的优化: 被2010年 science 文章引用: Atomic-Level Characterization of the Structural Dynamics of Proteins ( tip3p water )因此倾向于采用 amber99sb-ildn 。 最新的是amber ff14sb力场: https://pubs.acs.org/doi/pdf/10.1021/acs.jctc.5b00255 ff14sb 仍然是用的 tip3p water model 力场,说明 tip3p 还是可以继续用的。此外, ff14sb 目前已经推荐,但 gromacs 力场 package 里还没有默认 include 进来,而 ff99sb-ildn 也已经足够好 。 关于 ff99SB* 和 ff03* 力场, proposed in this paper: https://pubs.acs.org/doi/pdf/10.1021/jp901540t 该力场知名度一般,升级版 ff99SB*-ildn 在 D.E.Shaw 的文章里第一次提,详细讨论见: http://archive.ambermd.org/201409/0260.html 此外, Review article 说: the force fields that are capable of reversibly folding globular proteins such as AMBER99sb with * and/or ILDN 165c modifications and CHARMM22* 166 in the studies described in refs 170b and 170d may not necessarily be the most suitable for characterizing the metastable states of disordered ensembles. (其实这里不包括 amber99sb-ildn,170b 是 amber99sb*-ildn , 170d 是 charmm22*。 )
this method and its associated software, FUERZA, were described in 1996. The essential concept is to analyze the interactions between all pairs of atoms in the molecule and to determine which ones are pairwise stable . In related applications, only the bond and angle parameters based on the Seminario method are used, while the dihedral and improper torsion terms involving metal ions are set to zero. 【redundantinternal coordinates】 A given molecule can be described by various sets of internal coordinates , and these may give different values for a particular internal coordinate force constant. 这一点,什么意思? 文章 里举了一个很好的例子: \0 A quick HF/STO-3G calculation yields a force constant for the CO bond of 0.54 a.u . if the chosen internals are the CO bond, CN bond, and NCO angle; 0.57 a.u . if we choose the CO bond, NO bond, and CON angle; and 0.26 a.u . if we choose redundant internal coordinates CO, NO, CN, CNO, NOC, and OCN. But why? How to understand? This is due to fact that physical force constants are tensors of rank 2, but we want to use them in practical MD simulations as scalars. so 该文章 实现了 cartesian坐标下的高效率。进一步redundantinternal coordinates可以实现更高的效率,详见该 博文 。
about 1-4 interactions two implement methods: 1. a list of in ffnonbond.itp, gen_pair=no , only useful for force field like Gromos. 2. For OPLS and Amber, 1-4 interactions are scaled by a factor, so theirs is no in ffnonbond.itp, and gen_pair=yes . It should be noted that the coulomb interaction is also scaled. exclusion=3 is default for all force fields, because are included in all force fields.
水模型: spc: 1. H. J. C. Berendsen, J. P. M. Postma, W. F. van Gunsteren, and J. Hermans, In Intermolecular Forces , edited by B. Pullman (Reidel, Dordrecht, 1981 ), p. 331. 很好。 spce: H 的电荷,由 0.41 增大到 0.4238 。 Berendsen, H. J. C. ; Grigera, J. R.; Straatsma, T. P. (1987). The missing term in effective pair potentials. J. Phys. Chem . 91 : 6269–6271. doi : 10.1021/j100308a038 . The SPC/E model results in a better density and diffusion constant than the SPC model. tip3p 模型 :比 spce 还早。 但一直在应用 ! H 的电荷取 0.417. Jorgensen, W. L.; Chandrasekhar, J.; Madura, J. D.; Impey, R. W.; Klein (1983). Comparison of simple potential functions for simulating liquid water. J. Chem. Phys . 79 : 926–935. Bibcode : 1983JChPh..79..926J . doi : 10.1063/1.445869 . Three-site model (TIP3P) has better performance in calculating specific heats . 也就是有温度变化的情况。 tip4p: same paper with tip3p H 电荷取 0.52 tip4p-ew: H 电荷修正为 0.52422. for use with Ewald summation methods Horn, H. W.; Swope, W. C.; Pitera, J. W.; Madura, J. D.; Dick, T. J.; Hura, G. L.; Head-Gordon, T. (2004). Development of an improved four-site water model for biomolecular simulations: TIP4P-Ew. J. Chem. Phys . 120 : 9665–9678. Bibcode : 2004JChPh.120.9665H . doi : 10.1063/1.1683075 . 看起来不错。 tip5p: H 电荷为 0.241 。 Mahoney, M. W.; Jorgensen (2000). A five-site model for liquid water and the reproduction of the density anomaly by rigid, nonpolarizable potential functions. J. Chem. Phys . 112 : 8910–8922. Bibcode : 2000JChPh.112.8910M . doi : 10.1063/1.481505 . The TIP5P-E model is a reparameterization of TIP5P for use with Ewald sums. but not included in Gromacs package. 关于 water 力场的很好的介绍: https://wenku.baidu.com/view/db4ee821336c1eb91a375da5.html?pn=51 acpye.py 生成 gromacs 格式 top 文件的时候,离子的 atomtype 类型不统一,是 bug ,我自己注意修改就好了。
amber1984 第一个版本, parameters are obtained for implicit water models. vdw 参数, OPLS ( Jorgensen, W. L., et al, JACS, 1988 )最成功, and is for explicit water models. It inspires the amber 1994 版本 (该版本,去掉了 Hydrogen-bond 项)。 * more about the hydrogen-bond项:(from here ) Hagler et al. observed that an explicit hydrogen-bond term was not needed in their model. 591 In subsequent work of Lifson et al., they also observed that hydrogen bonds were well described using a model consisting of an electrostatic term and a VDW interaction term. 595 (591) Hagler, A.; Huler, E.; Lifson, S. Energy Functions for Peptides and Proteins. I. Derivation of a Consistent Force Field Including the Hydrogen Bond from Amide Crystals. J. Am. Chem. Soc. 1974 , 96 , 5319 − 5327. (595) Lifson, S.; Hagler, A.; Dauber, P. Consistent Force Field Studies of Intermolecular Forces in Hydrogen-Bonded Crystals. 1. Carboxylic Acids, Amides, and the C: O. Cntdot.. Cntdot.. Cntdot. H- Hydrogen Bonds. J. Am. Chem. Soc. 1979 , 101 , 5111 − 5121. ***** 1994 版本里 vdw 用的是 6-12 势的 R*(r m =2 1/6 σ =1.122 σ ) 和 ε ,但由于 σ和ε是两粒子之间的,具体表示在单粒子上时,还要考虑combine-rule。 对 CA ,文献给出 R*=1.9080 Å , ε=0.086; gromacs 里 σ = 0.339967 nm , ε=0.359824 对于 σ, 是因为 combine rule 不同,前者是直接相加 ,后者还要乘 1/2. 对于 ε ,单位不同,前者是 Kcal/mol ,后者是 KJ/mol ,combine-rule一致。 gromacs package中amber1994 版本 Cu 和 Zn 就都有了,是哪里来的? 经分析,发现 Zn2+ 的参数与 1991 年 Merz 的 JACS 文章中一致。 而Cu的参数,考虑到 σ应该比Zn2+要小一点, 不应该是Cu2+离子,此外,ions.itp里面也没有。 amber package里的amber94力场也没有Zn2+离子,应该以amber package为准。 按这里,似乎 amber 软件包里是有 Cu2+ 的,只是 gromacs 包里没有。 http://www.quimica.urv.es/~bo/MOLMOD/General/Forcefields/AMBER.html#AAT amber package里的离子力场方面更新比较快。使用xleap工具,进行体系的构建,可以发现已经把PME下优化的一价和二价离子参数加进去了。 一价的工作: https://pubs.acs.org/doi/abs/10.1021/jp8001614 Cl 已经变掉了,而且变化很显著。 二价的工作: https://pubs.acs.org/doi/10.1021/ct400146w
大体详见: http://blog.sciencenet.cn/blog-485752-1119787.html 逻辑: 分 5 大步,全程是 amber 程序包的格式。 § prmtop - a file containing a description of the molecular topology and the necessary force field parameters. § inpcrd (or a restrt from a previous run) - a file containing a description of the atom coordinates and optionally velocities and current periodic box dimensions. § mdin - the sander input file consisting of a series of namelists and control variables that determine the options and type of simulation to be run. 所以,目标是得到 prmtop 文件。 tleap -s -f 1OKL_tleap.in 1OKL_tleap.out 这一步,可以得到 prmtop 和 inpcrd 文件。 要得到 1OKL_tleap.in , MCPB.py -i 1OKL.in -s 4 其中 1OKL.in 是用户自己创建的:内容一般为: original_pdb 1OKL_fixed_H.pdb group_name 1OKL cut_off 2.8 ion_ids 4018 ion_mol2files ZN.mol2 naa_mol2files MNS.mol2 frcmod_files MNS.frcmod large_opt 1 这些文件都是需要的。 MNS 指的是一个小分子,一般体系里没有,那么就不需要 MNS.mol2, MNS.frcmod 了。 ZN.mol2 怎么得到? antechamber -fi pdb -fo mol2 -i ZN.pdb -o ZN_pre.mol2 -at amber -pf y 1OKL_fixed_H.pdb 怎么得到? pdb4amber -i 1OKL_H.pdb -o 1OKL_fixed_H.pdb 其中, 1OKL_H.pdb 这样得到: cat 1OKL_Hpp_fixed.pdb ZN.pdb MNS_H_fixed.pdb 1OKL_H.pdb 而 1OKL_Hpp_fixed.pdb: ambpdb -p 0.15_80_10_pH6.5_1OKL.top -c 0.15_80_10_pH6.5_1OKL.crd 1OKL_Hpp.pdb 其中 top 和 crd 文件由 H++ server 得到。 不用 H++ server 也行,有时候死活过不去 ,直接跳过ambpdb和pdb4amber,手动得到1OKL_fixed_H.pdb。 手动修改掉 pdb2gmx 得到的 pdb 文件的时候,主要是 HD1,HD2 à HD2,HD3 等。 注意 CU.mol2 中 ATOM 那一行,第 8 , 9 列是必要的,第 8 列要和 pdb 文件对应。自动生成的可能不对。 距离怎么定?修改 trpz.in ,默认是 2.8 Å. 再回来考虑这一步: MCPB.py -i 1OKL.in -s 4 这一步生成的 1OKL_tleap.in ,需要调用 ZN1.mol2 相比 ZN.mol2 , ZN1.mol2 的电荷是 resp 电荷。这一步还生成 1OKL_mcpbpy.pdb ,在 tleap 中用到,注意检查 leap.log , 1OKL_mcpbpy.pdb 的原子标号 可能会出错。 ZN1.mol2 这样得到: MCPB.py -i 1OKL.in -s 3 此外,1OKL_tleap.in还调用了1OKL_mcpbpy.frcmod,该文件包含了residues的电荷、bond、angle、dihedral参数。 该步还需要1OKL_large_mk.log,这样得到: g09 1OKL_large_mk.com 1OKL_large_mk.log 而frcmod文件这样得到: MCPB.py -i 1OKL.in -s 2 注意这里有 2,2z 和 2e 可选, 2z=z-matrix ,支持 more than Zn, 2e=empirical method, 目前只支持 Zn 。 2= Seminario method , 默认的。 这一步是将之前的 gaussian 计算,转换成 force constant 等参数。 因此该步需要 1OKL_small_opt.fchk ,不然会报错。该文件这样得到: formchk 1OKL_small_opt.chk 1OKL_small_opt.fchk 其中 1OKL_small_opt.chk 文件需要: g09 1OKL_small_opt.com 1OKL_small_opt.log( 默认用双核 ) g09 1OKL_small_fc.com 1OKL_small_fc.log 参见 com 文件内容,两步中都对 1OKL_small_opt.chk 进行了修改。 而这里对 com 文件,需要: MCPB.py -i 1OKL.in -s 1 来生成。 注意: 只能在 N 和 Cu2+ 之间形成共价键? 我的情况不适合,因为: Generally the metal ion related parameters are follow these creteria: A) The bond force constants between an metal ion and its ligating atoms are less than 200 kcal/(mol*Angstrom^2), and the eqlibirum bond distances are less than 2.8 Angstrom ;