[gmx-users] GAFF and alkynes (Justin Lemkul)

Rebeca García Fandiño regafan at hotmail.com
Fri Feb 13 20:52:36 CET 2015


Hi,
sorry for insisting about this issue...but it is not clear for me yet. I have increased (in 10 times) the force constant of the two angles involved in my molecule (those involving the alkyne) and I obtain a correct geometry without no errors during the simulation, and perfectly linear.

Is there anything dangerous using this approximation? I don´t understand very well Justin´s answer to my previous suggestion:

> A restraint is just a biasing potential towards some desired outcome.  It's not 
> necessarily even different in form from a normal harmonic interaction, but I'd 
> be hesitant to try to use some angle definition and a restraint on top of that. 
>   The GROMACS code indicates there is a linear angle potential, but it doesn't 
> seem to be documented.  This was likely added to circumvent the previous 
> problems with these linear species.

You suggest that maybe it is not a good idea but I cannot understand very well the reasons. Do you see anything very 'weird'  increasing the force constant of an angle during a simulation to obtain a linear functionality?

Thanks a lot for your help.

Best wishes,

Rebeca.


> From: gromacs.org_gmx-users-request at maillist.sys.kth.se
> Subject: gromacs.org_gmx-users Digest, Vol 130, Issue 38
> To: gromacs.org_gmx-users at maillist.sys.kth.se
> Date: Fri, 13 Feb 2015 04:53:24 +0100
> 
> Send gromacs.org_gmx-users mailing list submissions to
> 	gromacs.org_gmx-users at maillist.sys.kth.se
> 
> To subscribe or unsubscribe via the World Wide Web, visit
> 	https://maillist.sys.kth.se/mailman/listinfo/gromacs.org_gmx-users
> or, via email, send a message with subject or body 'help' to
> 	gromacs.org_gmx-users-request at maillist.sys.kth.se
> 
> You can reach the person managing the list at
> 	gromacs.org_gmx-users-owner at maillist.sys.kth.se
> 
> When replying, please edit your Subject line so it is more specific
> than "Re: Contents of gromacs.org_gmx-users digest..."
> 
> 
> Today's Topics:
> 
>    1. Re: GAFF and alkynes (Justin Lemkul)
>    2. Re: Semiisotropic pressure coupling (Justin Lemkul)
>    3. Re: Semiisotropic pressure coupling (shivangi nangia)
>    4. Re: Semiisotropic pressure coupling (shivangi nangia)
>    5. g_cluster (gromos method) creates clusters with members
>       having RMSD greater than the cutoff to the cluster "middle"
>       (Christopher Neale)
> 
> 
> ----------------------------------------------------------------------
> 
> Message: 1
> Date: Thu, 12 Feb 2015 20:36:30 -0500
> From: Justin Lemkul <jalemkul at vt.edu>
> To: gmx-users at gromacs.org
> Subject: Re: [gmx-users] GAFF and alkynes
> Message-ID: <54DD551E.8040201 at vt.edu>
> Content-Type: text/plain; charset=windows-1252; format=flowed
> 
> 
> 
> On 2/12/15 9:43 AM, Rebeca Garc?a Fandi?o wrote:
> > Thanks a lot for the suggestions, Justin.
> > In your CO2 tutorial I have seen the comment "One cannot effectively build this
> > molecule in the traditional sense, as there are algorithmic reasons why
> > an angle of 180? is not stable during a simulation."
> >
> > I was wondering what happens if I use angle restraints, for example. Would it lead to an unstable simulation? I suppose that if people use virtual sites to treat with linear molecules instead of applying angle restraints, there should be problems with this approximation...is this correct?
> 
> A restraint is just a biasing potential towards some desired outcome.  It's not 
> necessarily even different in form from a normal harmonic interaction, but I'd 
> be hesitant to try to use some angle definition and a restraint on top of that. 
>   The GROMACS code indicates there is a linear angle potential, but it doesn't 
> seem to be documented.  This was likely added to circumvent the previous 
> problems with these linear species.
> 
> -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
> 
> jalemkul at outerbanks.umaryland.edu | (410) 706-7441
> http://mackerell.umaryland.edu/~jalemkul
> 
> ==================================================
> 
> 
> ------------------------------
> 
> Message: 2
> Date: Thu, 12 Feb 2015 20:41:11 -0500
> From: Justin Lemkul <jalemkul at vt.edu>
> To: gmx-users at gromacs.org
> Subject: Re: [gmx-users] Semiisotropic pressure coupling
> Message-ID: <54DD5637.50206 at vt.edu>
> Content-Type: text/plain; charset=windows-1252; format=flowed
> 
> 
> 
> On 2/12/15 8:15 PM, shivangi nangia wrote:
> > Hello All,
> >
> > I am doing all-atom simulation with charmm 36 ff (obtained from reverse CG
> > using the script provided by Tsjerk)
> > The system consists of popc, 21 AA protein, water and ions.
> >
> > After reverse CG I did short EM and 5 ns NVT simulation.
> > However, when I do NPT simulation, the area per lipid drops from 67 ang2 to
> > 57 within 1 ns ( x and y dimensions decrese and z increases).
> >
> > I am using semiisotropic pressure coupling.
> >
> > I have a question regarding that.
> > I noticed in example .mdps and the tutorial by Justin Lemkul
> >
> > http://www.bevanlab.biochem.vt.edu/Pages/Personal/justin/gmx-tutorials/membrane_protein/Files/npt.mdp
> >
> >
> > http://www.bevanlab.biochem.vt.edu/Pages/Personal/justin/gmx-tutorials/membrane_protein/Files/md.mdp
> >
> >
> > tau_p is specified, the ref_p and compresibility is same, does not that
> > mean essentially I am doing isotropic coupling?
> >
> 
> No, that does not indicate isotropic pressure coupling.  The full pressure 
> tensor is used to calculate pressures in the x-y and z dimensions separately. 
> The values of tau_p and compressibility are simply used for calculating the 
> response time of the barostat.
> 
> The problem is that your nonbonded settings are incorrect:
> 
> > ; nblist cut-off
> > rlist                    = 1.2
> >
> > ; OPTIONS FOR ELECTROSTATICS AND VDW
> > ; Method for doing electrostatics
> > coulombtype              = PME
> > rcoulomb-switch          = 0
> > rcoulomb                 = 1.2
> > fourierspacing           = 0.16
> > pme_order                = 4
> > optimize-fft            = yes
> > ; Method for doing Van der Waals
> > vdw-type                 = switch
> > ; cut-off lengths
> > rvdw-switch              = 1.0
> > rvdw                     = 1.2
> >
> 
> Note that potential switching leads to errors in the lipid force field.  You 
> need a buffered neighbor list with force switching.  The exact parameters that 
> you should use are listed on 
> http://www.gromacs.org/Documentation/Terminology/Force_Fields/CHARMM.  For 
> lipids, you should absolutely not deviate from these settings.  The protein and 
> nucleic acid force fields (and the remainder of CHARMM additive models) should 
> use these settings, though they are more forgiving (i.e. potential switching 
> does not lead to noticeable artifacts).  With lipids, the requirements are 
> pretty absolute.  People frequently report a drop in APL; every time they do, 
> it's because they're using the wrong settings.
> 
> -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
> 
> jalemkul at outerbanks.umaryland.edu | (410) 706-7441
> http://mackerell.umaryland.edu/~jalemkul
> 
> ==================================================
> 
> 
> ------------------------------
> 
> Message: 3
> Date: Thu, 12 Feb 2015 22:30:03 -0500
> From: shivangi nangia <shivangi.nangia at gmail.com>
> To: Discussion list for GROMACS users <gmx-users at gromacs.org>
> Subject: Re: [gmx-users] Semiisotropic pressure coupling
> Message-ID:
> 	<CAH137_LODQ5hzCnqSF-QBe1xuy+es_vG0hv5isrvJyike3Pw1Q at mail.gmail.com>
> Content-Type: text/plain; charset=UTF-8
> 
> Dear Justin,
> >
> 
> Thanks for the reply.
> 
> The link you have mentioned says the parameters are for GROMACS 5.0
> 
> I am using an older version, 4.6.1
> 
> Is there another link/suggestions for that version.
> 
> Thanks so much,
> 
> Best,
> sxn
> 
> >
> >
> > On 2/12/15 8:15 PM, shivangi nangia wrote:
> >
> >> Hello All,
> >>
> >> I am doing all-atom simulation with charmm 36 ff (obtained from reverse CG
> >> using the script provided by Tsjerk)
> >> The system consists of popc, 21 AA protein, water and ions.
> >>
> >> After reverse CG I did short EM and 5 ns NVT simulation.
> >> However, when I do NPT simulation, the area per lipid drops from 67 ang2
> >> to
> >> 57 within 1 ns ( x and y dimensions decrese and z increases).
> >>
> >> I am using semiisotropic pressure coupling.
> >>
> >> I have a question regarding that.
> >> I noticed in example .mdps and the tutorial by Justin Lemkul
> >>
> >> http://www.bevanlab.biochem.vt.edu/Pages/Personal/justin/
> >> gmx-tutorials/membrane_protein/Files/npt.mdp
> >>
> >>
> >> http://www.bevanlab.biochem.vt.edu/Pages/Personal/justin/
> >> gmx-tutorials/membrane_protein/Files/md.mdp
> >>
> >>
> >> tau_p is specified, the ref_p and compresibility is same, does not that
> >> mean essentially I am doing isotropic coupling?
> >>
> >>
> > No, that does not indicate isotropic pressure coupling.  The full pressure
> > tensor is used to calculate pressures in the x-y and z dimensions
> > separately. The values of tau_p and compressibility are simply used for
> > calculating the response time of the barostat.
> >
> > The problem is that your nonbonded settings are incorrect:
> >
> >  ; nblist cut-off
> >> rlist                    = 1.2
> >>
> >> ; OPTIONS FOR ELECTROSTATICS AND VDW
> >> ; Method for doing electrostatics
> >> coulombtype              = PME
> >> rcoulomb-switch          = 0
> >> rcoulomb                 = 1.2
> >> fourierspacing           = 0.16
> >> pme_order                = 4
> >> optimize-fft            = yes
> >> ; Method for doing Van der Waals
> >> vdw-type                 = switch
> >> ; cut-off lengths
> >> rvdw-switch              = 1.0
> >> rvdw                     = 1.2
> >>
> >>
> > Note that potential switching leads to errors in the lipid force field.
> > You need a buffered neighbor list with force switching.  The exact
> > parameters that you should use are listed on http://www.gromacs.org/
> > Documentation/Terminology/Force_Fields/CHARMM.  For lipids, you should
> > absolutely not deviate from these settings.  The protein and nucleic acid
> > force fields (and the remainder of CHARMM additive models) should use these
> > settings, though they are more forgiving (i.e. potential switching does not
> > lead to noticeable artifacts).  With lipids, the requirements are pretty
> > absolute.  People frequently report a drop in APL; every time they do, it's
> > because they're using the wrong settings.
> >
> > -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
> >
> > jalemkul at 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-request at gromacs.org.
> >
> >
> 
> 
> ------------------------------
> 
> Message: 4
> Date: Thu, 12 Feb 2015 22:34:31 -0500
> From: shivangi nangia <shivangi.nangia at gmail.com>
> To: Discussion list for GROMACS users <gmx-users at gromacs.org>
> Subject: Re: [gmx-users] Semiisotropic pressure coupling
> Message-ID:
> 	<CAH137_LOL_sKLWZA_JdesieOVb-FW4CD4nziKkchCeAm7WgA=Q at mail.gmail.com>
> Content-Type: text/plain; charset=UTF-8
> 
> >
> > Dear Justin,
> >
> 
> Thanks for the reply.
> 
> The link you have mentioned says the parameters are for GROMACS 5.0
> 
> I am using an older version, 4.6.1
> 
> Is there another link/suggestions for that version.
> 
> Sorry, I forgot to mention.
> 
> On trying to use those parameters with my version of GROMACS, I get the
> error:
> 
> ERROR 1 [file md.mdp, line 78]:
>   Invalid enum 'force-switch' for variable vdw-modifier, using
>   'Potential-shift-Verlet'
>   Next time use one of: 'Potential-shift-Verlet' 'Potential-shift' 'None'
>   'Potential-switch' 'Exact-cutoff'
> 
> 
> The mdp options manual says:
> 
> *Force-switch*Smoothly switches the forces to zero between *rvdw-switch*
>  and *rvdw*. This shifts the potential shift over the whole range and
> switches it to zero at the cut-off. Note that this is more expensive to
> calculate than a plain cut-off and it is not required for energy
> conservation, since *Potential-shift* conserves energy just as well
> 
> Can I use Potential-shift instead of my version?
> Kindly suggesr.
> Thanks in advance.
> 
> 
> 
> On Thu, Feb 12, 2015 at 10:30 PM, shivangi nangia <shivangi.nangia at gmail.com
> > wrote:
> 
> >
> >
> > Dear Justin,
> >>
> >
> > Thanks for the reply.
> >
> > The link you have mentioned says the parameters are for GROMACS 5.0
> >
> > I am using an older version, 4.6.1
> >
> > Is there another link/suggestions for that version.
> >
> > Thanks so much,
> >
> > Best,
> > sxn
> >
> >>
> >>
> >> On 2/12/15 8:15 PM, shivangi nangia wrote:
> >>
> >>> Hello All,
> >>>
> >>> I am doing all-atom simulation with charmm 36 ff (obtained from reverse
> >>> CG
> >>> using the script provided by Tsjerk)
> >>> The system consists of popc, 21 AA protein, water and ions.
> >>>
> >>> After reverse CG I did short EM and 5 ns NVT simulation.
> >>> However, when I do NPT simulation, the area per lipid drops from 67 ang2
> >>> to
> >>> 57 within 1 ns ( x and y dimensions decrese and z increases).
> >>>
> >>> I am using semiisotropic pressure coupling.
> >>>
> >>> I have a question regarding that.
> >>> I noticed in example .mdps and the tutorial by Justin Lemkul
> >>>
> >>> http://www.bevanlab.biochem.vt.edu/Pages/Personal/justin/
> >>> gmx-tutorials/membrane_protein/Files/npt.mdp
> >>>
> >>>
> >>> http://www.bevanlab.biochem.vt.edu/Pages/Personal/justin/
> >>> gmx-tutorials/membrane_protein/Files/md.mdp
> >>>
> >>>
> >>> tau_p is specified, the ref_p and compresibility is same, does not that
> >>> mean essentially I am doing isotropic coupling?
> >>>
> >>>
> >> No, that does not indicate isotropic pressure coupling.  The full
> >> pressure tensor is used to calculate pressures in the x-y and z dimensions
> >> separately. The values of tau_p and compressibility are simply used for
> >> calculating the response time of the barostat.
> >>
> >> The problem is that your nonbonded settings are incorrect:
> >>
> >>  ; nblist cut-off
> >>> rlist                    = 1.2
> >>>
> >>> ; OPTIONS FOR ELECTROSTATICS AND VDW
> >>> ; Method for doing electrostatics
> >>> coulombtype              = PME
> >>> rcoulomb-switch          = 0
> >>> rcoulomb                 = 1.2
> >>> fourierspacing           = 0.16
> >>> pme_order                = 4
> >>> optimize-fft            = yes
> >>> ; Method for doing Van der Waals
> >>> vdw-type                 = switch
> >>> ; cut-off lengths
> >>> rvdw-switch              = 1.0
> >>> rvdw                     = 1.2
> >>>
> >>>
> >> Note that potential switching leads to errors in the lipid force field.
> >> You need a buffered neighbor list with force switching.  The exact
> >> parameters that you should use are listed on http://www.gromacs.org/
> >> Documentation/Terminology/Force_Fields/CHARMM.  For lipids, you should
> >> absolutely not deviate from these settings.  The protein and nucleic acid
> >> force fields (and the remainder of CHARMM additive models) should use these
> >> settings, though they are more forgiving (i.e. potential switching does not
> >> lead to noticeable artifacts).  With lipids, the requirements are pretty
> >> absolute.  People frequently report a drop in APL; every time they do, it's
> >> because they're using the wrong settings.
> >>
> >> -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
> >>
> >> jalemkul at 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-request at gromacs.org.
> >>
> >>
> >
> 
> 
> ------------------------------
> 
> Message: 5
> Date: Fri, 13 Feb 2015 03:53:17 +0000
> From: Christopher Neale <chris.neale at alum.utoronto.ca>
> To: "gromacs.org_gmx-users at maillist.sys.kth.se"
> 	<gromacs.org_gmx-users at maillist.sys.kth.se>
> Subject: [gmx-users] g_cluster (gromos method) creates clusters with
> 	members having RMSD greater than the cutoff to the cluster "middle"
> Message-ID: <1423799598212.17124 at alum.utoronto.ca>
> Content-Type: text/plain; charset="iso-8859-1"
> 
> Dear Users:
> 
> I have run g_cluster from gromacs 4.6.5 as follows:
> 
> g_cluster -s my.tpr -f tmp.xtc -method gromos -nofit -o rmsd-clust_nofit.xpm -g cluster_nofit.log -sz clust-size_nofit.xvg -cl clusters_nofit.pdb -n index.ndx -cutoff 0.275 -wcl 10 -cl finally.xtc
> 
> The top of cluster_nofit.log looks like this:
> Using gromos method for clustering
> Using RMSD cutoff 0.275 nm
> The RMSD ranges from 0.0247008 to 5.03412 nm
> Average RMSD is 0.45568
> Number of structures for matrix 12391
> Energy of the matrix is 4177.11 nm
> 
> Found 874 clusters
> 
> Writing middle structure for each cluster to finally.xtc
> Writing all structures for the first 10 clusters with more than 1 structures to finally.xtc%03d.xtc
> 
> cl. | #st  rmsd | middle rmsd | cluster members
>   1 | 4219  0.202 | 784700 .154 |  92900 103300 116100 125900 126500 133000 156600
> ...
>     |           |             | 693900 694000 694200 694300 694600 694700 695700
> ...
> 
> So the frame at 694700 is in the first cluster.
> 
> However, when I use g_rms to look at the rmsd between frame at 784700 and 694700, I find that the RMSD is 0.29 nm.
> 
> I am confused because frame 784700 is listed as the "middle" (which I assume is the cluster centroid) so frame 694700 should not be included in this cluster. 
> 
> When I put these two frames into a single .xtc with only 2 frames, g_cluster correctly puts them into two different clusters when I use a cutoff of 0.275 nm.
> 
> I understand that some of the members in a cluster can have a pairwise RMSD greater than the cutoff, but I thought that the "middle" of the cluster should have an RMSD less than the cutoff to all cluster members.
> 
> Can anybody help me figure out what I am missing?
> 
> Thank you,
> Chris.
> 
> 
> 
> ------------------------------
> 
> -- 
> 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-request at gromacs.org.
> 
> 
> End of gromacs.org_gmx-users Digest, Vol 130, Issue 38
> ******************************************************
 		 	   		  


More information about the gromacs.org_gmx-users mailing list