RE: [gmx-users] Why not PBC for implicit solvent?
Dear Justin (and Yun), Since our last discussion on implicit solvent simulations, I have read a bit more about how cutoffs in GBSA (or variants) implicit solvent simulations are used in other softwares such as AMBER, CHARMM and NAMD. It appears to me that many use switch or shift cutoffs. Have you tried this in your tests? Such boundary conditions should remove any non-conservation of energy problem. As you said, protein stability should also be tested though. Moreover, many use a non-zero salt concentration giving rise to a Debye screening (exp(-r*kappa), kappa^-1 ~ 1 nm under physiological concentration). Such screening makes the electrostatic energy tends to zero faster. In GROMACS, it is possible to specify a salt concentration for the implicit solvent scheme in the .mdp file, but such variable has not been coded for in the source code (as explained in the manual). From what I saw, this is still on the to-do list of the developers. Is it still under development? Finally, there as been a post on the Gromacs mailing list late 2012 about some energy differences between the AMBER and GROMACS implementation of OBC GBSA. It was found that the parameters of the gbsa.itp file should be changed to obtain a better energy agreement for the GB part. The SA part also showed some energy difference that could not be corrected for and were attributed to different surface accessible calculation algorithm. See : http://lists.gromacs.org/pipermail/gmx-users/2012-November/076230.html. Hope it helps, Sébastien > Date: Tue, 26 Feb 2013 06:56:40 -0500 > From: jalem...@vt.edu > To: gmx-users@gromacs.org > Subject: Re: [gmx-users] Why not PBC for implicit solvent? > > > > On 2/25/13 10:30 PM, Yun Shi wrote: > > Hi everyone, > > > > Previous posts mentioned setting pbc = none for MD simulations with > > implicit solvent. But I am trying to see the behavior of certain > > concentration of ligands (small molecules, no big biomolecules) in > > solvent, so I wonder if setting pbc = xyz would cause any problem for > > my system? > > > > I see no theoretical problem with running NVT with "pbc = xyz" with implicit > solvent, but definitely not NPT since the box will shrink inwards and lead to > periodicity artifacts (if it even remains stable at all). I usually set a > nonperiodic box because it allows me to use the infinite cutoff approach, > which > is the only one I have found to give sensible results. > > > Should I also stay with the normal cutoff values for vdw and coul > > interactions? > > > > Maybe, but do some serious testing before relying on the results. Maybe for > small molecules a finite cutoff will work, but for proteins I have tried > cutoffs > of 1.0, 2.0, 4.0, and even 8.0 nm and all have unfolded or distorted. > > -Justin > > -- > > > Justin A. Lemkul, Ph.D. > Research Scientist > Department of Biochemistry > Virginia Tech > Blacksburg, VA > jalemkul[at]vt.edu | (540) 231-9080 > http://www.bevanlab.biochem.vt.edu/Pages/Personal/justin > > > -- > gmx-users mailing listgmx-users@gromacs.org > http://lists.gromacs.org/mailman/listinfo/gmx-users > * Please search the archive at > http://www.gromacs.org/Support/Mailing_Lists/Search before posting! > * Please don't post (un)subscribe requests to the list. Use the > www interface or send it to gmx-users-requ...@gromacs.org. > * Can't post? Read http://www.gromacs.org/Support/Mailing_Lists -- gmx-users mailing listgmx-users@gromacs.org http://lists.gromacs.org/mailman/listinfo/gmx-users * Please search the archive at http://www.gromacs.org/Support/Mailing_Lists/Search before posting! * Please don't post (un)subscribe requests to the list. Use the www interface or send it to gmx-users-requ...@gromacs.org. * Can't post? Read http://www.gromacs.org/Support/Mailing_Lists
RE: [gmx-users] Implicit solvent
Dear Justin, Thank you very much for your answer. I did not take this as a criticism of implicit solvent model. I was just wondering what was the origin of the limitations you were speaking of. I did some small implicit solvent simulations, but I must say that I have not thoroughly looked at the impact of cutoffs on energy conservation. I will be careful about that. About the bug for the implicit solvent code in GROMACS, is it only found in the new version 4.6 or it is also present in earlier versions such as in the 4.5 series? Thanks, Sebastien > Date: Thu, 14 Feb 2013 10:48:51 -0500 > From: jalem...@vt.edu > To: gmx-users@gromacs.org > Subject: Re: [gmx-users] Implicit solvent > > > > On 2/14/13 10:21 AM, Sebastien Cote wrote: > > > > Dear Justin, > > I am not sure to follow you. You essentially say that it is better to avoid > > using implicit solvent i.e. the generalized Born-formalism implemented in > > GROMACS? > > For the case of optimizing and relaxing a system (expecting short MD), I > > agree that it might be preferable to use explicit solvent. > > But for some systems, such as the oligomerization of a protein from a > > random state, explicit solvent can make this problem intractable in MD > > (inefficient sampling of conformational space due to too many degrees of > > freedom from the solvent) or REMD (too many processors to use due to too > > many degrees of freedom). > > I saw some tutorials and workshops suggesting larger cutoffs when using > > implicit solvent as the non-bonded function do tend to zero at long range. > > Then, domain decomposition can be used over particle decomposition. > > Could you clarify your point please? Especially, your statement saying that > > : "you can only run on 1-2 processors in implicit solvent". > > There is an existing limitation in Gromacs, related to constraints, that > causes > the implicit solvent code to fail when run on CPU if you try to use more than > 2 > processors. It is my understanding that someone (probably Berk) is working > on a > better fix than what is in place now. DD does not work with nonperiodic unit > cells and simple grid searches; you have to use particle decomposition. > Therein > lies another limitation. I have found that only infinite cutoffs (i.e. > rlist=rcoulomb=rvdw=0) lead to stable simulations, while any finite cutoff > (even > several nm in length) leads to poor energy conservation and unstable > trajectories. > > You can use GPU acceleration for implicit solvent simulations via OpenMM, but > support for this feature is minimal, at best, at the present time. > > Please do not take my comments as a criticism of implicit solvent methods; > they > are perfectly useful for a number of cases. > > -Justin > > -- > > > Justin A. Lemkul, Ph.D. > Research Scientist > Department of Biochemistry > Virginia Tech > Blacksburg, VA > jalemkul[at]vt.edu | (540) 231-9080 > http://www.bevanlab.biochem.vt.edu/Pages/Personal/justin > > > -- > gmx-users mailing listgmx-users@gromacs.org > http://lists.gromacs.org/mailman/listinfo/gmx-users > * Please search the archive at > http://www.gromacs.org/Support/Mailing_Lists/Search before posting! > * Please don't post (un)subscribe requests to the list. Use the > www interface or send it to gmx-users-requ...@gromacs.org. > * Can't post? Read http://www.gromacs.org/Support/Mailing_Lists -- gmx-users mailing listgmx-users@gromacs.org http://lists.gromacs.org/mailman/listinfo/gmx-users * Please search the archive at http://www.gromacs.org/Support/Mailing_Lists/Search before posting! * Please don't post (un)subscribe requests to the list. Use the www interface or send it to gmx-users-requ...@gromacs.org. * Can't post? Read http://www.gromacs.org/Support/Mailing_Lists
RE: [gmx-users] Implicit solvent
Dear Justin, I am not sure to follow you. You essentially say that it is better to avoid using implicit solvent i.e. the generalized Born-formalism implemented in GROMACS? For the case of optimizing and relaxing a system (expecting short MD), I agree that it might be preferable to use explicit solvent. But for some systems, such as the oligomerization of a protein from a random state, explicit solvent can make this problem intractable in MD (inefficient sampling of conformational space due to too many degrees of freedom from the solvent) or REMD (too many processors to use due to too many degrees of freedom). I saw some tutorials and workshops suggesting larger cutoffs when using implicit solvent as the non-bonded function do tend to zero at long range. Then, domain decomposition can be used over particle decomposition. Could you clarify your point please? Especially, your statement saying that : "you can only run on 1-2 processors in implicit solvent". Thank you, Sebastien > Date: Thu, 14 Feb 2013 08:03:28 -0500 > From: jalem...@vt.edu > To: gmx-users@gromacs.org > Subject: Re: [gmx-users] Implicit solvent > > > > On 2/14/13 1:01 AM, Алексей Раевский wrote: > > Hi, dear developers! > > > > I want to ask you abou dynamics in implicit solvent. I have a complex of > > animal protein - dimer/trimer. After modeling by homology I have built > > another one from the plant organism. Dimer/trimer was constructed with > > superposition method. The question is - can I use imlicit solvent model to > > optimize and relax this complex? Is it real and where I can find suitable > > parameters of mdp? > > > > Infinite cutoffs with a nonperiodic unit cell are generally advisable for > running implicit solvent calculations. They don't parallelize well due to a > bug > that I believe is still being fixed, so you can only run on 1-2 processors in > implicit solvent. Performance for large systems actually scales better in > explicit solvent, especially if you use GPU acceleration. > > -Justin > > > -- > > > Justin A. Lemkul, Ph.D. > Research Scientist > Department of Biochemistry > Virginia Tech > Blacksburg, VA > jalemkul[at]vt.edu | (540) 231-9080 > http://www.bevanlab.biochem.vt.edu/Pages/Personal/justin > > > -- > gmx-users mailing listgmx-users@gromacs.org > http://lists.gromacs.org/mailman/listinfo/gmx-users > * Please search the archive at > http://www.gromacs.org/Support/Mailing_Lists/Search before posting! > * Please don't post (un)subscribe requests to the list. Use the > www interface or send it to gmx-users-requ...@gromacs.org. > * Can't post? Read http://www.gromacs.org/Support/Mailing_Lists -- gmx-users mailing listgmx-users@gromacs.org http://lists.gromacs.org/mailman/listinfo/gmx-users * Please search the archive at http://www.gromacs.org/Support/Mailing_Lists/Search before posting! * Please don't post (un)subscribe requests to the list. Use the www interface or send it to gmx-users-requ...@gromacs.org. * Can't post? Read http://www.gromacs.org/Support/Mailing_Lists
RE: [gmx-users] Problem with equilibrated lipid bilayer structure
Hello, You might want to use TIPS3P (special CHARMM TIP3P) instead of regular TIP3P because some researchers observed that regular TIP3P can impact the area per lipid significantly when simulating CHARMM lipid bilayer (see reference below). Concerning the CHARMM36 lipid force field in Gromacs, there is also an issue about cutoffs. The switch cutoff scheme in the charmm software is different than in the gromacs software. It might impact your results for some bilayers, while it does not for others. You might find that previous post on the mailing list enlightening (and the cited reference) : " Hi everyone,This message is just to let people know that there is finally a reference that we ask you to cite if you publish work using the CHARMM36 lipid force field contribution available from the GROMACS website. The reference is provided in a new version of the forcefield.doc file, and is also given below: Thomas J. Piggot, Ángel Piñeiro and Syma Khalid Molecular Dynamics Simulations of Phosphatidylcholine Membranes: A Comparative Force Field Study Journal of Chemical Theory and Computation DOI: 10.1021/ct3003157 http://pubs.acs.org/doi/abs/10.1021/ct3003157 A validation of the conversion for DPPC and POPC membrane systems (in terms of single point energies calculated in both GROMACS and NAMD) is provided in the supporting information of this work. Cheers Tom" Hope that helps, Sebastien > Date: Mon, 15 Oct 2012 13:48:47 +0800 > From: jernej.zi...@gmail.com > To: gmx-users@gromacs.org > Subject: [gmx-users] Problem with equilibrated lipid bilayer structure > > Hi. > I used CHARMM to prepare a cholesterol/POPC lipid bilayer. I > equilibrated it in CHARMM and the imported the resulting PDB to > GROMACS. In CHARMM I was the TIP3P water model, so before importing > the PDB I changed the atom types from the TIP3P's ones to SPCE ones. > Then I used this command to import the PDB: > pdb2gmx -f lipid-wat.pdb -o bilayer.gro -p bilayer.top -i > bilayer-posres.itp -ff charmm36cgenff -water spce -noter -v -renum > > After the import I center the system in the unit cell with: editconf > -f bilayer.gro -o bilayer.gro -c > > Problem 1: The system size (unit cell) is not properly > computed/detected. The following numbers are for the same structure > after the CHARMM equilibration. >CHARMM reports the size as: 10.45820 x 10.45825 x > 6.864382 (in nm; converted from Angstroms) >GROMACS states the size is: 12.30940 x 11.81980 x > 7.138100 (in nm) > > Problem 2: While I'm able to run MD in CHARMM if I start from the > equilibrated structure, I am unable to do so in GROMACS. As the system > apparently blows up with the cryptic message: > Program mdrun, VERSION 4.5.5 > Source code file: /build/buildd/gromacs-4.5.5/src/mdlib/pme.c, line: 538 > > Fatal error: > 5 particles communicated to PME node 0 are more than 2/3 times the > cut-off out of the domain decomposition cell of their charge group in > dimension y. > This usually means that your system is not well equilibrated. > For more information and tips for troubleshooting, please check the GROMACS > website at http://www.gromacs.org/Documentation/Errors > > I can minimize the system, but any attempt to run MD results in the > aforementioned message. Here's the MDP file I would like to use: > ;define = -DPOSRES ; position restrain the protein > ; Run parameters > integrator= md; leap-frog integrator > nsteps= 50; 2 * 50 = 1000 ps (1 ns) > dt= 0.002 ; 2 fs > ; Output control > nstxout = 100 ; save coordinates every 0.2 ps > nstvout = 100 ; save velocities every 0.2 ps > nstenergy = 100 ; save energies every 0.2 ps > nstlog= 100 ; update log file every 0.2 ps > ; Bond parameters > continuation = no; Restarting after NVT > constraint_algorithm = lincs ; holonomic constraints > constraints = all-bonds ; all bonds (even heavy atom-H bonds) > constrained > lincs_iter= 1 ; accuracy of LINCS > lincs_order = 4 ; also related to accuracy > ; Neighborsearching > ns_type = grid ; search neighboring grid cels > nstlist = 5 ; 10 fs > rlist = 1.2 ; short-range neighborlist cutoff (in nm) > rcoulomb = 1.2 ; short-range electrostatic cutoff (in nm) > rvdw = 1.2 ; short-range van der Waals cutoff (in nm) > ; Electrostatics > coulombtype = PME ; Particle Mesh Ewald for long-range > electrostatics > pme_order = 4 ; cubic interpolation > fourierspacing= 0.16 ; grid spacing for FFT > ; Temperature coupling is on > tcoupl= Nose-Hoover ; More accurate thermostat > tc-grps
RE: [gmx-users] CHARMM36 - Smaller Area per lipid for POPE - Why?
Thanks for the advices Chris. My peptide is known to be more favorably to PE than PC membrane that is why I am using POPE. Experimentally, the liquid phase transition is at 298K for POPE (if I am not mistaken). Is your 323K refer to some simulations? At first I wanted to use the new CHARMM36 lipids parameters because they are supposed to solve the previous CHARMM27 issue with the area per lipid. However, I am consistently obtained smaller APL then experiment and I am not able to reproduce the published APL obtained for POPE, even if I am starting from their equilibrated 80-POPE membrane and use same simulation conditions. That was the reason for starting this thread on the mailing list. Unfortunately, my peptide conformational space in solution is only well-represented by CHARMM27 (equivalently in CHARMM36), so I can not use Berger's lipid parameters with OPLS or GROMOS even if it would be preferable as they do not have APL inconsistency and are united-atom. I will made some tests in the NPAT ensemble. Perhaps the NPAT effects can be made neglegible by using bigger membrane compared to my peptide's size (?). Sebastien > From: chris.ne...@mail.utoronto.ca > To: gmx-users@gromacs.org > Date: Wed, 15 Aug 2012 17:29:29 + > Subject: [gmx-users] CHARMM36 - Smaller Area per lipid for POPE - Why? > > The area per lipid (APL) will certainly affect the free energy of > peptide/protein binding to a lipid bilayer. > I have not used charmm lipids extensively, but from what I understand they > older charmm lipids required > NPAT to get the correct APL. The newer charmm lipids were supposed to solve > that problem, but I have heard > it said that, though the problem has been alleviated to some extend, it still > remains. > > If I were you, I'd use POPC in place of POPE. POPE is notorious for giving > too-small APL's in simulations and I think > it even requires temperatures of 323 K to enter the liquid phase. > > That said, I don't have a specific answer to your question of whether there > are other affects of NPAT vs. NPT. > It is plausible that NPAT-based fluctuations could affect the pathway or the > kinetics. > > PS: I was not referring to lipid rafts, but the separate diffusion of the > upper and lower leaflets. Once the peptide is > fully inserted, if it spans both leaflets, this will tend to reduce this > leaflet-specific diffusion and would represent an > entropic penalty for binding (not sure how large). > > Chris. > > > > > Dear Peter, > > > > I also used h-bonds and I also switch LJ interaction from 0.8 nm to 1.2 nm > > (as in Klauda's paper). I will retry with a more solvated membrane. > > > > Would you have any thought on how the NPAT ensemble might affect > > peptide-membrane interactions like I am studying i.e. peptide is totally > > solvated, then adsorb, and finally may insert? The paper on > > peptide-membrane interaction like this usually use united-atom lipid in the > > NPT ensemble. Most of the work I have seen on Charmm membrane in the NPAT > > ensemble were for embedded membrane protein. > > Sorry, but I only have experience with large pre-embedded membrane proteins, > and those are governed both by signal sequences and post-translational > modification. > > Chris's last email on the subject might lead to the hypothesis that lipid > raft translation as the leaflets "slide" past one another could be a > contributing factor to adsorbption of your species. > > > > > Thanks, > > > > Sebastien > -- > gmx-users mailing list gmx-users@gromacs.org > http://lists.gromacs.org/mailman/listinfo/gmx-users > * Only plain text messages are allowed! > * Please search the archive at > http://www.gromacs.org/Support/Mailing_Lists/Search before posting! > * Please don't post (un)subscribe requests to the list. Use the > www interface or send it to gmx-users-requ...@gromacs.org. > * Can't post? Read http://www.gromacs.org/Support/Mailing_Lists -- gmx-users mailing listgmx-users@gromacs.org http://lists.gromacs.org/mailman/listinfo/gmx-users * Only plain text messages are allowed! * Please search the archive at http://www.gromacs.org/Support/Mailing_Lists/Search before posting! * Please don't post (un)subscribe requests to the list. Use the www interface or send it to gmx-users-requ...@gromacs.org. * Can't post? Read http://www.gromacs.org/Support/Mailing_Lists
RE: [gmx-users] Parameters for bonded interactions
Thanks Mark! This is exactly what I need. > Date: Wed, 15 Aug 2012 22:18:10 +1000 > From: mark.abra...@anu.edu.au > To: gmx-users@gromacs.org > Subject: Re: [gmx-users] Parameters for bonded interactions > > On 15/08/2012 8:44 PM, Sebastien Cote wrote: > > After checking the post-processed topology, it does not contain the > > information that I need. I would like to know the torsion angle parameters > > of each torsion angle, and then compare with the ffbonded file to see the > > corresponding four-atom groups ' X X X X '. > > > > Is there a way to convert the .tpr file to 'human readable'? > > Sorry, meant to say that last email. gmxdump is all that is available. > You'll have to work backwards from the atom numbers (starting counting > from zero!) to get the interaction number to look up the interaction > parameters from their list. > > Mark > > > > > Thanks again, > > > > Sebastien > > > > > >> Date: Wed, 15 Aug 2012 11:48:00 +1000 > >> From: mark.abra...@anu.edu.au > >> To: gmx-users@gromacs.org > >> Subject: Re: [gmx-users] Parameters for bonded interactions > >> > >> On 15/08/2012 9:46 AM, Sebastien Cote wrote: > >>> Dear Gromacs users, > >>> > >>> In the topology file of the protein, we see every two atoms that share a > >>> bond, three atoms that share a bond angle, and four atoms that share a > >>> torsion angle. However, the parameters (equilibrium value, energy > >>> constant, phase) are not explicitly shown as Gromacs fetch them in the > >>> ffbonded table. My question: Is there a simple way to see (without having > >>> to look at the code) which parameters are associated to which bond > >>> length, bond angle and torsion angle? > >>> > >>> More specifically, I would like to know which torsion angle parameters is > >>> specifically associated to each torsion angle. In the ffbonded table, > >>> there is sometimes different groups of four atoms that can be associated > >>> to the same torsion angle as there is the 'X' atom which can be any atom. > >> grompp -pp writes out a post-processed topology, which may be what you > >> need. > >> > >> Otherwise you can inspect the contents of the .tpr file to see which > >> interaction gets which parameters. That is not easy for a large system, > >> but is definitive. > >> > >> Mark > >> -- > >> gmx-users mailing list gmx-users@gromacs.org > >> http://lists.gromacs.org/mailman/listinfo/gmx-users > >> * Only plain text messages are allowed! > >> * Please search the archive at > >> http://www.gromacs.org/Support/Mailing_Lists/Search before posting! > >> * Please don't post (un)subscribe requests to the list. Use the > >> www interface or send it to gmx-users-requ...@gromacs.org. > >> * Can't post? Read http://www.gromacs.org/Support/Mailing_Lists > > -- > > gmx-users mailing list gmx-users@gromacs.org > > http://lists.gromacs.org/mailman/listinfo/gmx-users > > * Only plain text messages are allowed! > > * Please search the archive at > > http://www.gromacs.org/Support/Mailing_Lists/Search before posting! > > * Please don't post (un)subscribe requests to the list. Use the > > www interface or send it to gmx-users-requ...@gromacs.org. > > * Can't post? Read http://www.gromacs.org/Support/Mailing_Lists > > -- > gmx-users mailing list gmx-users@gromacs.org > http://lists.gromacs.org/mailman/listinfo/gmx-users > * Only plain text messages are allowed! > * Please search the archive at > http://www.gromacs.org/Support/Mailing_Lists/Search before posting! > * Please don't post (un)subscribe requests to the list. Use the > www interface or send it to gmx-users-requ...@gromacs.org. > * Can't post? Read http://www.gromacs.org/Support/Mailing_Lists -- gmx-users mailing listgmx-users@gromacs.org http://lists.gromacs.org/mailman/listinfo/gmx-users * Only plain text messages are allowed! * Please search the archive at http://www.gromacs.org/Support/Mailing_Lists/Search before posting! * Please don't post (un)subscribe requests to the list. Use the www interface or send it to gmx-users-requ...@gromacs.org. * Can't post? Read http://www.gromacs.org/Support/Mailing_Lists
RE: [gmx-users] Parameters for bonded interactions
After checking the post-processed topology, it does not contain the information that I need. I would like to know the torsion angle parameters of each torsion angle, and then compare with the ffbonded file to see the corresponding four-atom groups ' X X X X '. Is there a way to convert the .tpr file to 'human readable'? Thanks again, Sebastien > Date: Wed, 15 Aug 2012 11:48:00 +1000 > From: mark.abra...@anu.edu.au > To: gmx-users@gromacs.org > Subject: Re: [gmx-users] Parameters for bonded interactions > > On 15/08/2012 9:46 AM, Sebastien Cote wrote: > > Dear Gromacs users, > > > > In the topology file of the protein, we see every two atoms that share a > > bond, three atoms that share a bond angle, and four atoms that share a > > torsion angle. However, the parameters (equilibrium value, energy constant, > > phase) are not explicitly shown as Gromacs fetch them in the ffbonded > > table. My question: Is there a simple way to see (without having to look at > > the code) which parameters are associated to which bond length, bond angle > > and torsion angle? > > > > More specifically, I would like to know which torsion angle parameters is > > specifically associated to each torsion angle. In the ffbonded table, there > > is sometimes different groups of four atoms that can be associated to the > > same torsion angle as there is the 'X' atom which can be any atom. > > grompp -pp writes out a post-processed topology, which may be what you need. > > Otherwise you can inspect the contents of the .tpr file to see which > interaction gets which parameters. That is not easy for a large system, > but is definitive. > > Mark > -- > gmx-users mailing list gmx-users@gromacs.org > http://lists.gromacs.org/mailman/listinfo/gmx-users > * Only plain text messages are allowed! > * Please search the archive at > http://www.gromacs.org/Support/Mailing_Lists/Search before posting! > * Please don't post (un)subscribe requests to the list. Use the > www interface or send it to gmx-users-requ...@gromacs.org. > * Can't post? Read http://www.gromacs.org/Support/Mailing_Lists -- gmx-users mailing listgmx-users@gromacs.org http://lists.gromacs.org/mailman/listinfo/gmx-users * Only plain text messages are allowed! * Please search the archive at http://www.gromacs.org/Support/Mailing_Lists/Search before posting! * Please don't post (un)subscribe requests to the list. Use the www interface or send it to gmx-users-requ...@gromacs.org. * Can't post? Read http://www.gromacs.org/Support/Mailing_Lists
RE: [gmx-users] Parameters for bonded interactions
Thanks a lot Mark! > Date: Wed, 15 Aug 2012 11:48:00 +1000 > From: mark.abra...@anu.edu.au > To: gmx-users@gromacs.org > Subject: Re: [gmx-users] Parameters for bonded interactions > > On 15/08/2012 9:46 AM, Sebastien Cote wrote: > > Dear Gromacs users, > > > > In the topology file of the protein, we see every two atoms that share a > > bond, three atoms that share a bond angle, and four atoms that share a > > torsion angle. However, the parameters (equilibrium value, energy constant, > > phase) are not explicitly shown as Gromacs fetch them in the ffbonded > > table. My question: Is there a simple way to see (without having to look at > > the code) which parameters are associated to which bond length, bond angle > > and torsion angle? > > > > More specifically, I would like to know which torsion angle parameters is > > specifically associated to each torsion angle. In the ffbonded table, there > > is sometimes different groups of four atoms that can be associated to the > > same torsion angle as there is the 'X' atom which can be any atom. > > grompp -pp writes out a post-processed topology, which may be what you need. > > Otherwise you can inspect the contents of the .tpr file to see which > interaction gets which parameters. That is not easy for a large system, > but is definitive. > > Mark > -- > gmx-users mailing list gmx-users@gromacs.org > http://lists.gromacs.org/mailman/listinfo/gmx-users > * Only plain text messages are allowed! > * Please search the archive at > http://www.gromacs.org/Support/Mailing_Lists/Search before posting! > * Please don't post (un)subscribe requests to the list. Use the > www interface or send it to gmx-users-requ...@gromacs.org. > * Can't post? Read http://www.gromacs.org/Support/Mailing_Lists -- gmx-users mailing listgmx-users@gromacs.org http://lists.gromacs.org/mailman/listinfo/gmx-users * Only plain text messages are allowed! * Please search the archive at http://www.gromacs.org/Support/Mailing_Lists/Search before posting! * Please don't post (un)subscribe requests to the list. Use the www interface or send it to gmx-users-requ...@gromacs.org. * Can't post? Read http://www.gromacs.org/Support/Mailing_Lists
RE: [gmx-users] CHARMM36 - Smaller Area per lipid for POPE - Why?
Dear Peter, I also used h-bonds and I also switch LJ interaction from 0.8 nm to 1.2 nm (as in Klauda's paper). I will retry with a more solvated membrane. Would you have any thought on how the NPAT ensemble might affect peptide-membrane interactions like I am studying i.e. peptide is totally solvated, then adsorb, and finally may insert? The paper on peptide-membrane interaction like this usually use united-atom lipid in the NPT ensemble. Most of the work I have seen on Charmm membrane in the NPAT ensemble were for embedded membrane protein. Thanks, Sebastien > Date: Mon, 13 Aug 2012 16:12:29 -0500 > From: p...@uab.edu > To: gmx-users@gromacs.org > Subject: Re: [gmx-users] CHARMM36 - Smaller Area per lipid for POPE - Why? > > Oh something I didn't mention: for bond constraints I used h-bonds instead > of all-bonds. This may or may not make a difference (although I switched to > h-bonds based on the suggestion of some charmm/lipid thread on here from > a couple of years ago). > > On 2012-08-09 12:34:19PM -0300, Sebastien Cote wrote: > > > > Dear Peter, > > > > Did you use any different simulation conditions for your POPC membrane? I > > tried many different ones for POPE, without never reproducing Klauda's > > results. I may try yours on my POPE membrane. > > > > In my simulations, I want to study peptide-membrane interactions. The > > peptide is not embedded in the membrane. It is initially completely > > solvated without any interactions with the membrane. Then, I want to look > > at its adsorption and degree of insertion in the membrane. For that system, > > I can not remove the CoM motion of the protein alone, otherwise it will not > > adsorb and insert in the membrane. > > > > I may try (as you suggested) to remove CoM of the bottom leaflet on one > > hand, and the peptide-upperleaflet on the other hand. My peptide is not > > very long (17 to 35 amino acids), so I believe that remove the CoM of the > > peptide-upperleaflet/bottomleaflet will not have any pernicious effect. > > What do you think? > > > > Thanks for the suggestion, > > > > Sébastien > > > > > > > Date: Wed, 8 Aug 2012 20:19:56 -0500 > > > From: p...@uab.edu > > > To: gmx-users@gromacs.org > > > Subject: Re: [gmx-users] CHARMM36 - Smaller Area per lipid for POPE - Why? > > > > > > Personally, I could remove the COM of each leaflet when equilibrating the > > > bilayer by itself (and as a side note I am not experiencing a similar > > > problem > > > with POPC that you're having with POPE...). However, after the protein is > > > embedded, I have gotten good results for my protein, which extends from > > > the > > > water through the entire membrane into more water, by using a whole System > > > COM removal. The introduction of my particular embedded protein acts as a > > > physical coupling between the water layers with the lipids (not to > > > mention if > > > I choose to model the lipid raft localization crosslink, it will have to > > > happen anyway). If your protein doesn't extend fully past both layers of > > > the > > > membrane you may want to stick with just coupling a Membrane+Protein+1 > > > layer > > > of water or Membrane+Protein and Water separately (like in Justin's KALP15 > > > tutorial). You will have to decide what you think is physically realistic > > > based on the interaction between the water, membrane, and protein when the > > > protein is embedded. (if your protein is assymetrically embedded you may > > > even > > > use the following COM groups: protein+involved leaflet, second leaflet, > > > water). > > > > > > On 2012-08-09 09:38:01AM +1000, Mark Abraham wrote: > > > > On 9/08/2012 3:28 AM, Sebastien Cote wrote: > > > > > Thanks for the suggestion. I tried it, but for my system the gain is > > > > > not significant. > > > > > > > > > > I was aware that it is preferable to remove the centre-of-mass for > > > > > each leaflet separately. However, in my tests, I removed the > > > > > center-of-mass of the membrane because I intent to simulate > > > > > peptide-membrane interactions. In such case, the center-of-mass of > > > > > the protein-membrane system is usually removed. Is their any way to > > > > > remove the CoM motion of each leaflet separately on one hand, and > > > > &g
RE: [gmx-users] CHARMM36 - Smaller Area per lipid for POPE - Why?
Dear Felipe, I will take a lot at this paper. Thanks! Sebastien > Date: Tue, 4 Aug 012 9::8::7 +200< > From: luis.pinedadecas...@lnu.se > To: gmx-users@gromacs.org > Subject: Re: [gmx-users] CHARMM6 - Smaller Area per lipid for POPE - Why? > > Hi Sébastien, > > I found the following paper very instructive about this issue (simulated > areas per lipid in bilayers): > > Jensen, M. et al. Simulations of a membrane anchored peptide: structure, > dynamics, and influence on bilayer properties. Biophys. J. (004))6,, 556--5< > > Take maybe a look at it, if you haven't done it already. > > Regards, > > Felipe > > On 8//3//012 1::2 PM, Peter C. Lai wrote: > > Oh something I didn't mention: for bond constraints I used h-bonds instead > > of all-bonds. This may or may not make a difference (although I switched to > > h-bonds based on the suggestion of some charmm/lipid thread on here from > > a couple of years ago). > > > > On 012--8--9 2::4::9PPM -300,, Sebastien Cote wrote: > >> Dear Peter, > >> > >> Did you use any different simulation conditions for your POPC membrane? I > >> tried many different ones for POPE, without never reproducing Klauda's > >> results. I may try yours on my POPE membrane. > >> > >> In my simulations, I want to study peptide-membrane interactions. The > >> peptide is not embedded in the membrane. It is initially completely > >> solvated without any interactions with the membrane. Then, I want to look > >> at its adsorption and degree of insertion in the membrane. For that > >> system, I can not remove the CoM motion of the protein alone, otherwise it > >> will not adsorb and insert in the membrane. > >> > >> I may try (as you suggested) to remove CoM of the bottom leaflet on one > >> hand, and the peptide-upperleaflet on the other hand. My peptide is not > >> very long (7 to 5 amino acids), so I believe that remove the CoM of the > >> peptide-upperleaflet/bottomleaflet will not have any pernicious effect. > >> What do you think? > >> > >> Thanks for the suggestion, > >> > >> Sébastien > >> > >> > >>> Date: Wed, Aug 012 0::9::6 -500< > >>> From: p...@uab.edu > >>> To: gmx-users@gromacs.org > >>> Subject: Re: [gmx-users] CHARMM6 - Smaller Area per lipid for POPE - Why? > >>> > >>> Personally, I could remove the COM of each leaflet when equilibrating the > >>> bilayer by itself (and as a side note I am not experiencing a similar > >>> problem > >>> with POPC that you're having with POPE...). However, after the protein is > >>> embedded, I have gotten good results for my protein, which extends from > >>> the > >>> water through the entire membrane into more water, by using a whole System > >>> COM removal. The introduction of my particular embedded protein acts as a > >>> physical coupling between the water layers with the lipids (not to > >>> mention if > >>> I choose to model the lipid raft localization crosslink, it will have to > >>> happen anyway). If your protein doesn't extend fully past both layers of > >>> the > >>> membrane you may want to stick with just coupling a Membrane+Protein+ > >>> layer > >>> of water or Membrane+Protein and Water separately (like in Justin's KALP5< > >>> tutorial). You will have to decide what you think is physically realistic > >>> based on the interaction between the water, membrane, and protein when the > >>> protein is embedded. (if your protein is assymetrically embedded you may > >>> even > >>> use the following COM groups: protein+involved leaflet, second leaflet, > >>> water). > >>> > >>> On 012--8--9 9::8::1AAM +000,, Mark Abraham wrote: > >>>> On //8//012 ::8 AM, Sebastien Cote wrote: > >>>>> Thanks for the suggestion. I tried it, but for my system the gain is > >>>>> not significant. > >>>>> > >>>>> I was aware that it is preferable to remove the centre-of-mass for each > >>>>> leaflet separately. However, in my tests, I removed the center-of-mass > >>>>> of the membrane because I intent to simulate peptide-membrane > >>>>> interactions. In such case, the center-of-mass of the protein-membrane > >>>
[gmx-users] Parameters for bonded interactions
Dear Gromacs users, In the topology file of the protein, we see every two atoms that share a bond, three atoms that share a bond angle, and four atoms that share a torsion angle. However, the parameters (equilibrium value, energy constant, phase) are not explicitly shown as Gromacs fetch them in the ffbonded table. My question: Is there a simple way to see (without having to look at the code) which parameters are associated to which bond length, bond angle and torsion angle? More specifically, I would like to know which torsion angle parameters is specifically associated to each torsion angle. In the ffbonded table, there is sometimes different groups of four atoms that can be associated to the same torsion angle as there is the 'X' atom which can be any atom. Thanks for your help, Sebastien -- gmx-users mailing listgmx-users@gromacs.org http://lists.gromacs.org/mailman/listinfo/gmx-users * Only plain text messages are allowed! * Please search the archive at http://www.gromacs.org/Support/Mailing_Lists/Search before posting! * Please don't post (un)subscribe requests to the list. Use the www interface or send it to gmx-users-requ...@gromacs.org. * Can't post? Read http://www.gromacs.org/Support/Mailing_Lists
RE: [gmx-users] CHARMM36 - Smaller Area per lipid for POPE - Why?
Dear Peter, Did you use any different simulation conditions for your POPC membrane? I tried many different ones for POPE, without never reproducing Klauda's results. I may try yours on my POPE membrane. In my simulations, I want to study peptide-membrane interactions. The peptide is not embedded in the membrane. It is initially completely solvated without any interactions with the membrane. Then, I want to look at its adsorption and degree of insertion in the membrane. For that system, I can not remove the CoM motion of the protein alone, otherwise it will not adsorb and insert in the membrane. I may try (as you suggested) to remove CoM of the bottom leaflet on one hand, and the peptide-upperleaflet on the other hand. My peptide is not very long (17 to 35 amino acids), so I believe that remove the CoM of the peptide-upperleaflet/bottomleaflet will not have any pernicious effect. What do you think? Thanks for the suggestion, Sébastien > Date: Wed, 8 Aug 2012 20:19:56 -0500 > From: p...@uab.edu > To: gmx-users@gromacs.org > Subject: Re: [gmx-users] CHARMM36 - Smaller Area per lipid for POPE - Why? > > Personally, I could remove the COM of each leaflet when equilibrating the > bilayer by itself (and as a side note I am not experiencing a similar problem > with POPC that you're having with POPE...). However, after the protein is > embedded, I have gotten good results for my protein, which extends from the > water through the entire membrane into more water, by using a whole System > COM removal. The introduction of my particular embedded protein acts as a > physical coupling between the water layers with the lipids (not to mention if > I choose to model the lipid raft localization crosslink, it will have to > happen anyway). If your protein doesn't extend fully past both layers of the > membrane you may want to stick with just coupling a Membrane+Protein+1 layer > of water or Membrane+Protein and Water separately (like in Justin's KALP15 > tutorial). You will have to decide what you think is physically realistic > based on the interaction between the water, membrane, and protein when the > protein is embedded. (if your protein is assymetrically embedded you may even > use the following COM groups: protein+involved leaflet, second leaflet, > water). > > On 2012-08-09 09:38:01AM +1000, Mark Abraham wrote: > > On 9/08/2012 3:28 AM, Sebastien Cote wrote: > > > Thanks for the suggestion. I tried it, but for my system the gain is not > > > significant. > > > > > > I was aware that it is preferable to remove the centre-of-mass for each > > > leaflet separately. However, in my tests, I removed the center-of-mass of > > > the membrane because I intent to simulate peptide-membrane interactions. > > > In such case, the center-of-mass of the protein-membrane system is > > > usually removed. Is their any way to remove the CoM motion of each > > > leaflet separately on one hand, and peptide-membrane system CoM motion on > > > the other? > > > > See 7.3.3 of manual. > > > > Mark > > > > > > > > Thanks, > > > > > > Sebastien > > > > > > > > >> Date: Fri, 3 Aug 2012 11:10:22 -0400 > > >> Subject: Re: [gmx-users] CHARMM36 - Smaller Area per lipid for POPE - > > >> Why? > > >> From: da...@cornell.edu > > >> To: gmx-users@gromacs.org > > >> > > >> Hello, > > >> > > >> I ran into similar issues for a DPPC bilayer. It might be possible > > >> that the two leaflets of the bilayer are moving with respect to > > >> eachother. If this is not taken into account, these artificial > > >> velocities will mean the simulation thinks it is at a higher > > >> temperature than it really is. If possible, you might want to try > > >> subtracting the center of mass motion of each leaflet, rather than the > > >> center of mass motion of the entire bilayer. This will allow the > > >> system to equillibrate to the correct (higher) temperature, and should > > >> increase the area per lipid of the bilayer. > > >> > > >> Hope this helps. > > >> -David > > >> > > >> On Thu, Aug 2, 2012 at 8:22 AM, Sebastien Cote > > >> wrote: > > >>> > > >>> Dear Gromacs users, > > >>> > > >>> I did new tests on the POPE membrane with CHARMM36 parameters, but I > > >>> still always get area per lipid values that are smaller than > >
RE: [gmx-users] CHARMM36 - Smaller Area per lipid for POPE - Why?
Thanks for the suggestion. I tried it, but for my system the gain is not significant. I was aware that it is preferable to remove the centre-of-mass for each leaflet separately. However, in my tests, I removed the center-of-mass of the membrane because I intent to simulate peptide-membrane interactions. In such case, the center-of-mass of the protein-membrane system is usually removed. Is their any way to remove the CoM motion of each leaflet separately on one hand, and peptide-membrane system CoM motion on the other? Thanks, Sebastien > Date: Fri, 3 Aug 2012 11:10:22 -0400 > Subject: Re: [gmx-users] CHARMM36 - Smaller Area per lipid for POPE - Why? > From: da...@cornell.edu > To: gmx-users@gromacs.org > > Hello, > > I ran into similar issues for a DPPC bilayer. It might be possible > that the two leaflets of the bilayer are moving with respect to > eachother. If this is not taken into account, these artificial > velocities will mean the simulation thinks it is at a higher > temperature than it really is. If possible, you might want to try > subtracting the center of mass motion of each leaflet, rather than the > center of mass motion of the entire bilayer. This will allow the > system to equillibrate to the correct (higher) temperature, and should > increase the area per lipid of the bilayer. > > Hope this helps. > -David > > On Thu, Aug 2, 2012 at 8:22 AM, Sebastien Cote > wrote: > > > > > > Dear Gromacs users, > > > > I did new tests on the POPE membrane with CHARMM36 parameters, but I still > > always get area per lipid values that are smaller than experimental value > > by 4 to 6 Angstrom2. Here are my new tests. > > > > My initial configuration is an equilibrated POPE membrane with 80 lipids at > > 1 atm and 310K in NPT. It was taken from Klauda's website and it was > > obtained from the study in which the POPE parameters were tested (Klauda, > > J. B. et al. 2010 J. Phys. Chem. B, 114, 7830-7843). > > > > I use TIPS3P (Charmm's special TIP3P). My simulations parameters are > > similar to those used in a previous tread on the Gromacs mailing list > > (http://lists.gromacs.org/pipermail/gmx-users/2010-October/055161.html for > > DMPC, POPC and DPPC of 128 lipids each) : > > > > dt = 0.002 ps; rlist = 1.0 nm; rlistlong = 1.4 nm; coulombtype = pme; > > rcoulomb = 1.4 nm; vdwtype = switch or cutoff (see below); DispCorr = No; > > fourierspacing = 0.15 nm; pme_order = 6; tcoupl = nose-hoover; tau_t = 1.0 > > ps; ref_t = 310K; pcoupl = Parrinello-Rahman; pcoupltype = semiisotropic; > > tau_p = 5.0 ps; compressibility = 4.5e-5; ref_p = 1.0 atm; constraints = > > h-bonds; constraint_algorithm = LINCS. Nochargegrps was used when executing > > pdb2gmx. > > > > The simulation time of each simulation is 100 ns. I tried different VdW > > cutoff values, since it was previously mentioned that cutoff values for VdW > > may influence the area per lipid. The average value and standard deviation > > are calculated on the 20 to 100 ns time interval. > > > > 1- For VdW switch from 0.8 to 1.2 nm : The area per lipid is 54.8 +/- 1.6 > > A2. > > 2- For VdW switch from 1.1 to 1.2 nm : The area per lipid is 54.6 +/- 1.8 > > A2. > > 3- For VdW cutoff at 1.4 nm : The area per lipid is 55.9 +/- 1.6 A2. > > > > I also checked the influence of DispCorr with VdW switch from 0.8 to 1.2 nm > > : > > > > 1- Without DispCorr : The area per lipid is 54.8 +/- 1.6 A2. > > 2- With DispCorr : The area per lipid is 54.4 +/- 1.9 A2. > > > > I also checked the influence of PME cutoff with VdW switch from 0.8 to 1.2 > > nm : > > > > 1- For PME cutoff at 1.4 nm : The area per lipid is 54.8 +/- 1.6 A2. > > 2- For PME cutoff at 1.0 nm : The area per lipid is 56.4 +/- 1.5 A2. > > > > These values are smaller than 4-6 A2 when compared against the experimental > > value (59.75-60.75 A2) and the value obtained in Klauda's simulation (59.2 > > +/- 0.3 A2). DispCorr and LJ cutoff weakly impact the results. Reducing the > > PME cutoff seems to have the greatest effect, but the value obtained is > > still smaller than experimental value by 3-4 A2. > > > > I also tried other initial configurations, but the results were either very > > similar or worst. > > > > Larger membrane gave similar results for the mean values and smaller > > standard deviations. > > > > --- > > > > Have anyone else tried to simulate a CHARMM36 POPE membrane in Gromacs? Do > > you get similar results? > > > > Is
RE: [gmx-users] CHARMM36 - Smaller Area per lipid for POPE - Why?
Dear Gromacs users, I did new tests on the POPE membrane with CHARMM36 parameters, but I still always get area per lipid values that are smaller than experimental value by 4 to 6 Angstrom2. Here are my new tests. My initial configuration is an equilibrated POPE membrane with 80 lipids at 1 atm and 310K in NPT. It was taken from Klauda's website and it was obtained from the study in which the POPE parameters were tested (Klauda, J. B. et al. 2010 J. Phys. Chem. B, 114, 7830-7843). I use TIPS3P (Charmm's special TIP3P). My simulations parameters are similar to those used in a previous tread on the Gromacs mailing list (http://lists.gromacs.org/pipermail/gmx-users/2010-October/055161.html for DMPC, POPC and DPPC of 128 lipids each) : dt = 0.002 ps; rlist = 1.0 nm; rlistlong = 1.4 nm; coulombtype = pme; rcoulomb = 1.4 nm; vdwtype = switch or cutoff (see below); DispCorr = No; fourierspacing = 0.15 nm; pme_order = 6; tcoupl = nose-hoover; tau_t = 1.0 ps; ref_t = 310K; pcoupl = Parrinello-Rahman; pcoupltype = semiisotropic; tau_p = 5.0 ps; compressibility = 4.5e-5; ref_p = 1.0 atm; constraints = h-bonds; constraint_algorithm = LINCS. Nochargegrps was used when executing pdb2gmx. The simulation time of each simulation is 100 ns. I tried different VdW cutoff values, since it was previously mentioned that cutoff values for VdW may influence the area per lipid. The average value and standard deviation are calculated on the 20 to 100 ns time interval. 1- For VdW switch from 0.8 to 1.2 nm : The area per lipid is 54.8 +/- 1.6 A2. 2- For VdW switch from 1.1 to 1.2 nm : The area per lipid is 54.6 +/- 1.8 A2. 3- For VdW cutoff at 1.4 nm : The area per lipid is 55.9 +/- 1.6 A2. I also checked the influence of DispCorr with VdW switch from 0.8 to 1.2 nm : 1- Without DispCorr : The area per lipid is 54.8 +/- 1.6 A2. 2- With DispCorr : The area per lipid is 54.4 +/- 1.9 A2. I also checked the influence of PME cutoff with VdW switch from 0.8 to 1.2 nm : 1- For PME cutoff at 1.4 nm : The area per lipid is 54.8 +/- 1.6 A2. 2- For PME cutoff at 1.0 nm : The area per lipid is 56.4 +/- 1.5 A2. These values are smaller than 4-6 A2 when compared against the experimental value (59.75-60.75 A2) and the value obtained in Klauda's simulation (59.2 +/- 0.3 A2). DispCorr and LJ cutoff weakly impact the results. Reducing the PME cutoff seems to have the greatest effect, but the value obtained is still smaller than experimental value by 3-4 A2. I also tried other initial configurations, but the results were either very similar or worst. Larger membrane gave similar results for the mean values and smaller standard deviations. --- Have anyone else tried to simulate a CHARMM36 POPE membrane in Gromacs? Do you get similar results? Is a 3-4 A2 deviation from experiment likely to influence my membrane/peptide simulations? Would it then be preferable to go with CHARMM27 in the NPAT ensemble? At this point, I have no clue of how to reproduce correctly Klauda's results for POPE. Any suggestion is welcomed. Thanks, Sebastien > Date: Mon, 23 Jul 2012 16:06:40 -0500 > From: p...@uab.edu > To: gmx-users@gromacs.org > Subject: Re: [gmx-users] CHARMM36 - Smaller Area per lipid for POPE - Why? > > On 2012-07-23 02:34:31PM -0300, Sebastien Cote wrote: > > > > There is not much difference when using DispCorr or not. At least on the > > same time scale as the simulation with switch cutoff from 0.8 to 1.2 nm and > > on the same time scale. > > > > Should DispCorr be used in all membrane simulations? I thought that we > > should always use this correction. > > I alwasy thought it was actually forcefield dependent. I never use it with > CHARMM since the mdp files I used as the basis for mine didn't with C27, and > I get acceptable APL with POPC when using the same mdp with C36. I haven't > compared the codes for CHARMM to see if dispcorr is builtin to the gromacs > implementation or not, but the reason I brought it up is that on past > mailing list discussions about TIPS3P, there were reports of significant > density differences with and without dispcorr. > > > > > > Thanks, > > > > Sebastien > > > > > > > Date: Fri, 20 Jul 2012 12:47:44 -0500 > > > From: p...@uab.edu > > > To: gmx-users@gromacs.org > > > Subject: Re: [gmx-users] CHARMM36 - Smaller Area per lipid for POPE - Why? > > > > > > Did you play with DispCorr? > > > > > > On 2012-07-20 09:46:13AM -0300, Sebastien Cote wrote: > > > > > > > > Dear Gromacs users, > > > > > > > > My simulations on a POPE membrane using the CHARMM36 parameters are > > > > giving
RE: [gmx-users] CHARMM36 - Smaller Area per lipid for POPE - Why?
There is not much difference when using DispCorr or not. At least on the same time scale as the simulation with switch cutoff from 0.8 to 1.2 nm and on the same time scale. Should DispCorr be used in all membrane simulations? I thought that we should always use this correction. Thanks, Sebastien > Date: Fri, 20 Jul 2012 12:47:44 -0500 > From: p...@uab.edu > To: gmx-users@gromacs.org > Subject: Re: [gmx-users] CHARMM36 - Smaller Area per lipid for POPE - Why? > > Did you play with DispCorr? > > On 2012-07-20 09:46:13AM -0300, Sebastien Cote wrote: > > > > Dear Gromacs users, > > > > My simulations on a POPE membrane using the CHARMM36 parameters are giving > > ''area per lipid'' values well below the experimental value (59.75-60.75 > > Angstroms2). Is their someone else experiencing a similar problem? If yes, > > how did you solved it? > > > > I did the following : > > > > I used the CHARMM36 parameters kindly provided by Thomas J. Piggot on the > > Users contribution section on Gromacs website. > > My starting configuration was taken from : > > http://terpconnect.umd.edu/~jbklauda/research/download.html > > It is a POPE membrane of 80 lipids equilibrated in NPT at T=310K and P=1atm > > for 40 ns. It is taken from the article Klauda, J. B. et al. 2010 J. Phys. > > Chem. B, 114, 7830-7843. > > > > At first, I tested normal TIP3P vs. CHARMM TIP3P and saw that normal TIP3P > > gives smaller Area per lipid of about 2-3 Angstroms. This was also observed > > by T.J. Piggot (personnal communication) and Tieleman (Sapay, N. et al. > > 2010 J. Comp. Chem. 32, 1400-1410). So, I will present only the simulations > > using CHARMM TIP3P. As in Klauda's paper, my simulations are at 310K and 1 > > atm. As them, I used a switch cutoff for vdw, and I used normal cutoff for > > PME. The simulations are 20 ns. I can send my .mdp file for more details. I > > varied the switch condition on vdw : > > > > 1- For a switch from 0.8 to 1.2 (as in Klauda's paper), I got Area per > > lipid of about 56.5 Angstroms2; whereas they got 59.2 in their paper, > > matching the experimental value of 59.75-60.75. > > 2- For a switch from 1.0 to 1.2, I got Area per lipid of about 53.5 > > Angstroms2, which is smaller than the previous cutoff. This is surprising > > since a previous thread on gromacs-users mailing lists said that increasing > > the lower cutoff, increased the Area per lipid or had not impact on POPC of > > DPPC. > > 3- For a switch from 1.1 to 1.2, I got Area per lipid of about 55 > > Angstroms2. > > 4- For a hard cutoff at 1.4, I got Area per lipid of about 52 Angstroms2. > > > > I also tried to re-equilibrate the membrane in the NPAT ensemble for 10 ns > > at 310K and 1 atm. Then, when I launched the simulation in NPT, I ended up > > with different results : > > > > 1- Switch from 0.8 to 1.2 gave a smaller area per lipid of 54 Angstroms2. > > 2- Switch from 1.0 to 1.2 gave a larger area per lipid of 55 Angstroms2. > > 4- Hard cutoff at 1.4 gave a similar area per lipid of 52.5 Angstroms2. > > > > I looked at the POPE paramaters for CHARMM36 in Gromacs, and they agree > > with the published parameters. > > > > Am I doing anything wrong? Is their someone else experiencing a similar > > problem for POPE? If yes, how did you solved it? > > > > Should I instead use CHARMM27 parameters in the NPAT ensemble? I want to > > study the interaction between a peptide and the POPE membrane. I am > > troubled that the NPAT ensemble might influence my results in a bad way. > > Also, I can not use OPLS AA nor GROMOS for the protein interactions because > > these force fields are not giving the correct structural ensemble for my > > peptide in solution. > > > > I am willing to send more information if you need. > > > > Thanks a lot, > > Sincerely, > > > > Sébastien -- > > gmx-users mailing list gmx-users@gromacs.org > > http://lists.gromacs.org/mailman/listinfo/gmx-users > > * Only plain text messages are allowed! > > * Please search the archive at > > http://www.gromacs.org/Support/Mailing_Lists/Search before posting! > > * Please don't post (un)subscribe requests to the list. Use the > > www interface or send it to gmx-users-requ...@gromacs.org. > > * Can't post? Read http://www.gromacs.org/Support/Mailing_Lists > > -- > == > Peter C. Lai | University of Alabama-Birmingham > Progra
[gmx-users] CHARMM36 - Smaller Area per lipid for POPE - Why?
Dear Gromacs users, My simulations on a POPE membrane using the CHARMM36 parameters are giving ''area per lipid'' values well below the experimental value (59.75-60.75 Angstroms2). Is their someone else experiencing a similar problem? If yes, how did you solved it? I did the following : I used the CHARMM36 parameters kindly provided by Thomas J. Piggot on the Users contribution section on Gromacs website. My starting configuration was taken from : http://terpconnect.umd.edu/~jbklauda/research/download.html It is a POPE membrane of 80 lipids equilibrated in NPT at T=310K and P=1atm for 40 ns. It is taken from the article Klauda, J. B. et al. 2010 J. Phys. Chem. B, 114, 7830-7843. At first, I tested normal TIP3P vs. CHARMM TIP3P and saw that normal TIP3P gives smaller Area per lipid of about 2-3 Angstroms. This was also observed by T.J. Piggot (personnal communication) and Tieleman (Sapay, N. et al. 2010 J. Comp. Chem. 32, 1400-1410). So, I will present only the simulations using CHARMM TIP3P. As in Klauda's paper, my simulations are at 310K and 1 atm. As them, I used a switch cutoff for vdw, and I used normal cutoff for PME. The simulations are 20 ns. I can send my .mdp file for more details. I varied the switch condition on vdw : 1- For a switch from 0.8 to 1.2 (as in Klauda's paper), I got Area per lipid of about 56.5 Angstroms2; whereas they got 59.2 in their paper, matching the experimental value of 59.75-60.75. 2- For a switch from 1.0 to 1.2, I got Area per lipid of about 53.5 Angstroms2, which is smaller than the previous cutoff. This is surprising since a previous thread on gromacs-users mailing lists said that increasing the lower cutoff, increased the Area per lipid or had not impact on POPC of DPPC. 3- For a switch from 1.1 to 1.2, I got Area per lipid of about 55 Angstroms2. 4- For a hard cutoff at 1.4, I got Area per lipid of about 52 Angstroms2. I also tried to re-equilibrate the membrane in the NPAT ensemble for 10 ns at 310K and 1 atm. Then, when I launched the simulation in NPT, I ended up with different results : 1- Switch from 0.8 to 1.2 gave a smaller area per lipid of 54 Angstroms2. 2- Switch from 1.0 to 1.2 gave a larger area per lipid of 55 Angstroms2. 4- Hard cutoff at 1.4 gave a similar area per lipid of 52.5 Angstroms2. I looked at the POPE paramaters for CHARMM36 in Gromacs, and they agree with the published parameters. Am I doing anything wrong? Is their someone else experiencing a similar problem for POPE? If yes, how did you solved it? Should I instead use CHARMM27 parameters in the NPAT ensemble? I want to study the interaction between a peptide and the POPE membrane. I am troubled that the NPAT ensemble might influence my results in a bad way. Also, I can not use OPLS AA nor GROMOS for the protein interactions because these force fields are not giving the correct structural ensemble for my peptide in solution. I am willing to send more information if you need. Thanks a lot, Sincerely, Sébastien -- gmx-users mailing listgmx-users@gromacs.org http://lists.gromacs.org/mailman/listinfo/gmx-users * Only plain text messages are allowed! * Please search the archive at http://www.gromacs.org/Support/Mailing_Lists/Search before posting! * Please don't post (un)subscribe requests to the list. Use the www interface or send it to gmx-users-requ...@gromacs.org. * Can't post? Read http://www.gromacs.org/Support/Mailing_Lists
