Re: [gmx-users] itp file for trehalose
You can refer to following papers: HALVOR S. HANSEN, PHILIPPE H. HÜNENBERGER , JComput Chem 32: 998?1032, 2011 and ? ROBERTO D. LINS, PHILIPPE H. HUNENBERGER. J Comput Chem 26: 1400 ?1412, 2005 Nisha Quoting ABEL Stephane 175950 stephane.a...@cea.fr: Dear All, I am looking for a topology file (*.itp) for trehalose for simulations with the GROMOS53A6 (or more recent). Does someboby know where I can find it? Thank you in advance for your help Stéphane -- Stéphane Abel, PhD CEA Saclay DSV/IBITEC-S/SB2SM CNRS URA2096 91191 Saclay, FRANCE URL: http://www.st-abel.com -- -- 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] sugar force fields
Here are few papers you can refer to: Gromos forcefield for hexopyronase-based carbohydrates - Roberto D.Lins and Philippe H. Hunenberger 2005 (In gromos forcefield there are parameters for galactose) An improved OPLS?AA force field for carbohydrates- D Kony 2002 Hope it helps! Nisha P Quoting Michael Brunsteiner mbx0...@yahoo.com: Dear all, I understand that most major force fields (charmm, opls, amber, gromos?...) come with parameters for sugars. I wonder: 1) if any work has been done/published comparing the accuracy of these different sugar force fields 2) Of the force fields that come with gromacs is there one for which making sugar topologies is as straight forward as making topologies for amino acids, or if there is no such thing for which FF is making a sugar topology the least tedious? BTW: I want to model simple sugar molecules, monomers, dimers (no polymers) thanks! Michael -- 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] oligoglycines
Thanks Tsjerk! Looks like it is working! Quoting Tsjerk Wassenaar tsje...@gmail.com: Hi Nisha, For building you can also use pymol if you have it installed. On the command line you can issue: pymol -qcd 'editor.build_peptide(GGG);cmd.save(triglycine.pdb,not hydro)' Hope it helps, Tsjerk On Mon, Mar 28, 2011 at 6:05 PM, nishap.pa...@utoronto.ca wrote: Hello, I want to simulate n-glycines (diglycine, triglycine..etc) I tried to get the structure from PRODRG, but the program adds H's on the N-terminal instead of on C- terminal. Is there another program I could use, or a site where I could get the structure of oligoglycines? For example, diglycine (NH2-CH2-CO-NH-CH2-COOH), but PRODRG gives me ((NH3-CH2-CO-NH-CH2-COO) Thanks! Nisha P -- gmx-users mailing list gmx-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 thewww interface or send it to gmx-users-requ...@gromacs.org. Can't post? Read http://www.gromacs.org/Support/Mailing_Lists -- Tsjerk A. Wassenaar, Ph.D. post-doctoral researcher Molecular Dynamics Group * Groningen Institute for Biomolecular Research and Biotechnology * Zernike Institute for Advanced Materials University of Groningen The Netherlands -- 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
[gmx-users] oligoglycines
Hello, I want to simulate n-glycines (diglycine, triglycine..etc) I tried to get the structure from PRODRG, but the program adds H's on the N-terminal instead of on C- terminal. Is there another program I could use, or a site where I could get the structure of oligoglycines? For example, diglycine (NH2-CH2-CO-NH-CH2-COOH), but PRODRG gives me ((NH3-CH2-CO-NH-CH2-COO) Thanks! Nisha P -- 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] oligoglycines
Thanks Justin! Nisha P Quoting Justin A. Lemkul jalem...@vt.edu: nishap.pa...@utoronto.ca wrote: Hello, I want to simulate n-glycines (diglycine, triglycine..etc) I tried to get the structure from PRODRG, but the program adds H's on the N-terminal instead of on C- terminal. Is there another program I could use, or a site where I could get the structure of oligoglycines? For example, diglycine (NH2-CH2-CO-NH-CH2-COOH), but PRODRG gives me ((NH3-CH2-CO-NH-CH2-COO) Protonation state is easy to change with pdb2gmx when writing your topology. Use -ignh and/or -ter to set the proper state. You can also tune this with PRODRG (see the FAQs regarding ADDHYD and DELHYD commands), but since you can do it all with Gromacs anyway, it shouldn't be necessary. -Justin Thanks! Nisha P -- Justin A. Lemkul Ph.D. Candidate ICTAS Doctoral Scholar MILES-IGERT Trainee Department of Biochemistry Virginia Tech Blacksburg, VA jalemkul[at]vt.edu | (540) 231-9080 http://www.bevanlab.biochem.vt.edu/Pages/Personal/justin -- 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] Carbohydrate simulation: problem with the terminal atoms
I am not sure if you have looked through this paper, but it gives the parameters for sugars. A new GROMOS force field for hexopyranose-based carbohydrates by RD Lins, 2005 A biomolecular force field based on the free enthalpy of hydration and solvation: The GROMOS force?field parameter sets 53A5 and 53A6- C Oostenbrink, 2004 Hope that helps. Nisha P Quoting Somaye Badieyan badie...@yahoo.com: Hi everyone, I am using g53a6 force filed for the simulation of cellulose. However it seems that parametrization works only for the sugar blocked linked together and no terminal atoms for terminal residues are defined there (I mean the start and end atom, H and O atom at the beginning linked to c4 and H atom at the last Sugar residue linked to O1). I tried to solve the problem by making the topology file removing the start OH and end H group from the initial PDB file to have the topology file made by g53a6 and later add the parameters of OH and H to the topology file. The problem is at the time of minimization almost no change in energy (not converged even after 7325 step) with the maximum force (about 9e+03) on the H of new added OH group (the start OH group). and when I checked the output of minimization step I found the HO position and C2 and C3, O3 atoms are deformed and the defined angles/bonds/dihedral of this new OH are not kept (the other H group the the other side, end, is fine). I know the problem is due to the way I defined the parameters but I do not know what the problem exactly is. Since the original topology file without the terminal atoms start from 1, the terminal atoms are added at the end of gro file and topology file (although OH group is part of first residue): the molecule is a cellotetetraose and here is the added part to top file (Atom number: H at end: 57 , O at start:58, H at start:59): ; residue 1001 GLCB rtp GLCB q 0.0 1 CH1 1001 GLCB C4 1 0.332 13.019 ; qtot 0.232 2 CH1 1001 GLCB C3 2 0.232 13.019 ; qtot 0.464 3 OA 1001 GLCB O3 2 -0.642 15.9994 ; qtot -0.178 4 H 1001 GLCB HO3 2 0.41 1.008 ; qtot 0.232 . . . 56 OA 1004 GLCB O1 20 -0.642 15.9994 ; qtot -0.20 57 H 1004 GLCB HO1 20 0.282 1.008 ; qtot 0 ;added later to 1000 GLCB 58 OA 1001 GLCB O4 1 -0.542 15.9994 ; 59 H 1001 GLCB HO4 1 0.442 1.008 ; [ bonds ] ; ai aj funct c0 c1 c2 c3 1 2 2 gb_26 1 11 2 gb_26 . . . 55 56 2 gb_20 56 57 2 gb_1 58 1 2 gb_20 58 59 2 gb_1 [ angles ] ; ai aj ak funct c0 c1 c2 c3 2 1 11 2 ga_8 1 2 3 2 ga_9 1 2 5 2 ga_8 . . . 55 56 57 2 ga_12 11 1 58 2 ga_9 58 1 2 2 ga_9 1 58 59 2 ga_12 [ dihedrals ] ; ai aj ak al funct c0 c1 c2 c3 c4 c5 11 1 2 3 1 gd_17 11 1 2 5 1 gd_34 2 1 11 8 1 gd_34 . . . 59 58 1 2 1 gd_30 58 1 2 3 1 gd_18 58 1 2 5 1 gd_17 58 1 11 8 1 gd_17 [ dihedrals ] ; ai aj ak al funct c0 c1 c2 c3 . . . 11 2 58 1 2 gi_2 I prepare this parameters based on the other OH group in the sugars and parameters of link oxygen (the shared oxygen at the position of glycosidic bond that is expected to be O4 group of the next residue). I did not use the PRODRG since i need 53a6 and I found for carbohydrate there is too much different between the parameters in 53a6 and 43a1. Somewhere it ahs been mentioned to change the terminal database file (.tdb), however I think it may make the problem more complicated. Your help is appreciated. -Somaye ... Somayesadat Badieyan PhD Candidate and Research Assistant Biological Syatems Engineering 201 Seitz Hall, Virginia Tech Blacksburg, VA 24060 ... -- 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
Re: [gmx-users] Nucleosides OPLS-AA
I apologize, but I am trying to simulate Cytidine (not cytosine). The parameters are given for cytosine. Nisha Quoting David van der Spoel sp...@xray.bmc.uu.se: On 2011-03-18 20.14, nishap.pa...@utoronto.ca wrote: Hello, I am trying to use OPLS-AA force field for simulating nucleosides eg. cytosine, adenosine etc. I found parameters for nucleotide bases (eg. 1-methylcytosine) but I haven't been able to find parameters for nucleosides. Does anyone know where I can find parameters for nucleotides for OPLS-AA (if they do exist?). A paper citation would be helpful. Thanks. Nisha Patel It's all there in the atomtypes.atp file, a little fragment: opls_336 12.01100 ; Cytosine C4 Nucleotide base opls_337 12.01100 ; Cytosine C5 parameters: opls_338 12.01100 ; Cytosine C6 JACS,113,2810(1991) opls_3391.00800 ; Cytosine H-N1 -- David van der Spoel, Ph.D., Professor of Biology Dept. of Cell Molec. Biol., Uppsala University. Box 596, 75124 Uppsala, Sweden. Phone: +46184714205. sp...@xray.bmc.uu.sehttp://folding.bmc.uu.se -- 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] Nucleosides OPLS-AA
I know and I apologize, I mistyped. It is an isolated molecule. I did try to combine the parameters of OPLSAA of carbohydrates(Jorgensen, 1997) with the nucleotide base paper (Jorgensen, 1991). However, I am not sure the partial charge that I should use on the Carbon (C) atom of the sugar molecule that is linked with the Nitrogen (N) of nucleoside, because it is not an alcohol anymore. Quoting David van der Spoel sp...@xray.bmc.uu.se: On 2011-03-18 20.41, nishap.pa...@utoronto.ca wrote: I apologize, but I am trying to simulate Cytidine (not cytosine). The parameters are given for cytosine. That's not what you asked for, but in that case you have to combine the cytosine parameters with some kind of sugar ring. Is this in a polymer through the sugar rings as well, or an isolated molecule? Nisha Quoting David van der Spoel sp...@xray.bmc.uu.se: On 2011-03-18 20.14, nishap.pa...@utoronto.ca wrote: Hello, I am trying to use OPLS-AA force field for simulating nucleosides eg. cytosine, adenosine etc. I found parameters for nucleotide bases (eg. 1-methylcytosine) but I haven't been able to find parameters for nucleosides. Does anyone know where I can find parameters for nucleotides for OPLS-AA (if they do exist?). A paper citation would be helpful. Thanks. Nisha Patel It's all there in the atomtypes.atp file, a little fragment: opls_336 12.01100 ; Cytosine C4 Nucleotide base opls_337 12.01100 ; Cytosine C5 parameters: opls_338 12.01100 ; Cytosine C6 JACS,113,2810(1991) opls_339 1.00800 ; Cytosine H-N1 -- David van der Spoel, Ph.D., Professor of Biology Dept. of Cell Molec. Biol., Uppsala University. Box 596, 75124 Uppsala, Sweden. Phone: +46184714205. sp...@xray.bmc.uu.se http://folding.bmc.uu.se -- gmx-users mailing list gmx-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 -- David van der Spoel, Ph.D., Professor of Biology Dept. of Cell Molec. Biol., Uppsala University. Box 596, 75124 Uppsala, Sweden. Phone: +46184714205. sp...@xray.bmc.uu.sehttp://folding.bmc.uu.se -- 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] oplsaa galactose
Hi, I am trying to simulate hexopyronase using OPLS-AA forcefield using parameters from the paper: An improved OPLS?AA force field for carbohydrates, 2002 by Gunsteren. The set of torsional angle parameters in the paper are given in Kcal/mol, but because I am using OPLS-AA I would need to convert the functions into RB-types, right? Is there a way to convert these functions to RB-types? Or do I convert them from Kcal/mol to KJ/mol? Thanks. Nisha P Quoting nishap.pa...@utoronto.ca: Thanks Justin! Quoting Justin A. Lemkul jalem...@vt.edu: nishap.pa...@utoronto.ca wrote: Hello, I am trying to simulate hexopyronase using OPLS-AA forcefield using parameters from the paper: An improved OPLS?AA force field for carbohydrates, 2002 by Gunsteren. I was looking through the ffoplsaabon.itp file and there are some dihedral parameters for hexopyronase. I am not sure what the comment below means. Below are extra dihedrals for some special organic molecules. ; Since the atom types are identical to other dihedrals you have to specify ; them explicitly with a define if you happen to simulate this type of molecule. How would I mention it in my rtp file? I don't understand ' specify them explicitly with a define' This is my .rtp for Galactose See examples in aminoacids.rtp, e.g. ARG. -Justin [ GLA ] [ atoms ] Oopls_180 -0.400 1 C1 opls_1950.365 2 H1opls_1960.100 2 O1opls_154 -0.683 2 HO1 opls_1550.418 2 C2 opls_1580.205 3 H2opls_1760.060 3 O2opls_154 -0.683 3 HO2 opls_1550.418 3 C3 opls_1580.205 4 H3opls_1760.060 4 O3opls_154 -0.683 4 HO3 opls_1550.418 4 C4 opls_1580.205 5 H4opls_1760.060 5 O4opls_154 -0.683 5 HO4 opls_1550.418 5 C5 opls_1830.170 6 H5opls_1850.030 6 C6opls_1570.145 7 H61 opls_1760.060 7 H62 opls_1760.060 7 O6opls_154 -0.683 7 HO6 opls_1550.418 7 [ bonds ] OC1 C1 H1 C1 O1 O1 HO1 C1 C2 C2 H2 C2 O2 O2 HO2 C2 C3 C3 H3 C3 O3 O3 HO3 C3 C4 C4 H4 C4 O4 O4 HO4 C4 C5 C5 C6 C5 O C6 H61 C6 H62 C6 O6 O6 HO6 -- Justin A. Lemkul Ph.D. Candidate ICTAS Doctoral Scholar MILES-IGERT Trainee Department of Biochemistry Virginia Tech Blacksburg, VA jalemkul[at]vt.edu | (540) 231-9080 http://www.bevanlab.biochem.vt.edu/Pages/Personal/justin -- 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
[gmx-users] oplsaa galactose
Hello, I am trying to simulate hexopyronase using OPLS-AA forcefield using parameters from the paper: An improved OPLS?AA force field for carbohydrates, 2002 by Gunsteren. I was looking through the ffoplsaabon.itp file and there are some dihedral parameters for hexopyronase. I am not sure what the comment below means. Below are extra dihedrals for some special organic molecules. ; Since the atom types are identical to other dihedrals you have to specify ; them explicitly with a define if you happen to simulate this type of molecule. How would I mention it in my rtp file? I don't understand ' specify them explicitly with a define' This is my .rtp for Galactose [ GLA ] [ atoms ] Oopls_180 -0.400 1 C1opls_1950.365 2 H1opls_1960.100 2 O1opls_154 -0.683 2 HO1 opls_1550.418 2 C2opls_1580.205 3 H2opls_1760.060 3 O2opls_154 -0.683 3 HO2 opls_1550.418 3 C3opls_1580.205 4 H3opls_1760.060 4 O3opls_154 -0.683 4 HO3 opls_1550.418 4 C4opls_1580.205 5 H4opls_1760.060 5 O4opls_154 -0.683 5 HO4 opls_1550.418 5 C5opls_1830.170 6 H5opls_1850.030 6 C6opls_1570.145 7 H61 opls_1760.060 7 H62 opls_1760.060 7 O6opls_154 -0.683 7 HO6 opls_1550.418 7 [ bonds ] OC1 C1 H1 C1 O1 O1 HO1 C1 C2 C2 H2 C2 O2 O2 HO2 C2 C3 C3 H3 C3 O3 O3 HO3 C3 C4 C4 H4 C4 O4 O4 HO4 C4 C5 C5 C6 C5 O C6 H61 C6 H62 C6 O6 O6 HO6 -- 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] oplsaa galactose
Thanks Justin! Quoting Justin A. Lemkul jalem...@vt.edu: nishap.pa...@utoronto.ca wrote: Hello, I am trying to simulate hexopyronase using OPLS-AA forcefield using parameters from the paper: An improved OPLS?AA force field for carbohydrates, 2002 by Gunsteren. I was looking through the ffoplsaabon.itp file and there are some dihedral parameters for hexopyronase. I am not sure what the comment below means. Below are extra dihedrals for some special organic molecules. ; Since the atom types are identical to other dihedrals you have to specify ; them explicitly with a define if you happen to simulate this type of molecule. How would I mention it in my rtp file? I don't understand ' specify them explicitly with a define' This is my .rtp for Galactose See examples in aminoacids.rtp, e.g. ARG. -Justin [ GLA ] [ atoms ] Oopls_180 -0.400 1 C1 opls_1950.365 2 H1opls_1960.100 2 O1opls_154 -0.683 2 HO1 opls_1550.418 2 C2 opls_1580.205 3 H2opls_1760.060 3 O2opls_154 -0.683 3 HO2 opls_1550.418 3 C3 opls_1580.205 4 H3opls_1760.060 4 O3opls_154 -0.683 4 HO3 opls_1550.418 4 C4 opls_1580.205 5 H4opls_1760.060 5 O4opls_154 -0.683 5 HO4 opls_1550.418 5 C5 opls_1830.170 6 H5opls_1850.030 6 C6opls_1570.145 7 H61 opls_1760.060 7 H62 opls_1760.060 7 O6opls_154 -0.683 7 HO6 opls_1550.418 7 [ bonds ] OC1 C1 H1 C1 O1 O1 HO1 C1 C2 C2 H2 C2 O2 O2 HO2 C2 C3 C3 H3 C3 O3 O3 HO3 C3 C4 C4 H4 C4 O4 O4 HO4 C4 C5 C5 C6 C5 O C6 H61 C6 H62 C6 O6 O6 HO6 -- Justin A. Lemkul Ph.D. Candidate ICTAS Doctoral Scholar MILES-IGERT Trainee Department of Biochemistry Virginia Tech Blacksburg, VA jalemkul[at]vt.edu | (540) 231-9080 http://www.bevanlab.biochem.vt.edu/Pages/Personal/justin -- 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
[gmx-users] galactose bond stretch
Hello, I ran a simulation of one molecule of galactose (cyclic) in water. After the simulation run, when I checked the trajectory file in VMD, the bonds in the galactose molecule stretched and during the run changed back to its original starting form. I am using GROMOS force field ffG53A6 and parameters as mentioned for carbohydrates in the forcefield, so I am not sure what went wrong.I generated topology using 'pdb2gmx' command. Also GROMACS did not give any warning during the run and my simulation did not crash. Any insights? Thanks. Nisha P -- 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] galactose bond stretch
Thanks Justin! I am using constraints, but like you said it could be just PBC. I did compare some of my calculations to experimental values and they are fairly similar. Nisha P Quoting Justin A. Lemkul jalem...@vt.edu: nishap.pa...@utoronto.ca wrote: Hello, I ran a simulation of one molecule of galactose (cyclic) in water. After the simulation run, when I checked the trajectory file in VMD, the bonds in the galactose molecule stretched and during the run changed back to its original starting form. I am using GROMOS force field ffG53A6 and parameters as mentioned for carbohydrates in the forcefield, so I am not sure what went wrong.I generated topology using 'pdb2gmx' command. Also GROMACS did not give any warning during the run and my simulation did not crash. Any insights? You can check bond length distributions with g_bond. Are you using constraints? Are you sure this isn't just a periodicity artifact during visualization? If the bonds deviated severely from equilibrium values or constraint lengths, the simulation would have crashed, so I suspect that what you're seeing is just PBC. -Justin Thanks. Nisha P -- Justin A. Lemkul Ph.D. Candidate ICTAS Doctoral Scholar MILES-IGERT Trainee Department of Biochemistry Virginia Tech Blacksburg, VA jalemkul[at]vt.edu | (540) 231-9080 http://www.bevanlab.biochem.vt.edu/Pages/Personal/justin -- 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
[gmx-users] galactose
Hello, I am trying to simulate galactose in water using Gromos FFG53A6 force field. The molecular formula for galactose is C6H12O6, however in the FFG53A6.rtp file, the molecular formula is C6H10O5 for galactose-A and B. Why is it different? I would really appreciate some help! Thanks -Nisha P -- 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] galactose
Hello, I am trying to simulate galactose in water using Gromos FFG53A6 force field. The molecular formula for galactose is C6H12O6, however in the FFG53A6.rtp file, the molecular formula is C6H10O5 for galactose-A and B. Why is it different? I would really appreciate some help! Thanks -Nisha P -- 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] galactose
Thanks. I am going to look over it. But I was wondering is it possible to just simulate one molecule of sucrose (glucose+fructose) in water using any of the force fields by Gromacs? I realize I would have to add the parameters to the .rtp files, but as you mentioned that the force fields recognizes them as part of the larger polymer and so it doesn't take into account the missing atoms. -Nisha P Quoting Justin A. Lemkul jalem...@vt.edu: nishap.pa...@utoronto.ca wrote: Hello, I am trying to simulate galactose in water using Gromos FFG53A6 force field. The molecular formula for galactose is C6H12O6, however in the FFG53A6.rtp file, the molecular formula is C6H10O5 for galactose-A and B. Why is it different? I would really appreciate some help! The Gromos96 sugar parameters assume the residues are part of a larger polymer, thus condensation of sugars eliminates the missing H2O. There is a newer version of these parameters that was recently published, which may have some relevant details: dx.doi.org/10.1002/jcc.21675 -Justin Thanks -Nisha P -- Justin A. Lemkul Ph.D. Candidate ICTAS Doctoral Scholar MILES-IGERT Trainee Department of Biochemistry Virginia Tech Blacksburg, VA jalemkul[at]vt.edu | (540) 231-9080 http://www.bevanlab.biochem.vt.edu/Pages/Personal/justin -- 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
[gmx-users] Energy-groups?
Hello, I want to plot the interaction potential energy between my solute and solvent. In my .mdp file I did not mention anything under energygrps,so I am thinking it calculates the energies for the whole system. But is there a way I can extract say for example LJ-14 term between my solute and solvent using the same .edr file? Or would I have to specify my energygrps and run the simulation again. Thanks. Nisha -- 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] Energy-groups?
Quoting Justin A. Lemkul jalem...@vt.edu: Thanks Justin! I do want the non-bonded potential between my solute and solvent. So in my .mdp file I put my solute and solvent as energygrps and ran mdrun using this command: mdrun -s md1.tpr(including energygrps) -rerun md.xtc (my trajectory) Is that correct? I don't understand how it would be faster or any different from running the simulation from scratch, because in my .mdp file it still has time of 100ns. I just modified the energygrps. nishap.pa...@utoronto.ca wrote: Hello, I want to plot the interaction potential energy between my solute and solvent. In my .mdp file I did not mention anything under energygrps,so I am thinking it calculates the energies for the whole system. But is there a way I can extract say for example LJ-14 term between my solute and solvent using the same .edr file? Or would I have to specify my energygrps and run the simulation again. 1-4 interactions are intramolecular, so there should be no solute-solvent 1-4 term. If you want nonbonded potentials decomposed, yes, you have to rerun your trajectory (mdrun -rerun, not from scratch) with a .tpr file specifying the desired groups. -Justin Thanks. Nisha -- Justin A. Lemkul Ph.D. Candidate ICTAS Doctoral Scholar MILES-IGERT Trainee Department of Biochemistry Virginia Tech Blacksburg, VA jalemkul[at]vt.edu | (540) 231-9080 http://www.bevanlab.biochem.vt.edu/Pages/Personal/justin -- 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] Energy-groups?
I see. It did work . Thanks a lot Justin! Quoting Justin A. Lemkul jalem...@vt.edu: nishap.pa...@utoronto.ca wrote: Quoting Justin A. Lemkul jalem...@vt.edu: Thanks Justin! I do want the non-bonded potential between my solute and solvent. So in my .mdp file I put my solute and solvent as energygrps and ran mdrun using this command: mdrun -s md1.tpr(including energygrps) -rerun md.xtc (my trajectory) Is that correct? I don't understand how it would be faster or any different from running the simulation from scratch, because in my .mdp file it still has time of 100ns. I just modified the energygrps. It should be significantly faster. In this case, mdrun is not doing the integration, it's simply using the known positions to re-calculate energies. -Justin nishap.pa...@utoronto.ca wrote: Hello, I want to plot the interaction potential energy between my solute and solvent. In my .mdp file I did not mention anything under energygrps,so I am thinking it calculates the energies for the whole system. But is there a way I can extract say for example LJ-14 term between my solute and solvent using the same .edr file? Or would I have to specify my energygrps and run the simulation again. 1-4 interactions are intramolecular, so there should be no solute-solvent 1-4 term. If you want nonbonded potentials decomposed, yes, you have to rerun your trajectory (mdrun -rerun, not from scratch) with a .tpr file specifying the desired groups. -Justin Thanks. Nisha -- Justin A. Lemkul Ph.D. Candidate ICTAS Doctoral Scholar MILES-IGERT Trainee Department of Biochemistry Virginia Tech Blacksburg, VA jalemkul[at]vt.edu | (540) 231-9080 http://www.bevanlab.biochem.vt.edu/Pages/Personal/justin -- 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 -- Justin A. Lemkul Ph.D. Candidate ICTAS Doctoral Scholar MILES-IGERT Trainee Department of Biochemistry Virginia Tech Blacksburg, VA jalemkul[at]vt.edu | (540) 231-9080 http://www.bevanlab.biochem.vt.edu/Pages/Personal/justin -- 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
[gmx-users] LJ Plot
Hello, Is there an option to plot Lennard Jones potential? I tried looking through the list and manual but I did not find any suggestions on how I could plot a LJ 6-12 potential plot. Nisha -- 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] LJ Plot
Quoting Justin A. Lemkul jalem...@vt.edu: How does sigeps read the input file? I tried using sigeps to plot LJ for different solutes in solvent and it gives me the same C6 and C12 value, when it runs, which seems a weird. c6= 1.0e-03, c12= 1.0e-06 sigma = 0.0, epsilon = 0.0 Van der Waals minimum at 0, V = nan Is it the default? nishap.pa...@utoronto.ca wrote: Hello, Is there an option to plot Lennard Jones potential? I tried looking through the list and manual but I did not find any suggestions on how I could plot a LJ 6-12 potential plot. g_sigeps -Justin Nisha -- Justin A. Lemkul Ph.D. Candidate ICTAS Doctoral Scholar MILES-IGERT Trainee Department of Biochemistry Virginia Tech Blacksburg, VA jalemkul[at]vt.edu | (540) 231-9080 http://www.bevanlab.biochem.vt.edu/Pages/Personal/justin -- 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] LJ Plot
Quoting Justin A. Lemkul jalem...@vt.edu: I actually want to plot my simulation with ad without the attractive term C6 of van der Waals and superimpose them to see the difference. nishap.pa...@utoronto.ca wrote: Quoting Justin A. Lemkul jalem...@vt.edu: How does sigeps read the input file? I tried using sigeps to plot LJ for different solutes in solvent and it gives me the same C6 and C12 value, when it runs, which seems a weird. c6= 1.0e-03, c12= 1.0e-06 sigma = 0.0, epsilon = 0.0 Van der Waals minimum at 0, V = nan g_sigeps does not take any input file. It produces an LJ curve based on command line parameters. I guess this is not what you are after, although that's exactly what it sounded like from your last post. Perhaps you need to re-phrase your question. If you just want to plot van der Waals energy terms, use g_energy to pull them out of the .edr file. -Justin Is it the default? nishap.pa...@utoronto.ca wrote: Hello, Is there an option to plot Lennard Jones potential? I tried looking through the list and manual but I did not find any suggestions on how I could plot a LJ 6-12 potential plot. g_sigeps -Justin Nisha -- Justin A. Lemkul Ph.D. Candidate ICTAS Doctoral Scholar MILES-IGERT Trainee Department of Biochemistry Virginia Tech Blacksburg, VA jalemkul[at]vt.edu | (540) 231-9080 http://www.bevanlab.biochem.vt.edu/Pages/Personal/justin -- 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 -- Justin A. Lemkul Ph.D. Candidate ICTAS Doctoral Scholar MILES-IGERT Trainee Department of Biochemistry Virginia Tech Blacksburg, VA jalemkul[at]vt.edu | (540) 231-9080 http://www.bevanlab.biochem.vt.edu/Pages/Personal/justin -- 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] LJ Plot
Quoting ms deviceran...@gmail.com: I did two different simulation. One with van der Waals attractive term and one without van der Waals attractive term. And so to see the difference between the LJ potential between the two simulation I wanted to plot the LJ potential for both the simulations and superimpose them. Because if there is no attraction LJ-12 (i.e. 1/r^-6 set to zero )LJ potential plot would be different than the one with LJ 6-12 potential. On 06/12/10 23:26, nishap.pa...@utoronto.ca wrote: Quoting Justin A. Lemkul jalem...@vt.edu: I actually want to plot my simulation with ad without the attractive term C6 of van der Waals and superimpose them to see the difference. The expression plot my simulation makes no sense. You don't plot a simulation, you plot variables calculated from a trajectory. Also, what does it mean with and without the attractive term? You perhaps want to *do* two different simulations? nishap.pa...@utoronto.ca wrote: Quoting Justin A. Lemkul jalem...@vt.edu: How does sigeps read the input file? I tried using sigeps to plot LJ for different solutes in solvent and it gives me the same C6 and C12 value, when it runs, which seems a weird. c6 = 1.0e-03, c12 = 1.0e-06 sigma = 0.0, epsilon = 0.0 Van der Waals minimum at 0, V = nan g_sigeps does not take any input file. It produces an LJ curve based on command line parameters. I guess this is not what you are after, although that's exactly what it sounded like from your last post. Perhaps you need to re-phrase your question. If you just want to plot van der Waals energy terms, use g_energy to pull them out of the .edr file. -Justin Is it the default? nishap.pa...@utoronto.ca wrote: Hello, Is there an option to plot Lennard Jones potential? I tried looking through the list and manual but I did not find any suggestions on how I could plot a LJ 6-12 potential plot. g_sigeps -Justin Nisha -- Justin A. Lemkul Ph.D. Candidate ICTAS Doctoral Scholar MILES-IGERT Trainee Department of Biochemistry Virginia Tech Blacksburg, VA jalemkul[at]vt.edu | (540) 231-9080 http://www.bevanlab.biochem.vt.edu/Pages/Personal/justin -- gmx-users mailing list gmx-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 -- Justin A. Lemkul Ph.D. Candidate ICTAS Doctoral Scholar MILES-IGERT Trainee Department of Biochemistry Virginia Tech Blacksburg, VA jalemkul[at]vt.edu | (540) 231-9080 http://www.bevanlab.biochem.vt.edu/Pages/Personal/justin -- gmx-users mailing list gmx-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 -- Massimo Sandal, Ph.D. http://devicerandom.org -- 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] C6 attractive term off OPLSAA
I ran the simulation using a different force field ffG53a6. I modified the ffG53a6nb.itp file by changing the term C6 to zero for nonbonded parameters, but this time for urea in water. The simulation ran fine without any warning or exploding. I don't understand why it would work with one force-field and not with OPLS-AA. Also using mdrun -rerun, I would basically use the nb.itp file with nonbonded parameters as mentioned earlier and everything else is the same i.e. my .mdp parameters and do another run using my previous trajectory file? I tried to look through some of the posts for using -rerun but I don't understand how that would not still blow up the system. I would like to give it a try for sure but I am not quite sure how I would use the command as such. I really appreciate the help. Thanks. Quoting Mark Abraham mark.abra...@anu.edu.au: On 1/12/2010 5:06 AM, Justin A. Lemkul wrote: nishap.pa...@utoronto.ca wrote: I tried using the nonbonded parameters as defined below in my ffoplsaanb.itp file for methanol in water and this is the syntax I used: [nonbond_params] ;ijfuncc6c12 opls_154opls_11110.00E+0002.43E-006 opls_154opls_11210.00E+0000.00E+000 opls_154opls_112 10.00E+0000.00E+000 opls_155opls_111 10.00E+0000.00E+000 opls_155opls_112 10.00E+0000.00E+000 opls_155opls_112 10.00E+0000.00E+000 opls_156opls_111 10.00E+0002.70E-007 opls_156opls_112 10.00E+0000.00E+000 opls_156opls_11210.00E+0000.00E+000 opls_157opls_111 10.00E+0003.01E-006 opls_157opls_112 10.00E+0000.00E+000 When I tried to do mdrun, I got an error saying my system is exploding. I tried doing the mdrun without nonbonded parameters and it runs fine. So I am not sure if I am using the nonbond_params concept correctly. i.e. I want C6 to be zero between my solute (methanol) and solvent (water). This is the error: Warning: 1-4 interaction between 2 and 6 at distance 3.786 which is larger than the 1-4 table size 2.500 nm These are ignored for the rest of the simulation This usually means your system is exploding, if not, you should increase table-extension in your mdp file or with user tables increase the table size Suggestions? Your implementation of your concept is technically correct, but perhaps physically unreasonable. You're turning off the attractive LJ component and leaving only a repulsive component. It sounds about right that everything is flying apart. Indeed. The way to do this is with mdrun -rerun on a trajectory generated with the normal version of the forcefield, as I think I said in another thread yesterday. Mark -Justin Thanks. -Nisha P Quoting nishap.pa...@utoronto.ca: Okay I am going to give it a try. I just wanted to make sure I was calculating C6 and C12 correctly as well using sigma and epsilon according to rule 3 C12 = Sigma^(6)*C6 C6 = 4*sigma^(6)*epsilon Thanks Quoting Justin A. Lemkul jalem...@vt.edu: nishap.pa...@utoronto.ca wrote: Hello, I have a concern regarding C6 attractive term in LJ potential. I am using OPLS-AA force field, and I wish to turn off attractive term C6 by setting the parameters to zero. One of the suggestion was to add a [nonbond_params] in my ffoplsaanb.itp file and set the C6 to zero between the non-bonded pair. In my system for example, which consists of one methanol in water, I wish to set C6 term to zero between my solute and solvent. Since OPLSAA is all atom force field it treats each atom individually and has sigma and epsilon for each atom, so I am not sure how I would actually set my nonbond_params in my nb.itp file. I realize I need to convert each sigma and epsilon to C6 and C12, so say for example for methanol in water my [nonbond_params] would look something like this? [ nonbond_params ] ; ij func c6 c12 CTOW 1 0.00 calculated value for C12 here? CTHW1 1 0.00 CTHW2 1 0.00 CT is the carbon in Methanol. OW, HW1 an HW2 correspond to atoms in TIP3P water model Is that correct? Would I have to do that for each atom in methanol? Sounds about right to me. -Justin Any suggestions would be appreciated. Thanks. Nisha P -- Justin A. Lemkul Ph.D. Candidate ICTAS Doctoral Scholar MILES-IGERT Trainee Department of Biochemistry Virginia Tech Blacksburg, VA jalemkul[at]vt.edu | (540) 231-9080 http://www.bevanlab.biochem.vt.edu/Pages/Personal/justin -- 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
Re: [gmx-users] C6 attractive term off OPLSAA
I see. I actually want to compare my RDFs with C6 term off. Earlier I tried using force.c code file and turned C6 = 0, but when i compared my RDFs, it didn't look any different so I am not sure if it even worked at all or made any difference to the simulation, but again I was using OPLS-AA and the nb.itp file itself does not indicate value for C6 or C12. If there any other way I could achieve this? Quoting Justin A. Lemkul jalem...@vt.edu: nishap.pa...@utoronto.ca wrote: I ran the simulation using a different force field ffG53a6. I modified the ffG53a6nb.itp file by changing the term C6 to zero for nonbonded parameters, but this time for urea in water. The simulation ran fine without any warning or exploding. I don't understand why it would work with one force-field and not with OPLS-AA. I don't understand why it worked at all, frankly, but more than likely you've been lucky enough to make 53a6 work, but OPLS-AA failed. You're making significant alterations to the nonbonded interactions of the system. Also using mdrun -rerun, I would basically use the nb.itp file with nonbonded parameters as mentioned earlier and everything else is the same i.e. my .mdp parameters and do another run using my previous trajectory file? I tried to look through some of the posts for using -rerun but I don't understand how that would not still blow up the system. I would like to give it a try for sure but I am not quite sure how I would use the command as such. Using mdrun -rerun recalculates energies from an existing trajectory. You would generate a trajectory with a normal force field model, then create a new .tpr file that has your modified potential, such that you can decompose the energy terms cleverly. If your goal is to generate trajectories to study the effects of using a modified potential (i.e., C6 = 0), then you can't use -rerun. -Justin I really appreciate the help. Thanks. Quoting Mark Abraham mark.abra...@anu.edu.au: On 1/12/2010 5:06 AM, Justin A. Lemkul wrote: nishap.pa...@utoronto.ca wrote: I tried using the nonbonded parameters as defined below in my ffoplsaanb.itp file for methanol in water and this is the syntax I used: [nonbond_params] ;ijfuncc6c12 opls_154opls_11110.00E+0002.43E-006 opls_154opls_11210.00E+0000.00E+000 opls_154opls_112 10.00E+0000.00E+000 opls_155opls_111 10.00E+0000.00E+000 opls_155opls_112 10.00E+0000.00E+000 opls_155opls_112 10.00E+0000.00E+000 opls_156opls_111 10.00E+0002.70E-007 opls_156opls_112 10.00E+0000.00E+000 opls_156opls_11210.00E+0000.00E+000 opls_157opls_111 10.00E+0003.01E-006 opls_157opls_112 10.00E+0000.00E+000 When I tried to do mdrun, I got an error saying my system is exploding. I tried doing the mdrun without nonbonded parameters and it runs fine. So I am not sure if I am using the nonbond_params concept correctly. i.e. I want C6 to be zero between my solute (methanol) and solvent (water). This is the error: Warning: 1-4 interaction between 2 and 6 at distance 3.786 which is larger than the 1-4 table size 2.500 nm These are ignored for the rest of the simulation This usually means your system is exploding, if not, you should increase table-extension in your mdp file or with user tables increase the table size Suggestions? Your implementation of your concept is technically correct, but perhaps physically unreasonable. You're turning off the attractive LJ component and leaving only a repulsive component. It sounds about right that everything is flying apart. Indeed. The way to do this is with mdrun -rerun on a trajectory generated with the normal version of the forcefield, as I think I said in another thread yesterday. Mark -Justin Thanks. -Nisha P Quoting nishap.pa...@utoronto.ca: Okay I am going to give it a try. I just wanted to make sure I was calculating C6 and C12 correctly as well using sigma and epsilon according to rule 3 C12 = Sigma^(6)*C6 C6 = 4*sigma^(6)*epsilon Thanks Quoting Justin A. Lemkul jalem...@vt.edu: nishap.pa...@utoronto.ca wrote: Hello, I have a concern regarding C6 attractive term in LJ potential. I am using OPLS-AA force field, and I wish to turn off attractive term C6 by setting the parameters to zero. One of the suggestion was to add a [nonbond_params] in my ffoplsaanb.itp file and set the C6 to zero between the non-bonded pair. In my system for example, which consists of one methanol in water, I wish to set C6 term to zero between my solute and solvent. Since OPLSAA is all atom force field it treats each atom individually and has sigma and epsilon for each atom, so I am not sure how I would actually set my nonbond_params in my nb.itp file. I realize I need to
RE: [gmx-users] C6 attractive term off OPLSAA
I am going to give that a try. Thanks. Quoting Berk Hess g...@hotmail.com: Hi, Maybe it is not so clear from the topology table in the manual, but when you supply sigma and epsilon as non-bonded parameters, as for OPLS, the nonbond_params section also expects sigma and epsilon. So it seems you could not set only C6 to zero. However there is an, undocumented, trick: if you use a negative sigma, C6 is set to zero. Berk Date: Wed, 1 Dec 2010 10:42:11 -0500 From: nishap.pa...@utoronto.ca To: gmx-users@gromacs.org Subject: Re: [gmx-users] C6 attractive term off OPLSAA I see. I actually want to compare my RDFs with C6 term off. Earlier I tried using force.c code file and turned C6 = 0, but when i compared my RDFs, it didn't look any different so I am not sure if it even worked at all or made any difference to the simulation, but again I was using OPLS-AA and the nb.itp file itself does not indicate value for C6 or C12. If there any other way I could achieve this? Quoting Justin A. Lemkul jalem...@vt.edu: nishap.pa...@utoronto.ca wrote: I ran the simulation using a different force field ffG53a6. I modified the ffG53a6nb.itp file by changing the term C6 to zero for nonbonded parameters, but this time for urea in water. The simulation ran fine without any warning or exploding. I don't understand why it would work with one force-field and not with OPLS-AA. I don't understand why it worked at all, frankly, but more than likely you've been lucky enough to make 53a6 work, but OPLS-AA failed. You're making significant alterations to the nonbonded interactions of the system. Also using mdrun -rerun, I would basically use the nb.itp file with nonbonded parameters as mentioned earlier and everything else is the same i.e. my .mdp parameters and do another run using my previous trajectory file? I tried to look through some of the posts for using -rerun but I don't understand how that would not still blow up the system. I would like to give it a try for sure but I am not quite sure how I would use the command as such. Using mdrun -rerun recalculates energies from an existing trajectory. You would generate a trajectory with a normal force field model, then create a new .tpr file that has your modified potential, such that you can decompose the energy terms cleverly. If your goal is to generate trajectories to study the effects of using a modified potential (i.e., C6 = 0), then you can't use -rerun. -Justin I really appreciate the help. Thanks. Quoting Mark Abraham mark.abra...@anu.edu.au: On 1/12/2010 5:06 AM, Justin A. Lemkul wrote: nishap.pa...@utoronto.ca wrote: I tried using the nonbonded parameters as defined below in my ffoplsaanb.itp file for methanol in water and this is the syntax I used: [nonbond_params] ;ijfuncc6c12 opls_154opls_11110.00E+0002.43E-006 opls_154opls_11210.00E+0000.00E+000 opls_154opls_112 10.00E+0000.00E+000 opls_155opls_111 10.00E+0000.00E+000 opls_155opls_112 10.00E+0000.00E+000 opls_155opls_112 10.00E+0000.00E+000 opls_156opls_111 10.00E+0002.70E-007 opls_156opls_112 10.00E+0000.00E+000 opls_156opls_11210.00E+0000.00E+000 opls_157opls_111 10.00E+0003.01E-006 opls_157opls_112 10.00E+0000.00E+000 When I tried to do mdrun, I got an error saying my system is exploding. I tried doing the mdrun without nonbonded parameters and it runs fine. So I am not sure if I am using the nonbond_params concept correctly. i.e. I want C6 to be zero between my solute (methanol) and solvent (water). This is the error: Warning: 1-4 interaction between 2 and 6 at distance 3.786 which is larger than the 1-4 table size 2.500 nm These are ignored for the rest of the simulation This usually means your system is exploding, if not, you should increase table-extension in your mdp file or with user tables increase the table size Suggestions? Your implementation of your concept is technically correct, but perhaps physically unreasonable. You're turning off the attractive LJ component and leaving only a repulsive component. It sounds about right that everything is flying apart. Indeed. The way to do this is with mdrun -rerun on a trajectory generated with the normal version of the forcefield, as I think I said in another thread yesterday. Mark -Justin Thanks. -Nisha P Quoting nishap.pa...@utoronto.ca: Okay I am going to give it a try. I just wanted to make sure I was calculating C6 and C12 correctly as well using sigma and epsilon according to rule 3 C12 = Sigma^(6)*C6 C6 = 4*sigma^(6)*epsilon Thanks Quoting Justin A. Lemkul jalem...@vt.edu: nishap.pa...@utoronto.ca wrote: Hello, I have a concern regarding C6 attractive term in LJ potential. I am using OPLS-AA force field, and I wish to turn off attractive
[gmx-users] RE:C6 attractive term off OPLSAA
Would I change the sigma to negative in my nb.itp file the way it is set or do I need to include a [nonbond_params] section and make sigma negative in that specific section? because I realize from the manual OPLS-AA uses combination rule 3 but in the manual it shows to combine sigma for i and j using rule 2 and rule 1 and 3 correspond to C6 and C12 for i and j. Is that correct? Quoting nishap.pa...@utoronto.ca: I am going to give that a try. Thanks. Quoting Berk Hess g...@hotmail.com: Hi, Maybe it is not so clear from the topology table in the manual, but when you supply sigma and epsilon as non-bonded parameters, as for OPLS, the nonbond_params section also expects sigma and epsilon. So it seems you could not set only C6 to zero. However there is an, undocumented, trick: if you use a negative sigma, C6 is set to zero. Berk Date: Wed, 1 Dec 2010 10:42:11 -0500 From: nishap.pa...@utoronto.ca To: gmx-users@gromacs.org Subject: Re: [gmx-users] C6 attractive term off OPLSAA I see. I actually want to compare my RDFs with C6 term off. Earlier I tried using force.c code file and turned C6 = 0, but when i compared my RDFs, it didn't look any different so I am not sure if it even worked at all or made any difference to the simulation, but again I was using OPLS-AA and the nb.itp file itself does not indicate value for C6 or C12. If there any other way I could achieve this? Quoting Justin A. Lemkul jalem...@vt.edu: nishap.pa...@utoronto.ca wrote: I ran the simulation using a different force field ffG53a6. I modified the ffG53a6nb.itp file by changing the term C6 to zero for nonbonded parameters, but this time for urea in water. The simulation ran fine without any warning or exploding. I don't understand why it would work with one force-field and not with OPLS-AA. I don't understand why it worked at all, frankly, but more than likely you've been lucky enough to make 53a6 work, but OPLS-AA failed. You're making significant alterations to the nonbonded interactions of the system. Also using mdrun -rerun, I would basically use the nb.itp file with nonbonded parameters as mentioned earlier and everything else is the same i.e. my .mdp parameters and do another run using my previous trajectory file? I tried to look through some of the posts for using -rerun but I don't understand how that would not still blow up the system. I would like to give it a try for sure but I am not quite sure how I would use the command as such. Using mdrun -rerun recalculates energies from an existing trajectory. You would generate a trajectory with a normal force field model, then create a new .tpr file that has your modified potential, such that you can decompose the energy terms cleverly. If your goal is to generate trajectories to study the effects of using a modified potential (i.e., C6 = 0), then you can't use -rerun. -Justin I really appreciate the help. Thanks. Quoting Mark Abraham mark.abra...@anu.edu.au: On 1/12/2010 5:06 AM, Justin A. Lemkul wrote: nishap.pa...@utoronto.ca wrote: I tried using the nonbonded parameters as defined below in my ffoplsaanb.itp file for methanol in water and this is the syntax I used: [nonbond_params] ;ijfuncc6c12 opls_154opls_11110.00E+0002.43E-006 opls_154opls_11210.00E+0000.00E+000 opls_154opls_112 10.00E+0000.00E+000 opls_155opls_111 10.00E+0000.00E+000 opls_155opls_112 10.00E+0000.00E+000 opls_155opls_112 10.00E+0000.00E+000 opls_156opls_111 10.00E+0002.70E-007 opls_156opls_112 10.00E+0000.00E+000 opls_156opls_11210.00E+0000.00E+000 opls_157opls_111 10.00E+0003.01E-006 opls_157opls_112 10.00E+0000.00E+000 When I tried to do mdrun, I got an error saying my system is exploding. I tried doing the mdrun without nonbonded parameters and it runs fine. So I am not sure if I am using the nonbond_params concept correctly. i.e. I want C6 to be zero between my solute (methanol) and solvent (water). This is the error: Warning: 1-4 interaction between 2 and 6 at distance 3.786 which is larger than the 1-4 table size 2.500 nm These are ignored for the rest of the simulation This usually means your system is exploding, if not, you should increase table-extension in your mdp file or with user tables increase the table size Suggestions? Your implementation of your concept is technically correct, but perhaps physically unreasonable. You're turning off the attractive LJ component and leaving only a repulsive component. It sounds about right that everything is flying apart. Indeed. The way to do this is with mdrun -rerun on a trajectory generated with the normal version of the forcefield, as I think I said in another thread yesterday. Mark -Justin Thanks. -Nisha P Quoting nishap.pa...@utoronto.ca: Okay I am going to give it a try. I just wanted to make sure I
Re: [gmx-users] C6 attractive term off OPLSAA
I tried using the nonbonded parameters as defined below in my ffoplsaanb.itp file for methanol in water and this is the syntax I used: [nonbond_params ] ;i j funcc6 c12 opls_154opls_111 1 0.00E+000 2.43E-006 opls_154opls_112 1 0.00E+000 0.00E+000 opls_154opls_112 1 0.00E+000 0.00E+000 opls_155opls_111 1 0.00E+000 0.00E+000 opls_155opls_112 1 0.00E+000 0.00E+000 opls_155opls_112 1 0.00E+000 0.00E+000 opls_156opls_111 1 0.00E+000 2.70E-007 opls_156opls_112 1 0.00E+000 0.00E+000 opls_156opls_112 1 0.00E+000 0.00E+000 opls_157opls_111 1 0.00E+000 3.01E-006 opls_157opls_112 1 0.00E+000 0.00E+000 When I tried to do mdrun, I got an error saying my system is exploding. I tried doing the mdrun without nonbonded parameters and it runs fine. So I am not sure if I am using the nonbond_params concept correctly. i.e. I want C6 to be zero between my solute (methanol) and solvent (water). This is the error: Warning: 1-4 interaction between 2 and 6 at distance 3.786 which is larger than the 1-4 table size 2.500 nm These are ignored for the rest of the simulation This usually means your system is exploding, if not, you should increase table-extension in your mdp file or with user tables increase the table size Suggestions? Thanks. -Nisha P Quoting nishap.pa...@utoronto.ca: Okay I am going to give it a try. I just wanted to make sure I was calculating C6 and C12 correctly as well using sigma and epsilon according to rule 3 C12 = Sigma^(6)*C6 C6 = 4*sigma^(6)*epsilon Thanks Quoting Justin A. Lemkul jalem...@vt.edu: nishap.pa...@utoronto.ca wrote: Hello, I have a concern regarding C6 attractive term in LJ potential. I am using OPLS-AA force field, and I wish to turn off attractive term C6 by setting the parameters to zero. One of the suggestion was to add a [nonbond_params] in my ffoplsaanb.itp file and set the C6 to zero between the non-bonded pair. In my system for example, which consists of one methanol in water, I wish to set C6 term to zero between my solute and solvent. Since OPLSAA is all atom force field it treats each atom individually and has sigma and epsilon for each atom, so I am not sure how I would actually set my nonbond_params in my nb.itp file. I realize I need to convert each sigma and epsilon to C6 and C12, so say for example for methanol in water my [nonbond_params] would look something like this? [ nonbond_params ] ; ij func c6 c12 CTOW 1 0.00 calculated value for C12 here? CTHW1 1 0.00 CTHW2 1 0.00 CT is the carbon in Methanol. OW, HW1 an HW2 correspond to atoms in TIP3P water model Is that correct? Would I have to do that for each atom in methanol? Sounds about right to me. -Justin Any suggestions would be appreciated. Thanks. Nisha P -- Justin A. Lemkul Ph.D. Candidate ICTAS Doctoral Scholar MILES-IGERT Trainee Department of Biochemistry Virginia Tech Blacksburg, VA jalemkul[at]vt.edu | (540) 231-9080 http://www.bevanlab.biochem.vt.edu/Pages/Personal/justin -- 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
[gmx-users] C6 attractive term off OPLSAA
Hello, I have a concern regarding C6 attractive term in LJ potential. I am using OPLS-AA force field, and I wish to turn off attractive term C6 by setting the parameters to zero. One of the suggestion was to add a [nonbond_params] in my ffoplsaanb.itp file and set the C6 to zero between the non-bonded pair. In my system for example, which consists of one methanol in water, I wish to set C6 term to zero between my solute and solvent. Since OPLSAA is all atom force field it treats each atom individually and has sigma and epsilon for each atom, so I am not sure how I would actually set my nonbond_params in my nb.itp file. I realize I need to convert each sigma and epsilon to C6 and C12, so say for example for methanol in water my [nonbond_params] would look something like this? [ nonbond_params ] ; ij func c6 c12 CTOW 1 0.00 calculated value for C12 here? CTHW1 1 0.00 CTHW2 1 0.00 CT is the carbon in Methanol. OW, HW1 an HW2 correspond to atoms in TIP3P water model Is that correct? Would I have to do that for each atom in methanol? Any suggestions would be appreciated. Thanks. Nisha P -- 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] C6 attractive term off OPLSAA
Okay I am going to give it a try. I just wanted to make sure I was calculating C6 and C12 correctly as well using sigma and epsilon according to rule 3 C12 = Sigma^(6)*C6 C6 = 4*sigma^(6)*epsilon Thanks Quoting Justin A. Lemkul jalem...@vt.edu: nishap.pa...@utoronto.ca wrote: Hello, I have a concern regarding C6 attractive term in LJ potential. I am using OPLS-AA force field, and I wish to turn off attractive term C6 by setting the parameters to zero. One of the suggestion was to add a [nonbond_params] in my ffoplsaanb.itp file and set the C6 to zero between the non-bonded pair. In my system for example, which consists of one methanol in water, I wish to set C6 term to zero between my solute and solvent. Since OPLSAA is all atom force field it treats each atom individually and has sigma and epsilon for each atom, so I am not sure how I would actually set my nonbond_params in my nb.itp file. I realize I need to convert each sigma and epsilon to C6 and C12, so say for example for methanol in water my [nonbond_params] would look something like this? [ nonbond_params ] ; ij func c6 c12 CTOW 1 0.00 calculated value for C12 here? CTHW1 1 0.00 CTHW2 1 0.00 CT is the carbon in Methanol. OW, HW1 an HW2 correspond to atoms in TIP3P water model Is that correct? Would I have to do that for each atom in methanol? Sounds about right to me. -Justin Any suggestions would be appreciated. Thanks. Nisha P -- Justin A. Lemkul Ph.D. Candidate ICTAS Doctoral Scholar MILES-IGERT Trainee Department of Biochemistry Virginia Tech Blacksburg, VA jalemkul[at]vt.edu | (540) 231-9080 http://www.bevanlab.biochem.vt.edu/Pages/Personal/justin -- 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] electrostatic interaction
Hello, It doesn't contribute as such. I was using a different force field before OPLSAA and it was all atom so I could make the charges zero in the topology file, I wanted to do the same. Quoting Justin A. Lemkul jalem...@vt.edu: Mark Abraham wrote: On 16/11/2010 7:14 AM, nishap.pa...@utoronto.ca wrote: Hello, I want to turn off electrostatic interactions between CH4 and SOL in my system. I am using ffG53a6 forcefield for CH4 and spc for my water model. CH4 is an united atom and so I can't make the charges zero in the topology. Is there any other way I can turn off electrostatic interactions? Thanks. -Nisha Use an energy group exclusion. See manual. Is that even necessary? Shouldn't a united-atom CH4 have a zero charge, anyway? If so, does it contribute in any way to the Coulombic potential? -Justin -- Justin A. Lemkul Ph.D. Candidate ICTAS Doctoral Scholar MILES-IGERT Trainee Department of Biochemistry Virginia Tech Blacksburg, VA jalemkul[at]vt.edu | (540) 231-9080 http://www.bevanlab.biochem.vt.edu/Pages/Personal/justin -- 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
[gmx-users] electrostatic interaction
Hello, I want to turn off electrostatic interactions between CH4 and SOL in my system. I am using ffG53a6 forcefield for CH4 and spc for my water model. CH4 is an united atom and so I can't make the charges zero in the topology. Is there any other way I can turn off electrostatic interactions? Thanks. -Nisha -- 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] LJ potential
Hello, Is there a way to turn off the 6-term in LJ 6-12 potential? -Nisha P. -- 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] LJ potential
How can I do that using OPLSAA, because in the nb.itp file the values are for sigma/epsilon? Can I change C6=0 in the code file? Quoting David van der Spoel sp...@xray.bmc.uu.se: On 2010-09-27 18.37, nishap.pa...@utoronto.ca wrote: Hello, Is there a way to turn off the 6-term in LJ 6-12 potential? -Nisha P. Set the parameters to zero? -- David van der Spoel, Ph.D., Professor of Biology Dept. of Cell Molec. Biol., Uppsala University. Box 596, 75124 Uppsala, Sweden. Phone: +46184714205. sp...@xray.bmc.uu.sehttp://folding.bmc.uu.se -- 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
[gmx-users] Urea-water
Hello, I am simulating 2M urea in water with single methane (77 urea + 1926 water). I used genbox -ci -nmol to add the number of urea molecules and then used genbox to add water molecules (as suggested in the manual). When I tried to do energy minimization this is what I got : Steepest Descents converged to Fmax 1000 in 234 steps Potential Energy = -1.0054066e+05 So I tried g_energy to plot the potential. But I got average Potential as a positive number 2.84087e+06. I am not sure why the potential is positive, it probably is because of repulsion but is that correct? I would appreciate some help. Thanks. -Nisha P -- 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] Urea-water
I did look at the plot, and it shows that the curve is smoothing right under zero at -1.00e+05 but then why does the average potential show a positive number? Shouldn't that number be negative as well? Quoting David van der Spoel sp...@xray.bmc.uu.se: On 2010-09-20 19.58, nishap.pa...@utoronto.ca wrote: Hello, I am simulating 2M urea in water with single methane (77 urea + 1926 water). I used genbox -ci -nmol to add the number of urea molecules and then used genbox to add water molecules (as suggested in the manual). When I tried to do energy minimization this is what I got : Steepest Descents converged to Fmax 1000 in 234 steps Potential Energy = -1.0054066e+05 So I tried g_energy to plot the potential. But I got average Potential as a positive number 2.84087e+06. I am not sure why the potential is positive, it probably is because of repulsion but is that correct? Did you look at the graph? Please do, then you will understand. I would appreciate some help. Thanks. -Nisha P -- David van der Spoel, Ph.D., Professor of Biology Dept. of Cell Molec. Biol., Uppsala University. Box 596, 75124 Uppsala, Sweden. Phone: +46184714205. sp...@xray.bmc.uu.sehttp://folding.bmc.uu.se -- 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
[gmx-users] NAGA
Hello, I am trying to simulate N-methylglycine amide. I selected 'none' for the termini using oplsaa forcefield. I got the topology but I got the error saying: ERROR 1 [file topol.top, line 126]: No default Ryckaert-Bell. types I checked the line 126 in my topol file and it is the dihedral between the atoms C-C-C-N (5 1 9 7). This is my PDB file: ATOM 1 CA NAGA1 2.010 -3.420 3.690 1.00 20.00 C ATOM 2 HA1 NAGA1 2.621 -3.652 4.447 1.00 20.00 H ATOM 3 HA2 NAGA1 1.081 -3.720 3.908 1.00 20.00 H ATOM 4 HA3 NAGA1 2.321 -3.880 2.858 1.00 20.00 H ATOM 5 CB NAGA1 2.020 -1.900 3.470 1.00 20.00 C ATOM 6 OB NAGA1 2.720 -1.200 4.200 1.00 20.00 O ATOM 7 NE NAGA1 1.270 -1.350 2.500 1.00 20.00 N ATOM 8 HE NAGA1 1.280 -0.370 2.350 1.00 20.00 H ATOM 9 CC NAGA1 0.360 -2.060 1.560 1.00 20.00 C ATOM 10 HC1 NAGA1 -0.410 -1.462 1.336 1.00 20.00 H ATOM 11 HC2 NAGA1 0.020 -2.891 2.000 1.00 20.00 H ATOM 12 CD NAGA1 1.090 -2.440 0.270 1.00 20.00 C ATOM 13 OD NAGA1 2.270 -2.810 0.300 1.00 20.00 O ATOM 14 NB NAGA1 0.370 -2.340 -0.840 1.00 20.00 N ATOM 15 HB1 NAGA1 0.780 -2.570 -1.720 1.00 20.00 H ATOM 16 HB2 NAGA1 -0.580 -2.040 -0.790 1.00 20.00 H Is there a way for me to determine that specific dihedral in .bon.itp file? Nisha P -- gmx-users mailing listgmx-users@gromacs.org http://lists.gromacs.org/mailman/listinfo/gmx-users Please search the archive at http://www.gromacs.org/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/mailing_lists/users.php
[gmx-users] VDWradii.dat
Hello, I have a question about how gromacs assigns VDW radii for ions. I checked the vdwradii.dat file but it does not mention any ions, so where do the parameters for ions come from? I got the .pqr file for sodium and I got the atomic radii of 1.6652A, I am not sure where that value came from. I would really appreciate some help! Thanks. -Nisha P -- gmx-users mailing listgmx-users@gromacs.org http://lists.gromacs.org/mailman/listinfo/gmx-users Please search the archive at http://www.gromacs.org/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/mailing_lists/users.php
Re: [gmx-users] VDWradii.dat
I did use this command: editconf -f md.tpr -mead na.pqr. If I understand it correctly for -sig56 option, I would use the min rvdw from my mdp and divide it by 2 and use that value in the command line? Quoting David van der Spoel sp...@xray.bmc.uu.se: On 2010-08-30 20.32, nishap.pa...@utoronto.ca wrote: Hello, I have a question about how gromacs assigns VDW radii for ions. I checked the vdwradii.dat file but it does not mention any ions, so where do the parameters for ions come from? I got the .pqr file for sodium and I got the atomic radii of 1.6652A, I am not sure where that value came from. I would really appreciate some help! Thanks. -Nisha P Did you get the pqr file from editconf? Then read editconf -h. Check the -sig56 option. -- David van der Spoel, Ph.D., Professor of Biology Dept. of Cell Molec. Biol., Uppsala University. Box 596, 75124 Uppsala, Sweden. Phone: +46184714205. sp...@xray.bmc.uu.sehttp://folding.bmc.uu.se -- gmx-users mailing listgmx-users@gromacs.org http://lists.gromacs.org/mailman/listinfo/gmx-users Please search the archive at http://www.gromacs.org/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/mailing_lists/users.php -- gmx-users mailing listgmx-users@gromacs.org http://lists.gromacs.org/mailman/listinfo/gmx-users Please search the archive at http://www.gromacs.org/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/mailing_lists/users.php
Re: [gmx-users] VDWradii.dat
Also, I am using OPLS-AA force field with TIP3P water Quoting David van der Spoel sp...@xray.bmc.uu.se: On 2010-08-30 20.32, nishap.pa...@utoronto.ca wrote: Hello, I have a question about how gromacs assigns VDW radii for ions. I checked the vdwradii.dat file but it does not mention any ions, so where do the parameters for ions come from? I got the .pqr file for sodium and I got the atomic radii of 1.6652A, I am not sure where that value came from. I would really appreciate some help! Thanks. -Nisha P Did you get the pqr file from editconf? Then read editconf -h. Check the -sig56 option. -- David van der Spoel, Ph.D., Professor of Biology Dept. of Cell Molec. Biol., Uppsala University. Box 596, 75124 Uppsala, Sweden. Phone: +46184714205. sp...@xray.bmc.uu.sehttp://folding.bmc.uu.se -- gmx-users mailing listgmx-users@gromacs.org http://lists.gromacs.org/mailman/listinfo/gmx-users Please search the archive at http://www.gromacs.org/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/mailing_lists/users.php -- gmx-users mailing listgmx-users@gromacs.org http://lists.gromacs.org/mailman/listinfo/gmx-users Please search the archive at http://www.gromacs.org/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/mailing_lists/users.php
Re: [gmx-users] VDWradii.dat
Okay I see what you mean. Thanks alot! Quoting Justin A. Lemkul jalem...@vt.edu: nishap.pa...@utoronto.ca wrote: I did use this command: editconf -f md.tpr -mead na.pqr. If I understand it correctly for -sig56 option, I would use the min rvdw from my mdp and divide it by 2 and use that value in the command line? The value of rvdw is taken from the -rvdw flag, per the description in editconf -h and the source code. If you take the code block (line 608 in gmx_editconf.c): if (bSig56) sig6 = 2*c12/c6; else sig6 = c12/c6; vdw = 0.5*pow(sig6,1.0/6.0); You can easily calculate the radius value you obtained (after converting sigma/epsilon for the atom type to C6/C12). -Justin Quoting David van der Spoel sp...@xray.bmc.uu.se: On 2010-08-30 20.32, nishap.pa...@utoronto.ca wrote: Hello, I have a question about how gromacs assigns VDW radii for ions. I checked the vdwradii.dat file but it does not mention any ions, so where do the parameters for ions come from? I got the .pqr file for sodium and I got the atomic radii of 1.6652A, I am not sure where that value came from. I would really appreciate some help! Thanks. -Nisha P Did you get the pqr file from editconf? Then read editconf -h. Check the -sig56 option. -- David van der Spoel, Ph.D., Professor of Biology Dept. of Cell Molec. Biol., Uppsala University. Box 596, 75124 Uppsala, Sweden. Phone:+46184714205. sp...@xray.bmc.uu.sehttp://folding.bmc.uu.se -- gmx-users mailing listgmx-users@gromacs.org http://lists.gromacs.org/mailman/listinfo/gmx-users Please search the archive at http://www.gromacs.org/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/mailing_lists/users.php -- Justin A. Lemkul Ph.D. Candidate ICTAS Doctoral Scholar MILES-IGERT Trainee Department of Biochemistry Virginia Tech Blacksburg, VA jalemkul[at]vt.edu | (540) 231-9080 http://www.bevanlab.biochem.vt.edu/Pages/Personal/justin -- gmx-users mailing listgmx-users@gromacs.org http://lists.gromacs.org/mailman/listinfo/gmx-users Please search the archive at http://www.gromacs.org/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/mailing_lists/users.php -- gmx-users mailing listgmx-users@gromacs.org http://lists.gromacs.org/mailman/listinfo/gmx-users Please search the archive at http://www.gromacs.org/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/mailing_lists/users.php
Re: [gmx-users] VDW radii of Br and K+
Thanks! I did overlook that. Quoting David van der Spoel sp...@xray.bmc.uu.se: On 2010-08-24 22.53, nishap.pa...@utoronto.ca wrote: Hello, I converted a .tpr file for Bromine ion and Potassium ion using -mead option and editconf. It gives out the Van der Waals radius of the atom, but I don't understand why i got ~2.46A for Potassium and ~2.31A for Bromide. Isn't the bromide ion bigger than Potassium? And I used genion to add both the atoms. This is a force field quirk. From the oplsaa.ff/ffnonbonded.itp [ atomtypes ] ; full atom descriptions are available in ffoplsaa.atp ; name bond_typemasscharge ptype sigma epsilon opls_402 Br- 35 79.90400-1.000 A4.62376e-01 3.76560e-01 opls_408 K+ 19 39.09830 1.000 A4.93463e-01 1.37235e-03 You see that K+ has larger sigma than Br-. Thanks. -Nisha P -- David van der Spoel, Ph.D., Professor of Biology Dept. of Cell Molec. Biol., Uppsala University. Box 596, 75124 Uppsala, Sweden. Phone: +46184714205. sp...@xray.bmc.uu.sehttp://folding.bmc.uu.se -- gmx-users mailing listgmx-users@gromacs.org http://lists.gromacs.org/mailman/listinfo/gmx-users Please search the archive at http://www.gromacs.org/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/mailing_lists/users.php -- gmx-users mailing listgmx-users@gromacs.org http://lists.gromacs.org/mailman/listinfo/gmx-users Please search the archive at http://www.gromacs.org/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/mailing_lists/users.php
[gmx-users] VDW radii of Br and K+
Hello, I converted a .tpr file for Bromine ion and Potassium ion using -mead option and editconf. It gives out the Van der Waals radius of the atom, but I don't understand why i got ~2.46A for Potassium and ~2.31A for Bromide. Isn't the bromide ion bigger than Potassium? And I used genion to add both the atoms. Thanks. -Nisha P -- gmx-users mailing listgmx-users@gromacs.org http://lists.gromacs.org/mailman/listinfo/gmx-users Please search the archive at http://www.gromacs.org/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/mailing_lists/users.php
[gmx-users] NaCl
Hello, I want to simulate a simple NaCl ion in water. I know the method using genion which adds individual Na+ and Cl- ion, but I wish to simulate Na-Cl connected rather than free ions floating in water. When I ran the grompp command I got the error: No default Bond types So I am thinking, I need to add b0 and Kb values for NaCl? Is that correct? I would really appreciate some help! Thanks. -Nisha P -- gmx-users mailing listgmx-users@gromacs.org http://lists.gromacs.org/mailman/listinfo/gmx-users Please search the archive at http://www.gromacs.org/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/mailing_lists/users.php
Re: [gmx-users] NaCl
I believe the bonds are not defined in the OPLS-AA force field for NaCl. I know this might sound like a stupid question, but is there a way to determine the bond force constant (kb), I tried to look it up, but I couldn't find anything. Quoting Justin A. Lemkul jalem...@vt.edu: nishap.pa...@utoronto.ca wrote: Hello, I want to simulate a simple NaCl ion in water. I know the method using genion which adds individual Na+ and Cl- ion, but I wish to simulate Na-Cl connected rather than free ions floating in water. When I ran the grompp command I got the error: No default Bond types So I am thinking, I need to add b0 and Kb values for NaCl? Is that correct? I would really appreciate some help! This really is a case where you can try it and see much faster than it will take you to type the email asking for help :) Presumably, if you want a bond, and there are not parameters for it already built in, then yes, you need to actually define the bond. -Justin Thanks. -Nisha P -- Justin A. Lemkul Ph.D. Candidate ICTAS Doctoral Scholar MILES-IGERT Trainee Department of Biochemistry Virginia Tech Blacksburg, VA jalemkul[at]vt.edu | (540) 231-9080 http://www.bevanlab.biochem.vt.edu/Pages/Personal/justin -- gmx-users mailing listgmx-users@gromacs.org http://lists.gromacs.org/mailman/listinfo/gmx-users Please search the archive at http://www.gromacs.org/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/mailing_lists/users.php -- gmx-users mailing listgmx-users@gromacs.org http://lists.gromacs.org/mailman/listinfo/gmx-users Please search the archive at http://www.gromacs.org/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/mailing_lists/users.php
Re: [gmx-users] Conformational sampling
Thanks. I am definitely going to look into that. Quoting chris.ne...@utoronto.ca: Nisha: Simply applying available tools and seeing if the distribution looks okay is not a good plan. You should have a very well-defined idea of what you are trying to test and then pick a tool to get that done. For example, perhaps you have a reference that provides the converged distribution of a dihedral and you want to see if your run reproduces this value. You would then use some tool, or combination of tools, to calculate the dihedral values and create a histogram of them. Without an outside standard, it is certainly impossible to use a single run to prove that your sampling has converged. It is only possible to prove that it either (a) remains unconverged, or (b) that there is no evidence that it is unconverged. In this case, if you really care about it, the best thing to do is to run a few separate simulations starting from different starting conformations (be sure to start from different dihedral basins). Then plot the time-average of some observables (e.g. the distribution of sampling of a given dihedral) and the values from different runs should converge to the same distribution if you run long enough -- this is the true meaning of convergence, what many people do by way of things like block averaging is just a hack since one usually only has the resources to run a single trajectory. And don't forget that there are more things to converge than just dihedrals. You could also look at pairwise combinations of dihedrals, etc. Note that I'm trying to keep the advice general here because I assume that you're eventually going to tackle something more complex than glycine. Chris. Well I guess, what I am trying to get at is, for my system I want to make sure that 100ns has covered all the conformational changes within the molecule, although I know there is not that much conformational changes for glycine molecule, but I just wanted to confirm. I did run g_angle command and I got the theta values for different groups and the distribution looks okay according to me (I compared the average angle values that I obtained from g_angle to the actual angle values. Then that's probably the best you can do for a molecule as simple as glycine. I would certainly think that rotation about such few bonds would happen within even less time than 100 ns. -Justin Quoting Justin A. Lemkul jalemkul at vt.edu: nishap.patel at utoronto.ca wrote: So is there a way I can test for convergence for my zwitterion for 100ns run? I'm not yet clear what you're assessing or how you define convergence. In addition to what Chris said, you can look at dihedral transitions with g_angle. Surely there are a few dihedrals aside from standard phi/psi, but I don't know what that's going to tell you. -Justin Quoting Justin A. Lemkul jalemkul at vt.edu: nishap.patel at utoronto.ca wrote: Okay so I tried to analyze the torsion using g_chi and g_rama for one glycine zwitterion in water. for g_rama I didn't see anything in xmgrace, and same for g_chi. I used this command for g_chi g_chi -f traj.xtc -s gly.gro -phi -psi When I run the command it says 1 residue with dihedrals found , 2 dihedrals found. But when I open the log file its empty and so are the histo-phi/psiGLY.xvg plots. I tried using g_dih but it says: Found 0 phi-psi combinations For a zwitterion, these torsions don't exist. To measure phi, you need at least C-N-CA-N, and for psi N-CA-C-N. For a single zwitterion, you have only N-CA-C. You need at least a dipeptide. -Justin I am not sure how to check for all torsion convergence for my glycine zwitterion molecule. Am I missing something in the command line? -Nisha P Quoting chris.neale at utoronto.ca: Nisha, The approach is dictated by the goal. What do you want from this? and why are you doing it? e.g. if you want to test the FF, then it is a good idea to at least include US as a part of your strategy; if you want to determine if 100 ns of equilibrium sampling is sufficient to converge all torsions, then obviously you should not be doing US. In any event, would start by running some equilbrium simulations and analyzing the torsions with g_rama, g_chi, g_dih, etc. Be sure to read through the -h output of each program as there are some nuances. Chris. -- original message -- Hello, I would like to do conformational sampling for my simulation of one glycine in its zwitterionic form in water and obtain a PMF curve to see if the system is equilibrated and that all possible torsions are covered for my 100ns run. I am not sure how to approach this issue. Is there a tutorial I can follow? Do I need to do umbrella sampling and use WHAM to extract PMF? I would appreciate some help! Thanks. Nisha P. -- gmx-users mailing listgmx-users@gromacs.org http://lists.gromacs.org/mailman/listinfo/gmx-users Please search the archive at http://www.gromacs.org/search before posting! Please
Re: [gmx-users] Conformational sampling
Okay so I tried to analyze the torsion using g_chi and g_rama for one glycine zwitterion in water. for g_rama I didn't see anything in xmgrace, and same for g_chi. I used this command for g_chi g_chi -f traj.xtc -s gly.gro -phi -psi When I run the command it says 1 residue with dihedrals found , 2 dihedrals found. But when I open the log file its empty and so are the histo-phi/psiGLY.xvg plots. I tried using g_dih but it says: Found 0 phi-psi combinations I am not sure how to check for all torsion convergence for my glycine zwitterion molecule. Am I missing something in the command line? -Nisha P Quoting chris.ne...@utoronto.ca: Nisha, The approach is dictated by the goal. What do you want from this? and why are you doing it? e.g. if you want to test the FF, then it is a good idea to at least include US as a part of your strategy; if you want to determine if 100 ns of equilibrium sampling is sufficient to converge all torsions, then obviously you should not be doing US. In any event, would start by running some equilbrium simulations and analyzing the torsions with g_rama, g_chi, g_dih, etc. Be sure to read through the -h output of each program as there are some nuances. Chris. -- original message -- Hello, I would like to do conformational sampling for my simulation of one glycine in its zwitterionic form in water and obtain a PMF curve to see if the system is equilibrated and that all possible torsions are covered for my 100ns run. I am not sure how to approach this issue. Is there a tutorial I can follow? Do I need to do umbrella sampling and use WHAM to extract PMF? I would appreciate some help! Thanks. Nisha P. -- gmx-users mailing listgmx-users@gromacs.org http://lists.gromacs.org/mailman/listinfo/gmx-users Please search the archive at http://www.gromacs.org/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/mailing_lists/users.php -- gmx-users mailing listgmx-users@gromacs.org http://lists.gromacs.org/mailman/listinfo/gmx-users Please search the archive at http://www.gromacs.org/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/mailing_lists/users.php
Re: [gmx-users] Conformational sampling
So is there a way I can test for convergence for my zwitterion for 100ns run? Quoting Justin A. Lemkul jalem...@vt.edu: nishap.pa...@utoronto.ca wrote: Okay so I tried to analyze the torsion using g_chi and g_rama for one glycine zwitterion in water. for g_rama I didn't see anything in xmgrace, and same for g_chi. I used this command for g_chi g_chi -f traj.xtc -s gly.gro -phi -psi When I run the command it says 1 residue with dihedrals found , 2 dihedrals found. But when I open the log file its empty and so are the histo-phi/psiGLY.xvg plots. I tried using g_dih but it says: Found 0 phi-psi combinations For a zwitterion, these torsions don't exist. To measure phi, you need at least C-N-CA-N, and for psi N-CA-C-N. For a single zwitterion, you have only N-CA-C. You need at least a dipeptide. -Justin I am not sure how to check for all torsion convergence for my glycine zwitterion molecule. Am I missing something in the command line? -Nisha P Quoting chris.ne...@utoronto.ca: Nisha, The approach is dictated by the goal. What do you want from this? and why are you doing it? e.g. if you want to test the FF, then it is a good idea to at least include US as a part of your strategy; if you want to determine if 100 ns of equilibrium sampling is sufficient to converge all torsions, then obviously you should not be doing US. In any event, would start by running some equilbrium simulations and analyzing the torsions with g_rama, g_chi, g_dih, etc. Be sure to read through the -h output of each program as there are some nuances. Chris. -- original message -- Hello, I would like to do conformational sampling for my simulation of one glycine in its zwitterionic form in water and obtain a PMF curve to see if the system is equilibrated and that all possible torsions are covered for my 100ns run. I am not sure how to approach this issue. Is there a tutorial I can follow? Do I need to do umbrella sampling and use WHAM to extract PMF? I would appreciate some help! Thanks. Nisha P. -- gmx-users mailing listgmx-users@gromacs.org http://lists.gromacs.org/mailman/listinfo/gmx-users Please search the archive at http://www.gromacs.org/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/mailing_lists/users.php -- Justin A. Lemkul Ph.D. Candidate ICTAS Doctoral Scholar MILES-IGERT Trainee Department of Biochemistry Virginia Tech Blacksburg, VA jalemkul[at]vt.edu | (540) 231-9080 http://www.bevanlab.biochem.vt.edu/Pages/Personal/justin -- gmx-users mailing listgmx-users@gromacs.org http://lists.gromacs.org/mailman/listinfo/gmx-users Please search the archive at http://www.gromacs.org/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/mailing_lists/users.php -- gmx-users mailing listgmx-users@gromacs.org http://lists.gromacs.org/mailman/listinfo/gmx-users Please search the archive at http://www.gromacs.org/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/mailing_lists/users.php
Re: [gmx-users] Conformational sampling
Well I guess, what I am trying to get at is, for my system I want to make sure that 100ns has covered all the conformational changes within the molecule, although I know there is not that much conformational changes for glycine molecule, but I just wanted to confirm. I did run g_angle command and I got the theta values for different groups and the distribution looks okay according to me (I compared the average angle values that I obtained from g_angle to the actual angle values. Quoting Justin A. Lemkul jalem...@vt.edu: nishap.pa...@utoronto.ca wrote: So is there a way I can test for convergence for my zwitterion for 100ns run? I'm not yet clear what you're assessing or how you define convergence. In addition to what Chris said, you can look at dihedral transitions with g_angle. Surely there are a few dihedrals aside from standard phi/psi, but I don't know what that's going to tell you. -Justin Quoting Justin A. Lemkul jalem...@vt.edu: nishap.pa...@utoronto.ca wrote: Okay so I tried to analyze the torsion using g_chi and g_rama for one glycine zwitterion in water. for g_rama I didn't see anything in xmgrace, and same for g_chi. I used this command for g_chi g_chi -f traj.xtc -s gly.gro -phi -psi When I run the command it says 1 residue with dihedrals found , 2 dihedrals found. But when I open the log file its empty and so are the histo-phi/psiGLY.xvg plots. I tried using g_dih but it says: Found 0 phi-psi combinations For a zwitterion, these torsions don't exist. To measure phi, you need at least C-N-CA-N, and for psi N-CA-C-N. For a single zwitterion, you have only N-CA-C. You need at least a dipeptide. -Justin I am not sure how to check for all torsion convergence for my glycine zwitterion molecule. Am I missing something in the command line? -Nisha P Quoting chris.ne...@utoronto.ca: Nisha, The approach is dictated by the goal. What do you want from this? and why are you doing it? e.g. if you want to test the FF, then it is a good idea to at least include US as a part of your strategy; if you want to determine if 100 ns of equilibrium sampling is sufficient to converge all torsions, then obviously you should not be doing US. In any event, would start by running some equilbrium simulations and analyzing the torsions with g_rama, g_chi, g_dih, etc. Be sure to read through the -h output of each program as there are some nuances. Chris. -- original message -- Hello, I would like to do conformational sampling for my simulation of one glycine in its zwitterionic form in water and obtain a PMF curve to see if the system is equilibrated and that all possible torsions are covered for my 100ns run. I am not sure how to approach this issue. Is there a tutorial I can follow? Do I need to do umbrella sampling and use WHAM to extract PMF? I would appreciate some help! Thanks. Nisha P. -- gmx-users mailing listgmx-users@gromacs.org http://lists.gromacs.org/mailman/listinfo/gmx-users Please search the archive at http://www.gromacs.org/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/mailing_lists/users.php -- Justin A. Lemkul Ph.D. Candidate ICTAS Doctoral Scholar MILES-IGERT Trainee Department of Biochemistry Virginia Tech Blacksburg, VA jalemkul[at]vt.edu | (540) 231-9080 http://www.bevanlab.biochem.vt.edu/Pages/Personal/justin -- gmx-users mailing listgmx-users@gromacs.org http://lists.gromacs.org/mailman/listinfo/gmx-users Please search the archive at http://www.gromacs.org/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/mailing_lists/users.php -- Justin A. Lemkul Ph.D. Candidate ICTAS Doctoral Scholar MILES-IGERT Trainee Department of Biochemistry Virginia Tech Blacksburg, VA jalemkul[at]vt.edu | (540) 231-9080 http://www.bevanlab.biochem.vt.edu/Pages/Personal/justin -- gmx-users mailing listgmx-users@gromacs.org http://lists.gromacs.org/mailman/listinfo/gmx-users Please search the archive at http://www.gromacs.org/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/mailing_lists/users.php -- gmx-users mailing listgmx-users@gromacs.org http://lists.gromacs.org/mailman/listinfo/gmx-users Please search the archive at http://www.gromacs.org/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
[gmx-users] Conformational sampling
Hello, I would like to do conformational sampling for my simulation of one glycine in its zwitterionic form in water and obtain a PMF curve to see if the system is equilibrated and that all possible torsions are covered for my 100ns run. I am not sure how to approach this issue. Is there a tutorial I can follow? Do I need to do umbrella sampling and use WHAM to extract PMF? I would appreciate some help! Thanks. Nisha P. -- gmx-users mailing listgmx-users@gromacs.org http://lists.gromacs.org/mailman/listinfo/gmx-users Please search the archive at http://www.gromacs.org/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/mailing_lists/users.php
[gmx-users] combining replicate xtc files
Hello, I ran a simulation for 100ns on 2 nodes (16 cores). I tried to generate replicates with different starting trajectories and ran it on one core, so basically I have 64 replicates with different starting trajectories and I am using 8 nodes. I wish to combine all the .xtc files from each of the replicates and analyze it as a whole.xtc file. So I used this command : trjcat -f *.xtc -settime -o whole.xtc (and thereafter entering 'C' to continue, to combine all 64 xtc files with different starting time/trajectory) and I got a total of 350ns combining all the trajectories. I tried to get the rdf for 350ns, but the values stay lower than 1. For 100ns that I ran on 2 nodes I got values for rdf that fluctuate around 1. I tried to break-down 350ns trajectory in 100ns and I got values fluctuating around 1, for all set of 100ns (i.e. beginning 100ns, middle 100ns and last 100ns, similar to 100ns ran on 2 nodes). Shouldn't the rdf be constant after 100ns then? Why do my values stay lower than 1 for 350ns? I am not sure if I am combining the files incorrectly. Could anyone provide some insights? Thanks. -Nisha P -- gmx-users mailing listgmx-users@gromacs.org http://lists.gromacs.org/mailman/listinfo/gmx-users Please search the archive at http://www.gromacs.org/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/mailing_lists/users.php
Re: [gmx-users] Domain decomposition error
Thanks Mark, I did play around with the -npme option a little and turns out using -npme 1 works for the system this is mdrun command I have in the script file mdrun_openmpi -nosum -dlb yes -npme 1 -cpt 40 -maxh 48 -deffnm md Thanks. Quoting Mark Abraham mark.abra...@anu.edu.au: - Original Message - From: nishap.pa...@utoronto.ca Date: Monday, May 17, 2010 10:28 Subject: Re: [gmx-users] Domain decomposition error To: gmx-users@gromacs.org Thanks Justin. But how come it worked for methanol. The system is of the same size , and all the parameters are same, so I don't understand why it won't work for ethanol. I'd guess that the greater size of ethanol in combination with constraints is running afoul of the DD minimum cell-size requirements. As you will see in reading the .log file, for your ethanol, P-LINCS requires at least 0.497nm given the initial condition of your system. DD fudges that up by a factor of 1.25 to give some flexibility. Given that minimimum size requirement, only 6 cells can result from a partition in any dimension of a 4nm x 4nm x 4nm box. You've probably artificially required DD to use 14 processors with -npme 2. That requires a 14x1x1 or 7x2x1 DD, neither of which can be consistent with the combination of your box size and constraint usage. Methanol will have a smaller constraint (see its .log file to compare), so the DD will be legal. If you'd given your full mdrun command line in your post, then I wouldn't be guessing as much... The .log file recommends a range of useful solutions for you to consider - but it makes no suggests regarding your .mdp file. Roughly speaking, the .mdp file normally describes the model physics and controls the output conditions, and the (sometimes implicit) arguments of mdrun describe the implemention of the resulting algorithm. You're requiring an impossible implementation. Simplest is to not specify -npme unless you know you need to. For efficiency, both -npme and -np less -npme need to be sufficiently composite and preferably with a high common factor of two of their factors. If you leave mdrun alone, it will guess reasonably. For example, -npme=4 giving 12 DD nodes decomposed 4x3x1 will work admirably in your case, and I bet that's what mdrun picks. Mark Quoting Justin A. Lemkul jalem...@vt.edu: nishap.pa...@utoronto.ca wrote: Hello, I got this following error when I was trying to run a simulation of ethanol-water box size 4*4*4 nm (~6530 atoms). Fatal error: There is no domain decomposition for 14 nodes that is compatible with the given box and a minimum cell size of 0.62175 nm Change the number of nodes or mdrun option -rcon or -dds or your LINCS settings Look in the log file for details on the domain decomposition. I looked into my log file and this is what I got: Initializing Domain Decomposition on 16 nodes Dynamic load balancing: yes Will sort the charge groups at every domain (re)decomposition Initial maximum inter charge-group distances: two-body bonded interactions: 0.234 nm, LJ-14, atoms 1 9 multi-body bonded interactions: 0.234 nm, Angle, atoms 2 5 Minimum cell size due to bonded interactions: 0.257 nm Maximum distance for 5 constraints, at 120 deg. angles, all- trans: 0.497 nm Estimated maximum distance required for P-LINCS: 0.497 nm This distance will limit the DD cell size, you can override this with -rcon Guess for relative PME load: 0.13 Will use 14 particle-particle and 2 PME only nodes This is a guess, check the performance at the end of the log file Using 2 separate PME nodes Scaling the initial minimum size with 1/0.8 (option -dds) = 1.25 Optimizing the DD grid for 14 cells with a minimum initial size of 0.622 nm The maximum allowed number of cells is: X 6 Y 6 Z 6 I am using the same .mdp file that I used for methanol, and it is working fine, I don't understand why it's giving me problem for ethanol. I am running my simulation on 2 nodes. Not according to the log file messages. You're running on 16 nodes, with 14 for PP and 2 for PME. Your system is of insufficient size to be divided over this many PP nodes. Check the list archive for tips, or simply use less nodes for the calculation. There is also information in the manual about all of the settings that are mentionedin the log file. -Justin Suggestions? Thanks -Nisha P -- Justin A. Lemkul Ph.D. Candidate ICTAS Doctoral Scholar MILES-IGERT Trainee Department of Biochemistry Virginia Tech Blacksburg, VA jalemkul[at]vt.edu | (540) 231-9080 http://www.bevanlab.biochem.vt.edu/Pages/Personal/justin -- gmx-users mailing list gmx-users@gromacs.org http://lists.gromacs.org/mailman/listinfo/gmx-users Please search the archive at http://www.gromacs.org/search before posting! Please don't post (un)subscribe requests to the list. Use the www interface or send it to
[gmx-users] Domain decomposition error
Hello, I got this following error when I was trying to run a simulation of ethanol-water box size 4*4*4 nm (~6530 atoms). Fatal error: There is no domain decomposition for 14 nodes that is compatible with the given box and a minimum cell size of 0.62175 nm Change the number of nodes or mdrun option -rcon or -dds or your LINCS settings Look in the log file for details on the domain decomposition. I looked into my log file and this is what I got: Initializing Domain Decomposition on 16 nodes Dynamic load balancing: yes Will sort the charge groups at every domain (re)decomposition Initial maximum inter charge-group distances: two-body bonded interactions: 0.234 nm, LJ-14, atoms 1 9 multi-body bonded interactions: 0.234 nm, Angle, atoms 2 5 Minimum cell size due to bonded interactions: 0.257 nm Maximum distance for 5 constraints, at 120 deg. angles, all-trans: 0.497 nm Estimated maximum distance required for P-LINCS: 0.497 nm This distance will limit the DD cell size, you can override this with -rcon Guess for relative PME load: 0.13 Will use 14 particle-particle and 2 PME only nodes This is a guess, check the performance at the end of the log file Using 2 separate PME nodes Scaling the initial minimum size with 1/0.8 (option -dds) = 1.25 Optimizing the DD grid for 14 cells with a minimum initial size of 0.622 nm The maximum allowed number of cells is: X 6 Y 6 Z 6 I am using the same .mdp file that I used for methanol, and it is working fine, I don't understand why it's giving me problem for ethanol. I am running my simulation on 2 nodes. Suggestions? Thanks -Nisha P -- gmx-users mailing listgmx-users@gromacs.org http://lists.gromacs.org/mailman/listinfo/gmx-users Please search the archive at http://www.gromacs.org/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/mailing_lists/users.php
Re: [gmx-users] Domain decomposition error
Thanks Justin. But how come it worked for methanol. The system is of the same size , and all the parameters are same, so I don't understand why it won't work for ethanol. Quoting Justin A. Lemkul jalem...@vt.edu: nishap.pa...@utoronto.ca wrote: Hello, I got this following error when I was trying to run a simulation of ethanol-water box size 4*4*4 nm (~6530 atoms). Fatal error: There is no domain decomposition for 14 nodes that is compatible with the given box and a minimum cell size of 0.62175 nm Change the number of nodes or mdrun option -rcon or -dds or your LINCS settings Look in the log file for details on the domain decomposition. I looked into my log file and this is what I got: Initializing Domain Decomposition on 16 nodes Dynamic load balancing: yes Will sort the charge groups at every domain (re)decomposition Initial maximum inter charge-group distances: two-body bonded interactions: 0.234 nm, LJ-14, atoms 1 9 multi-body bonded interactions: 0.234 nm, Angle, atoms 2 5 Minimum cell size due to bonded interactions: 0.257 nm Maximum distance for 5 constraints, at 120 deg. angles, all-trans: 0.497 nm Estimated maximum distance required for P-LINCS: 0.497 nm This distance will limit the DD cell size, you can override this with -rcon Guess for relative PME load: 0.13 Will use 14 particle-particle and 2 PME only nodes This is a guess, check the performance at the end of the log file Using 2 separate PME nodes Scaling the initial minimum size with 1/0.8 (option -dds) = 1.25 Optimizing the DD grid for 14 cells with a minimum initial size of 0.622 nm The maximum allowed number of cells is: X 6 Y 6 Z 6 I am using the same .mdp file that I used for methanol, and it is working fine, I don't understand why it's giving me problem for ethanol. I am running my simulation on 2 nodes. Not according to the log file messages. You're running on 16 nodes, with 14 for PP and 2 for PME. Your system is of insufficient size to be divided over this many PP nodes. Check the list archive for tips, or simply use less nodes for the calculation. There is also information in the manual about all of the settings that are mentioned in the log file. -Justin Suggestions? Thanks -Nisha P -- Justin A. Lemkul Ph.D. Candidate ICTAS Doctoral Scholar MILES-IGERT Trainee Department of Biochemistry Virginia Tech Blacksburg, VA jalemkul[at]vt.edu | (540) 231-9080 http://www.bevanlab.biochem.vt.edu/Pages/Personal/justin -- gmx-users mailing listgmx-users@gromacs.org http://lists.gromacs.org/mailman/listinfo/gmx-users Please search the archive at http://www.gromacs.org/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/mailing_lists/users.php -- gmx-users mailing listgmx-users@gromacs.org http://lists.gromacs.org/mailman/listinfo/gmx-users Please search the archive at http://www.gromacs.org/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/mailing_lists/users.php
[gmx-users] trjcat settime
Hello, I have a concern regarding what -settime actually does. I ran replicates of same simulation with different starting trajectories, basically to get more sampling. Now I wish to combine all the .xtc files from each replicates so I can analyze it as a whole. Now what's the difference if I combine the .xtc files using just trjcat -f md.xtc md1.xtcmdN.xtc -o whole.xtc or if I use -settime option. Any insights would be appreciated. Thanks Nisha P -- gmx-users mailing listgmx-users@gromacs.org http://lists.gromacs.org/mailman/listinfo/gmx-users Please search the archive at http://www.gromacs.org/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/mailing_lists/users.php
[gmx-users] Hard Sphere?
Hello, I am not sure if anyone else has some experience with this, but I want to simulate a hard sphere in a solvent (water and heptane). I believe for the hard spheres, the attractive terms for LJ potential is very negligible so the only term to take into account is repulsion, from what I understand. Is there a way I can define hard sphere in gromacs? I would really appreciate some kind of direction. Thanks. -Nisha P -- gmx-users mailing listgmx-users@gromacs.org http://lists.gromacs.org/mailman/listinfo/gmx-users Please search the archive at http://www.gromacs.org/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/mailing_lists/users.php
[gmx-users] Magic error xtc?
Hello, I am running simulations and I got Magic Number error: Fatal error: Magic Number Error in XTC file (read 0, should be 1995) Checking file md.part0002.xtc Reading frame 0 time 48190.004 # Atoms 6533 Precision 0.001 (nm) Reading frame 178000 time 83790.008 Part one of the simulation ran fine, but it exceeded the run time so I ran it again, and got part0002.xtc file. Now it seems like my part0002.xtc file is corrupted, is there a way I can just run my simulation again from my first file without running the whole thing again? Thanks. N.P -- gmx-users mailing listgmx-users@gromacs.org http://lists.gromacs.org/mailman/listinfo/gmx-users Please search the archive at http://www.gromacs.org/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/mailing_lists/users.php
Re: [gmx-users] Magic error xtc?
Hi Justin, Well the way I continued my second part of the simulation was using my md.cpt file. I checked my .cpt file and this is what I got: Last frame -1 time 83786.961 , so I guess my second run updated the file, so does this mean I have to do the whole simulation again? or can I use any of the other files for example .trr file? I thought it might be that I don't have enough disk space? Quoting Justin A. Lemkul jalem...@vt.edu: nishap.pa...@utoronto.ca wrote: Hello, I am running simulations and I got Magic Number error: Fatal error: Magic Number Error in XTC file (read 0, should be 1995) Checking file md.part0002.xtc Reading frame 0 time 48190.004 # Atoms 6533 Precision 0.001 (nm) Reading frame 178000 time 83790.008 Part one of the simulation ran fine, but it exceeded the run time so I ran it again, and got part0002.xtc file. Now it seems like my part0002.xtc file is corrupted, is there a way I can just run my simulation again from my first file without running the whole thing again? That's the nice part about checkpointing. As long as you have the .cpt file from wherever you started the continuation, you can just re-run from that same point. That's really what you did in the first place. If, for some reason, you've over-written or deleted that .cpt file, then you can only start from some later .cpt file, or in the worst case, all over again. -Justin Thanks. N.P -- Justin A. Lemkul Ph.D. Candidate ICTAS Doctoral Scholar MILES-IGERT Trainee Department of Biochemistry Virginia Tech Blacksburg, VA jalemkul[at]vt.edu | (540) 231-9080 http://www.bevanlab.biochem.vt.edu/Pages/Personal/justin -- gmx-users mailing listgmx-users@gromacs.org http://lists.gromacs.org/mailman/listinfo/gmx-users Please search the archive at http://www.gromacs.org/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/mailing_lists/users.php -- gmx-users mailing listgmx-users@gromacs.org http://lists.gromacs.org/mailman/listinfo/gmx-users Please search the archive at http://www.gromacs.org/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/mailing_lists/users.php
Re: [gmx-users] Magic error xtc?
I am trying using prev_md.cpt file. I checked and it looks okay, so hopefully it will work. If not I will try using .trr files. Thanks Quoting Tsjerk Wassenaar tsje...@gmail.com: Hey, Otherwise it's also still possible to do old style continuation with a frame from the .trr/.edr file, provided you have those... Cheers, Tsjerk On Thu, Apr 22, 2010 at 7:55 PM, Justin A. Lemkul jalem...@vt.edu wrote: nishap.pa...@utoronto.ca wrote: Hi Justin, Well the way I continued my second part of the simulation was using my md.cpt file. I checked my .cpt file and this is what I got: Last frame -1 time 83786.961 , so I guess my second run updated the file, so does this mean I have to do the whole simulation again? or can I use any of the other files for example .trr file? I thought it might be that I don't have enough disk space? Checkpoint files will over-write each other, but a backup is kept as (whatever)_prev.cpt; you may be able to use this file to re-start your run. If md.cpt corresponds to a time where there is a viable frame in your trajectory, then there should be no problem re-starting it from this .cpt file. That is, if the frame at (or after) 83790 is corrupt, re-starting from time 83786 should present no problem. Maybe you're out of disk space, but only you know if that's the problem :) -Justin Quoting Justin A. Lemkul jalem...@vt.edu: nishap.pa...@utoronto.ca wrote: Hello, I am running simulations and I got Magic Number error: Fatal error: Magic Number Error in XTC file (read 0, should be 1995) Checking file md.part0002.xtc Reading frame 0 time 48190.004 # Atoms 6533 Precision 0.001 (nm) Reading frame 178000 time 83790.008 Part one of the simulation ran fine, but it exceeded the run time so I ran it again, and got part0002.xtc file. Now it seems like my part0002.xtc file is corrupted, is there a way I can just run my simulation again from my first file without running the whole thing again? That's the nice part about checkpointing. As long as you have the .cpt file from wherever you started the continuation, you can just re-run from that same point. That's really what you did in the first place. If, for some reason, you've over-written or deleted that .cpt file, then you can only start from some later .cpt file, or in the worst case, all over again. -Justin Thanks. N.P -- Justin A. Lemkul Ph.D. Candidate ICTAS Doctoral Scholar MILES-IGERT Trainee Department of Biochemistry Virginia Tech Blacksburg, VA jalemkul[at]vt.edu | (540) 231-9080 http://www.bevanlab.biochem.vt.edu/Pages/Personal/justin -- gmx-users mailing list gmx-us...@gromacs.org http://lists.gromacs.org/mailman/listinfo/gmx-users Please search the archive at http://www.gromacs.org/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/mailing_lists/users.php -- Justin A. Lemkul Ph.D. Candidate ICTAS Doctoral Scholar MILES-IGERT Trainee Department of Biochemistry Virginia Tech Blacksburg, VA jalemkul[at]vt.edu | (540) 231-9080 http://www.bevanlab.biochem.vt.edu/Pages/Personal/justin -- gmx-users mailing list gmx-us...@gromacs.org http://lists.gromacs.org/mailman/listinfo/gmx-users Please search the archive at http://www.gromacs.org/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/mailing_lists/users.php -- Tsjerk A. Wassenaar, Ph.D. post-doctoral researcher Molecular Dynamics Group Groningen Institute for Biomolecular Research and Biotechnology University of Groningen The Netherlands -- gmx-users mailing listgmx-users@gromacs.org http://lists.gromacs.org/mailman/listinfo/gmx-users Please search the archive at http://www.gromacs.org/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/mailing_lists/users.php -- gmx-users mailing listgmx-users@gromacs.org http://lists.gromacs.org/mailman/listinfo/gmx-users Please search the archive at http://www.gromacs.org/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/mailing_lists/users.php
[gmx-users] combining energy files
Hello, I am a running a simulation of 100ns, so I have two edr files md.edr and md.part0002.edr. When I calculate some of the energy values using g_energy for each of the .edr files individually I get the following result: For md.edr: Statistics over 27120401 steps [ 0. thru 54240.8008 ps ], 4 data sets All averages are over 271205 frames Energy Average RMSD Fluct. Drift Tot-Drift --- Temperature 297.9963.759733.75972 6.88223e-07 0.0373298 Pressure (bar) 1.66199262.287262.287 -2.24459e-05 -1.21748 Volume 66.16170.397280.39728 0 0 Density (SI)983.9245.902925.90292 -1.44217e-07 -0.00782243 Heat Capacity Cv: 12.4747 J/mol K (factor = 0.000159181) Isothermal Compressibility: 5.79816e-05 /bar Adiabatic bulk modulus:17246.9 bar For md.part0002.edr Statistics over 22986601 steps [ 54026.8008 thru 10.0078 ps ], 4 data sets All averages are over 229867 frames Energy Average RMSD Fluct. Drift Tot-Drift --- Temperature 297.9933.753983.75397 -5.66217e-07 -0.0260308 Pressure (bar) 1.51511263.692263.692 4.88234e-05 2.24457 Volume 66.1617 0.3984940.39849 1.29794e-07 0.00596703 Density (SI)983.924 5.93295.93284 -1.92684e-06 -0.088583 Heat Capacity Cv: 12.4747 J/mol K (factor = 0.000158697) Isothermal Compressibility: 5.8337e-05 /bar Adiabatic bulk modulus:17141.8 bar But then when I combined my .edr files using eneconv -f md.edr md.part002.edr. the output is saved into 'fixed.edr' and when I ran fixed.edr file through g_energy, I get negative values as follows : Statistics over 5001 steps [ 0. thru 10.0078 ps ], 4 data sets All averages are exact over 5001 steps Energy Average RMSD Fluct. Drift Tot-Drift --- Temperature -5.90251e-06 0 0 0 0 Pressure (bar) 4.46175e-06 0 0 -9.66961e-07 -0.0966961 Volume -1.32585e-06 0 0 0 0 Density (SI) -1.96389e-05 0 0 -2.119e-07 -0.02119 Heat Capacity Cv: 12.4718 J/mol K (factor = 0) Isothermal Compressibility: 0 /bar Adiabatic bulk modulus:inf bar I don't understand what happened there. Could someone please help me understand if I am combining the files correctly, or it is something I am missing. Thanks Nisha Patel -- gmx-users mailing listgmx-users@gromacs.org http://lists.gromacs.org/mailman/listinfo/gmx-users Please search the archive at http://www.gromacs.org/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/mailing_lists/users.php
Re: [gmx-users] combining energy files
I am using gromacs-4.0.5. Thanks. Quoting XAvier Periole x.peri...@rug.nl: That looks like a nice bug! You should specify which gmx version you are using. That would help getting it fixed. On Apr 13, 2010, at 21:34, nishap.pa...@utoronto.ca wrote: Hello, I am a running a simulation of 100ns, so I have two edr files md.edr and md.part0002.edr. When I calculate some of the energy values using g_energy for each of the .edr files individually I get the following result: For md.edr: Statistics over 27120401 steps [ 0. thru 54240.8008 ps ], 4 data sets All averages are over 271205 frames Energy Average RMSD Fluct. Drift Tot-Drift --- Temperature 297.9963.759733.75972 6.88223e-07 0.0373298 Pressure (bar) 1.66199262.287262.287 -2.24459e-05 -1.21748 Volume 66.16170.397280.39728 0 0 Density (SI)983.9245.902925.90292 -1.44217e-07 -0.00782243 Heat Capacity Cv: 12.4747 J/mol K (factor = 0.000159181) Isothermal Compressibility: 5.79816e-05 /bar Adiabatic bulk modulus:17246.9 bar For md.part0002.edr Statistics over 22986601 steps [ 54026.8008 thru 10.0078 ps ], 4 data sets All averages are over 229867 frames Energy Average RMSD Fluct. Drift Tot-Drift --- Temperature 297.9933.753983.75397 -5.66217e-07 -0.0260308 Pressure (bar) 1.51511263.692263.692 4.88234e-05 2.24457 Volume 66.1617 0.3984940.39849 1.29794e-07 0.00596703 Density (SI)983.924 5.93295.93284 -1.92684e-06 -0.088583 Heat Capacity Cv: 12.4747 J/mol K (factor = 0.000158697) Isothermal Compressibility: 5.8337e-05 /bar Adiabatic bulk modulus:17141.8 bar But then when I combined my .edr files using eneconv -f md.edr md.part002.edr. the output is saved into 'fixed.edr' and when I ran fixed.edr file through g_energy, I get negative values as follows : Statistics over 5001 steps [ 0. thru 10.0078 ps ], 4 data sets All averages are exact over 5001 steps Energy Average RMSD Fluct. Drift Tot-Drift --- Temperature -5.90251e-06 0 0 0 0 Pressure (bar) 4.46175e-06 0 0 -9.66961e-07 -0.0966961 Volume -1.32585e-06 0 0 0 0 Density (SI) -1.96389e-05 0 0 -2.119e-07 -0.02119 Heat Capacity Cv: 12.4718 J/mol K (factor = 0) Isothermal Compressibility: 0 /bar Adiabatic bulk modulus:inf bar I don't understand what happened there. Could someone please help me understand if I am combining the files correctly, or it is something I am missing. Thanks Nisha Patel -- gmx-users mailing listgmx-users@gromacs.org http://lists.gromacs.org/mailman/listinfo/gmx-users Please search the archive at http://www.gromacs.org/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/mailing_lists/users.php -- gmx-users mailing listgmx-users@gromacs.org http://lists.gromacs.org/mailman/listinfo/gmx-users Please search the archive at http://www.gromacs.org/search before posting! Please don't post (un)subscribe requests to the list. Use thewww interface or send it to gmx-users-requ...@gromacs.org. Can't post? Read http://www.gromacs.org/mailing_lists/users.php -- gmx-users mailing listgmx-users@gromacs.org http://lists.gromacs.org/mailman/listinfo/gmx-users Please search the archive at http://www.gromacs.org/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/mailing_lists/users.php
Re: [gmx-users] density tip3p
Hi Zuzana, You are right. I did get my density to ~0.986g/cm3, and I actually came across this other paper by Price and Brooks, and initially they get a density of Tip3p using pme at around 0.979g/cm3, but they tweak some of the parameters in oplsaa to get the correct density. I am not sure if you have tried that and if that helps at all, but also they are using a different simulation package. I am getting good rdf values for the water even with slightly low density, so I guess it is okay. Thanks Quoting Zuzana Benkova zuzana.benk...@savba.sk: Hi Nisha, I have performed a benchmark simulations of TIP3P water and I still obtained density 0.984g.cm-3. It looks that this model gives such a density. If you look at JCP 79, 926, 1983 you will see that the authors got 0.982g.cm-3 with MC. (Table 3) Try to evaluate alsoo other properties, each model of water has some deficiencies. I tried large number of cut-offs combinations LJ and QQ tyoes of interactions, but the density did not gange by more than 0.001g.cm-3. Greetings Zuzana Dňa 03/31/10, nishap.pa...@utoronto.ca napísal:Hello, So I looked at the manual and I changed my mdp parameters as follows to equilibrate at a constant pressure. RUN CONTROL PARAMETERS = integrator = sd tinit = 0 dt = 0.002 nsteps = 50 nstcomm = 1 nstxout = 0 nstvout = 0 nstfout = 0 nstlog = 500 nstenergy = 500 nstxtcout = 0 xtc-precision = 1000 nstlist = 5 ns_type = grid pbc = xyz rlist = 1.3 coulombtype = pme rcoulomb = 1.3 epsilon-r = 1 vdw-type = switch rvdw-switch = 1.0 rvdw = 1.2 DispCorr = EnerPres fourierspacing = 0.1 fourier_nx = 0 fourier_ny = 0 fourier_nz = 0 pme_order = 6 ewald_rtol = 1e-06 epsilon_surface = 0 optimize_fft = no tc_grps = system tau_t = 0.1 ref_t = 298 Pcoupl = berendsen tau_p = 1.0 compressibility = 4.5e-05 ref_p = 1.0 constraints = all-bonds constraint-algorithm = Lincs unconstrained-start = no lincs-order = 4 lincs-warnangle = 30 lincs-iter = 1 gen_vel = yes gen_temp = 298.0 gen_seed = 173529 I also tried using mdp parameters as mentioned in some of the gromacs tutorials. I ran this simulation for 1ns and I am still getting density around 985g/l, which is low for tip3p water (27nm^3, 900 water molecules). I checked some of the publications and the values they obtained were ~998g/l. I checked for vacuum in my system, but I didn't see anything out of the ordinary. Any insights? Thanks Nisha P Quoting nishap.pa...@utoronto.ca: Thanks. I will look up the manual again. Quoting Sander Pronk pr...@cbr.su.se: Hi Nisha, Looking at your .mdp, there are some issues that might lead to the behavior that you describe: First: you should try to look up the published densities for tip3p water at 300K - they might actually be close to what you get. Second: your neighbor list cut-off (rlist) might be too small to fully contain the charge groups (check the manual, section 7.3.11). Third: You haven't enabled long-range mean field correction for the pressure or energy. Expect the pressure to be strongly dependent on cut-off (see the same section). Sander On Mar 30, 2010, at 22:47 , nishap.pa...@utoronto.ca wrote: I have a box (3x3x3nm) of Tip3p water molecules ~900 and the density when I create the box using genbox is 997.177g/l. I did energy minimization run and the potential energy did converge smoothly, so I did NPT equilibration run of 100ps and I got the density value of 975g/l. Why does the density decrease after the run? these are the parameters that I used : RUN CONTROL PARAMETERS integrator = md tinit = 0 dt = 0.002 nsteps = 5 nstcomm = 0 nstxout = 0 nstvout = 0 nstfout = 0 nstlog = 100 nstenergy = 100 nstxtcout = 0 xtc_precision = 1000 nstlist = 5 ns-type = Grid pbc = xyz rlist = 1.1 coulombtype = pme rcoulomb = 1.1 epsilon-r = 1 vdw-type = switch rvdw-switch = 0.8 rvdw = 1.0 Tcoupl = V-rescale tc-grps = System
Re: [gmx-users] density tip3p
Hello, So I looked at the manual and I changed my mdp parameters as follows to equilibrate at a constant pressure. RUN CONTROL PARAMETERS = integrator = sd tinit= 0 dt = 0.002 nsteps = 50 nstcomm = 1 nstxout = 0 nstvout = 0 nstfout = 0 nstlog = 500 nstenergy= 500 nstxtcout= 0 xtc-precision= 1000 nstlist = 5 ns_type = grid pbc = xyz rlist= 1.3 coulombtype = pme rcoulomb = 1.3 epsilon-r= 1 vdw-type = switch rvdw-switch = 1.0 rvdw = 1.2 DispCorr = EnerPres fourierspacing = 0.1 fourier_nx = 0 fourier_ny = 0 fourier_nz = 0 pme_order= 6 ewald_rtol = 1e-06 epsilon_surface = 0 optimize_fft = no tc_grps = system tau_t= 0.1 ref_t= 298 Pcoupl = berendsen tau_p= 1.0 compressibility = 4.5e-05 ref_p= 1.0 constraints = all-bonds constraint-algorithm = Lincs unconstrained-start = no lincs-order = 4 lincs-warnangle = 30 lincs-iter = 1 gen_vel = yes gen_temp = 298.0 gen_seed = 173529 I also tried using mdp parameters as mentioned in some of the gromacs tutorials. I ran this simulation for 1ns and I am still getting density around 985g/l, which is low for tip3p water (27nm^3, 900 water molecules). I checked some of the publications and the values they obtained were ~998g/l. I checked for vacuum in my system, but I didn't see anything out of the ordinary. Any insights? Thanks Nisha P Quoting nishap.pa...@utoronto.ca: Thanks. I will look up the manual again. Quoting Sander Pronk pr...@cbr.su.se: Hi Nisha, Looking at your .mdp, there are some issues that might lead to the behavior that you describe: First: you should try to look up the published densities for tip3p water at 300K - they might actually be close to what you get. Second: your neighbor list cut-off (rlist) might be too small to fully contain the charge groups (check the manual, section 7.3.11). Third: You haven't enabled long-range mean field correction for the pressure or energy. Expect the pressure to be strongly dependent on cut-off (see the same section). Sander On Mar 30, 2010, at 22:47 , nishap.pa...@utoronto.ca wrote: I have a box (3x3x3nm) of Tip3p water molecules ~900 and the density when I create the box using genbox is 997.177g/l. I did energy minimization run and the potential energy did converge smoothly, so I did NPT equilibration run of 100ps and I got the density value of 975g/l. Why does the density decrease after the run? these are the parameters that I used : RUN CONTROL PARAMETERS integrator = md tinit= 0 dt = 0.002 nsteps = 5 nstcomm = 0 nstxout = 0 nstvout = 0 nstfout = 0 nstlog = 100 nstenergy= 100 nstxtcout= 0 xtc_precision= 1000 nstlist = 5 ns-type = Grid pbc = xyz rlist= 1.1 coulombtype = pme rcoulomb = 1.1 epsilon-r= 1 vdw-type = switch rvdw-switch = 0.8 rvdw = 1.0 Tcoupl = V-rescale tc-grps = System tau_t= 0.1 ref_t= 300 Pcoupl = berendsen tau_p= 0.5 compressibility = 4e-05 ref_p= 1.0 constraints = all-bonds constraint-algorithm = Lincs unconstrained-start = no lincs-order = 4 lincs-iter = 1 lincs-warnangle = 30 Am I missing something? Thanks Nisha Patel -- gmx-users mailing listgmx-users@gromacs.org http://lists.gromacs.org/mailman/listinfo/gmx-users Please search the archive at http://www.gromacs.org/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/mailing_lists/users.php -- gmx-users mailing listgmx-users@gromacs.org http://lists.gromacs.org/mailman/listinfo/gmx-users Please search the archive at http://www.gromacs.org/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
[gmx-users] Pressure average problem
Hello, I am simulating a box of water (TIP3P ) in a box of 3nm*3*3nm.My simulations are 100ns each. I am using different cut-off. i.e vdw-type = switch rvdw-switch = 0.8 rvdw = 1.0, 1.2 and 1.3 and these are my parameters for Pressure coupling: Pcoupl = Parrinello-Rahman Pcoupltype = isotropic tau_p= 2.0 compressibility = 4.5e-5 ref_p= 1.0 Now after my simulations, when I use g_energy to extract average pressure I get different values for all three cut-offs. For the cut-off of 1.0nm, the average value doesn't converge to 1.0atm but it is rather ~2.1atm, for 1.2cut off the avg. pressure is ~1.7 and for 1.3cut off avg. pressure is ~1.03. I don't understand why the pressure is not converging to 1.0atm for the cut-off of 1.0nm.. All the parameters are exactly the same, the only difference between these simulations are the cut-off values. Could anyone please give me some insight on how the pressure coupling works and if this behaviour seems normal? Thanks. Nisha Patel -- gmx-users mailing listgmx-users@gromacs.org http://lists.gromacs.org/mailman/listinfo/gmx-users Please search the archive at http://www.gromacs.org/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/mailing_lists/users.php
Re: [gmx-users] Pressure average problem
Quoting Mark Abraham mark.abra...@anu.edu.au: On 31/03/2010 5:19 AM, nishap.pa...@utoronto.ca wrote: Hello, I am simulating a box of water (TIP3P ) in a box of 3nm*3*3nm.My simulations are 100ns each. I am using different cut-off. i.e vdw-type = switch rvdw-switch = 0.8 rvdw = 1.0, 1.2 and 1.3 and these are my parameters for Pressure coupling: Pcoupl = Parrinello-Rahman Pcoupltype = isotropic tau_p = 2.0 compressibility = 4.5e-5 ref_p = 1.0 Now after my simulations, when I use g_energy to extract average pressure I get different values for all three cut-offs. For the cut-off of 1.0nm, the average value doesn't converge to 1.0atm but it is rather ~2.1atm, for 1.2cut off the avg. pressure is ~1.7 and for 1.3cut off avg. pressure is ~1.03. I don't understand why the pressure is not converging to 1.0atm for the cut-off of 1.0nm.. All the parameters are exactly the same, the only difference between these simulations are the cut-off values. Could anyone please give me some insight on how the pressure coupling works and if this behaviour seems normal? You've not told us critical information, like the rest of your .mdp file, or how you equilibrated before running your simulation from which you collected data. You could be doing all sorts of weird things and we simply couldn't know. Your electrostatic treatment is far more significant than your use of rvdw. Haphazardly changing simulation parameters is not a desirable thing to do. Find a suitable model of physics in a published calculation, evaluate its applicability to your situation, test its performance suitably, and then apply it. Mark Thanks Mark, I will look into it and hopefully it will make more sense then. -- gmx-users mailing listgmx-users@gromacs.org http://lists.gromacs.org/mailman/listinfo/gmx-users Please search the archive at http://www.gromacs.org/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/mailing_lists/users.php Nisha Patel -- gmx-users mailing listgmx-users@gromacs.org http://lists.gromacs.org/mailman/listinfo/gmx-users Please search the archive at http://www.gromacs.org/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/mailing_lists/users.php
[gmx-users] density tip3p
I have a box (3x3x3nm) of Tip3p water molecules ~900 and the density when I create the box using genbox is 997.177g/l. I did energy minimization run and the potential energy did converge smoothly, so I did NPT equilibration run of 100ps and I got the density value of 975g/l. Why does the density decrease after the run? these are the parameters that I used : RUN CONTROL PARAMETERS integrator = md tinit= 0 dt = 0.002 nsteps = 5 nstcomm = 0 nstxout = 0 nstvout = 0 nstfout = 0 nstlog = 100 nstenergy= 100 nstxtcout= 0 xtc_precision= 1000 nstlist = 5 ns-type = Grid pbc = xyz rlist= 1.1 coulombtype = pme rcoulomb = 1.1 epsilon-r= 1 vdw-type = switch rvdw-switch = 0.8 rvdw = 1.0 Tcoupl = V-rescale tc-grps = System tau_t= 0.1 ref_t= 300 Pcoupl = berendsen tau_p= 0.5 compressibility = 4e-05 ref_p= 1.0 constraints = all-bonds constraint-algorithm = Lincs unconstrained-start = no lincs-order = 4 lincs-iter = 1 lincs-warnangle = 30 Am I missing something? Thanks Nisha Patel -- gmx-users mailing listgmx-users@gromacs.org http://lists.gromacs.org/mailman/listinfo/gmx-users Please search the archive at http://www.gromacs.org/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/mailing_lists/users.php
Re: [gmx-users] density tip3p
Thanks. I will look up the manual again. Quoting Sander Pronk pr...@cbr.su.se: Hi Nisha, Looking at your .mdp, there are some issues that might lead to the behavior that you describe: First: you should try to look up the published densities for tip3p water at 300K - they might actually be close to what you get. Second: your neighbor list cut-off (rlist) might be too small to fully contain the charge groups (check the manual, section 7.3.11). Third: You haven't enabled long-range mean field correction for the pressure or energy. Expect the pressure to be strongly dependent on cut-off (see the same section). Sander On Mar 30, 2010, at 22:47 , nishap.pa...@utoronto.ca wrote: I have a box (3x3x3nm) of Tip3p water molecules ~900 and the density when I create the box using genbox is 997.177g/l. I did energy minimization run and the potential energy did converge smoothly, so I did NPT equilibration run of 100ps and I got the density value of 975g/l. Why does the density decrease after the run? these are the parameters that I used : RUN CONTROL PARAMETERS integrator = md tinit= 0 dt = 0.002 nsteps = 5 nstcomm = 0 nstxout = 0 nstvout = 0 nstfout = 0 nstlog = 100 nstenergy= 100 nstxtcout= 0 xtc_precision= 1000 nstlist = 5 ns-type = Grid pbc = xyz rlist= 1.1 coulombtype = pme rcoulomb = 1.1 epsilon-r= 1 vdw-type = switch rvdw-switch = 0.8 rvdw = 1.0 Tcoupl = V-rescale tc-grps = System tau_t= 0.1 ref_t= 300 Pcoupl = berendsen tau_p= 0.5 compressibility = 4e-05 ref_p= 1.0 constraints = all-bonds constraint-algorithm = Lincs unconstrained-start = no lincs-order = 4 lincs-iter = 1 lincs-warnangle = 30 Am I missing something? Thanks Nisha Patel -- gmx-users mailing listgmx-users@gromacs.org http://lists.gromacs.org/mailman/listinfo/gmx-users Please search the archive at http://www.gromacs.org/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/mailing_lists/users.php -- gmx-users mailing listgmx-users@gromacs.org http://lists.gromacs.org/mailman/listinfo/gmx-users Please search the archive at http://www.gromacs.org/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/mailing_lists/users.php Nisha Patel MSc Candidate Leslie Dan Faculty of Pharmacy Department of Pharmaceutical Sciences 144 College Street Room 1172 Toronto, ON M5S 3M2 Canada Telephone: 416-978-1536 -- gmx-users mailing listgmx-users@gromacs.org http://lists.gromacs.org/mailman/listinfo/gmx-users Please search the archive at http://www.gromacs.org/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/mailing_lists/users.php
Re: [gmx-users] Re: density vs time
What do you mean the value for methane? As in the charges? This is my topology file and I am using OPLS-AA for methane. Include forcefield parameters #include ffoplsaa.itp [ moleculetype ] ; Namenrexcl Methane3 [ atoms ] ; nr type resnr residue atom cgnr charge mass typeBchargeB massB 1 opls_138 1 METH C 1 -0.24 12.011 ; qtot -0.24 2 opls_140 1 METH H1 1 0.06 1.008 ; qtot -0.18 3 opls_140 1 METH H2 1 0.06 1.008 ; qtot -0.12 4 opls_140 1 METH H3 1 0.06 1.008 ; qtot -0.06 5 opls_140 1 METH H4 1 0.06 1.008 ; qtot 0.54 [ bonds ] ; aiaj functc0c1c2c3 1 2 1 1 3 1 1 4 1 1 5 1 [ angles ] ; aiajak functc0c1c2 c3 2 1 3 1 2 1 4 1 2 1 5 1 3 1 4 1 3 1 5 1 4 1 5 1 My RDF does look reasonable and the peaks are okay with the values published earlier, but because of what I am calculating (i.e. volume), it seems to be very sensitive and even 0.01 makes a difference. Thanks for suggestions. -Nisha P Quoting Vitaly V. Chaban vvcha...@gmail.com: I am trying to get an rdf graph actually, and my values are very close to one, but not exactly one, and I was wondering if there is some normalization issue with g_rdf? These are my values for rdf: 0 0 0.002 0 0.004 0 0.006 0 0.008 0 0.01 0 0.012 0 0.014 0 0.016 0 0.018 0 0.02 0 0.022 0 0.024 0 0.026 0 0.028 0 0.03 0 0.032 0 0.034 0 0.036 0 0.038 0 0.04 0 0.042 0 0.044 0 0.046 0 0.048 0 0.05 0 0.052 0 0.054 0 0.056 0 0.058 0 0.06 0 0.062 0 0.064 0 0.066 0 0.068 0 0.07 0 0.072 0 0.074 0 0.076 0 0.078 0 0.08 0 0.082 0 0.084 0 0.086 0 0.088 0 0.09 0 0.092 0 0.094 0 0.096 0 0.098 0 0.1 0 0.102 0 0.104 0 0.106 0 0.108 0 0.11 0 0.112 0 0.114 0 0.116 0 0.118 0 0.12 0 0.122 0 0.124 0 0.126 0 0.128 0 0.13 0 0.132 0 0.134 0 0.136 0 0.138 0 0.14 0 0.142 0 0.144 0 0.146 0 0.148 0 0.15 0 0.152 0 0.154 0 0.156 0 0.158 0 0.16 0 0.162 0 0.164 0 0.166 0 0.168 0 0.17 0 0.172 0 0.174 0 0.176 0 0.178 0 0.18 0 0.182 0 0.184 0 0.186 0 0.188 0 0.19 0 0.192 0 0.194 0 0.196 0 0.198 0 0.2 0 0.202 0 0.204 0 0.206 0 0.208 0 0.21 0 0.212 0 0.214 0 0.216 0 0.218 0 0.22 0 0.222 0 0.224 0 0.226 0 0.228 0 0.23 0 0.232 0 0.234 0 0.236 0 0.238 0 0.24 0 0.242 0 0.244 0 0.246 0 0.248 0 0.25 0 0.252 0 0.254 0 0.256 0 0.258 0 0.26 0 0.262 0 0.264 0 0.266 0 0.268 3.32195e-05 0.27 3.27293e-05 0.272 0.000193501 0.274 0.000317806 0.276 0.000689078 0.278 0.00126578 0.28 0.00197816 0.282 0.00348035 0.284 0.00653762 0.286 0.038 0.288 0.016713 0.29 0.0238596 0.292 0.0345034 0.294 0.0528889 0.296 0.0698775 0.298 0.0933922 0.3 0.118079 0.302
[gmx-users] density vs time
Hello, I am simulating one methane molecule in 899 water molecues in the box size of 3 3 3 nm (27nm^3). I would like to determine density vs time. Is there a way I can do that? I am running my simulation at constant Volume i.e. no pressure coupling. I tried using g_density but it gives me density vs box(nm). Any suggestions would be useful. Thanks. Nisha P -- gmx-users mailing listgmx-users@gromacs.org http://lists.gromacs.org/mailman/listinfo/gmx-users Please search the archive at http://www.gromacs.org/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/mailing_lists/users.php
Re: [gmx-users] density vs time
I am trying to get an rdf graph actually, and my values are very close to one, but not exactly one, and I was wondering if there is some normalization issue with g_rdf? These are my values for rdf: 0 0 0.002 0 0.004 0 0.006 0 0.008 0 0.01 0 0.012 0 0.014 0 0.016 0 0.018 0 0.02 0 0.022 0 0.024 0 0.026 0 0.028 0 0.03 0 0.032 0 0.034 0 0.036 0 0.038 0 0.04 0 0.042 0 0.044 0 0.046 0 0.048 0 0.05 0 0.052 0 0.054 0 0.056 0 0.058 0 0.06 0 0.062 0 0.064 0 0.066 0 0.068 0 0.07 0 0.072 0 0.074 0 0.076 0 0.078 0 0.08 0 0.082 0 0.084 0 0.086 0 0.088 0 0.09 0 0.092 0 0.094 0 0.096 0 0.098 0 0.1 0 0.102 0 0.104 0 0.106 0 0.108 0 0.11 0 0.112 0 0.114 0 0.116 0 0.118 0 0.12 0 0.122 0 0.124 0 0.126 0 0.128 0 0.13 0 0.132 0 0.134 0 0.136 0 0.138 0 0.14 0 0.142 0 0.144 0 0.146 0 0.148 0 0.15 0 0.152 0 0.154 0 0.156 0 0.158 0 0.16 0 0.162 0 0.164 0 0.166 0 0.168 0 0.17 0 0.172 0 0.174 0 0.176 0 0.178 0 0.18 0 0.182 0 0.184 0 0.186 0 0.188 0 0.19 0 0.192 0 0.194 0 0.196 0 0.198 0 0.2 0 0.202 0 0.204 0 0.206 0 0.208 0 0.21 0 0.212 0 0.214 0 0.216 0 0.218 0 0.22 0 0.222 0 0.224 0 0.226 0 0.228 0 0.23 0 0.232 0 0.234 0 0.236 0 0.238 0 0.24 0 0.242 0 0.244 0 0.246 0 0.248 0 0.25 0 0.252 0 0.254 0 0.256 0 0.258 0 0.26 0 0.262 0 0.264 0 0.266 0 0.268 3.32195e-05 0.27 3.27293e-05 0.272 0.000193501 0.274 0.000317806 0.276 0.000689078 0.278 0.00126578 0.28 0.00197816 0.282 0.00348035 0.284 0.00653762 0.286 0.038 0.288 0.016713 0.29 0.0238596 0.292 0.0345034 0.294 0.0528889 0.296 0.0698775 0.298 0.0933922 0.3 0.118079 0.302 0.156417 0.304 0.195336 0.306 0.242149 0.308 0.287833 0.31 0.351093 0.312 0.411434 0.314 0.473243 0.316 0.538749 0.3180.61093 0.32 0.688108 0.322 0.756748 0.324 0.837504 0.326 0.914596 0.328 0.993652 0.33 1.0663 0.3321.13958 0.3341.21522 0.3361.29358 0.3381.35543 0.341.42607 0.3421.49061 0.3441.54496 0.3461.59492 0.3481.64578 0.35 1.6887 0.3521.73861 0.354 1.7676 0.3561.78996 0.3581.82775 0.361.83798 0.3621.84435 0.3641.85567 0.3661.86664 0.3681.87024 0.371.86585 0.3721.85472 0.3741.84964 0.3761.84624 0.3781.82577 0.38 1.8069 0.3821.78844 0.3841.77257 0.3861.75025 0.3881.72797 0.391.69894 0.3921.67908 0.394 1.6537 0.3961.62812 0.3981.60552 0.41.57759 0.4021.55621 0.4041.53355 0.4061.49578 0.4081.48583 0.411.45567 0.4121.43467 0.4141.41543 0.416 1.3816 0.4181.36369 0.421.34185 0.4221.31562 0.4241.28986 0.4261.27127 0.4281.24739 0.431.23036 0.4321.21598 0.4341.20053
[gmx-users] waterbox
Hello, I am trying to use genbox to add water molecules to my system consisting of one molecule of ethane in the box size of 4 4 4nm( 64nm^3). If I understand it correctly, with the system of that size, I should have ~2080 water molecules in my system, but instead genbox adds ~2177 molecules of water, with higher density (1018.63 (g/l)). Is that correct? Thanks. Nisha Patel -- gmx-users mailing listgmx-users@gromacs.org http://lists.gromacs.org/mailman/listinfo/gmx-users Please search the archive at http://www.gromacs.org/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/mailing_lists/users.php
Re: [gmx-users] waterbox
Hi Justin, Well from the previous post: The number density for water should be about 32.5 waters per cubic nanometer. So if my box size is 4 4 4 =64nm^3, then shouldn't number of water molecules be 64*32.5 = 2080? May be I am not understanding it correctly. I am using Tip3p water model for my simulations. Thanks. Nisha Quoting Justin A. Lemkul jalem...@vt.edu: nishap.pa...@utoronto.ca wrote: Hello, I am trying to use genbox to add water molecules to my system consisting of one molecule of ethane in the box size of 4 4 4nm( 64nm^3). If I understand it correctly, with the system of that size, I should have ~2080 water molecules in my system, but instead genbox adds ~2177 molecules of water, with higher density (1018.63 (g/l)). Is that correct? Should be, but you haven't said what model water you're using. Assuming SPC, the box you've specified is 10x the size of spc216.gro, so I would expect that genbox would add roughly 216 x 10 = 2160 molecules, give or take a few. How did you come up with 2080 as the expected value? -Justin Thanks. Nisha Patel -- Justin A. Lemkul Ph.D. Candidate ICTAS Doctoral Scholar MILES-IGERT Trainee Department of Biochemistry Virginia Tech Blacksburg, VA jalemkul[at]vt.edu | (540) 231-9080 http://www.bevanlab.biochem.vt.edu/Pages/Personal/justin -- gmx-users mailing listgmx-users@gromacs.org http://lists.gromacs.org/mailman/listinfo/gmx-users Please search the archive at http://www.gromacs.org/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/mailing_lists/users.php Nisha Patel MSc Candidate Leslie Dan Faculty of Pharmacy Department of Pharmaceutical Sciences 144 College Street Room 1172 Toronto, ON M5S 3M2 Canada Telephone: 416-978-1536 -- gmx-users mailing listgmx-users@gromacs.org http://lists.gromacs.org/mailman/listinfo/gmx-users Please search the archive at http://www.gromacs.org/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/mailing_lists/users.php
Re: [gmx-users] ethanol bond types
Quoting David van der Spoel sp...@xray.bmc.uu.se: On 3/3/10 8:14 PM, nishap.pa...@utoronto.ca wrote: Hello, I am trying to simulate Ethanol in water using OPLSAA. I already tried it with all atoms, but now I want to try it with United atoms for ethanol, i.e. CH3,CH2,OH, and HO in my topology file. I created the topology but I got an error saying 'No bond types', so I checked ffoplsaabon.itp file, and as the error indicated I could not find any bond types. Is there a way for me to determine the bonds? This is my topology file for ethanol: [ moleculetype ] ; Namenrexcl Ethanol 3 [ atoms ] ; nr type resnr residue atom cgnr charge mass typeBchargeB massB 1 opls_068 1EOH CB 1 0 15.035 ; qtot 0 2 opls_081 1EOH CA 2 0.265 14.027 ; qtot 0.265 3 opls_078 1EOH OH 2 -0.715.9994 ; qtot -0.435 4 opls_079 1EOH HO 2 0.435 1.008 ; qtot 0 [ bonds ] ; aiaj functc0c1c2c3 1 2 1 2 3 1 3 4 1 [ pairs ] ; aiaj functc0c1c2c3 1 4 1 [ angles ] ; aiajak functc0c1c2 c3 1 2 3 1 2 3 4 1 [ dihedrals ] ; aiajakal functc0c1 c2c3c4c5 1 2 3 4 3 I checked one of the paper that has been published using the same parameters as follows: Table 1. Potential Parameters and Molecular Geometries of OPLS-UA and SPC/E qa atom or group ?, Å , kJ/mol OPLS R13.905 0.73220.000 R23.905 0.49370.265 -0.700 O 3.070 0.7113 H 0.000 0.0.435 I would really appreciate some suggestions, on how I should tackle the error. Thanks Nisha P The united atom parameters are really leftovers from long past. Jorgensen published his first all-atom alcohol simulations in 1988 IIRC. There is a methanol paper from 1983. Just search literature for Jorgensen, methanol ethanol and it will show up. Then you will to type in the parameters yourself, see chapter 5. -- David. Thanks David. These are the parameters I added. Do you think they seem reasonable? [ bondtypes ] ; ij func b0 kb C3C2 10.15300 259400.0 C2OH 10.14300 267800.0 [ angletypes ] ; ijk func th0 cth C3 C2 OH 1 108.00 100.000 ; C2 OH HO 1 108.50 110.000 [ dihedraltypes ] ; ijkl func coefficients C3 C2 OH HO 3 0.0 3.49000 0.0 3.13000 0.0 0.0 And this is the topology file I created using the parameters and pdb2gmx command: ; Include forcefield parameters #include ffoplsaa.itp [ moleculetype ] ; Namenrexcl Ethanol 3 [ atoms ] ; nr type resnr residue atom cgnr charge mass typeBchargeB massB 1 opls_068 1EOH CB 1 0 15.035 ; qtot 0 2 opls_081 1EOH CA 2 0.265 14.027 ; qtot 0.265 3 opls_078 1EOH OH 2 -0.715.9994 ; qtot -0.435 4 opls_079 1EOH HO 2 0.435 1.008 ; qtot 0 [ bonds ] ; aiaj functc0c1c2c3 1 2 1 2 3 1 3 4 1 [ pairs ] ; aiaj functc0c1c2c3 1 4 1 [ angles ] ; aiajak functc0c1c2 c3 1 2 3 1 2 3 4 1 [ dihedrals ] ; aiajakal functc0c1 c2c3c4c5 1 2 3 4 3 Thanks -Nisha David van der Spoel, PhD, Professor of Biology Dept. of Cell and Molecular Biology, Uppsala University. Husargatan 3, Box 596, 75124 Uppsala, Sweden phone: 46 18 471 4205 fax: 46 18 511 755 sp...@xray.bmc.uu.sesp...@gromacs.org http://xray.bmc.uu.se/~spoel -- gmx-users mailing listgmx-users@gromacs.org http://lists.gromacs.org/mailman/listinfo/gmx-users Please search the archive at http://www.gromacs.org/search before posting! Please don't post (un)subscribe requests to the list. Use the www interface or send it to
[gmx-users] ethanol bond types
Hello, I am trying to simulate Ethanol in water using OPLSAA. I already tried it with all atoms, but now I want to try it with United atoms for ethanol, i.e. CH3,CH2,OH, and HO in my topology file. I created the topology but I got an error saying 'No bond types', so I checked ffoplsaabon.itp file, and as the error indicated I could not find any bond types. Is there a way for me to determine the bonds? This is my topology file for ethanol: [ moleculetype ] ; Namenrexcl Ethanol 3 [ atoms ] ; nr type resnr residue atom cgnr charge mass typeBchargeB massB 1 opls_068 1EOH CB 1 0 15.035 ; qtot 0 2 opls_081 1EOH CA 2 0.265 14.027 ; qtot 0.265 3 opls_078 1EOH OH 2 -0.715.9994 ; qtot -0.435 4 opls_079 1EOH HO 2 0.435 1.008 ; qtot 0 [ bonds ] ; aiaj functc0c1c2c3 1 2 1 2 3 1 3 4 1 [ pairs ] ; aiaj functc0c1c2c3 1 4 1 [ angles ] ; aiajak functc0c1c2 c3 1 2 3 1 2 3 4 1 [ dihedrals ] ; aiajakal functc0c1 c2c3c4c5 1 2 3 4 3 I checked one of the paper that has been published using the same parameters as follows: Table 1. Potential Parameters and Molecular Geometries of OPLS-UA and SPC/E qa atom or group ?, Å , kJ/mol OPLS R13.905 0.73220.000 R23.905 0.49370.265 -0.700 O 3.070 0.7113 H 0.000 0.0.435 I would really appreciate some suggestions, on how I should tackle the error. Thanks Nisha P -- gmx-users mailing listgmx-users@gromacs.org http://lists.gromacs.org/mailman/listinfo/gmx-users Please search the archive at http://www.gromacs.org/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/mailing_lists/users.php
[gmx-users] extending simulation confusion?
Hi, I am trying to extend my simulation time from 100ns - 300ns. So I used the command : tpbconv -s old.tpr -extend 10 -o new.tpr but I got this: Extending remaining runtime of by 20 ps (now 14992 steps) Writing statusfile with starting step 0 and length 14992 steps... time 0.000 and length 30.000 ps So this means, basically the simulation is starting all over again and now it will run for 300ns? Is there a way I can specify what step to run from? Thanks -Nisha P -- gmx-users mailing listgmx-users@gromacs.org http://lists.gromacs.org/mailman/listinfo/gmx-users Please search the archive at http://www.gromacs.org/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/mailing_lists/users.php
Re: [gmx-users] rdf
Hi Will, I have been using pme from the very start and it didn't help with the volume. Could you please elaborate on what you mean by average g(r) over every molecule? And because my rdf-s are not smooth enough even with 100ns, I get fluctuating data using rdf-s. -Nisha Quoting Will Glover will_glo...@yahoo.com: Hi, 100 ns is overkill for those liquids. It should be smooth enough at 1 ns, provided you average g(r) over every molecule. What do you mean the volume decreases? With time? As for the cut-off, it's well known that switching Coulomb terms leads to artifacts, and increasing the cut-off distance doesn't necessarily make things converge. See http://pubs.acs.org/doi/abs/10.1021/ct0502256 for example. Use PME for electrostatics. Regards, -- Will --- On Sat, 1/30/10, nishap.pa...@utoronto.ca nishap.pa...@utoronto.ca wrote: From: nishap.pa...@utoronto.ca nishap.pa...@utoronto.ca Subject: [gmx-users] rdf To: Discussion list for GROMACS users gmx-users@gromacs.org Date: Saturday, January 30, 2010, 9:45 PM Hi, I am doing rdf's of simple molecules. I ran my simulation of water and methane for 100ns to get a smoother curve for the rdf. I am trying to determine the volume, and after a certain cut-off I would assume my values to be constant (i.e. volume), but the values fluctuate alot (i.e. decreasing). If I am using 'switch', from say 0.8-0.9, shouldn't the values smooth off after 0.9nm? and stay constant after that? Any insights would be helpful. Thanks. Nisha P --gmx-users mailing list gmx-users@gromacs.org http://lists.gromacs.org/mailman/listinfo/gmx-users Please search the archive at http://www.gromacs.org/search before posting! Please don't post (un)subscribe requests to the list. Use thewww interface or send it to gmx-users-requ...@gromacs.org. Can't post? Read http://www.gromacs.org/mailing_lists/users.php -- gmx-users mailing listgmx-users@gromacs.org http://lists.gromacs.org/mailman/listinfo/gmx-users Please search the archive at http://www.gromacs.org/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/mailing_lists/users.php Nisha Patel MSc Candidate Leslie Dan Faculty of Pharmacy Department of Pharmaceutical Sciences 144 College Street Room 1172 Toronto, ON M5S 3M2 Canada Telephone: 416-978-1536 -- gmx-users mailing listgmx-users@gromacs.org http://lists.gromacs.org/mailman/listinfo/gmx-users Please search the archive at http://www.gromacs.org/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/mailing_lists/users.php
Re: [gmx-users] rdf
Hi, I saved my rdf plot as a pdf file and I have attached it to this email. Hopefully it will work. My system is one methane with 893 molecules of tip3p water in a cubic box of size 3 3 3 (nm). My topology file is : ; Include forcefield parameters #include ffoplsaa.itp [ moleculetype ] ; Namenrexcl Methane3 [ atoms ] ; nr type resnr residue atom cgnr charge mass typeBchargeB massB 1 opls_138 1 METH C 1 -0.24 12.011 ; qtot -0.24 2 opls_140 1 METH H1 1 0.06 1.008 ; qtot -0.18 3 opls_140 1 METH H2 1 0.06 1.008 ; qtot -0.12 4 opls_140 1 METH H3 1 0.06 1.008 ; qtot -0.06 5 opls_140 1 METH H4 1 0.06 1.008 ; qtot 0.54 [ bonds ] ; aiaj functc0c1c2c3 1 2 1 1 3 1 1 4 1 1 5 1 [ angles ] ; aiajak functc0c1c2 c3 2 1 3 1 2 1 4 1 2 1 5 1 3 1 4 1 3 1 5 1 4 1 5 1 ; Include water topology #include tip3p.itp [ system ] ; Name methane pair in water [ molecules ] ; Compound#mols Methane 1 SOL 893 and my grompp file at constant volume is: ; VARIOUS PREPROCESSING OPTIONS title= Methane in water cpp = /lib/cpp ; RUN CONTROL PARAMETERS integrator = md tinit= 0 dt = 0.002 nsteps = 5000 comm_mode= Linear nstcomm = 1 nstxout = 1 nstvout = 1 nstfout = 0 nstlog = 1000 nstenergy= 1 nstxtcout= 1000 xtc_precision= 1000 xtc-grps = System energygrps = nstlist = 5 ns-type = Grid pbc = xyz rlist= 1.4 coulombtype = pme rcoulomb = 1.4 epsilon_rf = 1 epsilon_r= 1 vdw-type = cut-off rvdw = 1.4 Tcoupl = V-rescale tc-grps = System tau_t= 1.0 ref_t= 298 Pcoupl = No Pcoupltype = Isotropic tau_p= 1.0 compressibility = 4.5e-5 ref_p= 1.0 gen_vel = no gen_temp = 298 gen_seed = 173529 constraints = all-bonds constraint-algorithm = Lincs unconstrained-start = yes lincs-order = 4 lincs-iter = 2 lincs-warnangle = 30 Does this makes sense? And Chris, I don't understand what you mean by defining tail of my rdf. -Nisha Quoting chris.ne...@utoronto.ca: make a figure of your rdf, post it online somewhere (I use photobucket) and reply to the list with a link to your figure. Make sure your figure is well labeled and give us a thorough description of what you see that you don't like. Also tell us what your unit-cell box vectors are. Also describe your system thoroughly, e.g. one methane in a box of 100 water molecules... I know that the meaning of values fluctuate alot (i.e. decreasing). is immediately obvious to you, but for those of us who haven't seen the graph it is pretty tough to figure out what is actually going on. And off the bat, I'd say that the tail of your RDF is just poorly defined (data coming from few waters in the figurative corners of the unit cell for a central methane) and you are seeing the noise that results from this -- Try regenerating your data by block averaging into 0-33 ns, 33-66 ns, and 66-100 ns, then plot the average and standard deviation of your results... I bet this will show you what's going on. Chris. Hi Will, I have been using pme from the very start and it didn't help with the volume. Could you please elaborate on what you mean by average g(r) over every molecule? And because my rdf-s are not smooth enough even with 100ns, I get fluctuating data using rdf-s. I think this is a case where it is important to refer back to previous posts. Is this related to your system that has one molecule of methane solvated in water? If I recall, Chris already told you why the values you see aren't smooth: http://lists.gromacs.org/pipermail/gmx-users/2010-January/048297.html Related to the volume issue, what type of pressure coupling are you applying? Is your treatment of van der Waals appropriate (have you applied the switching function properly)? Is switch appropriate for your chosen force field? Are the cutoff's you are using those
Re: [gmx-users] rdf
okay this should work http://docs.google.com/fileview?id=0B_QNmYARywUiNzE0MDA0NWMtYjFiMy00MTRjLThkNDMtNTk3NDgyNjcxZTgxhl=en -Nisha Quoting chris.ne...@utoronto.ca: I didn't get the attachment. Please upload the image to some website and post the link. --original message -- I saved my rdf plot as a pdf file and I have attached it to this email. Hopefully it will work. My system is one methane with 893 molecules of tip3p water in a cubic box of size 3 3 3 (nm). My topology file is : ... -- gmx-users mailing listgmx-users@gromacs.org http://lists.gromacs.org/mailman/listinfo/gmx-users Please search the archive at http://www.gromacs.org/search before posting! Please don't post (un)subscribe requests to the list. Use thewww interface or send it to gmx-users-requ...@gromacs.org. Can't post? Read http://www.gromacs.org/mailing_lists/users.php Nisha Patel MSc Candidate Leslie Dan Faculty of Pharmacy Department of Pharmaceutical Sciences 144 College Street Room 1172 Toronto, ON M5S 3M2 Canada Telephone: 416-978-1536 -- gmx-users mailing listgmx-users@gromacs.org http://lists.gromacs.org/mailman/listinfo/gmx-users Please search the archive at http://www.gromacs.org/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/mailing_lists/users.php
[gmx-users] rdf
Hi, I am doing rdf's of simple molecules. I ran my simulation of water and methane for 100ns to get a smoother curve for the rdf. I am trying to determine the volume, and after a certain cut-off I would assume my values to be constant (i.e. volume), but the values fluctuate alot (i.e. decreasing). If I am using 'switch', from say 0.8-0.9, shouldn't the values smooth off after 0.9nm? and stay constant after that? Any insights would be helpful. Thanks. Nisha P -- gmx-users mailing listgmx-users@gromacs.org http://lists.gromacs.org/mailman/listinfo/gmx-users Please search the archive at http://www.gromacs.org/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/mailing_lists/users.php
[gmx-users] methane in water
Hello, I am running rdf simulations, and I tried running simple methane using OPLS-AA force-field for 100ns. I was hoping to get a smoother curve for such a small molecule in 100ns. Could someone please look over my mdp parameters and suggest if I am missing something. Thanks. Nisha Patel ; RUN CONTROL PARAMETERS = integrator = sd ; start time and timestep in ps = tinit= 0 dt = 0.002 nsteps = 5000 ; number of steps for center of mass motion removal = nstcomm = 100 ; OUTPUT CONTROL OPTIONS = ; Output frequency for coords (x), velocities (v) and forces (f) = nstxout = 10 nstvout = 10 nstfout = 0 ; Output frequency for energies to log file and energy file = nstlog = 10 nstenergy= 500 ; Output frequency and precision for xtc file = nstxtcout= 1 xtc-precision= 1000 ; This selects the subset of atoms for the xtc file. You can = ; select multiple groups. By default all atoms will be written. = ;xtc_grps = ; NEIGHBORSEARCHING PARAMETERS = ; nblist update frequency = nstlist = 10 ; ns algorithm (simple or grid) = ns_type = grid ; Periodic boundary conditions: xyz or none = pbc = xyz ; nblist cut-off = rlist= 1.0 domain-decomposition = no ; OPTIONS FOR ELECTROSTATICS AND VDW = ; Method for doing electrostatics = coulombtype = pme-switch ;rcoulomb-switch = 0.9 rcoulomb = 0.9 ; Dielectric constant (DC) for cut-off or DC of reaction field = epsilon-r= 1 ; Method for doing Van der Waals = vdw-type = switch ; cut-off lengths= rvdw-switch = 0.8 rvdw = 0.9 ; Apply long range dispersion corrections for Energy and Pressure = DispCorr = EnerPres ; Spacing for the PME/PPPM FFT grid = fourierspacing = 0.1 ; FFT grid size, when a value is 0 fourierspacing will be used = fourier_nx = 0 fourier_ny = 0 fourier_nz = 0 ; EWALD/PME/PPPM parameters = pme_order= 6 ewald_rtol = 1e-06 epsilon_surface = 0 optimize_fft = no ;restraints dihre= yes dihre-fc=1 nstdihreout=1000 disre= Simple disre_fc=1 ;OPTIONS FOR TEMPERATURE COUPLING tc_grps = system tau_t= 2.0 ref_t= 300 ;OPTIONS FOR PRESSURE COUPLING Pcoupl = Isotropic tau_p= 1.0 compressibility = 4.5e-05 ref_p= 1.0 ; Free energy control stuff free_energy = yes init_lambda = 0.0 delta_lambda = 0 sc_alpha =0.5 sc-power =1.0 sc-sigma = 0.3 ; GENERATE VELOCITIES FOR STARTUP RUN = gen_vel = yes gen_temp = 300 gen_seed = 171533 ; OPTIONS FOR BONDS = constraints = all-bonds ; Type of constraint algorithm = constraint-algorithm = Lincs ; Do not constrain the start configuration = unconstrained-start = no ; Relative tolerance of shake = shake-tol= 0.0001 ; Highest order in the expansion of the constraint coupling matrix = lincs-order = 12 ; Lincs will write a warning to the stderr if in one step a bond = ; rotates over more degrees than = lincs-warnangle = 30 ; ; File 'topol.top' was generated ; By user: user (1000) ; On host: Tigrans-Lab ; At date: Fri Jan 22 13:54:29 2010 ; ; This is your topology file ; methane pair in water ; ; Include forcefield parameters #include ffoplsaa.itp [ moleculetype ] ; Namenrexcl Methane3 [ atoms ] ; nr type resnr residue atom cgnr charge mass typeB chargeB massB 1 opls_138 1 METH C 1 -0.24 12.011 ; qtot -0.24 2 opls_140 1 METHHC1 1 0.06 1.008 ; qtot -0.18 3 opls_140 1 METHHC2 1 0.06 1.008 ; qtot -0.12 4 opls_140 1 METHHC3 1 0.06 1.008 ; qtot -0.06 5 opls_140 1 METHHC4 1 0.06 1.008 ; qtot 0 [ bonds ] ; aiaj functc0c1c2c3 1 2 1 1 3 1 1 4 1 1 5 1 [ angles ] ; aiajak functc0c1c2c3 2 1 3 1 2 1 4 1 2 1 5 1 3 1 4 1 3 1 5 1 4 1 5 1 ; Include Position restraint file #ifdef POSRES #include posre.itp #endif ; Include water topology #include tip3p.itp
[gmx-users] cut-off
Hello, I wanted to know, if say my box is 30A (cubic) and if I set my cut-off lengths i.e. rvdw-switch = 0.8 and rvdw =0.9 , does that mean, the simulation will cut-0ff at 9A? Should I be changing that to 1.5nm (to get half the box cut-off)? I am a little confused about that. Thanks. Nisha Patel -- gmx-users mailing listgmx-users@gromacs.org http://lists.gromacs.org/mailman/listinfo/gmx-users Please search the archive at http://www.gromacs.org/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/mailing_lists/users.php
Re: [gmx-users] cut-off
I see what you mean, so say if I am doing RDF, the values beyond 9A would be bogus? Quoting Justin A. Lemkul jalem...@vt.edu: nishap.pa...@utoronto.ca wrote: Hello, I wanted to know, if say my box is 30A (cubic) and if I set my cut-off lengths i.e. rvdw-switch = 0.8 and rvdw =0.9 , does that mean, the simulation will cut-0ff at 9A? Should I be changing that to 1.5nm (to get half the box cut-off)? I am a little confused about that. Setting rvdw = 0.9 means that van der Waals interactions beyond 0.9 nm are zero. It has no effect on, for example, Coulombic interactions. Most force fields have defined cutoffs that should be used in order to be consistent with the original derivation of the parameter set, so it is generally inadvisable to make ad hoc changes to the cutoffs in order to meet some arbitrary criterion. The minimum image convention specifies that your longest cutoff must always be less than half the smallest box vector. If you set up a system with a cutoff equal to exactly one half of a box vector, if that box vector decreases even slightly (i.e., under the influence of pressure coupling), then you will be calculating spurious forces. -Justin Thanks. Nisha Patel -- Justin A. Lemkul Ph.D. Candidate ICTAS Doctoral Scholar MILES-IGERT Trainee Department of Biochemistry Virginia Tech Blacksburg, VA jalemkul[at]vt.edu | (540) 231-9080 http://www.bevanlab.biochem.vt.edu/Pages/Personal/justin -- gmx-users mailing listgmx-users@gromacs.org http://lists.gromacs.org/mailman/listinfo/gmx-users Please search the archive at http://www.gromacs.org/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/mailing_lists/users.php Nisha Patel MSc Candidate Leslie Dan Faculty of Pharmacy Department of Pharmaceutical Sciences 144 College Street Room 1172 Toronto, ON M5S 3M2 Canada Telephone: 416-978-1536 -- gmx-users mailing listgmx-users@gromacs.org http://lists.gromacs.org/mailman/listinfo/gmx-users Please search the archive at http://www.gromacs.org/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/mailing_lists/users.php
[gmx-users] Ammonia in Water
Hello, I created a topology file for 1 ammonia in water. Could someone please have a look at the file and let me know if it looks okay. Thanks Nisha Patel ; Include forcefield parameters #include ffoplsaa.itp [ moleculetype ] ; Namenrexcl Ammonia 3 [ atoms ] ; nr type resnr residue atom cgnr charge mass typeBchargeB massB 1 opls_237 1NH3 N 1 -0.7614.0067 ; qtot -0.76 2 opls_240 1NH3 H1 1 0.38 1.008 ; qtot -0.38 3 opls_240 1NH3 H2 1 0.38 1.008 ; qtot 0 4 opls_240 1NH3 H3 1 0.38 1.008 ; qtot 0.38 [ bonds ] ; aiaj functc0c1c2c3 1 2 1 1 3 1 1 4 1 [ angles ] ; aiajak functc0c1c2 c3 2 1 3 1 2 1 4 1 3 1 4 1 ; Include Position restraint file #ifdef POSRES #include posre.itp #endif ; Include water topology #include spc.itp [ system ] ; Name ammonia in water [ molecules ] ; Compound#mols Ammonia 1 SOL 345 -- gmx-users mailing listgmx-users@gromacs.org http://lists.gromacs.org/mailman/listinfo/gmx-users Please search the archive at http://www.gromacs.org/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/mailing_lists/users.php
[gmx-users] Error during minimization
I have a question regarding the simulations that I am running. I tried creating a pdb file with one methane and water molecules, but the system keeps showing 'bad contact' error when I run the minimization step and I don't understand what I am doing wrong. This is my methane file : TITLE methane in water CRYST1 33.010 33.010 33.010 90.00 90.00 90.00 P 1 1 MODEL1 ATOM 1 CH4 CH4 1 19.520 25.090 16.410 1.00 0.00 So I use this file and run editconf -f methane1.pdb -o box.pdb -bt cubic. This creates a box with volume of around 35 nm^3 Then I add water molecules using : genbox -cp box.pdb -cs spc216 -o methwa.pdb -p topol.top This adds around 1191 water molecules to the existing methane file. So my box is now (33.010 33.010 33.010 with total 1192 molecules). Then I run minimization on this system using minim.mdp as follows: integrator= steep emtol= 1.0 nsteps= 100 nstenergy= 10 nstxtcout= 10 xtc_grps= CH4 SOL energygrps= CH4 SOL nstlist= 5 ns_type= simple rlist= 1.0 coulombtype= cut-off rcoulomb= 1.0 rvdw= 1.0 constraints= none pbc= no grompp -f minim.mdp -c methwa.pdb -p topol.top -v mdrun -v And it runs minimization with bad contact errors and this is the result I get: Steepest Descents converged to machine precision in 32 steps, but did not reach the requested Fmax 100. Potential Energy = -4.9481789e+04 Maximum force = 3.2827865e+06 on atom 500 Norm of force = 6.9951227e+04 Any suggestions? Nisha.P -- gmx-users mailing listgmx-users@gromacs.org http://lists.gromacs.org/mailman/listinfo/gmx-users Please search the archive at http://www.gromacs.org/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/mailing_lists/users.php
[gmx-users] Helpful Information Energy Minimization errors
Hello gmx users, After running energy minimization steps, I got the following error: teepest Descents: Tolerance (Fmax) = 1.0e+00 Number of steps= 100 Step=0, Dmax= 1.0e-02 nm, Epot= -2.84750e+04 Fmax= 1.53401e+04, atom= 2501 Step=1, Dmax= 1.0e-02 nm, Epot= -3.17123e+04 Fmax= 7.71044e+03, atom= 2501 Step=2, Dmax= 1.2e-02 nm, Epot= -3.49418e+04 Fmax= 3.42505e+03, atom= 2822 Step=3, Dmax= 1.4e-02 nm, Epot= -3.82660e+04 Fmax= 1.91814e+03, atom= 2822 Step=4, Dmax= 1.7e-02 nm, Epot= -4.10025e+04 Fmax= 1.19099e+03, atom= 2300 Step=5, Dmax= 2.1e-02 nm, Epot= -4.34223e+04 Fmax= 7.57338e+02, atom= 2300 Step=6, Dmax= 2.5e-02 nm, Epot= -4.56573e+04 Fmax= 9.04904e+02, atom= 470 Step=7, Dmax= 3.0e-02 nm, Epot= -4.68373e+04 Fmax= 1.87844e+03, atom= 470 Step=8, Dmax= 3.6e-02 nm, Epot= -4.73705e+04 Fmax= 8.91642e+02, atom= 470 t = 0.009 ps: Water molecule starting at atom 470 can not be settled. Check for bad contacts and/or reduce the timestep. Back Off! I just backed up step9b.pdb to ./#step9b.pdb.1# Back Off! I just backed up step9c.pdb to ./#step9c.pdb.1# Wrote pdb files with previous and current coordinates Step= 10, Dmax= 2.1e-02 nm, Epot= -4.79050e+04 Fmax= 9.54474e+02, atom= 140 Step= 11, Dmax= 2.6e-02 nm, Epot= -4.83740e+04 Fmax= 1.1e+03, atom= 14 Step= 12, Dmax= 3.1e-02 nm, Epot= -4.87475e+04 Fmax= 1.27989e+03, atom= 14 Step= 13, Dmax= 3.7e-02 nm, Epot= -4.90791e+04 Fmax= 1.53201e+03, atom= 500 Step= 14, Dmax= 4.5e-02 nm, Epot= -4.92883e+04 Fmax= 2.07694e+03, atom= 500 t = 0.015 ps: Water molecule starting at atom 500 can not be settled. Check for bad contacts and/or reduce the timestep. Back Off! I just backed up step15b.pdb to ./#step15b.pdb.1# Back Off! I just backed up step15c.pdb to ./#step15c.pdb.1# Wrote pdb files with previous and current coordinates Step= 15, Dmax= 5.3e-02 nm, Epot= -4.94818e+04 Fmax= 3.28279e+06, atom= 500 Step= 31, Dmax= 2.0e-06 nm, Epot= -4.93789e+04 Fmax= 6.50838e+03, atom= 500 Stepsize too small, or no change in energy. Converged to machine precision, but not to the requested precision Fmax 1 Double precision normally gives you higher accuracy. You might need to increase your constraint accuracy, or turn off constraints alltogether (set constraints = none in mdp file) writing lowest energy coordinates. Back Off! I just backed up confout.gro to ./#confout.gro.3# Steepest Descents converged to machine precision in 32 steps, but did not reach the requested Fmax 1. Potential Energy = -4.9481789e+04 Maximum force = 3.2827865e+06 on atom 500 Norm of force = 6.9951227e+04 so I added define = -DFLEXIBLE to my minim.mdp file and the minimization step ran smoothly. If you are having problems and getting bad contacts, try adding -Dflexible and hopefully it work fine. Thanks Chris for your help. Nisha P -- gmx-users mailing listgmx-users@gromacs.org http://lists.gromacs.org/mailman/listinfo/gmx-users Please search the archive at http://www.gromacs.org/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/mailing_lists/users.php
[gmx-users] pdb file for ammonia
Does anyone know where I can get a pdb file for ammonia? I tried searching but I am not sure if there is a website for the pdb files I can use in gromac. Any suggestions would be helpful. Thanks. Nisha -- gmx-users mailing listgmx-users@gromacs.org http://lists.gromacs.org/mailman/listinfo/gmx-users Please search the archive at http://www.gromacs.org/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/mailing_lists/users.php
Re: [gmx-users] Check for bad contacts and/or reduce the timestep
Hi, I tried energy minimization and I am getting the error while I am running the minimization step. But then I still tried to carry on with normal md-simulation, and it was fine, but I am sure it won't give me the result that it should. Is there a reason for that? I have tried very low energy steps as well, but I am still getting the error. -Nisha Quoting Thomas Schlesier schl...@uni-mainz.de: You can try a normal md-simulation after the energy minimisation, with a very low timestep (around two orders of magnitude lower, for 5ps). For this use no T-couple, or if with, then with a very low temperature. I had a case where i had a box of well ordered benzene molecules and only energyminimization didn't helpt much (In normal md simulation i got many lincs warning form the beginning). After the short NVE-md-simulation the benzene molecules weren't ordered and after a second energy minimization it was possible to do a normal md simulation at 300 K without any lincs warnings. Could be possible that the same protocol works also in your case, else follow Justin's advice. Greetings Thomas Message: 3 Date: Mon, 14 Dec 2009 14:09:47 -0500 From: Justin A. Lemkul jalem...@vt.edu Subject: Re: [gmx-users] Check for bad contacts and/or reduce the timestep To: Gromacs Users' List gmx-users@gromacs.org Message-ID: 4b268d7b.1060...@vt.edu Content-Type: text/plain; charset=ISO-8859-1; format=flowed nishap.pa...@utoronto.ca wrote: snip Steepest Descents converged to machine precision in 34 steps, but did not reach the requested Fmax 10. Potential Energy = -3.7206969e+04 Maximum force = 2.0111207e+04 on atom 1217 Norm of force = 6.9283209e+02 I don't understand constraints I need to turn off, since I haven't mentioned any in my grompp file. This is just generic advice that mdrun provides. The fact is your starting structure contains atomic overlap or clashes that cannot be resolved by energy minimization. Have a look at the trajectory and see if you can get some insight into where things are going wrong (since problematic atom numbers are being printed, and perhaps re-consider how you built your system. -Justin -- Justin A. Lemkul Ph.D. Candidate ICTAS Doctoral Scholar MILES-IGERT Trainee Department of Biochemistry Virginia Tech Blacksburg, VA jalemkul[at]vt.edu | (540) 231-9080 http://www.bevanlab.biochem.vt.edu/Pages/Personal/justin -- gmx-users mailing listgmx-users@gromacs.org http://lists.gromacs.org/mailman/listinfo/gmx-users Please search the archive at http://www.gromacs.org/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/mailing_lists/users.php -- gmx-users mailing listgmx-users@gromacs.org http://lists.gromacs.org/mailman/listinfo/gmx-users Please search the archive at http://www.gromacs.org/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/mailing_lists/users.php
Re: [gmx-users] Check for bad contacts and/or reduce the timestep
Hi Justin, I am trying to plot rdf of one methane molecule in water (spc water, 1095(water molecule in my file). I did create the file using genbox (31.010 31.010 31.010). This is the energy minimization file I am using: title = Methane in water cpp = /lib/cpp include = -I../top integrator = steep emtol = 1.0 nsteps = 100 nstenergy = 10 nstxtcout = 10 xtc_grps= CH4 SOL energygrps = CH4 SOL nstlist = 5 ns_type = simple rlist = 1.0 coulombtype = cut-off rcoulomb= 1.0 rvdw= 1.0 constraints = none pbc = no I tried reducing the steps (to 9 instead of 100) and it works but it doesn't make sense. I am going to try to build the system again. Thanks alot of for help. -Nisha Quoting Justin A. Lemkul jalem...@vt.edu: nishap.pa...@utoronto.ca wrote: Hi, I tried energy minimization and I am getting the error while I am running the minimization step. But then I still tried to carry on with normal md-simulation, and it was fine, but I am sure it won't give me the result that it should. Is there a reason for that? I have tried very low energy steps as well, but I am still getting the error. If minimization still reports the same problems, then what I said before is still true - you have bad contacts somewhere. Have you looked at the trajectory to see where these problematic atoms are and what might be causing the clashes? Sometimes running sufficiently gentle equilibration can relax away some of these bad contacts, but that is rarely the case and often problems will manifest later on in the course of MD. I would seriously consider re-evaluating however you built your system (which you haven't told us) and try to get the energy minimization to actually work before just plowing ahead with fingers crossed. -Justin -- Justin A. Lemkul Ph.D. Candidate ICTAS Doctoral Scholar MILES-IGERT Trainee Department of Biochemistry Virginia Tech Blacksburg, VA jalemkul[at]vt.edu | (540) 231-9080 http://www.bevanlab.biochem.vt.edu/Pages/Personal/justin -- gmx-users mailing listgmx-users@gromacs.org http://lists.gromacs.org/mailman/listinfo/gmx-users Please search the archive at http://www.gromacs.org/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/mailing_lists/users.php -- gmx-users mailing listgmx-users@gromacs.org http://lists.gromacs.org/mailman/listinfo/gmx-users Please search the archive at http://www.gromacs.org/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/mailing_lists/users.php
Re: [gmx-users] Check for bad contacts and/or reduce the timestep
Hi, You are right, I had to delete some water molecules, it is working fine now. Thanks. -Nisha P Quoting Justin A. Lemkul jalem...@vt.edu: nishap.pa...@utoronto.ca wrote: Hi Justin, I am trying to plot rdf of one methane molecule in water (spc water, 1095(water molecule in my file). I did create the file using genbox (31.010 31.010 31.010). This is the energy minimization file I am using: Surely a 31-nm cubic box is excessive for a system of 1095 water molecules. You are aware that Gromacs specifies all units in nm, not Angstrom, right? As for the parameters in the .mdp file, I don't see anything that would necessarily cause the instability you're seeing. A better coulombtype method (in general) is PME, but this may not have a significant effect during the EM procedure. title= Methane in water cpp= /lib/cpp include= -I../top integrator= steep emtol= 1.0 nsteps= 100 nstenergy= 10 nstxtcout= 10 xtc_grps= CH4 SOL energygrps= CH4 SOL nstlist= 5 ns_type= simple rlist= 1.0 coulombtype= cut-off rcoulomb= 1.0 rvdw= 1.0 constraints= none pbc= no I tried reducing the steps (to 9 instead of 100) and it works but it doesn't make sense. I am going to try to build the system again. Reducing the number of steps will actually make the energy minimization worse! If anything, use *more* steps to try to make the structure better. -Justin Thanks alot of for help. -Nisha Quoting Justin A. Lemkul jalem...@vt.edu: nishap.pa...@utoronto.ca wrote: Hi, I tried energy minimization and I am getting the error while I am running the minimization step. But then I still tried to carry on with normal md-simulation, and it was fine, but I am sure it won't give me the result that it should. Is there a reason for that? I have tried very low energy steps as well, but I am still getting the error. If minimization still reports the same problems, then what I said before is still true - you have bad contacts somewhere. Have you looked at the trajectory to see where these problematic atoms are and what might be causing the clashes? Sometimes running sufficiently gentle equilibration can relax away some of these bad contacts, but that is rarely the case and often problems will manifest later on in the course of MD. I would seriously consider re-evaluating however you built your system (which you haven't told us) and try to get the energy minimization to actually work before just plowing ahead with fingers crossed. -Justin -- Justin A. Lemkul Ph.D. Candidate ICTAS Doctoral Scholar MILES-IGERT Trainee Department of Biochemistry Virginia Tech Blacksburg, VA jalemkul[at]vt.edu | (540) 231-9080 http://www.bevanlab.biochem.vt.edu/Pages/Personal/justin -- gmx-users mailing listgmx-users@gromacs.org http://lists.gromacs.org/mailman/listinfo/gmx-users Please search the archive at http://www.gromacs.org/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/mailing_lists/users.php -- Justin A. Lemkul Ph.D. Candidate ICTAS Doctoral Scholar MILES-IGERT Trainee Department of Biochemistry Virginia Tech Blacksburg, VA jalemkul[at]vt.edu | (540) 231-9080 http://www.bevanlab.biochem.vt.edu/Pages/Personal/justin -- gmx-users mailing listgmx-users@gromacs.org http://lists.gromacs.org/mailman/listinfo/gmx-users Please search the archive at http://www.gromacs.org/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/mailing_lists/users.php -- gmx-users mailing listgmx-users@gromacs.org http://lists.gromacs.org/mailman/listinfo/gmx-users Please search the archive at http://www.gromacs.org/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/mailing_lists/users.php
[gmx-users] RDF methanol and water
HI, Okay, so I am trying to plot rdf of one methanol this time in water. Is there a way to select methanol molecule as a whole (i.e. get rdf from the center of the molecule), instead of Me, O and H. Even though I put Me, O and H in same group, I am getting a weird rdf plot (Attached to this message). -Thanks Nisha P # Grace project file # @version 50122 @page size 792, 612 @page scroll 5% @page inout 5% @link page off @map font 8 to Courier, Courier @map font 10 to Courier-Bold, Courier-Bold @map font 11 to Courier-BoldOblique, Courier-BoldOblique @map font 9 to Courier-Oblique, Courier-Oblique @map font 4 to Helvetica, Helvetica @map font 6 to Helvetica-Bold, Helvetica-Bold @map font 7 to Helvetica-BoldOblique, Helvetica-BoldOblique @map font 15 to Helvetica-Narrow, Helvetica-Narrow @map font 16 to Helvetica-Narrow-Bold, Helvetica-Narrow-Bold @map font 17 to Helvetica-Narrow-BoldOblique, Helvetica-Narrow-BoldOblique @map font 18 to Helvetica-Narrow-Oblique, Helvetica-Narrow-Oblique @map font 5 to Helvetica-Oblique, Helvetica-Oblique @map font 20 to NewCenturySchlbk-Bold, NewCenturySchlbk-Bold @map font 21 to NewCenturySchlbk-BoldItalic, NewCenturySchlbk-BoldItalic @map font 22 to NewCenturySchlbk-Italic, NewCenturySchlbk-Italic @map font 23 to NewCenturySchlbk-Roman, NewCenturySchlbk-Roman @map font 24 to Palatino-Bold, Palatino-Bold @map font 25 to Palatino-BoldItalic, Palatino-BoldItalic @map font 26 to Palatino-Italic, Palatino-Italic @map font 27 to Palatino-Roman, Palatino-Roman @map font 12 to Symbol, Symbol @map font 2 to Times-Bold, Times-Bold @map font 3 to Times-BoldItalic, Times-BoldItalic @map font 1 to Times-Italic, Times-Italic @map font 0 to Times-Roman, Times-Roman @map font 33 to ZapfChancery-MediumItalic, ZapfChancery-MediumItalic @map font 13 to ZapfDingbats, ZapfDingbats @map color 0 to (255, 255, 255), white @map color 1 to (0, 0, 0), black @map color 2 to (255, 0, 0), red @map color 3 to (0, 255, 0), green @map color 4 to (0, 0, 255), blue @map color 5 to (255, 255, 0), yellow @map color 6 to (188, 143, 143), brown @map color 7 to (220, 220, 220), grey @map color 8 to (148, 0, 211), violet @map color 9 to (0, 255, 255), cyan @map color 10 to (255, 0, 255), magenta @map color 11 to (255, 165, 0), orange @map color 12 to (114, 33, 188), indigo @map color 13 to (103, 7, 72), maroon @map color 14 to (64, 224, 208), turquoise @map color 15 to (0, 139, 0), green4 @reference date 0 @date wrap off @date wrap year 1950 @default linewidth 1.0 @default linestyle 1 @default color 1 @default pattern 1 @default font 0 @default char size 1.00 @default symbol size 1.00 @default sformat %.8g @background color 0 @page background fill on @timestamp off @timestamp 0.03, 0.03 @timestamp color 1 @timestamp rot 0 @timestamp font 0 @timestamp char size 1.00 @timestamp def Wed Dec 16 15:44:45 2009 @r0 off @link r0 to g0 @r0 type above @r0 linestyle 1 @r0 linewidth 1.0 @r0 color 1 @r0 line 0, 0, 0, 0 @r1 off @link r1 to g0 @r1 type above @r1 linestyle 1 @r1 linewidth 1.0 @r1 color 1 @r1 line 0, 0, 0, 0 @r2 off @link r2 to g0 @r2 type above @r2 linestyle 1 @r2 linewidth 1.0 @r2 color 1 @r2 line 0, 0, 0, 0 @r3 off @link r3 to g0 @r3 type above @r3 linestyle 1 @r3 linewidth 1.0 @r3 color 1 @r3 line 0, 0, 0, 0 @r4 off @link r4 to g0 @r4 type above @r4 linestyle 1 @r4 linewidth 1.0 @r4 color 1 @r4 line 0, 0, 0, 0 @g0 on @g0 hidden false @g0 type XY @g0 stacked false @g0 bar hgap 0.00 @g0 fixedpoint off @g0 fixedpoint type 0 @g0 fixedpoint xy 0.00, 0.00 @g0 fixedpoint format general general @g0 fixedpoint prec 6, 6 @with g0 @world 0, 0, 2, 2 @stack world 0, 0, 0, 0 @znorm 1 @view 0.15, 0.15, 1.15, 0.85 @title Radial distribution @title font 0 @title size 1.50 @title color 1 @subtitle MeO - SOL @subtitle font 0 @subtitle size 1.00 @subtitle color 1 @xaxes scale Normal @yaxes scale Normal @xaxes invert off @yaxes invert off @xaxis on @xaxis type zero false @xaxis offset 0.00 , 0.00 @xaxis bar on @xaxis bar color 1 @xaxis bar linestyle 1 @xaxis bar linewidth 1.0 @xaxis label r @xaxis label layout para @xaxis label place auto @xaxis label char size 1.00 @xaxis label font 0 @xaxis label color 1 @xaxis label place normal @xaxis tick on @xaxis tick major 0.5 @xaxis tick minor ticks 1 @xaxis tick default 6 @xaxis tick place rounded true @xaxis tick in @xaxis tick major size 1.00 @xaxis tick major color 1 @xaxis tick major linewidth 1.0 @xaxis tick major linestyle 1 @xaxis tick major grid off @xaxis tick minor color 1 @xaxis tick minor linewidth 1.0 @xaxis tick minor linestyle 1 @xaxis tick minor grid off @xaxis tick minor size 0.50 @xaxis ticklabel on @xaxis ticklabel format general @xaxis ticklabel prec 5 @xaxis
[gmx-users] Check for bad contacts and/or reduce the timestep
Hi, I am trying to run a simulation with 1360 Tip4p water molecules and one methane molecule (Box size: 35.500 35.3371 35.500 and Total 5541 atoms) but I am getting the error: t = 0.000 ps: Water molecule starting at atom 5433 can not be settled. Check for bad contacts and/or reduce the timestep. Here is my grompp. file : ; RUN CONTROL PARAMETERS = integrator = md ; start time and timestep in ps = tinit= 0.0 dt = 0.002 nsteps = 25000 ; number of steps for center of mass motion removal = nstcomm = 1 ENERGY MINIMIZATION OPTIONS emtol= 1 emstep = 0.01 ; OUTPUT CONTROL OPTIONS = ; Output frequency for coords (x), velocities (v) and forces (f) = nstxout = 1000 nstvout = 1000 nstfout = 0 ; Output frequency for energies to log file and energy file = nstlog = 0 nstenergy= 100 ; Output frequency and precision for xtc file = nstxtcout= 100 xtc_precision= 1000 ; This selects the subset of atoms for the xtc file. = ; Only the first group gets written out, it does not make sense = ; to have multiple groups. By default all atoms will be written = xtc_grps = ; Selection of energy groups = energygrps = CH4 SOL ; NEIGHBORSEARCHING PARAMETERS = ; nblist update frequency = nstlist = 1.0 ; ns algorithm (simple or grid) = ns_type = grid ; nblist cut-off = rlist= 1.0 ; OPTIONS FOR ELECTROSTATICS AND VDW = ; Method for doing electrostatics = coulomb_type = pme rcoulomb_switch = 0 rcoulomb = 1.0 ; Dielectric constant (DC) for cut-off or DC of reaction field = epsilon_r= 1.0 epsilon_rf = 0 ; Method for doing Van der Waals = vdw_type = Cut-off ; cut-off lengths= rvdw_switch = 0 rvdw = 1.0 dispcorr = enerpres ; OPTIONS FOR WEAK COUPLING ALGORITHMS = ; Temperature coupling = tcoupl = v-rescale ; Groups to couple separately = tc-grps = system ; Time constant (ps) and reference temperature (K) = tau_t= 0.1 ref_t= 300 ; Pressure coupling = Pcoupl = Berendsen ; Time constant (ps), compressibility (1/bar) and reference P (bar) = tau_p= 1.0 compressibility = 4.6e-5 ref_p= 1 ; GENERATE VELOCITIES FOR STARTUP RUN = gen_vel = yes gen_temp = 300 gen_seed = 171533 I am new to this, so any inputs would be appreciated. Thanks Nisha -- gmx-users mailing listgmx-users@gromacs.org http://lists.gromacs.org/mailman/listinfo/gmx-users Please search the archive at http://www.gromacs.org/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/mailing_lists/users.php
Re: [gmx-users] Check for bad contacts and/or reduce the timestep
- Forwarded message from nishap.pa...@utoronto.ca - Date: Mon, 14 Dec 2009 14:06:30 -0500 From: nishap.pa...@utoronto.ca Reply-To: nishap.pa...@utoronto.ca Subject: Re: [gmx-users] Check for bad contacts and/or reduce the timestep To: Justin A. Lemkul jalem...@vt.edu I tried to run minimization step and I got the following error: Steepest Descents: Tolerance (Fmax) = 1.0e+01 Number of steps=25000 Step=1, Dmax= 1.0e-02 nm, Epot= 1.57470e+09 Fmax= 2.12566e+11, atom= 4973 Step=2, Dmax= 1.2e-02 nm, Epot= 2.42785e+08 Fmax= 2.02821e+10, atom= 4985 Step=3, Dmax= 1.4e-02 nm, Epot= 6.83525e+07 Fmax= 2.24629e+09, atom= 4973 Step=4, Dmax= 1.7e-02 nm, Epot= 1.36192e+07 Fmax= 2.34481e+08, atom= 4985 Step=5, Dmax= 2.1e-02 nm, Epot= 3.37948e+06 Fmax= 2.50839e+07, atom= 4973 Step=6, Dmax= 2.5e-02 nm, Epot= 9.58580e+05 Fmax= 2.79175e+06, atom= 1321 Step=7, Dmax= 3.0e-02 nm, Epot= 2.98861e+05 Fmax= 3.83032e+06, atom= 4973 Step=8, Dmax= 3.6e-02 nm, Epot= 2.15478e+05 Fmax= 5.79850e+05, atom= 4973 Step=9, Dmax= 4.3e-02 nm, Epot= 7.16581e+04 Fmax= 1.12895e+06, atom= 4973 Step= 10, Dmax= 5.2e-02 nm, Epot= 4.39860e+04 Fmax= 7.57820e+04, atom= 101 t = 0.022 ps: Water molecule starting at atom 101 can not be settled. Check for bad contacts and/or reduce the timestep. Back Off! I just backed up step11b.pdb to ./#step11b.pdb.6# Back Off! I just backed up step11c.pdb to ./#step11c.pdb.6# Wrote pdb files with previous and current coordinates Step= 11, Dmax= 6.2e-02 nm, Epot= -1.94430e+04 Fmax= 7.60249e+04, atom= 102 Step= 12, Dmax= 7.4e-02 nm, Epot= -2.29926e+04 Fmax= 3.21409e+05, atom= 4917 t = 0.026 ps: Water molecule starting at atom 4917 can not be settled. Check for bad contacts and/or reduce the timestep. Back Off! I just backed up step13b.pdb to ./#step13b.pdb.6# Back Off! I just backed up step13c.pdb to ./#step13c.pdb.6# Wrote pdb files with previous and current coordinates Step= 14, Dmax= 4.5e-02 nm, Epot= -2.80812e+04 Fmax= 1.42056e+04, atom= 9577 Step= 15, Dmax= 5.3e-02 nm, Epot= -3.28009e+04 Fmax= 1.37135e+05, atom= 4893 t = 0.032 ps: Water molecule starting at atom 4893 can not be settled. Check for bad contacts and/or reduce the timestep. Back Off! I just backed up step16b.pdb to ./#step16b.pdb.6# Back Off! I just backed up step16c.pdb to ./#step16c.pdb.6# Wrote pdb files with previous and current coordinates Step= 16, Dmax= 6.4e-02 nm, Epot= -3.72070e+04 Fmax= 2.01112e+04, atom= 1217 t = 0.034 ps: Water molecule starting at atom 4917 can not be settled. Check for bad contacts and/or reduce the timestep. Back Off! I just backed up step17b.pdb to ./#step17b.pdb.2# Back Off! I just backed up step17c.pdb to ./#step17c.pdb.2# Wrote pdb files with previous and current coordinates Step= 33, Dmax= 1.2e-06 nm, Epot= -1.10597e+04 Fmax= 1.29469e+06, atom= 4893 Stepsize too small, or no change in energy. Converged to machine precision, but not to the requested precision Fmax 10 Double precision normally gives you higher accuracy. You might need to increase your constraint accuracy, or turn off constraints alltogether (set constraints = none in mdp file) writing lowest energy coordinates. Steepest Descents converged to machine precision in 34 steps, but did not reach the requested Fmax 10. Potential Energy = -3.7206969e+04 Maximum force = 2.0111207e+04 on atom 1217 Norm of force = 6.9283209e+02 I don't understand constraints I need to turn off, since I haven't mentioned any in my grompp file. -Nisha nishap.pa...@utoronto.ca wrote: Hi, I am trying to run a simulation with 1360 Tip4p water molecules and one methane molecule (Box size: 35.500 35.3371 35.500 and Total 5541 atoms) but I am getting the error: t = 0.000 ps: Water molecule starting at atom 5433 can not be settled. Check for bad contacts and/or reduce the timestep. Have you done sufficient energy minimization? General information can be found within the wiki pages of the Gromacs site, i.e.: http://www.gromacs.org/Documentation/Errors#LINCS.2fSETTLE.2fSHAKE_warnings -Justin Here is my grompp. file : ; RUN CONTROL PARAMETERS = integrator = md ; start time and timestep in ps = tinit= 0.0 dt = 0.002 nsteps = 25000 ; number of steps for center of mass motion removal = nstcomm = 1 ENERGY MINIMIZATION OPTIONS emtol= 1 emstep = 0.01 ; OUTPUT CONTROL OPTIONS = ; Output frequency for coords (x), velocities (v) and forces (f) = nstxout = 1000 nstvout = 1000 nstfout = 0 ; Output frequency for energies to log file and energy file = nstlog = 0 nstenergy= 100 ; Output frequency and precision for xtc file = nstxtcout= 100 xtc_precision= 1000 ; This selects the subset of atoms for