1.       convert the pdb-file to a gromacs structure file(.gro / .pdb) and a gromacs topology file(.top)

这一步需要gromacs的预处理程序pdb2gmx

标准的输入命令应该是：

＃ pdb2gmx -ignh -f speptide.pdb -p speptide.top -o speptide.gro

这个命令之后有一个force field 的选择 系统给出了11种力场，具体以后实际计算中如何选择力场还是应该多参阅文献和gromacs说明书的。

这步之后可以在屏幕看到：

                ——— PLEASE NOTE ————

You have succesfully generated a topology from: speptide.pdb.

The G43a1 force field and the spc water model are used.

Note that the default mechanism for selecting a force fields has

changed, starting from GROMACS version 3.2.0

                ——— ETON ESAELP ————

至此，完成了从pdb到gro和top的转变。另外，我也看过有的作这一步命令时候并不进行top文件的转换，因为在genbox（对盒子加水溶剂化时还要对top文件进行修改）

2.       solvate the peptide in water

这一步需要用到editconf和genbox两个程序，命令的标准书写应该是：

# editconf -bt cubic -f speptide.gro -o speptide_box.gro -c -d 0.5

可以在屏幕上看到

Read 191 atoms

Volume: 8.17315 nm^3, corresponds to roughly 3600 electrons

No velocities found

    system size :  1.774  3.372  1.367 (nm)

    diameter    :  3.517               (nm)

    center      :  2.650  1.453  2.417 (nm)

    box vectors :  1.774  3.372  1.366 (nm)

    box angles  :  90.00  90.00  90.00 (degrees)

    box volume  :   8.17               (nm^3)

    shift       : -0.391  0.806 -0.158 (nm)

new center      :  2.259  2.259  2.259 (nm)

new box vectors :  4.517  4.517  4.517 (nm)

new box angles  :  90.00  90.00  90.00 (degrees)

new box volume  :  92.18               (nm^3)

 

因为程序的默认，所以这一步命令可以写作：

# editconf –f speptide –o –d 0.5

这里－f之后的sp文件是默认gro文件或者pdb文件，－o后面生成的output file是放入盒子的体系的gro文件 按照这样的格式，输出的文件被默认命名为－out.gro

 

体系放入盒子之后在盒子中注入水

# genbox -cp speptide_box.gro -cs -p speptide.top -o speptide_water.gro

也可简写为 ：

# genbox –cp out –cs –p sprptide –o b4em

Screen–

Output configuration contains 9035 atoms in 2967 residues

Volume                 :     92.1813 (nm^3)

Density                :     994.449 (g/l)

Number of SOL molecules:   2948  

 

之后看加水之后的top文件 用命令 #  tail speptide.top

可以看到

 [ system ]

; Name

Protein in water

 

[ molecules ]

; Compound        #mols

Protein             1

SOL              2948

其中SOL后面的数字2948就是genbox加入盒子的水分子的数目。

Editconf程序的另一个用途是讲gro文件转化回pdb 这时可以讲speptide_water.gro转化回pdb观察

＃ editconf –f speptide_water.gro –o speptide.pdb

拖回本机 用spbdv或者vmd观察 。

The next step is to generate index file . in the tutor of  gromacs , we are told that there are a set of index groups to select, unfortunately, I didi not find them ,so I have to use make_ndx to generate one .

# make_ndx –f b4em

 

3.       perform an enery  minimization of the peptide in solvent

Now the simulation system is almost ready. Before we can start the dynamics, we must perform an energy minimization, to alleviate any bad contacts (atoms overlapping such that a significant repulsion would result, causing numerical problems in the simulation) that might be present in the system.

# grompp –v –f em.mdp –c speptide_water.gro –p speptide.top –o speptide_em.tpr

Or use this command for short:

# grompp –v –f em –c b4em –p speptide –o em

After this command ,the bad contacs have been removed .so we can do enery minimization now.

# mdrun –v –s speptide_em.tpr –o speptide_em.trr –c after_em.gro –g emlog.log

Or

# mdrun –v –s em –o em –c after_em –g emlog

From screen ,we can see these ——

 Steepest Descents:

   Tolerance (Fmax)   =  2.00000e+03

   Number of steps    =          100

Step=    0, Dmax= 1.0e-02 nm, Epot= -8.02562e+04 Fmax= 2.70468e+04, atom= 6792

Step=    1, Dmax= 1.0e-02 nm, Epot= -8.69641e+04 Fmax= 1.24531e+04, atom= 6657

Step=    2, Dmax= 1.2e-02 nm, Epot= -9.48067e+04 Fmax= 5.54072e+03, atom= 9009

Step=    3, Dmax= 1.4e-02 nm, Epot= -1.02743e+05 Fmax= 2.86957e+03, atom= 8766

Step=    4, Dmax= 1.7e-02 nm, Epot= -1.08941e+05 Fmax= 1.53075e+04, atom= 150

Step=    5, Dmax= 2.1e-02 nm, Epot= -1.10083e+05 Fmax= 1.31653e+04, atom= 150

Step=    6, Dmax= 2.5e-02 nm, Epot= -1.10564e+05 Fmax= 2.16454e+04, atom= 150

Step=    7, Dmax= 3.0e-02 nm, Epot= -1.11528e+05 Fmax= 1.82398e+04, atom= 150

Step=    9, Dmax= 1.8e-02 nm, Epot= -1.12767e+05 Fmax= 3.87548e+03, atom= 120

Step=   11, Dmax= 1.1e-02 nm, Epot= -1.13598e+05 Fmax= 1.23013e+04, atom= 120

Step=   12, Dmax= 1.3e-02 nm, Epot= -1.14464e+05 Fmax= 5.37656e+03, atom= 120

Step=   13, Dmax= 1.5e-02 nm, Epot= -1.14787e+05 Fmax= 1.74494e+04, atom= 120

Step=   14, Dmax= 1.9e-02 nm, Epot= -1.15874e+05 Fmax= 8.35820e+03, atom= 120

Step=   16, Dmax= 1.1e-02 nm, Epot= -1.16420e+05 Fmax= 6.65421e+03, atom= 120

Step=   17, Dmax= 1.3e-02 nm, Epot= -1.16807e+05 Fmax= 1.08106e+04, atom= 120

Step=   18, Dmax= 1.6e-02 nm, Epot= -1.17248e+05 Fmax= 1.12834e+04, atom= 120

Step=   19, Dmax= 1.9e-02 nm, Epot= -1.17425e+05 Fmax= 1.45145e+04, atom= 120

Step=   20, Dmax= 2.3e-02 nm, Epot= -1.17590e+05 Fmax= 1.67089e+04, atom= 120

Step=   22, Dmax= 1.4e-02 nm, Epot= -1.18845e+05 Fmax= 2.65542e+03, atom= 121

Step=   23, Dmax= 1.7e-02 nm, Epot= -1.19198e+05 Fmax= 2.39072e+04, atom= 120

Step=   24, Dmax= 2.0e-02 nm, Epot= -1.20588e+05 Fmax= 6.17310e+03, atom= 120

Step=   26, Dmax= 1.2e-02 nm, Epot= -1.20901e+05 Fmax= 1.04224e+04, atom= 120

Step=   27, Dmax= 1.4e-02 nm, Epot= -1.21229e+05 Fmax= 8.89909e+03, atom= 120

Step=   28, Dmax= 1.7e-02 nm, Epot= -1.21346e+05 Fmax= 1.51125e+04, atom= 120

Step=   29, Dmax= 2.1e-02 nm, Epot= -1.21686e+05 Fmax= 1.29390e+04, atom= 120

Step=   31, Dmax= 1.2e-02 nm, Epot= -1.22216e+05 Fmax= 3.21419e+03, atom= 120

Step=   33, Dmax= 7.5e-03 nm, Epot= -1.22523e+05 Fmax= 6.52631e+03, atom= 120

Step=   34, Dmax= 8.9e-03 nm, Epot= -1.22796e+05 Fmax= 6.31916e+03, atom= 120

Step=   35, Dmax= 1.1e-02 nm, Epot= -1.22998e+05 Fmax= 8.04922e+03, atom= 120

Step=   36, Dmax= 1.3e-02 nm, Epot= -1.23167e+05 Fmax= 9.93619e+03, atom= 120

Step=   37, Dmax= 1.5e-02 nm, Epot= -1.23317e+05 Fmax= 1.10029e+04, atom= 120

Step=   39, Dmax= 9.3e-03 nm, Epot= -1.23857e+05 Fmax= 1.88770e+03, atom= 810

 

writing lowest energy coordinates.

 

Steepest Descents converged to Fmax < 2000 in 40 steps

Potential Energy  = -1.2385705e+05

Maximum force     =  1.8876959e+03 on atom 81

Norm of force     =  1.3663420e+04

NOTICE!!!!!!!!!

If the potential enery after minimization is lower than -1.1e+05kJ/mol, it is acceptable and the structure can be used for MD caculations.