[gmx-users] Dihedral calculations in topology
Dear GMX users, I have a few questions about something I don’t fully understand. I am trying to work with polymers, where the basics are written by Justin Lemkul via this link: https://mailman-1.sys.kth.se/pipermail/gromacs.org_gmx-users/2009-March/040125.html This works up to the point where I want to use grompp, where it gives errors regarding the Proper Dihedral types; No default Proper Dih. types. When consulting the lines where these errors refers to, it is caused by the three dihedrals defined in my topology: [ dihedrals ] ; aiajakal functc0c1 c2c3c4c5 2 1 5 8 1 1 5 811 1 5 81112 1 Now, according to the manual table 5.5, page 138 (v. 2016.1), these three dihedrals are defined as funct=1, which corresponds with a proper dihedral, while by default the OPLS ff only knows R-B dihedrals, which are funct=3. Up to this point my story is correct right? Then, why does pdb2gmx select this dihedral type? In the above topology [ai aj ak al] = [1 5 8 11] refers to a [opls_135 opls_136 opls_136 opls_135] dihedral, which I believe is a [CT CT CT CT] dihedral. This dihedral is already defined by ffbonded.itp as follows: CT CT CT CT 3 2.92880 -1.46440 0.20920 -1.67360 0.0 0.0 ; hydrocarbon all-atom Am I correct about this part? Is it valid to manually change the funct=1 in the above topology to funct=3 to use this dihedral type? And if it is valid, how can I prevent pdb2gmx to assign a funct=1 to the topology? Mark -- Gromacs Users mailing list * Please search the archive at http://www.gromacs.org/Support/Mailing_Lists/GMX-Users_List before posting! * Can't post? Read http://www.gromacs.org/Support/Mailing_Lists * For (un)subscribe requests visit https://maillist.sys.kth.se/mailman/listinfo/gromacs.org_gmx-users or send a mail to gmx-users-requ...@gromacs.org.
Re: [gmx-users] Dihedral calculations in topology
Hi, On Wed, Dec 21, 2016 at 9:17 PM Kamps, M. wrote: > Dear GMX users, > > I have a few questions about something I don’t fully understand. I am > trying to work with polymers, where the basics are written by Justin > Lemkul via this link: > > https://mailman-1.sys.kth.se/pipermail/gromacs.org_gmx-users/2009-March/040125.html > > This works up to the point where I want to use grompp, where it gives > errors regarding the Proper Dihedral types; No default Proper Dih. > types. When consulting the lines where these errors refers to, it is > caused by the three dihedrals defined in my topology: > > [ dihedrals ] > ; aiajakal functc0c1 > c2c3c4c5 > 2 1 5 8 1 > 1 5 811 1 > 5 81112 1 > > Now, according to the manual table 5.5, page 138 (v. 2016.1), these > three dihedrals are defined as funct=1, which corresponds with a > proper dihedral, while by default the OPLS ff only knows R-B > dihedrals, which are funct=3. Up to this point my story is correct > right? > The normal [bondedtypes] defaults in oplsaa.ff/aminoacids.rtp uses dihedral type 3, yes. Whether that's applicable to your case we can't tell. Which file are your .rtp entries coming from? The terminal output of pdb2gmx and grompp say a lot of useful things about which it would be nice not to guess. Mark Then, why does pdb2gmx select this dihedral type? In the above > topology [ai aj ak al] = [1 5 8 11] refers to a [opls_135 opls_136 > opls_136 opls_135] dihedral, which I believe is a [CT CT CT CT] > dihedral. This dihedral is already defined by ffbonded.itp as follows: > CT CT CT CT 3 2.92880 -1.46440 0.20920 > -1.67360 0.0 0.0 ; hydrocarbon all-atom > Am I correct about this part? > > Is it valid to manually change the funct=1 in the above topology to > funct=3 to use this dihedral type? And if it is valid, how can I > prevent pdb2gmx to assign a funct=1 to the topology? > > Mark > -- > Gromacs Users mailing list > > * Please search the archive at > http://www.gromacs.org/Support/Mailing_Lists/GMX-Users_List before > posting! > > * Can't post? Read http://www.gromacs.org/Support/Mailing_Lists > > * For (un)subscribe requests visit > https://maillist.sys.kth.se/mailman/listinfo/gromacs.org_gmx-users or > send a mail to gmx-users-requ...@gromacs.org. -- Gromacs Users mailing list * Please search the archive at http://www.gromacs.org/Support/Mailing_Lists/GMX-Users_List before posting! * Can't post? Read http://www.gromacs.org/Support/Mailing_Lists * For (un)subscribe requests visit https://maillist.sys.kth.se/mailman/listinfo/gromacs.org_gmx-users or send a mail to gmx-users-requ...@gromacs.org.
Re: [gmx-users] Dihedral calculations in topology
Dear Mark, Thanks for your reply! I am replicating Polyethylene, of which the [polymer.rtp] and the [polymer.hdb] are exactly the same as described by Justin Lemkul in this link: https://mailman-1.sys.kth.se/pipermail/gromacs.org_gmx-users/2009-March/040125.html So, my rtp entries are: ; Terminal PE residue (chain "begin") ;C1 is -CH3 [ EthB ] [ atoms ] C1opls_135 -0.180 1 H11opls_1400.060 1 H12opls_1400.060 1 H13opls_1400.0601 C2opls_136 -0.120 2 H21opls_1400.060 2 H22opls_1400.060 2 [ bonds ] C1 H11 C1 H12 C1 H13 C1 C2 C2 H21 C2 H22 C2 +C1 ; Terminal PE residue (chain "end") ;C2 is -CH3 [ EthE ] [ atoms ] C1opls_136 -0.120 1 H11opls_1400.060 1 H12opls_1400.060 1 C2opls_135 -0.180 2 H21opls_1400.060 2 H22opls_1400.060 2 H23opls_1400.0602 [ bonds ] C1 -C2 C1 H11 C1 H12 C1 C2 C2 H21 C2 H22 C2 H23 And the used hdb entries are: EthB 2 3 4 H1 C1 C2 +C1 2 6 H2 C2 C1 +C1 EthE 2 2 6 H1 C1 C2 -C2 3 4 H2 C2 C1 -C2 The used pdb file is: ATOM 1 C1 EthB1 0.000 0.000 0.000 ATOM 2 C2 EthB1 1.273 0.847 0.000 ATOM 3 C1 EthE1 2.546 0.000 0.000 ATOM 4 C2 EthE1 3.818 0.847 0.000 Which is a simple C4H10 molecule, composed out of two residues, with a total of three dihedrals. The output given by grompp is: Setting the LD random seed to -962199660 Generated 337431 of the 337431 non-bonded parameter combinations Generating 1-4 interactions: fudge = 0.5 Generated 337431 of the 337431 1-4 parameter combinations ERROR 1 [file topol_polymer.itp, line 119]: No default Proper Dih. types ERROR 2 [file topol_polymer.itp, line 120]: No default Proper Dih. types ERROR 3 [file topol_polymer.itp, line 121]: No default Proper Dih. types Where lines 119 till 121 correspond with the three dihedrals in my topology (of the polymer). The output of pdb2gmx is as follows: Opening force field file ./oplsaa_polymer.ff/aminoacids.r2b Reading C4H10.pdb... Read 4 atoms Analyzing pdb file Splitting chemical chains based on TER records or chain id changing. There are 1 chains and 0 blocks of water and 2 residues with 4 atoms chain #res #atoms 1 ' ' 2 4 All occupancy fields zero. This is probably not an X-Ray structure Opening force field file ./oplsaa_polymer.ff/atomtypes.atp Atomtype 823 Reading residue database... (oplsaa_polymer) Opening force field file ./oplsaa_polymer.ff/Polymer.rtp Reading .rtp file without '[ bondedtypes ]' directive, Will proceed as if the entry was: [ bondedtypes ] ; bonds angles dihedrals impropers all_dihedrals nr_exclusions HH14 remove_dih 1 1 1 2 0 3 1 1 Residue 4 Sorting it all out... Opening force field file ./oplsaa_polymer.ff/aminoacids.rtp Residue 55 Sorting it all out... Opening force field file ./oplsaa_polymer.ff/surface.rtp Using default: not generating all possible dihedrals Using default: excluding 3 bonded neighbors Using default: generating 1,4 H--H interactions Using default: removing proper dihedrals found on the same bond as a proper dihedral Residue 62Warning: file does not end with a newline, last line: PTPT 0.000 0 Residue 63 Sorting it all out... Opening force field file ./oplsaa_polymer.ff/Polymer.hdb Opening force field file ./oplsaa_polymer.ff/aminoacids.hdb Opening force field file ./oplsaa_polymer.ff/aminoacids.n.tdb Opening force field file ./oplsaa_polymer.ff/aminoacids.c.tdb Processing chain 1 (4 atoms, 2 residues) Warning: Starting residue EthB1 in chain not identified as Protein/RNA/DNA. Warning: Starting residue EthE1 in chain not identified as Protein/RNA/DNA. Problem with chain definition, or missing terminal residues. This chain does not appear to contain a recognized chain molecule. If this is incorrect, you can edit residuetypes.dat to modify the behavior. 8 out of 8 lines of specbond.dat converted successfully Checking for duplicate atoms Generating any missing hydrogen atoms and/or adding termini. Now there are 2 residues with 14 atoms Making bonds... Number of bonds was 14, now 13 Generating angles, dihedrals and pairs... Before cleaning: 27 pairs Before cleaning: 27 dihedrals Making cmap torsions... There are3 dihedrals,0 impropers, 24 angles 27 pairs, 13 bonds and 0 virtual sites Total mass 58.124 a.m.u. Total charge -0.000 e Writing topology Writing coordinate file... -- Gromacs Users mailing list * Please search the archive at http://www.gromacs.org/Support/Mailing_Lists/GMX-Users_List before posting! * Can't post? Read http://www.gromacs.org/Support/Mailing
Re: [gmx-users] Dihedral calculations in topology
Hi, pdb2gmx told you that it saw that your rtp file didn't have default bondedtypes, so it was going to use certain values, and that means default dihedral type 1. If you want to use the same defaults as aminoacids.rtp, then copy them over to the rtp file you're using :-) Mark On Wed, Dec 21, 2016 at 10:12 PM Kamps, M. wrote: > Dear Mark, > > Thanks for your reply! > I am replicating Polyethylene, of which the [polymer.rtp] and the > [polymer.hdb] are exactly the same as described by Justin Lemkul in > this link: > https://mailman-1.sys.kth.se/pipermail/gromacs.org_gmx-users/2009-March/040125.html > > So, my rtp entries are: > ; Terminal PE residue (chain "begin") > ;C1 is -CH3 > [ EthB ] > [ atoms ] > C1opls_135 -0.180 1 >H11opls_1400.060 1 >H12opls_1400.060 1 >H13opls_1400.0601 > C2opls_136 -0.120 2 >H21opls_1400.060 2 >H22opls_1400.060 2 > [ bonds ] > C1 H11 > C1 H12 > C1 H13 > C1 C2 > C2 H21 > C2 H22 > C2 +C1 > ; Terminal PE residue (chain "end") > ;C2 is -CH3 > [ EthE ] > [ atoms ] > C1opls_136 -0.120 1 >H11opls_1400.060 1 >H12opls_1400.060 1 > C2opls_135 -0.180 2 >H21opls_1400.060 2 >H22opls_1400.060 2 >H23opls_1400.0602 > [ bonds ] > C1 -C2 > C1 H11 > C1 H12 > C1 C2 > C2 H21 > C2 H22 > C2 H23 > > And the used hdb entries are: > EthB 2 > 3 4 H1 C1 C2 +C1 > 2 6 H2 C2 C1 +C1 > EthE 2 > 2 6 H1 C1 C2 -C2 > 3 4 H2 C2 C1 -C2 > > The used pdb file is: > ATOM 1 C1 EthB1 0.000 0.000 0.000 > ATOM 2 C2 EthB1 1.273 0.847 0.000 > ATOM 3 C1 EthE1 2.546 0.000 0.000 > ATOM 4 C2 EthE1 3.818 0.847 0.000 > > Which is a simple C4H10 molecule, composed out of two residues, with a > total of three dihedrals. The output given by grompp is: > Setting the LD random seed to -962199660 > Generated 337431 of the 337431 non-bonded parameter combinations > Generating 1-4 interactions: fudge = 0.5 > Generated 337431 of the 337431 1-4 parameter combinations > ERROR 1 [file topol_polymer.itp, line 119]: > No default Proper Dih. types > ERROR 2 [file topol_polymer.itp, line 120]: > No default Proper Dih. types > ERROR 3 [file topol_polymer.itp, line 121]: > No default Proper Dih. types > Where lines 119 till 121 correspond with the three dihedrals in my > topology (of the polymer). > > The output of pdb2gmx is as follows: > Opening force field file ./oplsaa_polymer.ff/aminoacids.r2b > Reading C4H10.pdb... > Read 4 atoms > Analyzing pdb file > Splitting chemical chains based on TER records or chain id changing. > There are 1 chains and 0 blocks of water and 2 residues with 4 atoms > > chain #res #atoms > 1 ' ' 2 4 > > All occupancy fields zero. This is probably not an X-Ray structure > Opening force field file ./oplsaa_polymer.ff/atomtypes.atp > Atomtype 823 > Reading residue database... (oplsaa_polymer) > Opening force field file ./oplsaa_polymer.ff/Polymer.rtp > Reading .rtp file without '[ bondedtypes ]' directive, > Will proceed as if the entry was: > [ bondedtypes ] > ; bonds angles dihedrals impropers all_dihedrals nr_exclusions > HH14 remove_dih > 1 1 1 2 0 3 >1 1 > > Residue 4 > Sorting it all out... > Opening force field file ./oplsaa_polymer.ff/aminoacids.rtp > Residue 55 > Sorting it all out... > Opening force field file ./oplsaa_polymer.ff/surface.rtp > Using default: not generating all possible dihedrals > Using default: excluding 3 bonded neighbors > Using default: generating 1,4 H--H interactions > Using default: removing proper dihedrals found on the same bond as a > proper dihedral > Residue 62Warning: file does not end with a newline, last line: > PTPT 0.000 0 > Residue 63 > Sorting it all out... > Opening force field file ./oplsaa_polymer.ff/Polymer.hdb > Opening force field file ./oplsaa_polymer.ff/aminoacids.hdb > Opening force field file ./oplsaa_polymer.ff/aminoacids.n.tdb > Opening force field file ./oplsaa_polymer.ff/aminoacids.c.tdb > Processing chain 1 (4 atoms, 2 residues) > Warning: Starting residue EthB1 in chain not identified as Protein/RNA/DNA. > Warning: Starting residue EthE1 in chain not identified as Protein/RNA/DNA. > Problem with chain definition, or missing terminal residues. > This chain does not appear to contain a recognized chain molecule. > If this is incorrect, you can edit residuetypes.dat to modify the behavior. > 8 out of 8 lines of specbond.dat converted successfully > Checking for duplicate atoms > Generating any missing hydrogen atoms and/or adding termini. > Now there are 2 residues with
Re: [gmx-users] Dihedral calculations in topology
Mark, Thanks again for the swift reply! Stupid of me to miss that message! I noticed some other code that was generated as output of pdb2gmx, which I do not fully understand. Can you maybe explain some of it? Processing chain 1 (4 atoms, 2 residues) Warning: Starting residue EthB1 in chain not identified as Protein/RNA/DNA. Warning: Starting residue EthE1 in chain not identified as Protein/RNA/DNA. Problem with chain definition, or missing terminal residues. This chain does not appear to contain a recognized chain molecule. If this is incorrect, you can edit residuetypes.dat to modify the behavior. I understand the two warnings. To my directory I added a file residuetypes.dat, where I specificy that (among others) EthB and EthE (which are the ending and beginning residues of polyethylene) are Polymer. What does the rest of this message mean? I don't understand how I could fix this message? Or even what it wants to do, since i've defined my own terminal residues. Is it not possible to define my own 'group'? residuetypes.dat: EthBPolymer EthEPolymer Mark -- Gromacs Users mailing list * Please search the archive at http://www.gromacs.org/Support/Mailing_Lists/GMX-Users_List before posting! * Can't post? Read http://www.gromacs.org/Support/Mailing_Lists * For (un)subscribe requests visit https://maillist.sys.kth.se/mailman/listinfo/gromacs.org_gmx-users or send a mail to gmx-users-requ...@gromacs.org.
Re: [gmx-users] Dihedral calculations in topology
On 12/21/16 6:20 AM, Mark Abraham wrote: Hi, pdb2gmx told you that it saw that your rtp file didn't have default bondedtypes, so it was going to use certain values, and that means default dihedral type 1. If you want to use the same defaults as aminoacids.rtp, then copy them over to the rtp file you're using :-) Which is (just for the sake of being pedantic) exactly what the linked post says to do. My .rtp entries there are just the residues that are needed. That is not intended or stated to be a complete .rtp file. -Justin Mark On Wed, Dec 21, 2016 at 10:12 PM Kamps, M. wrote: Dear Mark, Thanks for your reply! I am replicating Polyethylene, of which the [polymer.rtp] and the [polymer.hdb] are exactly the same as described by Justin Lemkul in this link: https://mailman-1.sys.kth.se/pipermail/gromacs.org_gmx-users/2009-March/040125.html So, my rtp entries are: ; Terminal PE residue (chain "begin") ;C1 is -CH3 [ EthB ] [ atoms ] C1opls_135 -0.180 1 H11opls_1400.060 1 H12opls_1400.060 1 H13opls_1400.0601 C2opls_136 -0.120 2 H21opls_1400.060 2 H22opls_1400.060 2 [ bonds ] C1 H11 C1 H12 C1 H13 C1 C2 C2 H21 C2 H22 C2 +C1 ; Terminal PE residue (chain "end") ;C2 is -CH3 [ EthE ] [ atoms ] C1opls_136 -0.120 1 H11opls_1400.060 1 H12opls_1400.060 1 C2opls_135 -0.180 2 H21opls_1400.060 2 H22opls_1400.060 2 H23opls_1400.0602 [ bonds ] C1 -C2 C1 H11 C1 H12 C1 C2 C2 H21 C2 H22 C2 H23 And the used hdb entries are: EthB 2 3 4 H1 C1 C2 +C1 2 6 H2 C2 C1 +C1 EthE 2 2 6 H1 C1 C2 -C2 3 4 H2 C2 C1 -C2 The used pdb file is: ATOM 1 C1 EthB1 0.000 0.000 0.000 ATOM 2 C2 EthB1 1.273 0.847 0.000 ATOM 3 C1 EthE1 2.546 0.000 0.000 ATOM 4 C2 EthE1 3.818 0.847 0.000 Which is a simple C4H10 molecule, composed out of two residues, with a total of three dihedrals. The output given by grompp is: Setting the LD random seed to -962199660 Generated 337431 of the 337431 non-bonded parameter combinations Generating 1-4 interactions: fudge = 0.5 Generated 337431 of the 337431 1-4 parameter combinations ERROR 1 [file topol_polymer.itp, line 119]: No default Proper Dih. types ERROR 2 [file topol_polymer.itp, line 120]: No default Proper Dih. types ERROR 3 [file topol_polymer.itp, line 121]: No default Proper Dih. types Where lines 119 till 121 correspond with the three dihedrals in my topology (of the polymer). The output of pdb2gmx is as follows: Opening force field file ./oplsaa_polymer.ff/aminoacids.r2b Reading C4H10.pdb... Read 4 atoms Analyzing pdb file Splitting chemical chains based on TER records or chain id changing. There are 1 chains and 0 blocks of water and 2 residues with 4 atoms chain #res #atoms 1 ' ' 2 4 All occupancy fields zero. This is probably not an X-Ray structure Opening force field file ./oplsaa_polymer.ff/atomtypes.atp Atomtype 823 Reading residue database... (oplsaa_polymer) Opening force field file ./oplsaa_polymer.ff/Polymer.rtp Reading .rtp file without '[ bondedtypes ]' directive, Will proceed as if the entry was: [ bondedtypes ] ; bonds angles dihedrals impropers all_dihedrals nr_exclusions HH14 remove_dih 1 1 1 2 0 3 1 1 Residue 4 Sorting it all out... Opening force field file ./oplsaa_polymer.ff/aminoacids.rtp Residue 55 Sorting it all out... Opening force field file ./oplsaa_polymer.ff/surface.rtp Using default: not generating all possible dihedrals Using default: excluding 3 bonded neighbors Using default: generating 1,4 H--H interactions Using default: removing proper dihedrals found on the same bond as a proper dihedral Residue 62Warning: file does not end with a newline, last line: PTPT 0.000 0 Residue 63 Sorting it all out... Opening force field file ./oplsaa_polymer.ff/Polymer.hdb Opening force field file ./oplsaa_polymer.ff/aminoacids.hdb Opening force field file ./oplsaa_polymer.ff/aminoacids.n.tdb Opening force field file ./oplsaa_polymer.ff/aminoacids.c.tdb Processing chain 1 (4 atoms, 2 residues) Warning: Starting residue EthB1 in chain not identified as Protein/RNA/DNA. Warning: Starting residue EthE1 in chain not identified as Protein/RNA/DNA. Problem with chain definition, or missing terminal residues. This chain does not appear to contain a recognized chain molecule. If this is incorrect, you can edit residuetypes.dat to modify the behavior. 8 out of 8 lines of specbond.dat converted successfully Checking for duplicate atoms Generating any missing hydrogen atoms and/or adding termini. Now ther
Re: [gmx-users] Dihedral calculations in topology
On 12/21/16 9:21 AM, Kamps, M. wrote: Mark, Thanks again for the swift reply! Stupid of me to miss that message! I noticed some other code that was generated as output of pdb2gmx, which I do not fully understand. Can you maybe explain some of it? Processing chain 1 (4 atoms, 2 residues) Warning: Starting residue EthB1 in chain not identified as Protein/RNA/DNA. Warning: Starting residue EthE1 in chain not identified as Protein/RNA/DNA. Problem with chain definition, or missing terminal residues. This chain does not appear to contain a recognized chain molecule. If this is incorrect, you can edit residuetypes.dat to modify the behavior. I understand the two warnings. To my directory I added a file residuetypes.dat, where I specificy that (among others) EthB and EthE (which are the ending and beginning residues of polyethylene) are Polymer. What does the rest of this message mean? I don't understand how I could fix this message? Or even what it wants to do, since i've defined my own terminal residues. Is it not possible to define my own 'group'? residuetypes.dat: EthBPolymer EthEPolymer There are certain types of molecules that pdb2gmx understands. If you're not doing something that is Protein/DNA/RNA, then just call it "Other" and it will work fine. The "warning" probably shouldn't be called that. It's merely information that is used to point out instances when the user may not be doing something consistent, e.g. a modified amino acid hasn't been added as protein, so there's a discontinuity in the chain. If all you have is polymer, that's fine. -Justin -- == Justin A. Lemkul, Ph.D. Ruth L. Kirschstein NRSA Postdoctoral Fellow Department of Pharmaceutical Sciences School of Pharmacy Health Sciences Facility II, Room 629 University of Maryland, Baltimore 20 Penn St. Baltimore, MD 21201 jalem...@outerbanks.umaryland.edu | (410) 706-7441 http://mackerell.umaryland.edu/~jalemkul == -- Gromacs Users mailing list * Please search the archive at http://www.gromacs.org/Support/Mailing_Lists/GMX-Users_List before posting! * Can't post? Read http://www.gromacs.org/Support/Mailing_Lists * For (un)subscribe requests visit https://maillist.sys.kth.se/mailman/listinfo/gromacs.org_gmx-users or send a mail to gmx-users-requ...@gromacs.org.