[gmx-users] Conflicting atom names when using g_density
Hi Everyone, I am having an awkward situation. I want to use g_density to get the density profile of a lipid bilayer. At first glance, it looks like a trivial task. But for my case, it's not because I have conflicting atom names from two different molecules. I am using GROMOS 53a6. In DPPC topology, I have C1 C2 C3 N4...O7... In another lipid topology, say lipid B, I have C1 C2 C3 N4... C7... The C1 in DPPC and C1 in lipid B could have different atomic mass as GROMOS is united atom representation. C1 of DPPC could be a CH2, C2 of lipid B could be a CH3. Then I don't know how to specify the entry in the .dat file: C1 = ? Is it possible to calculate the density profile of the whole bilayer by just invoking g_density once? Or in my awkward situation, I have to calculate the density profile of each molecule type, then add them together? Then I will be concerned with the numerical errors involved. Thanks. Bin -- 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] Molecule removal by g_membed
Hi All, I am simulating a system of lipopeptides embedded into a lipid bilayer. I have a question about how g_membed works. My bilayer has three kinds of lipids. Let's call them lipid_1, lipid_2, lipid_3. Sodium ions are used to balance the charges. When I invoke g_membed, it will ask for which group to embed in the membrane. In my case, obviously it is the group Lipopeptide. Then g_membed asks for a group to embed Lipopeptide into. Should it be the joint group of lipid_1, lipid_2 and lipid_3? Or something else? I used make_ndx to generate an index group called Membrane by merging group Lipid_1, Lipid_2 and Lipid_3. It seems it works with g_membed. It gives the following information: Will remove 2 Lipid_1 molecules Will remove 1 Lipid_2 molecules Will remove 3 Lipid_3 molecules Will remove 7 SOL molecules Will remove 0 NA molecules Will remove 0 Lipopeptide molecules Apparently the behaviour of g_membed is what I want. But what concerns me is I don't know whether the removal of 0 Sodium ions is just a coincidence or it's a deterministic action by g_membed. Of course I don't want the removal of ions as it would cause a charged system. So my question is under what conditions will g_membed remove ions? Thank you in advance. Regards Bin -- 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: Use pull code to restrain the COM
Hi Justin, Thank you for your quick reply. In Setting A, grompp reports Pull group natoms pbc atom distance at start reference at t=0 0 0 0 1 14623185-0.000 -0.000 -3.268 6.668 6.693 6.492 The values of reference at t=0 are from the .mdp file. They are the COM of pull group 1 calculated by g_traj -ox -com. But I can not understand the values of distance at start. Where does -3.268 come from? If pull group 0, the reference group is set to be 6.668 6.693 6.492, and grompp can calculate the COM of pull group 1, the distance should be 0.0 0.0 0.0 instead of -0.000 -0.000 -3.268. In Setting B, grompp reports Pull group natoms pbc atomdistance at start reference at t=0 0 0 0 1 14623185 -0.000 -0.000 -3.268 0.000 0.000 -3.268 It seems Setting B doesn't agree with Setting A. Then pull-start=yes doesn't have the effect as I guessed. My situation is a bit complicated and challenging. I am simulating the oligomerization of lipo-peptides inside a lipid bilayer. There are several lipo-peptides in the system. Initially I placed them in the water, i.e., outside the bilayer. Then I found it takes too long (at the scale of microseconds even seconds) for them to spontaneously insert into the bilayer and form an oligomer. Then I artificially placed the lipopeptides into the bilayer, and put them together to create an artificial oligomer. I want to see if they can stay there. But obviously the artificial insertion and oligomerization is unfavourable from the viewpoint of energy. So the lipopeptides have a strong tendency to leave the bilayer. Now I want to restrain the COM of each lipopeptide (at least in z direction) for some time, but give them conformational and orientational degrees of freedom as much as possible, to give rise to the possibility that the artificial oligomer can find a favourable conformation and stay there after the removal of restraints. Unfortunately GROMACS can not restrain the COMs of several groups simultaneously. Then in one short simulation (say 100ps), I can only restrain the COM of one lipopeptide and freeze the others (in z direction). The whole process proceeds as: MD (restrain COM of peptide A, freeze BCD) - CG EM - MD (restrain COM of B, freeze ACD) - CG EM - MD (restrain COM of C, freeze ABD) - CG EM - MD (restrain COM of D, freeze ABC) - Next cyle - From my current experience, the EM step is necessary as in MD steps the freezed lipopeptides developed a large amount of strain. If the strain can not get released, the MD runs will explode in several hundreds of picoseconds (GROMACS reports failures of LINCS and eventually segmentation fault, an obvious indication that the system has exploded.) And pull-k1= 50 is indeed very small. The restrained lipopeptides can displace in z direction by a quite large amount in only hundreds of picoseconds. Could you suggest a different approach, or some suggestions to refine this process? Thank you very much. Best Regards Bin 2013/8/9 gmx-users-requ...@gromacs.org Send gmx-users mailing list submissions to gmx-users@gromacs.org To subscribe or unsubscribe via the World Wide Web, visit http://lists.gromacs.org/mailman/listinfo/gmx-users or, via email, send a message with subject or body 'help' to gmx-users-requ...@gromacs.org You can reach the person managing the list at gmx-users-ow...@gromacs.org When replying, please edit your Subject line so it is more specific than Re: Contents of gmx-users digest... Today's Topics: 1. Aw: [gmx-users] Umbrella sampling - position restraints (lloyd riggs) 2. Use pull code to restrain the COM (Bin Liu) 3. Re: Use pull code to restrain the COM (Justin Lemkul) 4. Re: Trying to replicate Aqvist's results (solvation free energy). (Andr? Farias de Moura) 5. Re: Trying to replicate Aqvist's results (solvation free energy). (Heymman) 6. g_wham -sym (Shima Arasteh) -- Message: 1 Date: Thu, 8 Aug 2013 23:19:59 +0200 (CEST) From: lloyd riggs lloyd.ri...@gmx.ch Subject: Aw: [gmx-users] Umbrella sampling - position restraints To: Discussion list for GROMACS users gmx-users@gromacs.org Message-ID: trinity-790b128a-a41d-4a6e-87ac-908747401fc1-1375996798996@3capp-gmx-bs33 Content-Type: text/plain; charset=us-ascii An HTML attachment was scrubbed... URL: http://lists.gromacs.org/pipermail/gmx-users/attachments/20130808/fe9e3133/attachment-0001.html -- Message: 2 Date: Thu, 8 Aug 2013 18:43:09 -0400 From: Bin Liu fdusuperstr...@gmail.com Subject: [gmx-users] Use pull code to restrain the COM To: gmx-users@gromacs.org Message-ID: caha8nljjpyduptamu50tzxkhtwd+gkpur87qkcnzwj3v37s...@mail.gmail.com Content-Type: text
[gmx-users] Use pull code to restrain the COM
Hi All, I would like to restrain the COM of a molecule, say a protein, in my simulation. I found the pull code can do the job. However, I am not sure about several parameters in the .mdp file. For example, if I want to restrain the COM of the protein in only z direction, but not rigidly, I can specify the following parameters: *; Setting A* *pull= umbrella* *pull-geometry = position* *pull-dim= N N Y ; only applied to Z direction* *pull-ngroups= 1 * *pull-group1 = Protein* *pull-start = no* *pull-init1 = 3.0 3.0 3.0 ; suppose the initial position of the COM is 3.0 3.0 3.0* *pull-k1 = 50* I am not confident about my understanding of the settings. Here I didn't specify *pull-group0*. Therefore the position of the reference point is 0.0 0.0 0.0. But *pull-init1* moves the reference point to *3.0 3.0 3.0. *I don't know the reasonable range for *pull-k1 *too. I want the COM can only have a small deviation from its initial position. Is 50 too large or too small? And I figured out another approach. *; Setting B* *pull= umbrella* *pull-geometry = position* *pull-dim= N N Y ; only applied to Z direction* *pull-ngroups= 1 * *pull-group1 = Protein* *pull-start = yes* *pull-k1 = 50* * * The GROMACS manual says *pull-start=yes *adds the COM distance of the starting conformation to *pull-init. *Does it mean Setting B is equivalent to Setting A? If I set *pull-geometry= position, *does *pull-start=yes *add the initial position of the COM to *pull-init* as a vector? If it's true, then I don't need to calculate the initial position of the COM. Thanks a lot. Best Bin -- 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: gmx-users Digest, Vol 111, Issue 53
Hi Trayder, Thank you for the suggestion. I will try the plugin soon. Yes, I can play the two trajectories simultaneously, one with protein disabled. In this way, I indeed can get a very good presentation of what is happening. However, it's quite inconvenient to set Graphical Representations in VMD to two trajectories every time. And having two trajectories to set complicates the Tcl script to automate the visualization. That's why I ask people how to apply -nojump to a part of a system. Cheers Bin Message: 2 Date: Thu, 11 Jul 2013 11:51:51 +1000 From: Trayder Thomas trayder.tho...@monash.edu Subject: Re: [gmx-users] How to apply trjconv -nojump to a part of a system To: Discussion list for GROMACS users gmx-users@gromacs.org Message-ID: CAHSz8K7Tc02KXqr6yHQ=c= rssocu5gdgxbkrpbk-+-h93yw...@mail.gmail.com Content-Type: text/plain; charset=ISO-8859-1 VMD might do what you want with the PBC tools plugin (installed by default). http://www.ks.uiuc.edu/Research/vmd/plugins/pbctools/ unwrap being the equivalent of -nojump Otherwise, couldn't you just view your 2 trajectories simultaneously, one with protein the other not? -Trayder On Wed, Jul 10, 2013 at 9:13 AM, Bin Liu fdusuperstr...@gmail.com wrote: Hi All, For the convenience of visualization, I need to remove the jump of one component (say a protein) of the system at the boundary. I don't need to, or say I need not to remove the jump of the other components (say a lipid bilayer), since otherwise the system will look falling apart. I noticed I can cluster a part of a system, then output all the atoms in the system in which only the part is clustered, and the other components unchanged. Does GROMACS have similar function when *-nojump* is used? If this can not be accomplished directly, is there a way to circumvent it? I figured out a way, but haven't implemented it. I can plug the coordinates of the protein treated with *-nojump* into the trajectory of the whole system which is not treated with *-nojump*. It is kind of substituting the coordinates of one component in one trajectory with the coordinates of the same component in another trajectory. Is there anyone aware of a tool or a script to do this job? Thank you very much. Regards Bin -- 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] AVX2 SIMD intrinsics speed boost
Hi all, If my understanding is correct, GROMACS parallelization and acceleration page indicates AVX2 SIMD intrinsics can offer a speed boost on a Haswell CPU. I was wondering how much performance gain we can expect from it. In another word, what's the approximate speed increase if we run a simulation with AVX2 SIMD intrinsics on a Haswell CPU (say i7 4770K) than on an Ivy Bridge CPU of the same clock (say i7 3770K) with the current AVX SIMD intrinsics? And is there a timeline for the release of AVX2 SIMD intrinsics? This information is crucial if we want to assemble a machine with balanced CPU and GPU performance. My current machine has i7 3770K (3.5GHz, stock frequency) and Geforce 650 Ti (768 CUDA cores, 1032MHz). When I ran simulations with rcoulomb=1.0 and rvdw=1.0, I got this at the end of the log file: *Force evaluation time GPU/CPU: 1.762 ms/1.150 ms = 1.531* * * It seems I need a GPU with 50% more CUDA cores. In the best scenario, If AVX2 can give 30% speed boost, and I can successfully overclock 4770K to 4.5GHz, I need 1900 CUDA cores( 130%*(4.5GHz/3.5GHz)*1.531*768 cores) at the same frequency to get balanced CPU and GPU performance. Then I will need a GeForce GTX 780 (2304 CUDA cores at 863MHz, equivalent to 1925 CUDA cores at 1032MHz). Since GROMACS is highly insensitive to memory clock and latency, I hope this naive arithmetic can give a good estimation which graphic card I should purchase. Best Bin -- 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 apply trjconv -nojump to a part of a system
Hi All, For the convenience of visualization, I need to remove the jump of one component (say a protein) of the system at the boundary. I don't need to, or say I need not to remove the jump of the other components (say a lipid bilayer), since otherwise the system will look falling apart. I noticed I can cluster a part of a system, then output all the atoms in the system in which only the part is clustered, and the other components unchanged. Does GROMACS have similar function when *-nojump* is used? If this can not be accomplished directly, is there a way to circumvent it? I figured out a way, but haven't implemented it. I can plug the coordinates of the protein treated with *-nojump* into the trajectory of the whole system which is not treated with *-nojump*. It is kind of substituting the coordinates of one component in one trajectory with the coordinates of the same component in another trajectory. Is there anyone aware of a tool or a script to do this job? Thank you very much. Regards Bin -- 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] The effect of coulomb-modifier and vdw-modifier in Verlet cutoff scheme
Hi All, In GROMACS 4.6.x series, Verlet cutoff scheme is introduced to enable OpenMP parallelization and GPU acceleration. Then some new run parameters are introduced to control the use of Verlet cutoff scheme. However, I noticed the GROMACS manual doesn't give in-depth knowledge on some parameters. For example, I am sure whether and how to use the coulomb-modifier and vdw-modifier parameters with Verlet cutoff scheme. My shallow experience told me these two parameters won't affect the simulation much in term of speed and accuracy. I did a rough benchmark on a DMPC128 bilayer (323K) system to check the area per lipid against (traditional) group cutoff scheme results. Group cutoff scheme: 0.656 (0.008) nm^2 (Uncertainty) Verlet cutoff scheme (rcoulomb=rvdw=1.0, coulomb-modifier = Potential-shift-Verlet, vdw-modifier = Potential-shift-Verlet) 0.643nm^2 Verlet cutoff scheme (rcoulomb=rvdw=1.0, coulomb-modifier = None, vdw-modifier = None) 0.645nm^2 Basically the results from two Verlet cutoff parameter sets are indistinguishable. I am looking forward to your help to give me some insight into this question. Thank you very much. Regards Bin -- 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:The effect of coulomb-modifier and vdw-modifier in Verlet cutoff scheme
Dear Mark, Could you elaborate on your answer? In my group cutoff scheme, I used ns_type = grid ; search neighboring grid cels nstlist = 5 ; 10 fs rlist = 1.3 ; short-range neighborlist cutoff (in nm) rcoulomb= 1.3 ; short-range electrostatic cutoff (in nm) rvdw= 1.0 ; short-range van der Waals cutoff (in nm) vdwtype = Shift rvdw_switch = 0.9 What is the advantage of turning on coulomb-modifier and vdw-modifier in terms of physical or chemical accuracy of simulations? Thanks. Bin Those modifiers shift only the potential, as manual 7.3 points out. So the forces and sampling are unaffected, so it does not surprise me that APL is unaffected by the use of such a shift. If your group cutoff scheme was unbuffered (rlist = max(rcoulcomb,rvdw) and nstlist 1), then if the observed difference is significant, then that could be the reason. -- 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
