[gmx-users] Error in free energy calculation

[Tasneem Kausar]

Dear All

I am trying to perform free energy calculations of a protein drug system.
I am using Gromacs-5.1.4 for this purpose. decoulpe molecule type i.e. drug
molecule is defined in separate drug.itp file. From the discussions of
mailing list I have found that this is not a problem. To obtain the
equilibrium values for protein-ligand restraint 15 ns simulation was
performed. Here is the last section of topology file.
 ; Include Position restraint file
#ifdef POSRES
#include "posre.itp"
#endif
#include "drug.itp"
[ intermolecular_interactions ]
[ bonds ]
;    i     j  type     r0A     r1A     r2A    fcA    r0B     r1B     r2B
fcB
  1444  3118    10     0.591   0.591   10.0   0.0    0.591   0.591   10.0
4184.000

[ angle_restraints ]
;   ai    aj    ak    al  type    thA      fcA    multA  thB      fcB
multB
  1431  1444  3118  1444     1    52.35    0.0    1    52.35    41.840    1
  1444  3118  3120  3118     1   118.14    0.0    1   118.14    41.840    1

[ dihedral_restraints ]
;   ai    aj    ak    al  type    phiA     dphiA  fcA    phiB      dphiB
fcB
  1444  1431  1444  3118     1    -99.41   0.0    0.0    -99.41    0.0
41.840
  1431  1444  3118  3120     1     17.49    0.0    0.0    17.49    0.0
41.840
  1444  3118  3120  3119     1     -6.49   0.0    0.0     -6.49    0.0
41.840
; Include water topology
#include "gromos54a7.ff/spc.itp"

#ifdef POSRES_WATER
; Position restraint for each water oxygen
[ position_restraints ]
;  i funct       fcx        fcy        fcz
   1    1       1000       1000       1000
#endif

; Include topology for ions
#include "gromos54a7.ff/ions.itp"

[ system ]
; Name
DUMM in water t= 10000.00000
[ molecules ]
; Compound        #mols
Protein_chain_A     1
DRG                 1
SOL             17685

Is there any default position for [ intermolecular_interaction ] in
topology file.
Since I am using a drug.itp file for drug molecule. Numbering in itp file
starts from 1. Is atom number in itp file needed to be modified?


grompp gives many warning  and 1 error
Back Off! I just backed up mdout.mdp to ./#mdout.mdp.9#

WARNING 1 [file em_steep_0.mdp, line 78]:
  Unknown left-hand 'pull_geometry' in parameter file



WARNING 2 [file em_steep_0.mdp, line 78]:
  Unknown left-hand 'pull_dim' in parameter file



WARNING 3 [file em_steep_0.mdp, line 78]:
  Unknown left-hand 'pull_start' in parameter file



WARNING 4 [file em_steep_0.mdp, line 78]:
  Unknown left-hand 'pull_init1' in parameter file



WARNING 5 [file em_steep_0.mdp, line 78]:
  Unknown left-hand 'pull_ngroups' in parameter file



WARNING 6 [file em_steep_0.mdp, line 78]:
  Unknown left-hand 'pull_group0' in parameter file



WARNING 7 [file em_steep_0.mdp, line 78]:
  Unknown left-hand 'pull_group1' in parameter file



WARNING 8 [file em_steep_0.mdp, line 78]:
  Unknown left-hand 'pull_k1' in parameter file



WARNING 9 [file em_steep_0.mdp, line 78]:
  Unknown left-hand 'pull_kB1' in parameter file



NOTE 1 [file em_steep_0.mdp]:
  With Verlet lists the optimal nstlist is >= 10, with GPUs >= 20. Note
  that with the Verlet scheme, nstlist has no effect on the accuracy of
  your simulation.

Setting the LD random seed to 1259182625
Generated 226 of the 1711 non-bonded parameter combinations

-------------------------------------------------------
Program gmx grompp, VERSION 5.1.4
Source code file:
/home/shahid/Software/gromacs-5.1.4/src/gromacs/gmxpreprocess/topio.c,
line: 755

Fatal error:
Syntax error - File complex.top, line 19956
Last line read:
'[ intermolecular_interactions ]'
Invalid order for directive intermolecular_interactions
For more information and tips for troubleshooting, please check the GROMACS
website at http://www.gromacs.org/Documentation/Errors
-------------------------------------------------------

Here are mdp parameters

; Options for the decoupling
sc-alpha                 = 0.5
sc-coul                  = no
sc-power                 = 1
sc-sigma                 = 0.3
couple-moltype           = DRG
couple-lambda0           = none
couple-lambda1           = vdw-q
couple-intramol          = no
nstdhdl                  = 10
; No velocities during EM
gen_vel                  = no
; options for bonds
constraints              = all_bonds
constraint-algorithm     = lincs
continuation             = no
lincs-order              = 12
pull           = no
pull_geometry  = distance
pull_dim       = N N Y
pull_start     = yes
pull_init1     = 0.654
pull_ngroups   = 1
pull_group0    = a_1444
pull_group1    = a_3118
pull_k1        = 0.0   ; kJ*mol^(-1)*nm^(-2)
pull_kB1       = 4184  ; kJ*mol^(-1)*nm^(-2)

since I have defined restraint parameters for protein and drug. then what
is need of pulling code and why?

Kindly resolve the issue.
Thanks in Advance.
Regards