[gmx-users] Error in free energy calculation
Hannes Loeffler
Hannes.Loeffler at stfc.ac.uk
Thu Mar 16 12:19:44 CET 2017
I can't help you much because I have not really paid attention how
restraints in alchemical free energy simulations are handled in newer
versions of Gromacs. I understand that [intermolecular_interactions]
is a new section that allows using restraints based on global indices.
The error message suggest that you have put this section at the wrong
place in the topology file. I guess the manual will have information
on where grompp expects this to be.
The warnings come from the fact to you have set pull=no but other pull
control parameters are found (and ignored as not needed).
I guess from your input that you are trying to use something like
Stefan Boresch' VBA method. The bond restraint seems excessively
high. I also note that you attempt to recouple the DRG molecule in one
step (couple-lambda1 = vdw-q). This will probably mean that you should
also use a softcore potential for the Coulomb interactions.
Alternatively, you may want to recouple electrostatic and vdW
interactions separately. There is literature that suggests the latter
approach to be more efficient.
On Thu, 16 Mar 2017 15:19:48 +0530
Tasneem Kausar <tasneemkausar12 at gmail.com> wrote:
> 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
More information about the gromacs.org_gmx-users
mailing list