[gmx-users] problem in topology with AMBER99SB

Justin A. Lemkul jalemkul at vt.edu
Mon Sep 5 19:31:46 CEST 2011



Anna Marabotti wrote:
> Dear gmx-users,
> I'm trying to create the topology for a ligand using Amber99SB force 
> field, but I'm experimenting several problems. Here what I did:
>  - I recovered the .mol2 file of my ligand (comp1.mol2)
> - I checked and added H with Amber tool reduce: reduce comp1.mol2 > 
> comp1-H.mol2
> - I added charges using antechamber: antechamber -i comp1-H.mol2 -fi 
> mol2 -o comp1-H_C.mol2 -fo mol2 -c bcc -s2
> - I created the frcmod parameter file: parmcheck -i comp1-H_C.mol2 -f 
> mol2 -o comp1-H_C.frcmod 
> - then I used: xleap -s -f /opt/amber11/dat/leap/cmd/leaprc.ff99SB
>                             source leaprc.gaff
>                             lig = loadmol2 comp1-H_C.mol2
>                             check lig
>                             loadamberparams comp1-H_C.frcmod
>                             saveamberparm lig comp1.prmtop comp1.inpcrd
> - I created the .top file with amb2gmx.pl: amb2gmx.pl --prmtop 
> comp1.prmtop --crd comp1.inpcrd --outname comp1
> then
> mv comp1.top comp1.itp
> 
> My comp1.itp file is as following:
> 
> ; comp1.top created by rdparm2gmx.pl Mon Sep 5 12:45:09 CEST 2011
> 
> [ defaults ]
> 
> ; nbfunc comb-rule gen-pairs fudgeLJ fudgeQQ
> 
> 1 2 yes 0.5 0.8333
> 
> [ atomtypes ]
> 
> ;name bond_type mass charge ptype sigma epsilon
> 
> c3 c3 0.0000 0.0000 A 3.39967e-01 4.57730e-01
> 
> n n 0.0000 0.0000 A 3.25000e-01 7.11280e-01
> 
> ha ha 0.0000 0.0000 A 2.59964e-01 6.27600e-02
> 
> ca ca 0.0000 0.0000 A 3.39967e-01 3.59824e-01
> 
> hn hn 0.0000 0.0000 A 1.06908e-01 6.56888e-02
> 
> hc hc 0.0000 0.0000 A 2.64953e-01 6.56888e-02
> 
> c c 0.0000 0.0000 A 3.39967e-01 3.59824e-01
> 
> o o 0.0000 0.0000 A 2.95992e-01 8.78640e-01
> 
> [ moleculetype ]
> 
> ; Name nrexcl
> 
> solute 3
> 
> [ atoms ]
> 
> ; nr type resnr residue atom cgnr charge mass typeB chargeB
> 
> 1 o 1 UNK O1 1 -0.58510 16.000000
> 
> 2 o 1 UNK O2 2 -0.47810 16.000000
> 
> 3 o 1 UNK O3 3 -0.47710 16.000000
> 
> 4 n 1 UNK N1 4 -0.47010 14.000000
> 
> ...............................................................................
> 
> [ bonds ]
> 
> ; ai aj funct r k
> 
> 4 33 1 1.0090e-01 3.4326e+05
> 
> 7 26 1 1.0920e-01 2.8225e+05
> 
> 7 25 1 1.0920e-01 2.8225e+05
> 
> 7 24 1 1.0920e-01 2.8225e+05
> 
> ...............................................................................
> 
> [ pairs ]
> 
> ; ai aj funct
> 
> 1 33 1
> 
> 4 36 1
> 
> 4 34 1
> 
> 5 33 1
> 
> ...............................................................................
> 
> [ angles ]
> 
> ; ai aj ak funct theta cth
> 
> 5 9 32 1 1.1005e+02 3.8802e+02
> 
> 5 9 31 1 1.1005e+02 3.8802e+02
> 
> 5 9 30 1 1.1005e+02 3.8802e+02
> 
> 5 8 29 1 1.1005e+02 3.8802e+02
> 
> ................................................................................
> 
> [ dihedrals ]
> 
> ;i j k l func C0 ... C5
> 
> 1 6 4 33 3 29.28800 -8.36800 -20.92000 0.00000 0.00000 0.00000 ;
> 
> 4 12 17 36 3 30.33400 0.00000 -30.33400 0.00000 0.00000 0.00000 ;
> 
> 4 12 14 34 3 30.33400 0.00000 -30.33400 0.00000 0.00000 0.00000 ;
> 
> 5 6 4 33 3 20.92000 0.00000 -20.92000 0.00000 0.00000 0.00000 ;
> 
> 6 5 9 32 3 0.65270 1.95811 0.00000 -2.61082 0.00000 0.00000 ;
> 
> ..............................................................................
> 
> [ system ]
> 
> 40 system
> 
> [ molecules ]
> 
> ; Compound nmols
> 
> solute 1
> 
>  
> 
> I pasted the coordinates of comp1.pdb into those of protein.pdb, II 
> removed directives [ system ] and [ molecules ] and I included 
> "comp1.itp" in the "topol.top" file obtained with pdb2gmx (Gromacs 
> 4.5.4) using ff Amber99SB, on the protein alone:
> 
> topol.top:
> 
> ; Include forcefield parameters
> 
> #include "amber99sb.ff/forcefield.itp"
> 
> #include "comp1.itp"
> 
> ; Include chain topologies
> 
> #include "topol_Protein_chain_C.itp"
> 
> #include "topol_Other.itp"
> 
> ; Include water topology
> 
> #include "amber99sb.ff/tip3p.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 "amber99sb.ff/ions.itp"
> 
> [ system ]
> 
> ; Name
> 
> Protein in water
> 
> [ molecules ]
> 
> ; Compound #mols
> 
> Protein_chain_C 1
> 
> Other 1
> 
> UNK 1
> 
>  
> 
> Then, I launched editconf-genbox-grompp to neutralize system:
> 
> /opt/gromacs45/bin/grompp -f em.mdp -c prot-lig_solv.gro 
> -o prot-lig_solv.tpr -p
> 
> The result was:
> 
> Fatal error:
> 
> Syntax error - File comp1.itp, line 3
> 
> Last line read:
> 
> '[ defaults ]'
> 
> Invalid order for directive defaults
> 
>  
> 
> I had a check on the manual and it seems to me that all directives are 
> in the correct order! However, I removed
> 

The [defaults] directive can appear only once in a topology and must be the very 
first entry.  There is one in forcefield.itp and then you added another in your 
ligand .itp file.  That's why grompp complains.

> directive [ defaults ] and re-launched grompp:
> 
> /opt/gromacs45/bin/grompp -f em.mdp -c prot-lig_solv.gro 
> -o prot-lig_solv.tpr -p
> 
> Result:
> 
> Fatal error:
> 
> Syntax error - File comp1.itp, line 7
> 
> Last line read:
> 
> '[ atomtypes ]'
> 
> Invalid order for directive atomtypes
> 
>  
> 
> I removed directive [ atomtypes ] and re-launched grompp:
> 
> /opt/gromacs45/bin/grompp -f em.mdp -c prot-lig_solv.gro 
> -o prot-lig_solv.tpr -p
> 
> Result:
> 
> Fatal error:
> 
> Atomtype n not found
> 
>  
> 
> What's wrong? Why Gromacs claims that directives are not in correct 
> order? Could you please help me?
> 

There is a very specific hierarchy that must be maintained in a topology.  Force 
field-level directives must appear before any molecules are declared, system 
directives must follow all molecule declarations, etc.

You need the [atomtypes] directive for the topology to be valid, so keep it in 
there and work with that topology.  It seems to me that the order should be 
correct, as shown, since the ligand's new [atomtypes] have to precede all 
[moleculetype] directives.  I have a similar case that I'm actually working on 
right now that works fine, i.e.:

; Include forcefield parameters
#include "amber03.ff/forcefield.itp"

; Include ligand
#include "ligand.itp"

...then the protein stuff from pdb2gmx.  The ligand.itp file contains an 
[atomtypes] directive, which, in this position, is legal.  Which version of 
Gromacs are you using?  If it's not 4.5.4, then I'll suggest that you upgrade 
and try again.  If it is 4.5.4, then check the contents of the ligand .itp file 
carefully to make sure there aren't any remnants of the [defaults] directive 
lying around.

-Justin

-- 
========================================

Justin A. Lemkul
Ph.D. Candidate
ICTAS Doctoral Scholar
MILES-IGERT Trainee
Department of Biochemistry
Virginia Tech
Blacksburg, VA
jalemkul[at]vt.edu | (540) 231-9080
http://www.bevanlab.biochem.vt.edu/Pages/Personal/justin

========================================



More information about the gromacs.org_gmx-users mailing list