[gmx-users] FW: trjconv -dump problem
Dear Gromacs users, I minimized a protein structure 1IARcompleted_WT.pdb, getting, among others, the files 1IARcompleted_WT_minimized.trr and 1IARcompleted_WT_minimized_potential_energy.xvg . Looking at the latter file showed the last frame to be at 217 ps: . . @title Gromacs Energies @xaxis label Time (ps) @yaxis label (kJ/mol) . . @ s0 legend Potential 0.00 -31997.519531 0.00 -33810.406250 200.00 -69850.609375 217.00 -69898.031250 I wanted to extract only this last frame from the .trr file, so used trjconv -f 1IARcompleted_WT_minimized.trr -o 1IARcompleted_WT_minimized_217.trr -dump 217 However, this seems to have produced a file with the t = 200 ps conformation, though the dump parameter was recorded, as trjconv output was: . . Option Type Value Description -- . . -dumptime 217 Dump frame nearest specified time (ps) . . Will write trr: Trajectory in portable xdr format trn version: GMX_trn_file (single precision) Reading frame 2 time 200.000 Dumping frame at t= 200 ps Reading frame 3 time 217.000 . . Also, using -dump 200 gave an identical file to the above. Any idea why the expected timeframe isn't reproduced? I'm using gromacs 4.5.3. Thanks, Ehud Schreiber. -- 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] FW: trjconv -dump problem
On 14/02/2012 7:59 PM, Ehud Schreiber wrote: Dear Gromacs users, I minimized a protein structure 1IARcompleted_WT.pdb, getting, among others, the files 1IARcompleted_WT_minimized.trr and 1IARcompleted_WT_minimized_potential_energy.xvg . Looking at the latter file showed the last frame to be at 217 ps: . . @title Gromacs Energies @xaxis label Time (ps) @yaxis label (kJ/mol) . . @ s0 legend Potential 0.00 -31997.519531 0.00 -33810.406250 200.00 -69850.609375 217.00 -69898.031250 I wanted to extract only this last frame from the .trr file, so used trjconv -f 1IARcompleted_WT_minimized.trr -o 1IARcompleted_WT_minimized_217.trr -dump 217 However, this seems to have produced a file with the t = 200 ps conformation, though the dump parameter was recorded, as trjconv output was: . . Option Type Value Description -- . . -dumptime 217 Dump frame nearest specified time (ps) . . Will write trr: Trajectory in portable xdr format trn version: GMX_trn_file (single precision) Reading frame 2 time 200.000 Dumping frame at t= 200 ps Reading frame 3 time 217.000 . . Also, using -dump 200 gave an identical file to the above. Any idea why the expected timeframe isn't reproduced? I'm using gromacs 4.5.3. Your nstxout parameter means not every frame is written. Prior to the implementation of checkpointing, the final frame was written to the .trr regardless of nstxout, but that no longer occurs. The final frame is in your checkpoint file, and you can use that anywhere you might use a coordinate file - including trjconv to get a simple coordinate file from it. Mark -- 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] Npt equilibration of the membrane-mimicking CCl4 layer
This also was solved by the some extra minimisation steps. I've forced with another problem :D During npt equilibration my system have slightly expanded so my desired volume and density were perturbed. I've noticed the below options in npt wich could help me ref_p= 1 1 compressibility = 4.5e-5 i'm using this compressibility value because I'm modelling the lipid-like environment so I think that I must increase pressure. Could you remind me the dependence of pressure from density and volume for liquids ? :) James 2012/2/14 James Starlight jmsstarli...@gmail.com It seems that I've fixed that problem by reduce vdv radii for Cl during defining of my box Eventually I've obtained box with the desired density than I've delete vdvradii.dat for my wor dir by when I've launched equilibration I've oibtained Fatal error: Too many LINCS warnings (1598) If you know what you are doing you can adjust the lincs warning threshold in your mdp file I've never seen this before I'm using 1.o cutoff for pme and 1.4 for vdv my LINKS parameters are ; Bond parameters continuation= no; first dynamics run constraint_algorithm = lincs; holonomic constraints constraints= all-bonds; all bonds (even heavy atom-H bonds) constrained lincs_iter= 1; accuracy of LINCS lincs_order= 4; also related to accuracy How I could solve it? James 2012/2/14 James Starlight jmsstarli...@gmail.com Mark, I've checked only density value with 500 molecules Ccl4 I have density that is twisely less that I need ( in accordance to the literature ). Also I've checked my box visually and found that the box is not properly tightly packed so I dont know why genbox didnt add some extra mollecules :( In other words I wounder to know if there is any way to add some extra molecules to the pre defined box to make my system more tighly packed ( to short distance between existing molecules and place new ones in the new space ) ? James 2012/2/14 Mark Abraham mark.abra...@anu.edu.au On 14/02/2012 4:57 PM, James Starlight wrote: Justin, Firstly I've created the box of desired size with only 500 molecules ( I need 1000) Than I've tried to add extra 200 molecules by means of Genbox genbox -cp super_box.gro -ci Ccl4.gro -nmol 200 -o new_solv.gro but no molecules have been added Added 0 molecules (out of 200 requested) of Cl4 ... then there are no gaps large enough to insert your molecules. Either make gaps, or check out genbox -h for advice on defining the radii. also I've tried genbox -cp super_box.gro -cp Ccl4.gro -nmol 200 -o new_solv.gro Two -cp options is not what you want, and -nmol probably only works with -ci. but system were crashed with message Reading solute configuration God Rules Over Mankind, Animals, Cosmos and Such Containing 2500 atoms in 500 residues Initialising van der waals distances... WARNING: masses and atomic (Van der Waals) radii will be determined based on residue and atom names. These numbers can deviate from the correct mass and radius of the atom type. Reading solvent configuration God Rules Over Mankind, Animals, Cosmos and Such solvent configuration contains 5 atoms in 1 residues Is there any ways to add extra mollecules to the pre defined box ? Yes - but there has to be room for them. Mark -- gmx-users mailing listgmx-users@gromacs.org http://lists.gromacs.org/**mailman/listinfo/gmx-usershttp://lists.gromacs.org/mailman/listinfo/gmx-users Please search the archive at http://www.gromacs.org/** Support/Mailing_Lists/Searchhttp://www.gromacs.org/Support/Mailing_Lists/Searchbefore 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_Listshttp://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] RE: trjconv -dump problem
Dear Mark (or anybody else interested), The .trr file does include the final (t = 217) frame - first, trjconv said: Reading frame 3 time 217.000 and second this is verified by converting the whole trajectory to .gro: trjconv -f 1IARcompleted_WT_minimized.trr -s 1IARcompleted_WT_minimized.tpr -o 1IARcompleted_WT_minimized_path.gro It therefore seems the behavior is a bug, as the last frame is there, and is needed, especially in minimizations. An output in a .gro format is not sufficient for a further minimization because of its limited accuracy format of coordinates. Also, my mdrun did not produce a checkpoint file (I'm not sure whether because I didn't ask to or because the run was shorter than 15 minutes). Thanks, Ehud Schreiber. -- Message: 2 Date: Tue, 14 Feb 2012 21:58:10 +1100 From: Mark Abraham mark.abra...@anu.edu.au Subject: Re: [gmx-users] FW: trjconv -dump problem To: Discussion list for GROMACS users gmx-users@gromacs.org Message-ID: 4f3a3e42.6000...@anu.edu.au Content-Type: text/plain; charset=iso-8859-1 On 14/02/2012 7:59 PM, Ehud Schreiber wrote: Dear Gromacs users, I minimized a protein structure 1IARcompleted_WT.pdb, getting, among others, the files 1IARcompleted_WT_minimized.trr and 1IARcompleted_WT_minimized_potential_energy.xvg . Looking at the latter file showed the last frame to be at 217 ps: . . @title Gromacs Energies @xaxis label Time (ps) @yaxis label (kJ/mol) . . @ s0 legend Potential 0.00 -31997.519531 0.00 -33810.406250 200.00 -69850.609375 217.00 -69898.031250 I wanted to extract only this last frame from the .trr file, so used trjconv -f 1IARcompleted_WT_minimized.trr -o 1IARcompleted_WT_minimized_217.trr -dump 217 However, this seems to have produced a file with the t = 200 ps conformation, though the dump parameter was recorded, as trjconv output was: . . Option Type Value Description -- . . -dumptime 217 Dump frame nearest specified time (ps) . . Will write trr: Trajectory in portable xdr format trn version: GMX_trn_file (single precision) Reading frame 2 time 200.000 Dumping frame at t= 200 ps Reading frame 3 time 217.000 . . Also, using -dump 200 gave an identical file to the above. Any idea why the expected timeframe isn't reproduced? I'm using gromacs 4.5.3. Your nstxout parameter means not every frame is written. Prior to the implementation of checkpointing, the final frame was written to the .trr regardless of nstxout, but that no longer occurs. The final frame is in your checkpoint file, and you can use that anywhere you might use a coordinate file - including trjconv to get a simple coordinate file from it. Mark -- 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] GROMOS96 53a6 and PRODRG topologies
Dear: My receptor is membrane proteins,the I get protein-ligand complex.I want to do MD simulations in DPPC.Can I use the GROMOS96 53a6 force field modified in order to include Berger’s parameters for lipids?The topology for the ligand was created employing the server PRODRG 2.5 Beta.How to change the charges? -- 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] GROMOS96 53a6 and PRODRG topologies
xiaojiong wrote: Dear: My receptor is membrane proteins,the I get protein-ligand complex.I want to do MD simulations in DPPC.Can I use the GROMOS96 53a6 force field modified in order to include Berger’s parameters for lipids? Sure, but there are other choices, as well. The topology for the ligand was created employing the server PRODRG 2.5 Beta. How to change the charges? http://pubs.acs.org/doi/abs/10.1021/ci100335w -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] Npt equilibration of the membrane-mimicking CCl4 layer
I assume that you energy minimisd the system, but still have atomic clashes? One thing which helped me in a similar case, was a short simulation at low temperature with a really small timestep (about 3-5 magnitudes smaller than the normal timestep). With this the atoms which clashes move away from each other, but due to the very small timestep, they are not fast as rockets and shouldn't lead to an exploding system. Then when most clashes are resolved you can continue with a normal equilibration. Greetings Thomas It seems that I've fixed that problem by reduce vdv radii for Cl during defining of my box Eventually I've obtained box with the desired density than I've delete vdvradii.dat for my wor dir by when I've launched equilibration I've oibtained Fatal error: Too many LINCS warnings (1598) If you know what you are doing you can adjust the lincs warning threshold in your mdp file I've never seen this before I'm using 1.o cutoff for pme and 1.4 for vdv my LINKS parameters are ; Bond parameters continuation= no; first dynamics run constraint_algorithm = lincs; holonomic constraints constraints= all-bonds; all bonds (even heavy atom-H bonds) constrained lincs_iter= 1; accuracy of LINCS lincs_order= 4; also related to accuracy How I could solve it? James 2012/2/14 James Starlight jmsstarli...@gmail.com Mark, I've checked only density value with 500 molecules Ccl4 I have density that is twisely less that I need ( in accordance to the literature ). Also I've checked my box visually and found that the box is not properly tightly packed so I dont know why genbox didnt add some extra mollecules :( In other words I wounder to know if there is any way to add some extra molecules to the pre defined box to make my system more tighly packed ( to short distance between existing molecules and place new ones in the new space ) ? James 2012/2/14 Mark Abraham mark.abra...@anu.edu.au On 14/02/2012 4:57 PM, James Starlight wrote: Justin, Firstly I've created the box of desired size with only 500 molecules ( I need 1000) Than I've tried to add extra 200 molecules by means of Genbox genbox -cp super_box.gro -ci Ccl4.gro -nmol 200 -o new_solv.gro but no molecules have been added Added 0 molecules (out of 200 requested) of Cl4 ... then there are no gaps large enough to insert your molecules. Either make gaps, or check out genbox -h for advice on defining the radii. also I've tried genbox -cp super_box.gro -cp Ccl4.gro -nmol 200 -o new_solv.gro Two -cp options is not what you want, and -nmol probably only works with -ci. but system were crashed with message Reading solute configuration God Rules Over Mankind, Animals, Cosmos and Such Containing 2500 atoms in 500 residues Initialising van der waals distances... WARNING: masses and atomic (Van der Waals) radii will be determined based on residue and atom names. These numbers can deviate from the correct mass and radius of the atom type. Reading solvent configuration God Rules Over Mankind, Animals, Cosmos and Such solvent configuration contains 5 atoms in 1 residues Is there any ways to add extra mollecules to the pre defined box ? Yes - but there has to be room for them. Mark -- gmx-users mailing listgmx-users@gromacs.org http://lists.gromacs.org/**mailman/listinfo/gmx-usershttp://lists.gromacs.org/mailman/listinfo/gmx-users Please search the archive at http://www.gromacs.org/** Support/Mailing_Lists/Searchhttp://www.gromacs.org/Support/Mailing_Lists/Searchbefore 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_Listshttp://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] Domain decomposition in FEP
Dear all, I'm trying to run FEP calculation of the ligand in the protein in water. Performing NVT-dynamics for Lennard-Johnes perturbation using soft core results in the following error: There is no domain decomposition for 20 nodes... I've checked out the log-file and found this: Initial maximum inter charge-group distances: two-body bonded interactions: 1.474 nm, Bond, atoms 4371 4391 multi-body bonded interactions: 0.522 nm, Ryckaert-Bell., atoms 3747 3756 Minimum cell size due to bonded interactions: 1.621 nm Guess for relative PME load: 0.29 Will use 20 particle-particle and 12 PME only nodes This is a guess, check the performance at the end of the log file Using 12 separate PME nodes Scaling the initial minimum size with 1/0.8 (option -dds) = 1.25 Optimizing the DD grid for 20 cells with a minimum initial size of 2.026 nm The maximum allowed number of cells is: X 3 Y 3 Z 3 --- Program mdrun_mpi, VERSION 4.5.1 Source code file: domdec.c, line: 6428 Fatal error: There is no domain decomposition for 20 nodes that is compatible with the given box and a minimum cell size of 2.02621 nm Change the number of nodes or mdrun option -rdd or -dds Look in the log file for details on the domain decomposition For more information and tips for troubleshooting, please check the GROMACS website at http://www.gromacs.org/Documentation/Errors --- I understand that distance of 1.474 nm is rather big for bond! BUT I can't find such a distance between atoms 4371 4391. The fragment of gro-file for such index interval looks like this: 640LIG C4 4371 2.897 3.416 3.072 640LIG N3 4372 3.296 3.109 3.434 640LIG N6 4373 3.416 2.999 3.595 640LIG C6 4374 2.696 2.762 3.601 640LIG C7 4375 3.025 3.232 3.150 640LIG C8 4376 3.208 3.142 3.338 640LIG C9 4377 2.758 2.644 3.649 640LIG C10 4378 2.767 2.884 3.595 640LIG N4 4379 3.319 2.826 3.731 640LIG C11 4380 2.968 2.768 3.677 640LIG C12 4381 2.894 2.648 3.686 640LIG C13 4382 2.901 2.887 3.638 640LIG N5 4383 3.100 2.766 3.707 640LIG C14 4384 3.420 2.908 3.693 640LIG O1 4385 3.018 3.627 3.068 640LIG C15 4386 3.028 3.095 3.182 640LIG C16 4387 3.123 3.049 3.275 640LIG O2 4388 2.566 2.759 3.559 640LIG O3 4389 2.699 2.993 3.545 640LIG O4 4390 2.680 2.530 3.655 640LIG C17 4391 2.791 3.426 3.184 640LIG C18 4392 2.834 3.453 2.937 and the actual distance is 1.55 A - quite OK for single C-C bond. And another interesting fact is that it doesn't occur during FEP calculation of the same system without soft-core and aimed at only charge perturbation. I will appreciate any help regarding this. Alexey Zeifman -- 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] Number of windows in umbrella sampling
shahid nayeem wrote: Thanks Justin. I will try again. But please refer to some protocol if you know and one last question that before doing umbrella sampling simulation how can one be sure that the pulling is good and one should go ahead with selecting window and doing umbrella sampling. In the end when you see histo.xvg and profile.xvg , it cost a lot computationally and if you don't get it right these computational resources are wasted. There are certainly simulations that have been published using umbrella sampling on a variety of systems. I'm sure some relatively quick literature searching will lead you to some suitable ideas. Spending a few hours doing some reading will save you months of wasted CPU time. -Justin Shahid nayeem On Mon, Feb 13, 2012 at 8:12 PM, Justin A. Lemkul jalem...@vt.edu mailto:jalem...@vt.edu wrote: shahid nayeem wrote: Thanks for quick reply. I have created mutant of a complex by changing interface residue in VMD. These mutant are experimentally known to show less binding affinity. I want to reproduce these results with umbrella sampling. Now I am sending profile and histo file for wt and mutant. Please suggest where i am wrong. The PMF curves look poorly converged. Your reaction coordinates are not the same for both the WT and mutant (you appear to have a far shorter reaction coordinate for the mutant). The energy minimum is also ill-defined for the mutant. As for the reason behind these phenomena, I cannot say, nor do I have time to sort through your data and try to work it out for you. Refer to the literature, find similar protocols, and proceed from there. -Justin Shahid Nayeem On Mon, Feb 13, 2012 at 7:15 PM, Justin A. Lemkul jalem...@vt.edu mailto:jalem...@vt.edu mailto:jalem...@vt.edu mailto:jalem...@vt.edu wrote: shahid nayeem wrote: Dear Justin I am doing umbrella pulling simulation of a protein complex wt and mutant. I expect mutant to give lower deltG value. I am attaching a tif file of energy vs time curve of wt and mutant protein on pulling simulation. These energies are obtained by g_energy and selecting 11 option which is COM pulling energy. In this curve the first peak decreases and again rises at longer time. How many windows should be selected. As expected the peak of COM pulling energy is lower in mutants. Please explain why the energy again rises at higher time. Should I use the windows upto 160ps only because thereafter in both curve there is rise in energy value. the pull code used is as follows. What you have obtained is a path-dependent energy that may or may not signify anything useful - it almost certainly does not. I cannot offer an explanation of the sharp increase towards the end of the simulation other than to speculate that your box is of insufficient size and you're encountering PBC issues. As for how many windows are necessary, it's also impossible to tell. You need enough windows to adequately sample the reaction coordinate. Thus, it is decided based on how far you need to separate the two species, how strong the force constant is during the US simulations, the nature of the interactions in the system and how fast they converge, etc. -Justin pull_geometry = distance pull_dim= Y N Y pull_vec1 = 0.75 0 1 pull_start = yes pull_ngroups= 1 pull_group0 = Chain_B pull_group1 = Chain_A pull_rate1 = 0.01 pull_k1 = 1000 -- ==== Justin A. Lemkul Ph.D. Candidate ICTAS Doctoral Scholar MILES-IGERT Trainee Department of Biochemistry Virginia Tech Blacksburg, VA jalemkul[at]vt.edu http://vt.edu http://vt.edu | (540) 231-9080 http://www.bevanlab.biochem.vt.edu/Pages/Personal/justin http://vt.edu/Pages/Personal/justin http://www.bevanlab.biochem.__vt.edu/Pages/Personal/justin http://www.bevanlab.biochem.vt.edu/Pages/Personal/justin ==== -- gmx-users mailing listgmx-users@gromacs.org mailto:gmx-users@gromacs.org mailto:gmx-users@gromacs.org mailto:gmx-users@gromacs.org http://lists.gromacs.org/mailman/listinfo/gmx-users
Re: [gmx-users] Number of windows in umbrella sampling
Justin A. Lemkul wrote: shahid nayeem wrote: Thanks Justin. I will try again. But please refer to some protocol if you know and one last question that before doing umbrella sampling simulation how can one be sure that the pulling is good and one should go ahead with selecting window and doing umbrella sampling. In the end when you see histo.xvg and profile.xvg , it cost a lot computationally and if you don't get it right these computational resources are wasted. There are certainly simulations that have been published using umbrella sampling on a variety of systems. I'm sure some relatively quick literature searching will lead you to some suitable ideas. Spending a few hours doing some reading will save you months of wasted CPU time. I should, of course, mention that what you have done is not necessarily wasted, it's simply insufficient. Simulations can always be extended and new umbrella sampling windows added after the fact and incorporated into the analysis. A better experimental plan will save these headaches in the future. -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
Re: [gmx-users] Domain decomposition in FEP
Dear Alexey Zeifman, From my experience it is just the nodes have to break down into some sort of even multiple by dividing the domain decomposition size or vise versi. I think somone elses answere was just to try it with a couple sizes for the decomposition, as multiples of the number of nodes. also your PME reserved nodes are probably too much or your wasting space, but I am not sure if this effects the domain decomposition or if all 20 nodes are thrown into that? Stephan Watkins Original-Nachricht Datum: Tue, 14 Feb 2012 18:01:24 +0400 Von: Alexey Zeifman az...@mail.ru An: gmx-users@gromacs.org Betreff: [gmx-users] Domain decomposition in FEP Dear all, I'm trying to run FEP calculation of the ligand in the protein in water. Performing NVT-dynamics for Lennard-Johnes perturbation using soft core results in the following error: There is no domain decomposition for 20 nodes... I've checked out the log-file and found this: Initial maximum inter charge-group distances: two-body bonded interactions: 1.474 nm, Bond, atoms 4371 4391 multi-body bonded interactions: 0.522 nm, Ryckaert-Bell., atoms 3747 3756 Minimum cell size due to bonded interactions: 1.621 nm Guess for relative PME load: 0.29 Will use 20 particle-particle and 12 PME only nodes This is a guess, check the performance at the end of the log file Using 12 separate PME nodes Scaling the initial minimum size with 1/0.8 (option -dds) = 1.25 Optimizing the DD grid for 20 cells with a minimum initial size of 2.026 nm The maximum allowed number of cells is: X 3 Y 3 Z 3 --- Program mdrun_mpi, VERSION 4.5.1 Source code file: domdec.c, line: 6428 Fatal error: There is no domain decomposition for 20 nodes that is compatible with the given box and a minimum cell size of 2.02621 nm Change the number of nodes or mdrun option -rdd or -dds Look in the log file for details on the domain decomposition For more information and tips for troubleshooting, please check the GROMACS website at http://www.gromacs.org/Documentation/Errors --- I understand that distance of 1.474 nm is rather big for bond! BUT I can't find such a distance between atoms 4371 4391. The fragment of gro-file for such index interval looks like this: 640LIG C4 4371 2.897 3.416 3.072 640LIG N3 4372 3.296 3.109 3.434 640LIG N6 4373 3.416 2.999 3.595 640LIG C6 4374 2.696 2.762 3.601 640LIG C7 4375 3.025 3.232 3.150 640LIG C8 4376 3.208 3.142 3.338 640LIG C9 4377 2.758 2.644 3.649 640LIG C10 4378 2.767 2.884 3.595 640LIG N4 4379 3.319 2.826 3.731 640LIG C11 4380 2.968 2.768 3.677 640LIG C12 4381 2.894 2.648 3.686 640LIG C13 4382 2.901 2.887 3.638 640LIG N5 4383 3.100 2.766 3.707 640LIG C14 4384 3.420 2.908 3.693 640LIG O1 4385 3.018 3.627 3.068 640LIG C15 4386 3.028 3.095 3.182 640LIG C16 4387 3.123 3.049 3.275 640LIG O2 4388 2.566 2.759 3.559 640LIG O3 4389 2.699 2.993 3.545 640LIG O4 4390 2.680 2.530 3.655 640LIG C17 4391 2.791 3.426 3.184 640LIG C18 4392 2.834 3.453 2.937 and the actual distance is 1.55 A - quite OK for single C-C bond. And another interesting fact is that it doesn't occur during FEP calculation of the same system without soft-core and aimed at only charge perturbation. I will appreciate any help regarding this. Alexey Zeifman -- NEU: FreePhone 3-fach-Flat mit kostenlosem Smartphone! Jetzt informieren: http://www.gmx.net/de/go/freephone/ -- 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] Domain decomposition in FEP
Alexey Zeifman wrote: Dear all, I'm trying to run FEP calculation of the ligand in the protein in water. Performing NVT-dynamics for Lennard-Johnes perturbation using soft core results in the following error: There is no domain decomposition for 20 nodes... I've checked out the log-file and found this: Initial maximum inter charge-group distances: two-body bonded interactions: 1.474 nm, Bond, atoms 4371 4391 multi-body bonded interactions: 0.522 nm, Ryckaert-Bell., atoms 3747 3756 Minimum cell size due to bonded interactions: 1.621 nm Guess for relative PME load: 0.29 Will use 20 particle-particle and 12 PME only nodes This is a guess, check the performance at the end of the log file Using 12 separate PME nodes Scaling the initial minimum size with 1/0.8 (option -dds) = 1.25 Optimizing the DD grid for 20 cells with a minimum initial size of 2.026 nm The maximum allowed number of cells is: X 3 Y 3 Z 3 --- Program mdrun_mpi, VERSION 4.5.1 Source code file: domdec.c, line: 6428 Fatal error: There is no domain decomposition for 20 nodes that is compatible with the given box and a minimum cell size of 2.02621 nm Change the number of nodes or mdrun option -rdd or -dds Look in the log file for details on the domain decomposition For more information and tips for troubleshooting, please check the GROMACS website at http://www.gromacs.org/Documentation/Errors --- I understand that distance of 1.474 nm is rather big for bond! BUT I can't find such a distance between atoms 4371 4391. The fragment of gro-file for such index interval looks like this: 640LIG C4 4371 2.897 3.416 3.072 640LIG N3 4372 3.296 3.109 3.434 640LIG N6 4373 3.416 2.999 3.595 640LIG C6 4374 2.696 2.762 3.601 640LIG C7 4375 3.025 3.232 3.150 640LIG C8 4376 3.208 3.142 3.338 640LIG C9 4377 2.758 2.644 3.649 640LIGC10 4378 2.767 2.884 3.595 640LIG N4 4379 3.319 2.826 3.731 640LIGC11 4380 2.968 2.768 3.677 640LIGC12 4381 2.894 2.648 3.686 640LIGC13 4382 2.901 2.887 3.638 640LIG N5 4383 3.100 2.766 3.707 640LIGC14 4384 3.420 2.908 3.693 640LIG O1 4385 3.018 3.627 3.068 640LIGC15 4386 3.028 3.095 3.182 640LIGC16 4387 3.123 3.049 3.275 640LIG O2 4388 2.566 2.759 3.559 640LIG O3 4389 2.699 2.993 3.545 640LIG O4 4390 2.680 2.530 3.655 640LIGC17 4391 2.791 3.426 3.184 640LIGC18 4392 2.834 3.453 2.937 and the actual distance is 1.55 A - quite OK for single C-C bond. And another interesting fact is that it doesn't occur during FEP calculation of the same system without soft-core and aimed at only charge perturbation. I will appreciate any help regarding this. The general information about this error can be found at: http://www.gromacs.org/Documentation/Errors#There_is_no_domain_decomposition_for_n_nodes_that_is_compatible_with_the_given_box_and_a_minimum_cell_size_of_x_nm The DD setup is based on a host of factors. Seeing a complete .mdp file would help. -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
Re: [gmx-users] Domain decomposition in FEP
On 15/02/2012 1:44 AM, lloyd riggs wrote: Dear Alexey Zeifman, From my experience it is just the nodes have to break down into some sort of even multiple by dividing the domain decomposition size or vise versi. I think somone elses answere was just to try it with a couple sizes for the decomposition, as multiples of the number of nodes. also your PME reserved nodes are probably too much or your wasting space, but I am not sure if this effects the domain decomposition or if all 20 nodes are thrown into that? This is not really the issue. You can read in the OP's output that there are 20 PP and 12 PME nodes, which is a normal enough ratio, and fine for a 0.29 PME load... If the OP is confident that his initial configuration really doesn't have such a long bond (no input checkpoint file messing things up?) he should try again with GROMACS 4.5.5. Mark Stephan Watkins Original-Nachricht Datum: Tue, 14 Feb 2012 18:01:24 +0400 Von: Alexey Zeifmanaz...@mail.ru An: gmx-users@gromacs.org Betreff: [gmx-users] Domain decomposition in FEP Dear all, I'm trying to run FEP calculation of the ligand in the protein in water. Performing NVT-dynamics for Lennard-Johnes perturbation using soft core results in the following error: There is no domain decomposition for 20 nodes... I've checked out the log-file and found this: Initial maximum inter charge-group distances: two-body bonded interactions: 1.474 nm, Bond, atoms 4371 4391 multi-body bonded interactions: 0.522 nm, Ryckaert-Bell., atoms 3747 3756 Minimum cell size due to bonded interactions: 1.621 nm Guess for relative PME load: 0.29 Will use 20 particle-particle and 12 PME only nodes This is a guess, check the performance at the end of the log file Using 12 separate PME nodes Scaling the initial minimum size with 1/0.8 (option -dds) = 1.25 Optimizing the DD grid for 20 cells with a minimum initial size of 2.026 nm The maximum allowed number of cells is: X 3 Y 3 Z 3 --- Program mdrun_mpi, VERSION 4.5.1 Source code file: domdec.c, line: 6428 Fatal error: There is no domain decomposition for 20 nodes that is compatible with the given box and a minimum cell size of 2.02621 nm Change the number of nodes or mdrun option -rdd or -dds Look in the log file for details on the domain decomposition For more information and tips for troubleshooting, please check the GROMACS website at http://www.gromacs.org/Documentation/Errors --- I understand that distance of 1.474 nm is rather big for bond! BUT I can't find such a distance between atoms 4371 4391. The fragment of gro-file for such index interval looks like this: 640LIG C4 4371 2.897 3.416 3.072 640LIG N3 4372 3.296 3.109 3.434 640LIG N6 4373 3.416 2.999 3.595 640LIG C6 4374 2.696 2.762 3.601 640LIG C7 4375 3.025 3.232 3.150 640LIG C8 4376 3.208 3.142 3.338 640LIG C9 4377 2.758 2.644 3.649 640LIGC10 4378 2.767 2.884 3.595 640LIG N4 4379 3.319 2.826 3.731 640LIGC11 4380 2.968 2.768 3.677 640LIGC12 4381 2.894 2.648 3.686 640LIGC13 4382 2.901 2.887 3.638 640LIG N5 4383 3.100 2.766 3.707 640LIGC14 4384 3.420 2.908 3.693 640LIG O1 4385 3.018 3.627 3.068 640LIGC15 4386 3.028 3.095 3.182 640LIGC16 4387 3.123 3.049 3.275 640LIG O2 4388 2.566 2.759 3.559 640LIG O3 4389 2.699 2.993 3.545 640LIG O4 4390 2.680 2.530 3.655 640LIGC17 4391 2.791 3.426 3.184 640LIGC18 4392 2.834 3.453 2.937 and the actual distance is 1.55 A - quite OK for single C-C bond. And another interesting fact is that it doesn't occur during FEP calculation of the same system without soft-core and aimed at only charge perturbation. I will appreciate any help regarding this. Alexey Zeifman -- 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] rotational autocorrelation function
Hello, I have done COM pulling simulation to pull a small molecule (5-Fluorouracil) through a path. I want to see the rotational freedom for the pulled over molecule during the SMD simulation. Is it meaningful to calculate the rotational autocorrelation for a molecule which observes external force? Regards, Vijay. -- 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] RE: trjconv -dump problem
On 14/02/2012 11:39 PM, Ehud Schreiber wrote: Dear Mark (or anybody else interested), The .trr file does include the final (t = 217) frame - first, trjconv said: Reading frame 3 time 217.000 and second this is verified by converting the whole trajectory to .gro: trjconv -f 1IARcompleted_WT_minimized.trr -s 1IARcompleted_WT_minimized.tpr -o 1IARcompleted_WT_minimized_path.gro It therefore seems the behavior is a bug, as the last frame is there, and is needed, especially in minimizations. An output in a .gro format is not sufficient for a further minimization because of its limited accuracy format of coordinates. I looked at the code. This is an artefact of trjconv having to try to work out the spacing of the frames in your .trr to then work out which frame is probably going to be the closest. It assumes the frame will be equally spaced, and accepts the first frame that is closer than half the spacing it computes from the first two frames. Since your frames are presumably 100, 200 and 217, the spacing is 100 and it accepts 200 because it's the first one it finds within 100/2 of 217. You can work around this in lots of ways (trjconv -b 217; using nstxout=0 in the first EM; using nstxout=1 in the first EM; using nstxout extremely large in the first EM; taking the EM output .trr and giving it to the next grompp -t without doing anything). Also, my mdrun did not produce a checkpoint file (I'm not sure whether because I didn't ask to or because the run was shorter than 15 minutes). OK, maybe EM has a different behaviour (no checkpoint, thus write last configuration always). Mark Thanks, Ehud Schreiber. -- Message: 2 Date: Tue, 14 Feb 2012 21:58:10 +1100 From: Mark Abrahammark.abra...@anu.edu.au Subject: Re: [gmx-users] FW: trjconv -dump problem To: Discussion list for GROMACS usersgmx-users@gromacs.org Message-ID:4f3a3e42.6000...@anu.edu.au Content-Type: text/plain; charset=iso-8859-1 On 14/02/2012 7:59 PM, Ehud Schreiber wrote: Dear Gromacs users, I minimized a protein structure 1IARcompleted_WT.pdb, getting, among others, the files 1IARcompleted_WT_minimized.trr and 1IARcompleted_WT_minimized_potential_energy.xvg . Looking at the latter file showed the last frame to be at 217 ps: . . @title Gromacs Energies @xaxis label Time (ps) @yaxis label (kJ/mol) . . @ s0 legend Potential 0.00 -31997.519531 0.00 -33810.406250 200.00 -69850.609375 217.00 -69898.031250 I wanted to extract only this last frame from the .trr file, so used trjconv -f 1IARcompleted_WT_minimized.trr -o 1IARcompleted_WT_minimized_217.trr -dump 217 However, this seems to have produced a file with the t = 200 ps conformation, though the dump parameter was recorded, as trjconv output was: . . Option Type Value Description -- . . -dumptime 217 Dump frame nearest specified time (ps) . . Will write trr: Trajectory in portable xdr format trn version: GMX_trn_file (single precision) Reading frame 2 time 200.000 Dumping frame at t= 200 ps Reading frame 3 time 217.000 . . Also, using -dump 200 gave an identical file to the above. Any idea why the expected timeframe isn't reproduced? I'm using gromacs 4.5.3. Your nstxout parameter means not every frame is written. Prior to the implementation of checkpointing, the final frame was written to the .trr regardless of nstxout, but that no longer occurs. The final frame is in your checkpoint file, and you can use that anywhere you might use a coordinate file - including trjconv to get a simple coordinate file from it. Mark -- 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] Npt equilibration of the membrane-mimicking CCl4 layer
On 14/02/2012 11:01 PM, James Starlight wrote: This also was solved by the some extra minimisation steps. I've forced with another problem :D During npt equilibration my system have slightly expanded so my desired volume and density were perturbed. I've noticed the below options in npt wich could help me ref_p= 1 1 compressibility = 4.5e-5 i'm using this compressibility value because I'm modelling the lipid-like environment so I think that I must increase pressure. Could you remind me the dependence of pressure from density and volume for liquids ? :) Your forcefield, simulation cell contents and .mdp settings will determine the equilibrium density. Whether you need to do anything depends on whether you've made a statistically significant post-equilibration measurement of your average density. Haphazardly increasing the reference pressure for the coupling will reduce the volume, but now you are simulating at that pressure. See http://www.gromacs.org/Documentation/Terminology/Pressure for background info. Mark James 2012/2/14 James Starlight jmsstarli...@gmail.com mailto:jmsstarli...@gmail.com It seems that I've fixed that problem by reduce vdv radii for Cl during defining of my box Eventually I've obtained box with the desired density than I've delete vdvradii.dat for my wor dir by when I've launched equilibration I've oibtained Fatal error: Too many LINCS warnings (1598) If you know what you are doing you can adjust the lincs warning threshold in your mdp file I've never seen this before I'm using 1.o cutoff for pme and 1.4 for vdv my LINKS parameters are ; Bond parameters continuation= no; first dynamics run constraint_algorithm = lincs; holonomic constraints constraints= all-bonds; all bonds (even heavy atom-H bonds) constrained lincs_iter= 1; accuracy of LINCS lincs_order= 4; also related to accuracy How I could solve it? James 2012/2/14 James Starlight jmsstarli...@gmail.com mailto:jmsstarli...@gmail.com Mark, I've checked only density value with 500 molecules Ccl4 I have density that is twisely less that I need ( in accordance to the literature ). Also I've checked my box visually and found that the box is not properly tightly packed so I dont know why genbox didnt add some extra mollecules :( In other words I wounder to know if there is any way to add some extra molecules to the pre defined box to make my system more tighly packed ( to short distance between existing molecules and place new ones in the new space ) ? James 2012/2/14 Mark Abraham mark.abra...@anu.edu.au mailto:mark.abra...@anu.edu.au On 14/02/2012 4:57 PM, James Starlight wrote: Justin, Firstly I've created the box of desired size with only 500 molecules ( I need 1000) Than I've tried to add extra 200 molecules by means of Genbox genbox -cp super_box.gro -ci Ccl4.gro -nmol 200 -o new_solv.gro but no molecules have been added Added 0 molecules (out of 200 requested) of Cl4 ... then there are no gaps large enough to insert your molecules. Either make gaps, or check out genbox -h for advice on defining the radii. also I've tried genbox -cp super_box.gro -cp Ccl4.gro -nmol 200 -o new_solv.gro Two -cp options is not what you want, and -nmol probably only works with -ci. but system were crashed with message Reading solute configuration God Rules Over Mankind, Animals, Cosmos and Such Containing 2500 atoms in 500 residues Initialising van der waals distances... WARNING: masses and atomic (Van der Waals) radii will be determined based on residue and atom names. These numbers can deviate from the correct mass and radius of the atom type. Reading solvent configuration God Rules Over Mankind, Animals, Cosmos and Such solvent configuration contains 5 atoms in 1 residues Is there any ways to add extra mollecules to the pre defined box ? Yes - but there has to be room for them. Mark -- gmx-users mailing list gmx-users@gromacs.org mailto: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
[gmx-users] what does the Coulomb (SR) mean under PME?
Hi all, Anybody can help me out! Thanks in advance? In my test simulation, there are one Na^+ and one Cl^- (the distance of 1nm) in vacuum in a very big simulation box (10*10*10 nm^3). I calculated the energies under two different conditions, coulombtype=PME vs coulombtype=cut-off. 1. coulombtype=PME Coul. recip. -139.191 -- 0 0 (kJ/mol) Coul-SR:NA-CL*-0.0322805* -- 0 0 (kJ/mol) LJ-SR:NA-CL -0.00100285 -- 0 0 (kJ/mol) 2. coulombtype=cut-off Coul-SR:NA-CL * -138.935* -- 0 0 (kJ/mol) LJ-SR:NA-CL -0.00100285 -- 0 0 (kJ/mol) I did some calculations. Under coulombtype=cut-off, E(Coul-SR) = f*e_i*e_j/r with f=138.935, which is what I expect. However under coulombtype=PME, how does Gromacs get the value of -0.0322805kJ/mol for *SR* term? If I understand it correctly, SR means short-range, which is within the cut-off distance. Since the cut-off distance (1.2nm here) is much smaller than the simulation box length (10nm), there is not short-range interaction with the image box. Some details: Gromacs 4.5.5 PBC is used. distance of Na^+ and Cl^-: 1nm Na^+ and Cl^- are frozen in XYZ dimensions. Only 1 step is run. best regards, Baofu -- 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] Umbrella sampling and WHAM along a curved pathway?
Dear all, I would like to calculate the PMF along a curved reaction pathway using umbrella sampling. I just wonder if it is appropriate to use g_wham to extract the PMF along the curved pathway? Any help would be appreciated. Hao Jiang -- 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] Question regarding freeze groups and their use
Hi all, I know there has been a *lot* of discussion on the mailing list on using freezegrps and potential pitfalls, but after having read much of this discussion the past couple of days and the manual I still find myself with questions and problems I'm hoping someone can help with. What I want to do, specifically, is have my protein frozen (*not* just constrained with some large force, I need actual freezing of every protein atom) but have the solvent still react to the frozen protein. The purpose behind this is to compare the water behavior near a particular part of a protein when it is dynamic and when it is not and thereby learn something about the communication between protein and solvent in this case. However, it seems that I'm not finding a rational set of mdp options to allow this, at least not a set for which dynamics will run with any speed. I have tried many of the things people have discussed on the mailing list, including turning off constraints, turning off pressure coupling, reducing the Protein heat bath temperature to 0 K, and using energygrp_excl. Here are the parameters I just tried: integrator = md nsteps = 100 dt = 0.001 nstxout = 1 nstvout = 1 nstxtcout = 1 nstenergy = 1 nstlog = 1 continuation= yes constraints = none ns_type = grid nstlist = 5 rlist = 1.0 rcoulomb= 1.0 rvdw= 1.0 coulombtype = PME pme_order = 4 fourierspacing = 0.16 tcoupl = V-rescale tc-grps = Protein Non-Protein tau_t = 0.1 0.1 ref_t = 0 300 pcoupl = no pbc = xyz DispCorr= EnerPres gen_vel = no freezegrps = Protein freezedim = Y Y Y ;;; This ran, but *incredibly* slowly (I'm using a single 32 processor node, and only got 255 steps in 10 wallclock minutes). From the discussions I've read on the mailing list it seems that PME is somewhat less than ideal for simulations involving frozen atoms, but I used PME for electrostatics in my non-frozen simulations and would like to keep that consistent if possible (I do get a grompp warning of course). I am also getting warnings in the log file like: DD step 254 vol min/aver 0.183! load imb.: force 2395.3% pme mesh/force 1.243 which I believe must be due to the frozen protein causing difficulties with domain decomposition. Can anyone offer any advice regarding these issues? I know this has been discussed often, but nothing I'm finding in the archives is particularly relevant. -- -- J. Nathan Scott, Ph.D. Postdoctoral Fellow Department of Chemistry and Biochemistry Montana State University -- 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] water channel
Hi Gmxers, Happy Valentine's Day! Sorry, I am still simulating a protein that has water channel. I was just wondering if there is a way to calculate the water density in the channel throughout my trajectory. I was lucky enough to run a NPT simulation. So I guess I can use g_energy to directly get density. But since the number of the channel water is quite dynamic, I do not think I can separate energy group for them in my mdp. Could someone help me with this? Thanks, Yao -- 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] water channel
You can look at gromacs tool g_flux and g_count at: https://github.com/orbeckst/g_count Jianguo From: Yao Yao ya...@ymail.com To: gmx-users@gromacs.org gmx-users@gromacs.org Sent: Wednesday, 15 February 2012, 10:38 Subject: [gmx-users] water channel Hi Gmxers, Happy Valentine's Day! Sorry, I am still simulating a protein that has water channel. I was just wondering if there is a way to calculate the water density in the channel throughout my trajectory. I was lucky enough to run a NPT simulation. So I guess I can use g_energy to directly get density. But since the number of the channel water is quite dynamic, I do not think I can separate energy group for them in my mdp. Could someone help me with this? Thanks, Yao -- 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-- 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] charmm27 in gromacs
Dear Gmx Developer or Users, Can anyone explain which section is for the parameters of improper dihedral angle on the file of ffbonded.itp ? On the file of ffbonded.itp, there is not any comment to differentiate the proper and improper dihedral angle. E.g. on thie file of ffbonded.itp in the directory of charmm27.ff There are two sections about dihedral angle: [ dihedraltypes ] ; i j k l funcphi0cp mult C CT1 NH1 C * 9* 180.00 0.8368 1 .. [ dihedraltypes ] ; i j k l funcq0 cq CPB CPA NPH CPA* 2* 0. 174.0544 .. There is also no explaination the explaination about the function type 2. It seems difficult for users to edit the parameters. Thanks for the explainations! -- 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] charmm27 in gromacs
Tom wrote: Dear Gmx Developer or Users, Can anyone explain which section is for the parameters of improper dihedral angle on the file of ffbonded.itp ? On the file of ffbonded.itp, there is not any comment to differentiate the proper and improper dihedral angle. E.g. on thie file of ffbonded.itp in the directory of charmm27.ff There are two sections about dihedral angle: [ dihedraltypes ] ; i j k l funcphi0cp mult C CT1 NH1 C * 9* 180.00 0.8368 1 .. [ dihedraltypes ] ; i j k l funcq0 cq CPB CPA NPH CPA* 2* 0. 174.0544 .. There is also no explaination the explaination about the function type 2. Table 5.5 of the manual lists all function types, what they are, and what the required parameters are. Functional forms are discussed in section 4.2. It seems difficult for users to edit the parameters. Comprehensive explanations of all nuts and bolts are beyond what can be included in the manual or force field files. A deep knowledge of the underlying force field (which may span several published papers, other manuals, etc) is required to understand any of the force fields implemented in Gromacs (or any other software, for that matter). -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] Water Shell Density
Hi Gmxers, Is there a way to calculate the density of water in a protein hydration layer, like from 5 A to 10 A (radius) from the protein surface? Thanks, Yao-- 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] about the umbrella sampling potential
Dear all, I got confused about the usage of pull-geometry and pull-dim (pull-vector) settings in the umbrella sampling simulations (not the pull / SMD simulations). In an umbrella sampling potential, U(r)=0.5*(r-r0)*(r-r0), are r and r0 referring to (A) the distance between the pulled group and the reference group, or (B) the distance between the pulled group and the reference group *projected* to a vector pre-defined by pull-dim or pull-vector? Note these two distances are not equal unless the pulled group stay strictly on the predefined vector, which I believe is not always the case. If A is true, in principle group1 should be able to move freely within a spherical shell centered at group0, right? What is the point to specify pull-geometry and pull-dim in this case? If B is true, group1 should in principle be able to move freely in plane perpendicular to the pre-defined vector. So if pull-dim= N N Y, the pulled group should be able to move freely in the X-Y plane, right? Any help would be greatly appreciated! Hao Jiang -- 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] Npt equilibration of the membrane-mimicking CCl4 layer
Mark, due to hight density the volume of my system have been slightly increased and during NPT phase I've obtained error Fatal error: One of the box vectors has become shorter than twice the cut-off length or box_yy-|box_zy| or box_zz has become smaller than the cut-off. I'm using 0.9 for electrostatic and 1.4 for vdw cutofs and the dimensions of my box was 6.5 3 3 on the initial step and 6.6 3 3.3 before crush :) I want prevent such expansion of my system by increasing of pressure and/ or compressibility but I have not found exact sollution yet. James 2012/2/14 Mark Abraham mark.abra...@anu.edu.au On 14/02/2012 11:01 PM, James Starlight wrote: This also was solved by the some extra minimisation steps. I've forced with another problem :D During npt equilibration my system have slightly expanded so my desired volume and density were perturbed. I've noticed the below options in npt wich could help me ref_p= 1 1 compressibility = 4.5e-5 i'm using this compressibility value because I'm modelling the lipid-like environment so I think that I must increase pressure. Could you remind me the dependence of pressure from density and volume for liquids ? :) Your forcefield, simulation cell contents and .mdp settings will determine the equilibrium density. Whether you need to do anything depends on whether you've made a statistically significant post-equilibration measurement of your average density. Haphazardly increasing the reference pressure for the coupling will reduce the volume, but now you are simulating at that pressure. See http://www.gromacs.org/Documentation/Terminology/Pressurefor background info. Mark James 2012/2/14 James Starlight jmsstarli...@gmail.com It seems that I've fixed that problem by reduce vdv radii for Cl during defining of my box Eventually I've obtained box with the desired density than I've delete vdvradii.dat for my wor dir by when I've launched equilibration I've oibtained Fatal error: Too many LINCS warnings (1598) If you know what you are doing you can adjust the lincs warning threshold in your mdp file I've never seen this before I'm using 1.o cutoff for pme and 1.4 for vdv my LINKS parameters are ; Bond parameters continuation= no; first dynamics run constraint_algorithm = lincs; holonomic constraints constraints= all-bonds; all bonds (even heavy atom-H bonds) constrained lincs_iter= 1; accuracy of LINCS lincs_order= 4; also related to accuracy How I could solve it? James 2012/2/14 James Starlight jmsstarli...@gmail.com Mark, I've checked only density value with 500 molecules Ccl4 I have density that is twisely less that I need ( in accordance to the literature ). Also I've checked my box visually and found that the box is not properly tightly packed so I dont know why genbox didnt add some extra mollecules :( In other words I wounder to know if there is any way to add some extra molecules to the pre defined box to make my system more tighly packed ( to short distance between existing molecules and place new ones in the new space ) ? James 2012/2/14 Mark Abraham mark.abra...@anu.edu.au On 14/02/2012 4:57 PM, James Starlight wrote: Justin, Firstly I've created the box of desired size with only 500 molecules ( I need 1000) Than I've tried to add extra 200 molecules by means of Genbox genbox -cp super_box.gro -ci Ccl4.gro -nmol 200 -o new_solv.gro but no molecules have been added Added 0 molecules (out of 200 requested) of Cl4 ... then there are no gaps large enough to insert your molecules. Either make gaps, or check out genbox -h for advice on defining the radii. also I've tried genbox -cp super_box.gro -cp Ccl4.gro -nmol 200 -o new_solv.gro Two -cp options is not what you want, and -nmol probably only works with -ci. but system were crashed with message Reading solute configuration God Rules Over Mankind, Animals, Cosmos and Such Containing 2500 atoms in 500 residues Initialising van der waals distances... WARNING: masses and atomic (Van der Waals) radii will be determined based on residue and atom names. These numbers can deviate from the correct mass and radius of the atom type. Reading solvent configuration God Rules Over Mankind, Animals, Cosmos and Such solvent configuration contains 5 atoms in 1 residues Is there any ways to add extra mollecules to the pre defined box ? Yes - but there has to be room for them. Mark -- 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 list
[gmx-users] Adding new residues and pdb2gmx
Hi. I've downloaded the charmm36 (gromacs-charmm36.ff_4.5.4.tgz) and added some residues to it by editing the lipids.rtp file. I also plan to use pdb2gmx to convert CHARMM-generated PDB files to GRO format. The Gromacs manual suggests one should create a residuetypes.dat file in the parent directory containing all the residues present in the RTP files. After doing so I used CHARMM (using CHARMM36 lipid forcefield) to generate the PDB of the lipid LPPC and run: pdb2gmx -f lppc.pdb -water none -noter -ff charmm36 -v. pdb2gmx appears to hang (100 % usage of one computer core) at the last line of the following messages: Option Filename Type Description -f lppc.pdb InputStructure file: gro g96 pdb tpr etc. -o conf.gro Output Structure file: gro g96 pdb etc. -p topol.top Output Topology file -i posre.itp Output Include file for topology -n clean.ndx Output, Opt. Index file -q clean.pdb Output, Opt. Structure file: gro g96 pdb etc. Option Type Value Description -- -[no]h bool no Print help info and quit -[no]version bool no Print version info and quit -niceint0 Set the nicelevel -chainsepenum id_or_ter Condition in PDB files when a new chain and molecule_type should be started: id_or_ter, id_and_ter, ter, id or interactive -ff string charmm36 Force field, interactive by default. Use -h for information. -water enum noneWater model to use: select, none, spc, spce, tip3p, tip4p or tip5p -[no]inter bool no Set the next 8 options to interactive -[no]ss bool no Interactive SS bridge selection -[no]ter bool no Interactive termini selection, iso charged -[no]lys bool no Interactive lysine selection, iso charged -[no]arg bool no Interactive arginine selection, iso charged -[no]asp bool no Interactive aspartic Acid selection, iso charged -[no]glu bool no Interactive glutamic Acid selection, iso charged -[no]gln bool no Interactive glutamine selection, iso neutral -[no]his bool no Interactive histidine selection, iso checking H-bonds -angle real 135 Minimum hydrogen-donor-acceptor angle for a H-bond (degrees) -distreal 0.3 Maximum donor-acceptor distance for a H-bond (nm) -[no]una bool no Select aromatic rings with united CH atoms on phenylalanine, tryptophane and tyrosine -[no]ignhbool no Ignore hydrogen atoms that are in the coordinate file -[no]missing bool no Continue when atoms are missing, dangerous -[no]v bool yes Be slightly more verbose in messages -posrefc real 1000Force constant for position restraints -vsite enum noneConvert atoms to virtual sites: none, hydrogens or aromatics -[no]heavyh bool no Make hydrogen atoms heavy -[no]deuterate bool no Change the mass of hydrogens to 2 amu -[no]chargegrp bool yes Use charge groups in the .rtp file -[no]cmapbool yes Use cmap torsions (if enabled in the .rtp file) -[no]renum bool no Renumber the residues consecutively in the output -[no]rtpres bool no Use .rtp entry names as residue names Using the Charmm36 force field in directory ./charmm36.ff Opening force field file ./charmm36.ff/aminoacids.r2b Opening force field file ./charmm36.ff/rna.r2b Reading lppc.pdb... Read 70 atoms Analyzing pdb file Splitting PDB chains based on TER records or changing chain id. There are 1 chains and 0 blocks of water and 1 residues with 70 atoms chain #res #atoms 1 ' ' 1 70 All occupancies are one Opening force field file ./charmm36.ff/atomtypes.atp Atomtype 1 (after that I abort the command with CTRL+C) - - - The LPPC residue is present in both the ./residuetypes.dat and ./charmm36/lipids.rtp files. What's wrong? I would really like to increase the verbosity of the pdb2gmx command but it seems the -v switch has no/little effect. Thanks in advance, Jernej Zidar -- 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] How to use Buckingham potentials ?
Dear Gromacs users, I am planing to use buckingham potential for the non-bonded interactions of my system. I know that by changing the nbfunc to 2 in [ defaults ] directive of topology will allow to use the Buckingham potential , But I don't know how to specify the A, B,C values in the ffnonbonded.itp file. what values should I specify in the [ atomtypes ] directive . As i know only A,B,C values for the pairs of atoms so i have tried by specifying A,B,C values both in [ atomtypes ] and [ pairtypes ] as [ atomtypes ] ;name at.num mass chargeptype A B C CC32A 6 12.011000.1030A 132277.5784 33.057851240.002710395 HCA2 1 1.00800 0.0355 A 29962.4608 41.580041580.000212547 OC30A 8 15.999400 -0.3480 A 243922.5976 40.2414486921 0.000803746 [ pairtypes ] ; i j func A BC OC30A OC30A 2 243922.597640.2414486921 0.000803746 OC30A CC32A 2 179625.814436.297640653 0.001476115 OC30A HCA22 85489.9984 40.899795501 0.000413379 CC32A CC32A 2 132277.578433.057851240.002710395 CC32A HCA22 62955.3928 36.832412523 0.000759396 HCA2HCA22 29962.4608 41.580041580.000212547 ;### I have given the A(in KJ/mol), B(in nm^-1), C(in Kj/mol*nm^6) ;D = A, B = 1/P, C = E Doing so , grompp is showing the following errors ERROR 5 [file ffnonbondedpeo.itp, line 13]: Too many parameters or not enough parameters for topology B ERROR 6 [file ffnonbondedpeo.itp, line 14]: Too many parameters or not enough parameters for topology B Generated 6 of the 6 non-bonded parameter combinations Excluding 3 bonded neighbours molecule type 'polymer' ERROR 7 [file topol.top, line 1652]: ERROR: The cut-off length is longer than half the shortest box vector or longer than the smallest box diagonal element. Increase the box size or decrease rlist. -- Program grompp, VERSION 4.5.5 Source code file: grompp.c, line: 1372 Fatal error: There were 7 errors in input file --- Any help will be highly appreciated. Regards, Ramesh Cheerla -- 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] Npt equilibration of the membrane-mimicking CCl4 layer
On 15/02/2012 4:45 PM, James Starlight wrote: Mark, due to hight density the volume of my system have been slightly increased and during NPT phase I've obtained error Fatal error: One of the box vectors has become shorter than twice the cut-off length or box_yy-|box_zy| or box_zz has become smaller than the cut-off. I'm using 0.9 for electrostatic and 1.4 for vdw cutofs and the dimensions of my box was 6.5 3 3 on the initial step and 6.6 3 3.3 before crush :) I want prevent such expansion of my system by increasing of pressure and/ or compressibility but I have not found exact sollution yet. Your system is dangerously small for those cut-offs if your initial density is not correct for your model physics. Your y and z dimensions only just contain a full cut-off sphere. You should also make sure you are following the advice about choice of P-coupling algorithm in manual 3.4.9, and consider using a very small integration time step. I remain unconvinced by this thread that you have generated a starting configuration that does not have atomic clashes. Mark James 2012/2/14 Mark Abraham mark.abra...@anu.edu.au mailto:mark.abra...@anu.edu.au On 14/02/2012 11:01 PM, James Starlight wrote: This also was solved by the some extra minimisation steps. I've forced with another problem :D During npt equilibration my system have slightly expanded so my desired volume and density were perturbed. I've noticed the below options in npt wich could help me ref_p= 1 1 compressibility = 4.5e-5 i'm using this compressibility value because I'm modelling the lipid-like environment so I think that I must increase pressure. Could you remind me the dependence of pressure from density and volume for liquids ? :) Your forcefield, simulation cell contents and .mdp settings will determine the equilibrium density. Whether you need to do anything depends on whether you've made a statistically significant post-equilibration measurement of your average density. Haphazardly increasing the reference pressure for the coupling will reduce the volume, but now you are simulating at that pressure. See http://www.gromacs.org/Documentation/Terminology/Pressure for background info. Mark James 2012/2/14 James Starlight jmsstarli...@gmail.com mailto:jmsstarli...@gmail.com It seems that I've fixed that problem by reduce vdv radii for Cl during defining of my box Eventually I've obtained box with the desired density than I've delete vdvradii.dat for my wor dir by when I've launched equilibration I've oibtained Fatal error: Too many LINCS warnings (1598) If you know what you are doing you can adjust the lincs warning threshold in your mdp file I've never seen this before I'm using 1.o cutoff for pme and 1.4 for vdv my LINKS parameters are ; Bond parameters continuation= no; first dynamics run constraint_algorithm = lincs; holonomic constraints constraints= all-bonds; all bonds (even heavy atom-H bonds) constrained lincs_iter= 1; accuracy of LINCS lincs_order= 4; also related to accuracy How I could solve it? James 2012/2/14 James Starlight jmsstarli...@gmail.com mailto:jmsstarli...@gmail.com Mark, I've checked only density value with 500 molecules Ccl4 I have density that is twisely less that I need ( in accordance to the literature ). Also I've checked my box visually and found that the box is not properly tightly packed so I dont know why genbox didnt add some extra mollecules :( In other words I wounder to know if there is any way to add some extra molecules to the pre defined box to make my system more tighly packed ( to short distance between existing molecules and place new ones in the new space ) ? James 2012/2/14 Mark Abraham mark.abra...@anu.edu.au mailto:mark.abra...@anu.edu.au On 14/02/2012 4:57 PM, James Starlight wrote: Justin, Firstly I've created the box of desired size with only 500 molecules ( I need 1000) Than I've tried to add extra 200 molecules by means of Genbox genbox -cp super_box.gro -ci Ccl4.gro -nmol 200 -o new_solv.gro but no molecules have been added Added 0 molecules (out of 200 requested) of Cl4 ... then there are no gaps large enough to insert your molecules. Either make gaps, or check out genbox -h for
Re: [gmx-users] How to use Buckingham potentials ?
On 15/02/2012 5:42 PM, ramesh cheerla wrote: Dear Gromacs users, I am planing to use buckingham potential for the non-bonded interactions of my system. I know that by changing the nbfunc to 2 in [ defaults ] directive of topology will allow to use the Buckingham potential , But I don't know how to specify the A, B,C values in the ffnonbonded.itp file. With [ nonbond_params], not [atomtypes] or [pairtypes]. See table 5.4 of section 5.7.1 of manual. Mark what values should I specify in the [ atomtypes ] directive . As i know only A,B,C values for the pairs of atoms so i have tried by specifying A,B,C values both in [ atomtypes ] and [ pairtypes ] as [ atomtypes ] ;name at.num mass chargeptype A B C CC32A 6 12.011000.1030A 132277.5784 33.057851240.002710395 HCA2 1 1.00800 0.0355 A 29962.4608 41.580041580.000212547 OC30A 8 15.999400 -0.3480 A 243922.5976 40.2414486921 0.000803746 [ pairtypes ] ; i j func A BC OC30A OC30A 2 243922.597640.2414486921 0.000803746 OC30A CC32A 2 179625.814436.297640653 0.001476115 OC30A HCA22 85489.9984 40.899795501 0.000413379 CC32A CC32A 2 132277.578433.057851240.002710395 CC32A HCA22 62955.3928 36.832412523 0.000759396 HCA2HCA22 29962.4608 41.580041580.000212547 ;### I have given the A(in KJ/mol), B(in nm^-1), C(in Kj/mol*nm^6) ;D = A, B = 1/P, C = E Doing so , grompp is showing the following errors ERROR 5 [file ffnonbondedpeo.itp, line 13]: Too many parameters or not enough parameters for topology B ERROR 6 [file ffnonbondedpeo.itp, line 14]: Too many parameters or not enough parameters for topology B Generated 6 of the 6 non-bonded parameter combinations Excluding 3 bonded neighbours molecule type 'polymer' ERROR 7 [file topol.top, line 1652]: ERROR: The cut-off length is longer than half the shortest box vector or longer than the smallest box diagonal element. Increase the box size or decrease rlist. -- Program grompp, VERSION 4.5.5 Source code file: grompp.c, line: 1372 Fatal error: There were 7 errors in input file --- Any help will be highly appreciated. Regards, Ramesh Cheerla -- 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] Adding new residues and pdb2gmx
On 15/02/2012 5:00 PM, Jernej Zidar wrote: Hi. I've downloaded the charmm36 (gromacs-charmm36.ff_4.5.4.tgz) and added some residues to it by editing the lipids.rtp file. I also plan to use pdb2gmx to convert CHARMM-generated PDB files to GRO format. The Gromacs manual suggests one should create a residuetypes.dat file in the parent directory containing all the residues present in the RTP files. After doing so I used CHARMM (using CHARMM36 lipid forcefield) to generate the PDB of the lipid LPPC and run: pdb2gmx -f lppc.pdb -water none -noter -ff charmm36 -v. pdb2gmx appears to hang (100 % usage of one computer core) at the last line of the following messages: Option Filename Type Description -f lppc.pdb InputStructure file: gro g96 pdb tpr etc. -o conf.gro Output Structure file: gro g96 pdb etc. -p topol.top Output Topology file -i posre.itp Output Include file for topology -n clean.ndx Output, Opt. Index file -q clean.pdb Output, Opt. Structure file: gro g96 pdb etc. Option Type Value Description -- -[no]h bool no Print help info and quit -[no]version bool no Print version info and quit -niceint0 Set the nicelevel -chainsepenum id_or_ter Condition in PDB files when a new chain and molecule_type should be started: id_or_ter, id_and_ter, ter, id or interactive -ff string charmm36 Force field, interactive by default. Use -h for information. -water enum noneWater model to use: select, none, spc, spce, tip3p, tip4p or tip5p -[no]inter bool no Set the next 8 options to interactive -[no]ss bool no Interactive SS bridge selection -[no]ter bool no Interactive termini selection, iso charged -[no]lys bool no Interactive lysine selection, iso charged -[no]arg bool no Interactive arginine selection, iso charged -[no]asp bool no Interactive aspartic Acid selection, iso charged -[no]glu bool no Interactive glutamic Acid selection, iso charged -[no]gln bool no Interactive glutamine selection, iso neutral -[no]his bool no Interactive histidine selection, iso checking H-bonds -angle real 135 Minimum hydrogen-donor-acceptor angle for a H-bond (degrees) -distreal 0.3 Maximum donor-acceptor distance for a H-bond (nm) -[no]una bool no Select aromatic rings with united CH atoms on phenylalanine, tryptophane and tyrosine -[no]ignhbool no Ignore hydrogen atoms that are in the coordinate file -[no]missing bool no Continue when atoms are missing, dangerous -[no]v bool yes Be slightly more verbose in messages -posrefc real 1000Force constant for position restraints -vsite enum noneConvert atoms to virtual sites: none, hydrogens or aromatics -[no]heavyh bool no Make hydrogen atoms heavy -[no]deuterate bool no Change the mass of hydrogens to 2 amu -[no]chargegrp bool yes Use charge groups in the .rtp file -[no]cmapbool yes Use cmap torsions (if enabled in the .rtp file) -[no]renum bool no Renumber the residues consecutively in the output -[no]rtpres bool no Use .rtp entry names as residue names Using the Charmm36 force field in directory ./charmm36.ff Opening force field file ./charmm36.ff/aminoacids.r2b Opening force field file ./charmm36.ff/rna.r2b Reading lppc.pdb... Read 70 atoms Analyzing pdb file Splitting PDB chains based on TER records or changing chain id. There are 1 chains and 0 blocks of water and 1 residues with 70 atoms chain #res #atoms 1 ' ' 1 70 All occupancies are one Opening force field file ./charmm36.ff/atomtypes.atp Atomtype 1 (after that I abort the command with CTRL+C) - - - The LPPC residue is present in both the ./residuetypes.dat and ./charmm36/lipids.rtp files. What's wrong? I would really like to increase the verbosity of the pdb2gmx command but it seems the -v switch has no/little effect. The most likely explanation is that you've not followed a file format somehow. A blank line at the end of a file you've modified might be required. If you've edited using some broken editor, dos2unix might be your friend. Check the original charmm36 forcefield works on some simple protein case. Then check that case works with your modified forcefield. Mark -- gmx-users mailing listgmx-users@gromacs.org http://lists.gromacs.org/mailman/listinfo/gmx-users Please search the archive at
Re: [gmx-users] Question regarding freeze groups and their use
On 15/02/2012 5:09 AM, J. Nathan Scott wrote: Hi all, I know there has been a *lot* of discussion on the mailing list on using freezegrps and potential pitfalls, but after having read much of this discussion the past couple of days and the manual I still find myself with questions and problems I'm hoping someone can help with. What I want to do, specifically, is have my protein frozen (*not* just constrained with some large force, I need actual freezing of every protein atom) but have the solvent still react to the frozen protein. The purpose behind this is to compare the water behavior near a particular part of a protein when it is dynamic and when it is not and thereby learn something about the communication between protein and solvent in this case. However, it seems that I'm not finding a rational set of mdp options to allow this, at least not a set for which dynamics will run with any speed. I have tried many of the things people have discussed on the mailing list, including turning off constraints, turning off pressure coupling, reducing the Protein heat bath temperature to 0 K, and using energygrp_excl. Here are the parameters I just tried: integrator = md nsteps = 100 dt = 0.001 nstxout = 1 nstvout = 1 nstxtcout = 1 nstenergy = 1 nstlog = 1 continuation= yes constraints = none dt = 0.0005 is normally used when there are no constraints. ns_type = grid nstlist = 5 rlist = 1.0 rcoulomb= 1.0 rvdw= 1.0 coulombtype = PME pme_order = 4 fourierspacing = 0.16 tcoupl = V-rescale tc-grps = Protein Non-Protein tau_t = 0.1 0.1 ref_t = 0 300 pcoupl = no pbc = xyz DispCorr= EnerPres gen_vel = no freezegrps = Protein freezedim = Y Y Y ;;; This ran, but *incredibly* slowly (I'm using a single 32 processor node, and only got 255 steps in 10 wallclock minutes). From the discussions I've read on the mailing list it seems that PME is somewhat less than ideal for simulations involving frozen atoms, but I used PME for electrostatics in my non-frozen simulations and would like to keep that consistent if possible (I do get a grompp warning of course). I am also getting warnings in the log file like: DD step 254 vol min/aver 0.183! load imb.: force 2395.3% pme mesh/force 1.243 which I believe must be due to the frozen protein causing difficulties with domain decomposition. Can anyone offer any advice regarding these issues? I know this has been discussed often, but nothing I'm finding in the archives is particularly relevant. Allowing a run to complete and inspecting the timing breakdown at the end of the .log file might be instructive. Manual 7.3.24 suggests the use of energy group exclusions. Mark -- 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] pdo making
Dear all I am caculating delta G by umbrella sampling with position pull in gromacs 4.5.after doing all of the steps when I want to use g_wham: g_wham -it tpr_files.dat -if pullf_files.dat -o -hist I have this error: found pull geometry position and more than 1 pull dimension(3) and it wants pdo files to do this analysis but I can not find any option -pd in this version of gromacs in mdrun command. What can I do to make pdo file in this version of gromacs? Thank you in progress Parto Haghighi -- 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] what does the Coulomb (SR) mean under PME?
On 15/02/2012 2:46 AM, Qiao Baofu wrote: Hi all, Anybody can help me out! Thanks in advance? In my test simulation, there are one Na^+ and one Cl^- (the distance of 1nm) in vacuum in a very big simulation box (10*10*10 nm^3). I calculated the energies under two different conditions, coulombtype=PME vs coulombtype=cut-off. 1. coulombtype=PME Coul. recip. -139.191 -- 0 0 (kJ/mol) Coul-SR:NA-CL *-0.0322805* -- 0 0 (kJ/mol) LJ-SR:NA-CL -0.00100285 -- 0 0 (kJ/mol) 2. coulombtype=cut-off Coul-SR:NA-CL *-138.935* -- 0 0 (kJ/mol) LJ-SR:NA-CL -0.00100285 -- 0 0 (kJ/mol) I did some calculations. Under coulombtype=cut-off, E(Coul-SR) = f*e_i*e_j/r with f=138.935, which is what I expect. However under coulombtype=PME, how does Gromacs get the value of -0.0322805kJ/mol for *SR* term? If I understand it correctly, SR means short-range, which is within the cut-off distance. Since the cut-off distance (1.2nm here) is much smaller than the simulation box length (10nm), there is not short-range interaction with the image box. Your two atoms are within the cut-off distance, so their interaction contributes to Coul-SR:NA-CL. The total energy for that interaction is distributed over the Coul-SR and Coul. recip terms (which is how PME works...), which you will note is approximately equal to Coul-SR for coulombtype cut-off. By construction, your periodic images are too far away to contribute to this term in either case (and you'd have been prevented by grompp from getting this far, if that were not true). Some details: Gromacs 4.5.5 PBC is used. distance of Na^+ and Cl^-: 1nm Na^+ and Cl^- are frozen in XYZ dimensions. Only 1 step is run. Using mdrun -rerun is often the best way to calculate single-point energies/forces. Mark -- 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
