[gmx-users] Re: Creating an atomtype where all nonbonded interactions are enumerated

Andrew DeYoung adeyoung at andrew.cmu.edu
Thu Mar 22 18:07:05 CET 2012


Hi Justin,

Thank you very much for all your time and help!

If you have time, can someone please help me understand how to read a .tpr
file?  I convertied my .tpr file to human-readable form using "gmxdump -s
topol.tpr".

Following an enumeration of all of the .mdp-like parameters, I see the
following section:

------------------------------------------------------------------
topology:
   name="IL"
   #atoms               = 6144
   molblock (0):
      moltype              = 0 "EMI"
      #molecules           = 256
      #atoms_mol           = 19
      #posres_xA           = 0
      #posres_xB           = 0
   molblock (1):
      moltype              = 1 "BF4"
      #molecules           = 256
      #atoms_mol           = 5
      #posres_xA           = 0
      #posres_xB           = 0
   ffparams:
      atnr=8
      ntypes=100
         functype[0]=LJ_SR, c6= 3.35274590e-03, c12= 3.95094276e-06
         functype[1]=LJ_SR, c6= 2.60915095e-03, c12= 3.84019631e-06
         functype[2]=LJ_SR, c6= 2.80388421e-03, c12= 4.30620821e-06
;...
         functype[63]=LJ_SR, c6= 0.00000000e+00, c12= 0.00000000e+00
         functype[64]=BONDS, b0A= 1.46600e-01, cbA= 2.82000e+05, b0B=
1.46600e-01, cbB= 2.82000e+05
;...
         functype[70]=BONDS, b0A= 1.08000e-01, cbA= 3.07106e+05, b0B=
1.08000e-01, cbB= 3.07106e+05
         functype[71]=ANGLES, thA= 1.26400e+02, ctA= 5.85200e+02, thB=
1.26400e+02, ctB= 5.85200e+02
;...
         functype[81]=ANGLES, thA= 1.12700e+02, ctA= 8.36800e+02, thB=
1.12700e+02, ctB= 8.36800e+02
         functype[82]=RBDIHS, rbcA[0]= 1.94599991e+01, rbcA[1]=
0.00000000e+00, rbcA[2]=-1.94599991e+01, rbcA[3]= 0.00000000e+00, rbcA[4]=
0.00000000e+00, rbcA[5]= 0.00000000e+00
rbcB[0]= 1.94599991e+01, rbcB[1]= 0.00000000e+00, rbcB[2]=-1.94599991e+01,
rbcB[3]= 0.00000000e+00, rbcB[4]= 0.00000000e+00, rbcB[5]= 0.00000000e+00
;...
         functype[88]=RBDIHS, rbcA[0]= 6.65499985e-01, rbcA[1]=
1.99650002e+00, rbcA[2]= 0.00000000e+00, rbcA[3]=-2.66199994e+00, rbcA[4]=
0.00000000e+00, rbcA[5]= 0.00000000e+00
rbcB[0]= 6.65499985e-01, rbcB[1]= 1.99650002e+00, rbcB[2]= 0.00000000e+00,
rbcB[3]=-2.66199994e+00, rbcB[4]= 0.00000000e+00, rbcB[5]= 0.00000000e+00
         functype[89]=LJ14, c6A= 1.30457547e-03, c12A= 1.92009816e-06, c6B=
1.30457547e-03, c12B= 1.92009816e-06
;...
         functype[97]=LJ14, c6A= 6.12890653e-05, c12A= 1.49631507e-08, c6B=
6.12890653e-05, c12B= 1.49631507e-08
         functype[98]=BONDS, b0A= 1.39300e-01, cbA= 2.42672e+05, b0B=
1.39300e-01, cbB= 2.42672e+05
         functype[99]=ANGLES, thA= 1.09500e+02, ctA= 4.18400e+02, thB=
1.09500e+02, ctB= 4.18400e+02
      reppow               = 12
      fudgeQQ              = 0.5
------------------------------------------------------------------

My question is, how can I tell which atomtypes or pairs (I am not sure if
the indices refer to atomtypes or pairs) each "functype[*]=LJ_SR, c6=*,
c12=*" refers to?  Each LJ_SR entry is indexed, but how do I determine what
the indices refer to?

In my .itp files which comprise my homemade force field, I define 22 unique
atomtypes, but here in the .tpr file, 63 LJ_SR entries are present.  Does
this mean that the LJ_SR entries refer to pairs of atoms?

Later in the .tpr file, the residue types are listed, along with the atom
names.  For example, 

------------------------------------------------------------------
   moltype (0):
      name="EMI"
      atoms:
         atom (19):
            atom[     0]={type=  0, typeB=  0, ptype=    Atom, m=
1.40067e+01, q= 1.50000e-01, mB= 1.40067e+01, qB= 1.50000e-01, resind=    0,
atomnumber=  7}
;...
            atom[    18]={type=  3, typeB=  3, ptype=    Atom, m=
1.00790e+00, q= 6.00000e-02, mB= 1.00790e+00, qB= 6.00000e-02, resind=    0,
atomnumber=  1}
         atom (19):
            atom[0]={name="NA"}
;...
            atom[18]={name="HC"}
         type (19):
            type[0]={name="NA",nameB="NA"}
;...
            type[18]={name="HC",nameB="HC"}
         residue (1):
            residue[0]={name="EMI", nr=1, ic=' '}
      cgs:
         nr=13
         cgs[0]={0..1}
;...
         cgs[12]={16..18}
      excls:
         nr=19
         nra=195
         excls[0][0..12]={0, 1, 2, 3, 4, 5, 6, 8, 9, 10, 11, 12, 13}
;...
         excls[18][187..194]={1, 6, 7, 14, 15, 16, 17, 18}
      Bond:
         nr: 57
         iatoms:
            0 type=64 (BONDS) 0 2
;...
            18 type=69 (BONDS) 7 18
      Angle:
         nr: 132
         iatoms:
            0 type=71 (ANGLES) 2 0 3
;...
            32 type=75 (ANGLES) 17 7 18
      Ryckaert-Bell.:
         nr: 165
         iatoms:
            0 type=82 (RBDIHS) 0 3 1 5
;...
            32 type=88 (RBDIHS) 15 6 7 18
      LJ-14:
         nr: 108
         iatoms:
            0 type=89 (LJ14) 0 6
;...
            35 type=97 (LJ14) 15 18
------------------------------------------------------------------

and similarly for the other residues.  But, I do not see LJ_SR entries in
these residue sections.  So, how can I determine what the LJ_SR entries are
referring to near the beginning of the file?  

I am sorry if this a silly question (which it might be).

Thank you for your time!

Andrew DeYoung
Carnegie Mellon University

-----Original Message-----
Andrew DeYoung wrote:
> Justin,
> 
> Thank you so much for your help!  That was really helpful.
> 
> It seems that I simply append my [ nonbond_params ] section to the end of
my
> ffnonbonded.itp file.  However, if I ONLY have this:
> 
> ---
> ; ffnonbonded.itp
> [ atomtypes ]
> ;type at.n	mass	charge  ptype  sig  eps
> ; all of the liquid atomtypes
> ; ...
> 
> [ nonbond_params ]
> ; i j func sig eps
> mywall N 1 0.29750 4.56912
> mywall C 1 0.31200 3.07524
> mywall H 1 0.26170 1.17675
> mywall B 1 0.31225 3.48667
> mywall F 1 0.28974 2.82420
> ---
> 
> then I get an error message saying that atomtype mywall does not exist.
> However, if I simply add a bogus entry for mywall in the [ atomtypes ]
> section of ffnonbonded.itp, then grompp runs fine:
> 

Right, you can't use an atomtype unless it is defined.

> ---
> ; ffnonbonded.itp
> [ atomtypes ]
> ;type at.n	mass	charge  ptype  sig  eps
> mywall 6  12.01100  0.000 A 0.00000e-1 0.00000e-1 ; bogus sig and eps
> ; all of the liquid atomtypes
> ; ...
> 
> [ nonbond_params ]
> ; i j func sig eps
> mywall N 1 0.29750 4.56912
> mywall C 1 0.31200 3.07524
> mywall H 1 0.26170 1.17675
> mywall B 1 0.31225 3.48667
> mywall F 1 0.28974 2.82420
> ---
> 
> Is there any way that I can verify that the code is actually using the
> parameters for mywall in the [ nonbond_params ] section, rather than the
> bogus parameters (0 for sig and eps) for mywall in the [ atomtypes ]
> section?
> 

Run gmxdump on your .tpr file.  All the relevant parameters should be in
there somewhere.

-Justin




More information about the gromacs.org_gmx-users mailing list