GROMACS里限制、冻结等参数的用法
在使用gromacs进行分子动力学模拟时候，我们常需要把某个原子、基团或者分子保持在固定位置，另一方面需要把分子中的某个键长、键角或者二面角限制住，使其保持某个长度或者角度。gromacs里有许多功能参数可以实现，下面进行简要的讲解。
1）冻结参数
需要冻结某个组分时候，可在mdp文件里面写入如下样式的参数：
freezegrps = protein
freezedim = Y Y Y
以上案例为冻结住蛋白（protein）组，后一项“Y”为在x，y，z方向冻结。若只想在x，y两个方向冻结，可以改为：freezedim = Y Y N，以此类推。
若需要冻结多个组，可在等号后面继续添加即可，如下所示，同时冻结蛋白和水：
freezegrps = protein SOL
freezedim = Y Y Y Y Y Y
有时候需要冻结特定组时候，需要用gmx make_ndx命令建立索引组，然后执行gmx grompp命令时添加 -n index.ndx 使得在mdp文件里面的设置组的名称能够被程序识别。如建立索引组的名称如果叫做：r_1，则冻结命令写为：
freezegrps = r_1
freezedim = Y Y Y
注：冻结之后可能会使体系能力非常高，导致模拟崩溃。建议先能力最小化之后再做冻结，或者检查出错的成键参数，把相应力常数改小。
2）成键项限制
使用这几项限制时候，必须放在需要限制的成键项所在分子的itp文件内，或者紧跟调用该分子的include命令后面。
2.1 距离限制：使用 [ distance_restraints ] 定义，案例如下：
[ distance_restraints ]
; ai aj type index type' low up1 up2 fac
10 16 1 0 1 0.0 0.3 0.4 1.0
10 28 1 1 1 0.0 0.3 0.4 1.0
10 46 1 1 1 0.0 0.3 0.4 1.0
ai 和 aj 列为要限制的粒子的原子编号（在该分子内的编号）。 type 列总是 1。在上述例子中，限制 10-28 和 10-46 的索引都是1，因此它们会被同时计算。对一起计算的那些限制的额外要求是，这几个限制必须位于连续的行中，中间不能有任何其他的限制干扰。列 low，up1 和 up2 的值分别对应下面方程中的r0，r1和r2的值：
对一些限制使用不同的力常数可能有益；这可通过fac列进行控制。对每个限制，参数文件中的力常数会乘以fac列中的值。
注：可以直接在[ bonds ]里用把两个原子的type改为6，也是添加一个谐振子势，但是该两个原子之间认为没有键连。
2.2 角度限制：使用 [ angle_restraints ] 定义，案例如下：
[ angle_restraints ]; i j k l type theta0 fc multiplicity23 677 34 677 1 112.0 3000 1
说明限制23-677-34三个原子形成的夹角在112.0度附近，力常数为3000。type没有用处，高版本会考虑删除。multiplicity多重性为1，默认1即可。
2.3 二面角限制：使用 [ dihedral_restraints ] 定义，案例如下：
[ dihedral_restraints ]; i j k l type label phi dphi kfac power4 12 8 11 1 1 180 0 1 2
上述限制的二面角是4-12-8-11组成的，type是限制类型，只能用1。label没有用处，高版本会考虑删除。phi是参考角，dphi是超过参考角多少度开始使用限制势，power没有用处，高版本会考虑删除。kfac乘上mdp中的dihre_fc将作为限制势力常数，所以可在mdp文件中做如下设置：
dihre = simple ; 简单（每个分子）的距离限制dihre_fc = 1000 ; 默认为1000，单位：kJ mol-1 nm-2，距离限制的力常数dihre_tau = 0.0 ; 单位：ps，距离限制运行平均的时间常数。值为零时关闭时间平均nstdihreout = 100 ; 单位：步，将限制所涉及的所有原子对之间的运行时间平均距离和瞬时距离，写入能量文件时的间隔步数。
3）位置限制
相对于最初位置进行限制，一般生成一个posre.itp文件，里面使用 [ position_restraints ]定义：
[ position_restraints ]
; ai funct fcx fcy fcz
1 1 1000 1000 1000 ; 限制在一个点
2 1 1000 0 1000 ; 限制在一条线（y轴）
2 1 1000 0 0 ; 限制在一个面（y-z平面）
一般该文件（posre.itp）在需要限制的整个分子itp后面写入，或者在调用该分子的include语句后面使用，使用方法为：
#ifdef POSRES
#include "posre.itp"
#endif
下列示意了错误用法和正确用法：
错误用法：
#include "topol_A.itp"
#include "topol_B.itp"
#include "ligand.itp"
#ifdef POSRES
#include "posre_A.itp"
#include "posre_B.itp"
#include "ligand_posre.itp"
#endif
正确用法：
#include "topol_A.itp"
#ifdef POSRES
#include "posre_A.itp"
#endif
#include "topol_B.itp"
#ifdef POSRES
#include "posre_B.itp"
#endif
#include "ligand.itp"
#ifdef POSRES
#include "ligand_posre.itp"
#endif
上述#ifdef名称都是使用的是：POSRES，所以可以在mdp文件里面添加如下一项，即可使限制生效：
define = -DPOSRES
其中等号后的“-D”是必须的，有时候只想单独限制某一个分子，那么可以使#ifdef后面的名称不一样，如改为：
#ifdef POSRES_LIG
#include "ligand_posre.itp"
#endif
此时，mdp文件里面需定义为：define = -DPOSRES_LIG，才使得该限制生效，若需要同时生效多个限制，则继续添加即可，如：
define = -DPOSRES_LIG -DPOSRES -DPOSRES_Agmx pdb2gmx命令生成蛋白等大分子itp时候，会自动生成、添加限制文件posre.itp在该分子itp里面，且名称默认为：POSRES 。使用下面命令可对任意分子进行限制文件制作：
gmx genrestr -f **.gro -o posre.itp -fc 1000 1000 1000**.gro里面应只包含生成限制文件分子的内容。或者使用**.pdb文件也行，gro/pdb文件务必和该分子的itp文件原子个数、顺序一致。生成posre.itp限制文件，力常数在x,y,z三个方向均为1000，可自行定义限制的方向。
同时，使用高版本时候，在运行gmx grompp命令时候，需要添加 -r **.gro 参数，**.gro文件可同-c参数后面的命令一致，若无限制则不加-r参数。如下案例：
gmx grompp -f nvt.mdp -c em.gro -p top.top -o nvt.tpr -r em.gro注：冻结和限制的区别在于，冻结是完完全全的不能运动或者振动，限制有时候会在其周围微小的变动，若限制的力常数设置得非常大，那么和冻结没有较大的区别了。
