Life is Worth Living

前一篇文章中讨论了键长键角的硬参数拟合，本文继续讨论力场参数拟合中的软参数拟合。软参数意义重大，软参数对分子间的相互作用和分子构象影响最大，硬参数影响较小，因此一般都先拟合软参数，再拟合硬参数。本章将对软参数的拟合方法进行简述，对软参数中的二面角参数进行详解。

前文提到了分子力场中的电荷模型的类型和历史沿革，接下来我们将讨论分子力场中的剩余部分的参数来源。除针对库仑势的电荷模型外，本文将讨论库仑势参数，键长，键角和二面角参数的拟合。力场的拟合可以大致分为硬参数拟合和软参数拟合，硬参数包含键长和键角参数，软参数包含范德华参数，库仑参数(电荷模型)和二面角参数。软参数对分子间的相互作用和分子构象影响最大，硬参数影响较小，因此一般都先拟合软参数，再拟合硬参数。即可按照库仑参数，范德华参数，二面角参数，键长和键角参数顺序拟合。本文以从简到难的顺序，先讨论键长和键角的拟合方法。

AMBER力场的御用电荷模型是限制性静电势电荷，Restrained Electro Static Potential (RESP)电荷。对于RESP电荷的思想来源推荐阅读以下两篇Sobereva的博文《RESP拟合静电势电荷的原理以及在Multiwfn中的计算》，以及《RESP2原子电荷的思想以及在Multiwfn中的计算》。针对GAFF小分子力场，Amber团队仍然首推RESP电荷，但是因为RESP电荷需要HF或B3LYP级别的量化计算，因此GAFF同时也提供一种快速的带键校正的半经验电荷模型AM1-BCC。关于各类电荷模型可阅读论文卢天,陈飞武等.原子电荷计算方法的对比[J].物理化学学报,2012(01):1-18.进行大致了解学习。本文将在以上文章基础上对RESP电荷模型和AM1-BCC电荷模型进行一定的补充讲解。

本文将以乙酰基噻唑分子为例，简单讲解GAFF力场小分子转换为Gromacs拓扑格式ITP/TOP应当包含的内容。

引言

在上一篇文章中我们了解了AMBER力场中的函数类型，根据我们所知，分子力场由函数形势和力场参数组成，本文将详细介绍GAFF小分子通用力场中如何描述不同化学环境中原子的力场参数。

引言

分子力场描述分子之间相互作用的势能与力。分子动力学模拟则通过分子力场逐步更新分子体系中所有原子的势能、速度和坐标，描述特定条件下的分子运动，从而获取我们关注的模拟体系微观或宏观信息。本博客将通过一系列文章以General Amber Force Field(GAFF)小分子通用力场为例，详细讲述小分子力场的的组成与拟合，对认识分子力场与分子力场优化提供一份基础教程。本篇内容作为开篇，将介绍分子力场的基础信息。

本文将以AMBER力场和Gromacs ITP文件格式为例，进行力场的介绍。之所以使用这样的组合，是因为AMBER力场几乎所有的力场搭建工具都打包在开源软件AmberTools中，因此配套的GAFF不需要借助在线工具即可离线构建小分子力场；选择Gromacs的力场格式ITP则是因为ITP文件比Amber的力场格式文件prmtop更加可读。

分子力场中，为描述复杂分子体系的电荷分布，常常会引入带有电荷的哑原子进行描述。最经典的哑原子分子比如TIP4P和TIP5P的水分子力场，4点水模型TIP4P将O的负电荷在H-O-H的夹角内的哑原子上，O的质量则留在O原子处；5点水模型TIP5P则将O的负电荷拆分在O原子的Lone Pair处，模拟O的孤电子对。

除去水分子外，AmberTools的MDGX模块提供了添加哑原子，或虚位点（Virtual Site）的方法，如教程35中对氯苯的苯环上哑原子，及Sigma Hole哑原子添加的方法。不过MDGX使用较为复杂，MDGX的电荷拟合教程也较复杂。本文将从Antechamber的RESP电荷拟合讲起，讨论添加哑原子后直接使用Antechamber进行哑原子电荷拟合的方法。

RESP电荷用于分子模拟，也是Amber/GAFF力场的御用电荷。Amber官方的的电荷教程以商业软件Gaussian和Jaguar为主，但目前使用开源的Psi4的Antechamber拟合教程并没有先例。除此以外，Multiwfn也是一种选择，Multiwfn支持多种量化格式的RESP和RESP2电荷的拟合。本文中，我们仅依靠python开发环境的Psi4和Antechamber进行有机小分子的RESP电荷拟合，以避开昂贵商用的Gaussian和Jaguar等软件.

My company are now using a docking engine which requires MOL2 as input. Since I use rdkit alot, I want to use rdkit to save mol2 file, while rdkit doesn’t support Tripos MOL2 format in current release. I don’t want to move to OpenBabel or other packages because I’m lazy ;P. So how about write one MOL2 writer?

I have found many email or issue of rdkit, and I found a pull request https://github.com/rdkit/rdkit/pull/415 after I write this script (Sad!). When I wrote these script, I haven’t notict this pull request, but I have read MOL2 scripts of UCSF Chimera and an introduction for MOL2 format. Chimera has many internal atom types which is perfectly matching MOL2 supported atom types, while rdkit don’t. For example, the amide bond. If I directly write bond type of amide bond, most MOL2 reader take this a rotatable single bond… Thus the only solution is use SMARTS pattern matching to assign MOL2 atom types.

Sometimes you want to see the 3D structure of molecule you are manipulating in jupyter, but you have to save it as a file with coordinates and open it in PyMol. This is very inconvenient. Although you can use PyMol to render image and show it in jupyter, but it is only a static image. Here I recommand a python package py3Dmol, which is a simple jupyter widget based on a JavaScript library 3Dmol.js. JavaScript can run in almost any browser, so you can use it in jupyter notebooks and in web pages.