R: Re: [gmx-users] Rotation Constraints - PMF - external potential
Dear Carsten, thank you very much for your support! At the beginning I did not follow your indications correctly but now, using the orientation restraints and the indication in the paper Keep it Flexible: Driving Macromolecular Rotary Motions in Atomistic Simulations with GROMACS Kutzner C, Czub J, and Grubmüller H J. Chem Theory and Comp. 7: 1381-1393 (2011) I believe I made it. Just few questions to be sure, if I can: how the vector v ( namely the axis rotation defined in the paper) is defined? giving a (large) number of atoms as rotation group, I suppose that the vector u identify the center of mass of the group, and in this point the v vector defines the direction that one gives with the rot_vec entry in mdp file. Is it correct? I suppose that in order to just keep the direction of the v axis during the dynamics, my choise have to be rote_rate=0, is it correct? Thanks again, Anna Messaggio originale Da: ckut...@gwdg.de Data: 25/07/2013 11.12 A: battis...@libero.itbattis...@libero.it Cc: Discussion list for GROMACS usersgmx-users@gromacs.org Ogg: Re: [gmx-users] Rotation Constraints - PMF - external potential On Jul 24, 2013, at 5:53 PM, battis...@libero.it wrote: Dear Carsten could you give me more information about your suggestions? I tried but probably I did not understand well what you meant. Hi Anna, I suggested to use the enforced rotation module of Gromacs 4.6 to restrain the orientation of your molecule(s). If you want to use the orientation restraints module instead, I am afraid I can not help you much with that, maybe someone else on this list? In order to avoid the rotation of the structure A and of the structure B, I have defined into the index file a group A_B that contains A+B and I have setted in the mdp file the following parameters: ; Orientation restraints: No or Yes orire= yes ; Orientation restraints force constant and tau for time averaging orire-fc = 500 orire-tau= 100 orire-fitgrp = A_B ; Output frequency for trace(SD) and S to energy file nstorireout = 100 As I have synthetically described in the first post , the structures A and B (characterized by a cylindrical shape) are defined by a number of 32 unit- structures that I call s. Into the itp is defined the topology for the s structure, and so in order to put an orientation restraints between atoms that are not included into the same itp file, I cannot put into the topology a section like that described into the manual 4.6.2 pag. 92 namely, [ orientation_restraints ], could I ? Could you tell me How I can fix the orientation of the systems A and B? Using the enforced rotation module you would choose an index group and an axis for each group that you want to fix the orientation, set the rotation angle to zero and choose an appropriate force constant. Appropriate potential functions would be the pivot-free ones if I understand your setting correctly. I don't understand the manual's explanation about the orire-fitgrp: (fit group for orientation restraining. This group of atoms is used to determine the rotation R of the system with respect to the reference orientation. The reference orientation is the starting conformation of the first subsystem. For a protein, backbone is a reasonable choice) How one have to give the group? using an index file or defining the group in the topology? This is the orire-fitgrp = A_B mdp file setting that you made. Best, Carsten Messaggio originale Da: ckut...@gwdg.de Data: 23/07/2013 13.09 A: battis...@libero.itbattis...@libero.it, Discussion list for GROMACS usersgmx-users@gromacs.org Ogg: Re: [gmx-users] Rotation Constraints - PMF Hi Anna, please have a look at the Enforced Rotation Section in the Gromacs 4.6 manual. You can restrain the angle of rotation about an axis by setting the rotation rate to zero. There is also a 4.5 add-on available with rotational restraints in the Gromacs git repository (branch rotation). For more info you may want to look at this page: http://www.mpibpc.mpg.de/grubmueller/rotation Best, Carsten On Jul 23, 2013, at 12:18 PM, battis...@libero.it wrote: Dear user and expert, I'd like ask you a suggestion about a problem that I will try present you schematically. I have got a structure s and I have generated the topolgy file itp for it. A number of separate s in turn generate a complex structure A, that is characterized by a cylindrical shape. Now, I constructed a system with two cylindrical structures, A and B (in total made by 64 s structures), and I'd like make an Umbrella Sampling calculation in order to study the PMF varying the distance between A and B. My problem is that I'd like fix the orientation of the axis of each structure A and B long the z axis, during the dynamics. So I need to put a force into the
RE: [gmx-users] Calculate interaction energy dynamically
From: davhak2...@hotmail.com To: gmx-users@gromacs.org Subject: RE: [gmx-users] Calculate interaction energy dynamically Date: Thu, 25 Jul 2013 09:13:14 + Date: Wed, 24 Jul 2013 09:35:26 -0400 From: jalem...@vt.edu To: gmx-users@gromacs.org Subject: Re: [gmx-users] Calculate interaction energy dynamically On 7/24/13 9:25 AM, Davit Hakobyan wrote: Thank you again for your time and help. Performing rerun on the original system passes without warnings: The subsequent g_energy also runs smoothly. One gets no warning also if passing the ordered *ordered.tpr and the new index (*ordered.ndx) files to the mdrun -rerun of the original system. This indeed leaves with the ordered trajectory file as the only problematic point which might cause the warnings. Good to know. Interestingly, the ordered *.trr file (product of trjorder command) is nearly twice as small as the original *.trr file. Since the coordinates in the ordered *.trr file look fine then probably it is the velocities which might be missing in it. Can the There's no might about it - your gmxcheck output from before shows that velocities are definitely absent, which should account for the difference in file size. warnings come from this issue ? In this case the g_energy should not be able to calculate the kinetic energy. If this is the case may one still rely on the numbers obtained by g_energy for the coulomb and vdw terms (despite the warning messages) ? This is a very weird issue. The .trr file clearly has step and time information (from your earlier post), but the resulting .edr file from mdrun -rerun does not have the correct step number (though it appears to have time). Your problem does not derive from anything to do with velocities, as far as I can tell. Proceed carefully and scrutinize the output to see if it indeed makes sense - you can verify energy values from snapshots that have not been ordered by creating index files that correspond to the same lipids, for instance. This is a very non-standard exercise, so it's probably never been tested and debugged. The energy calculation between the nearest neighbors at a given time frame based on the original trajectory file (with help of new index file) indeed gives the same values as for the ordered trajectory file. So one can probably rely on the energy output from the ordered trajectory despite the warnings. In case of more information about the warnings I will write in this thread. It seems that the problem appears when for the trjorder command one indicates *.trr instead of *.xtc as an output. When *.xtc is mentioned like trjorder ... -o output.xtc both the mdrun -rerun and g_energy subsequently pass without warnings. Thank you very much for all the help. -Justin -- == Justin A. Lemkul, Ph.D. Postdoctoral Fellow Department of Pharmaceutical Sciences School of Pharmacy Health Sciences Facility II, Room 601 University of Maryland, Baltimore 20 Penn St. Baltimore, MD 21201 jalem...@outerbanks.umaryland.edu | (410) 706-7441 == -- 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 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] Rotation Constraints - PMF - external potential
Hi Anna, On Jul 26, 2013, at 11:17 AM, battis...@libero.it wrote: Dear Carsten, thank you very much for your support! At the beginning I did not follow your indications correctly but now, using the orientation restraints and the indication in the paper I think you do not need orientation restraints at all, just enforced rotation. Keep it Flexible: Driving Macromolecular Rotary Motions in Atomistic Simulations with GROMACS Kutzner C, Czub J, and Grubmüller H J. Chem Theory and Comp. 7: 1381-1393 (2011) I believe I made it. Just few questions to be sure, if I can: how the vector v ( namely the axis rotation defined in the paper) is defined? giving a (large) number of atoms as rotation group, I suppose that the vector u identify the center of mass of the group, and in this point the v vector defines the direction that one gives with the rot_vec entry in mdp file. Is it correct? Yes, the rotation vector v you provide in the mdp file as rot_vec0, and u as rot_pivot0. It depends a bit on the rotation potential type you choose, for some of the potentials (the pivot-free ones) you do not have to provide a u vector, since it is calculated on the fly as the center of positions or as the center of mass (depending on whether you choose mass-weighting or not with the rot_massw0 parameter). If you have more than one protein that you want to constrain the rotation for, set rot_ngroups to a value larger than 1 and provide values for each group in rot_vec0, rot_vec1, and so on. I suppose that in order to just keep the direction of the v axis during the dynamics, my choise have to be rote_rate=0, is it correct? Yes, and you have to choose a spring constant (k, or rot_k0) that determines how close your protein stays at the given angle, i.e. how large you allow the fluctuations about the angle to be. Good luck! Carsten Thanks again, Anna Messaggio originale Da: ckut...@gwdg.de Data: 25/07/2013 11.12 A: battis...@libero.itbattis...@libero.it Cc: Discussion list for GROMACS usersgmx-users@gromacs.org Ogg: Re: [gmx-users] Rotation Constraints - PMF - external potential On Jul 24, 2013, at 5:53 PM, battis...@libero.it wrote: Dear Carsten could you give me more information about your suggestions? I tried but probably I did not understand well what you meant. Hi Anna, I suggested to use the enforced rotation module of Gromacs 4.6 to restrain the orientation of your molecule(s). If you want to use the orientation restraints module instead, I am afraid I can not help you much with that, maybe someone else on this list? In order to avoid the rotation of the structure A and of the structure B, I have defined into the index file a group A_B that contains A+B and I have setted in the mdp file the following parameters: ; Orientation restraints: No or Yes orire= yes ; Orientation restraints force constant and tau for time averaging orire-fc = 500 orire-tau= 100 orire-fitgrp = A_B ; Output frequency for trace(SD) and S to energy file nstorireout = 100 As I have synthetically described in the first post , the structures A and B (characterized by a cylindrical shape) are defined by a number of 32 unit- structures that I call s. Into the itp is defined the topology for the s structure, and so in order to put an orientation restraints between atoms that are not included into the same itp file, I cannot put into the topology a section like that described into the manual 4.6.2 pag. 92 namely, [ orientation_restraints ], could I ? Could you tell me How I can fix the orientation of the systems A and B? Using the enforced rotation module you would choose an index group and an axis for each group that you want to fix the orientation, set the rotation angle to zero and choose an appropriate force constant. Appropriate potential functions would be the pivot-free ones if I understand your setting correctly. I don't understand the manual's explanation about the orire-fitgrp: (fit group for orientation restraining. This group of atoms is used to determine the rotation R of the system with respect to the reference orientation. The reference orientation is the starting conformation of the first subsystem. For a protein, backbone is a reasonable choice) How one have to give the group? using an index file or defining the group in the topology? This is the orire-fitgrp = A_B mdp file setting that you made. Best, Carsten Messaggio originale Da: ckut...@gwdg.de Data: 23/07/2013 13.09 A: battis...@libero.itbattis...@libero.it, Discussion list for GROMACS usersgmx-users@gromacs.org Ogg: Re: [gmx-users] Rotation Constraints - PMF Hi Anna, please have a look at the Enforced Rotation Section in the Gromacs 4.6 manual. You can restrain the angle of rotation about an axis
[gmx-users] Reminder about US GROMACS workshop + soliciting presenters for talks and tutorials
Dear gmx-users list members- I'd like to remind you all about the 2013 GROMACS Workshop and Conference at the University of Virginia Sept 13th-15th; registration costs will rise from $60 to $95 on Aug 22nd. The original invitation, with link to the website, is at the end of this email. We are issuing a last call for contributed 20 min talks, with decisions on the schedule made next week, and would like to invite you to apply to present your GROMACS-related research. Email gromacsworkshop-registrat...@virginia.edu if you wish to submit a talk for consideration; see the preliminary schedule on the conference webpage for more information on the design of the talks. We still also want to give a opportunity to advanced users who may be interested in either presenting specific tutorial subjects at the workshop or alternatively having an online tutorial hosted in the GROMACS website. See the list of tutorial topics on the workshop schedule, or propose your own. We want to try to leverage the energy and experience of the GROMACS community to enrich the tutorial experience, rather than have the same people up there the whole weekend! Sincerely, The 2013 GROMACS USA Workshop and Conference Steering Committee Michael Shirts (chair) Angel Garcia Berk Hess Yu-Shan Lin Erik Lindahl Peter Kasson Original announcement: ~~~ We are pleased to announce the 2013 GROMACS USA Workshop and Conference at the University of Virginia in Charlottesville, Virginia, on September 13-15th. This will be the first full GROMACS workshop held here in the U.S. The workshop and conference will consist of two days of tutorials, discussions of future plans for GROMACS, face time with a large percentage of the developers, plenary software and application sessions, and a last day of development planning to which attendees are also invited to help influence the future directions of GROMACS. Please visit http://faculty.virginia.edu/gromacsworkshop for registration and for much more information about the program. For specific questions about registering or logistics after visiting the website, please write to gromacsworkshop-registrat...@virginia.edu. Michael Shirts Assistant Professor Department of Chemical Engineering University of Virginia michael.shi...@virginia.edu (434)-243-1821 -- 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] Calculate interaction energy dynamically
On 7/26/13 5:21 AM, Davit Hakobyan wrote: From: davhak2...@hotmail.com To: gmx-users@gromacs.org Subject: RE: [gmx-users] Calculate interaction energy dynamically Date: Thu, 25 Jul 2013 09:13:14 + Date: Wed, 24 Jul 2013 09:35:26 -0400 From: jalem...@vt.edu To: gmx-users@gromacs.org Subject: Re: [gmx-users] Calculate interaction energy dynamically On 7/24/13 9:25 AM, Davit Hakobyan wrote: Thank you again for your time and help. Performing rerun on the original system passes without warnings: The subsequent g_energy also runs smoothly. One gets no warning also if passing the ordered *ordered.tpr and the new index (*ordered.ndx) files to the mdrun -rerun of the original system. This indeed leaves with the ordered trajectory file as the only problematic point which might cause the warnings. Good to know. Interestingly, the ordered *.trr file (product of trjorder command) is nearly twice as small as the original *.trr file. Since the coordinates in the ordered *.trr file look fine then probably it is the velocities which might be missing in it. Can the There's no might about it - your gmxcheck output from before shows that velocities are definitely absent, which should account for the difference in file size. warnings come from this issue ? In this case the g_energy should not be able to calculate the kinetic energy. If this is the case may one still rely on the numbers obtained by g_energy for the coulomb and vdw terms (despite the warning messages) ? This is a very weird issue. The .trr file clearly has step and time information (from your earlier post), but the resulting .edr file from mdrun -rerun does not have the correct step number (though it appears to have time). Your problem does not derive from anything to do with velocities, as far as I can tell. Proceed carefully and scrutinize the output to see if it indeed makes sense - you can verify energy values from snapshots that have not been ordered by creating index files that correspond to the same lipids, for instance. This is a very non-standard exercise, so it's probably never been tested and debugged. The energy calculation between the nearest neighbors at a given time frame based on the original trajectory file (with help of new index file) indeed gives the same values as for the ordered trajectory file. So one can probably rely on the energy output from the ordered trajectory despite the warnings. In case of more information about the warnings I will write in this thread. It seems that the problem appears when for the trjorder command one indicates *.trr instead of *.xtc as an output. When *.xtc is mentioned like trjorder ... -o output.xtc both the mdrun -rerun and g_energy subsequently pass without warnings. Interesting. That leads me to suspect that there is something wrong with the way trjorder writes .trr frames. It writes .xtc files correctly, but incorrect step information in .trr files. Thanks for following up. I'll look into it and probably file a bug report on it. What Gromacs version are you using? Maybe you said that before, but I've forgotten. -Justin -- == Justin A. Lemkul, Ph.D. Postdoctoral Fellow Department of Pharmaceutical Sciences School of Pharmacy Health Sciences Facility II, Room 601 University of Maryland, Baltimore 20 Penn St. Baltimore, MD 21201 jalem...@outerbanks.umaryland.edu | (410) 706-7441 == -- 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
Fw: [gmx-users] Umbrella Sampling _ pulled ion
Thanks for your previous replies. All are always beneficial and I appreciate you. As I see in pullx.xvg file, the third column refers to the distance of pull group and reference group. All around -1.2. Is it sufficient to judge that the US has been performed properly? Thanks for your suggestions in advance. Sincerely, Shima - Forwarded Message - From: Justin Lemkul jalem...@vt.edu To: Discussion list for GROMACS users gmx-users@gromacs.org Sent: Thursday, July 25, 2013 4:25 PM Subject: Re: [gmx-users] Umbrella Sampling _ pulled ion On 7/25/13 7:52 AM, Shima Arasteh wrote: After running for 100 ps, I visualized the pullf_umbrella.xvg, in this plot I found the F value around 100 kJ/mol/ns. But force constant which I had set in md_US.mdp file was 4000 KJ/mol/ns. Does this difference show me that the US has not done properly? Force and force constant are different. The values in the pullf.xvg file indicate the magnitude of force required to maintain the position of the restrained ion. The pullx.xvg will probably be more useful in directly determining how effective the restraining potential was. -Justin ;define = -DPOSRES integrator = md ; leap-frog integrator nsteps =10 ; 1 * 10 = 100 ps dt = 0.001 ; 1 fs nstcomm = 10 tinit = 0 ; Output control nstxout = 5000 ; save coordinates every 100 ps nstvout = 5000 ; save velocities every 100 ps nstenergy = 1000 ; save energies every 2 ps nstfout = 5000 nstxtcout = 5000 ; every 10 ps continuation = yes ; first dynamics run constraint_algorithm = lincs ; holonomic constraints constraints = h-bonds ; all bonds (even heavy atom-H bonds) constrained lincs_iter = 1 ; accuracy of LINCS lincs_order = 4 ; also related to accuracy ; Neighborsearching ns_type = grid ; search neighboring grid cells nstlist = 5 ; 10 fs rlist = 1.2 ; short-range neighborlist cutoff (in nm) rlistlong = 1.4 rcoulomb = 1.2 ; short-range electrostatic cutoff (in nm) rvdw = 1.2 ; short-range van der Waals cutoff (in nm) vdwtype = switch rvdw_switch = 1.0 ; Electrostatics coulombtype = PME ; Particle Mesh Ewald for long-range electrostatics pme_order = 4 ; cubic interpolation fourierspacing = 0.12 ; grid spacing for FFT fourier_nx = 0 fourier_ny = 0 fourier_nz = 0 ewald_rtol = 1e-5 optimize_fft = yes ; Temperature coupling is on tcoupl = Nose-Hoover ; modified Berendsen thermostat tc-grps = Protein_POPC Water_Ces_CL ; two coupling groups - more accurate tau_t = 0.5 0.5 ; time constant, in ps ref_t = 310 310 ; reference temperature, one for each group, in K ; Pressure coupling is on pcoupl = Parrinello-Rahman ; no pressure coupling in NVT pcoupltype = semiisotropic tau_p = 4 ref_p = 1.01325 1.01325 compressibility = 4.5e-5 4.5e-5 ; Periodic boundary conditions pbc = xyz ; 3-D PBC ; Long-range dispersion correction DispCorr = EnerPres ; Pull code pull = umbrella pull_geometry = position pull_dim = N N Y pull_start = yes pull_ngroups = 1 pull_group0 = POPC pull_group1 = Ces_ion pull_init1 = 0.0 0.0 0.0 pull_rate1 = 0.0 0.0 0.00 pull_vec1 = 0 0 1 pull_k1 = 4000 ; kJ mol^-1 nm^-2 pull_nstxout = 1000 ; every 2 ps pull_nstfout = 1000 ; every 2 ps ; Velocity generation gen_vel = no ; assign velocities from Maxwell distribution refcoord_scaling = com Would you please let me know your suggestions ? Sincerely, Shima - Original Message - From: Justin Lemkul jalem...@vt.edu To: Shima Arasteh shima_arasteh2...@yahoo.com Cc: Sent: Wednesday, July 24, 2013 9:41 PM Subject: Re: [gmx-users] Umbrella Sampling _ pulled ion On 7/24/13 11:53 AM, Shima Arasteh wrote: Yes, Thanks. Would you give me a hint on this fact that how I would be sure that I am running a correct US ? with proper settings? Either the restraint distance is maintained with acceptable sampling around that distance (given the force constant) or it's not. To save time, I' d prefer to run the US.mdp just for one window. Do you agree with me that if I run an incorrect US for any of the windows, I would get an odd result for npt_US or md_US? It should be pretty easy to see. -Justin Many many thanks for your time and suggestions. Sincerely, Shima - Original Message - From: Justin Lemkul jalem...@vt.edu To: Shima Arasteh shima_arasteh2...@yahoo.com; Discussion list for GROMACS
Re: Fw: [gmx-users] Umbrella Sampling _ pulled ion
On 7/26/13 11:52 AM, Shima Arasteh wrote: Thanks for your previous replies. All are always beneficial and I appreciate you. As I see in pullx.xvg file, the third column refers to the distance of pull group and reference group. All around -1.2. Is it sufficient to judge that the US has been performed properly? I have no idea. Is that value correct based on your desired approach and anticipated outcomes? -Justin -- == Justin A. Lemkul, Ph.D. Postdoctoral Fellow Department of Pharmaceutical Sciences School of Pharmacy Health Sciences Facility II, Room 601 University of Maryland, Baltimore 20 Penn St. Baltimore, MD 21201 jalem...@outerbanks.umaryland.edu | (410) 706-7441 == -- 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: Fw: [gmx-users] Umbrella Sampling _ pulled ion
Please remember to keep the discussion on the list. On 7/26/13 12:04 PM, Shima Arasteh wrote: what I aim to see is that the pulled ion should be remained in center of the channel (peptide inserted in lipid bilayer). I didn't restrain the channel, it makes a little bit hard for me to judge exactly. Is position restraint of channels essential here? The data in pullx.xvg indicate that the distance between the reference group and the pulled group was maintained at -1.2 nm. Visual inspection of the trajectory will tell you whether or not that is the position you desire. Having noted the output of grompp (you did that, right?) you can match up the input specifications and target restrained distance and the actual outcome. Unless something has gone completely haywire, the outcome is very predictable. Learning to assess the validity of what you have done through careful record-keeping and cross-referencing with input files will prevent you from having to wait around for people on the Internet to check your work ;) -Justin -- == Justin A. Lemkul, Ph.D. Postdoctoral Fellow Department of Pharmaceutical Sciences School of Pharmacy Health Sciences Facility II, Room 601 University of Maryland, Baltimore 20 Penn St. Baltimore, MD 21201 jalem...@outerbanks.umaryland.edu | (410) 706-7441 == -- 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] Calculate interaction energy dynamically
Date: Fri, 26 Jul 2013 11:30:14 -0400 From: jalem...@vt.edu To: gmx-users@gromacs.org Subject: Re: [gmx-users] Calculate interaction energy dynamically On 7/26/13 5:21 AM, Davit Hakobyan wrote: From: davhak2...@hotmail.com To: gmx-users@gromacs.org Subject: RE: [gmx-users] Calculate interaction energy dynamically Date: Thu, 25 Jul 2013 09:13:14 + Date: Wed, 24 Jul 2013 09:35:26 -0400 From: jalem...@vt.edu To: gmx-users@gromacs.org Subject: Re: [gmx-users] Calculate interaction energy dynamically On 7/24/13 9:25 AM, Davit Hakobyan wrote: Thank you again for your time and help. Performing rerun on the original system passes without warnings: The subsequent g_energy also runs smoothly. One gets no warning also if passing the ordered *ordered.tpr and the new index (*ordered.ndx) files to the mdrun -rerun of the original system. This indeed leaves with the ordered trajectory file as the only problematic point which might cause the warnings. Good to know. Interestingly, the ordered *.trr file (product of trjorder command) is nearly twice as small as the original *.trr file. Since the coordinates in the ordered *.trr file look fine then probably it is the velocities which might be missing in it. Can the There's no might about it - your gmxcheck output from before shows that velocities are definitely absent, which should account for the difference in file size. warnings come from this issue ? In this case the g_energy should not be able to calculate the kinetic energy. If this is the case may one still rely on the numbers obtained by g_energy for the coulomb and vdw terms (despite the warning messages) ? This is a very weird issue. The .trr file clearly has step and time information (from your earlier post), but the resulting .edr file from mdrun -rerun does not have the correct step number (though it appears to have time). Your problem does not derive from anything to do with velocities, as far as I can tell. Proceed carefully and scrutinize the output to see if it indeed makes sense - you can verify energy values from snapshots that have not been ordered by creating index files that correspond to the same lipids, for instance. This is a very non-standard exercise, so it's probably never been tested and debugged. The energy calculation between the nearest neighbors at a given time frame based on the original trajectory file (with help of new index file) indeed gives the same values as for the ordered trajectory file. So one can probably rely on the energy output from the ordered trajectory despite the warnings. In case of more information about the warnings I will write in this thread. It seems that the problem appears when for the trjorder command one indicates *.trr instead of *.xtc as an output. When *.xtc is mentioned like trjorder ... -o output.xtc both the mdrun -rerun and g_energy subsequently pass without warnings. Interesting. That leads me to suspect that there is something wrong with the way trjorder writes .trr frames. It writes .xtc files correctly, but incorrect step information in .trr files. Thanks for following up. I'll look into it and probably file a bug report on it. What Gromacs version are you using? Maybe you said that before, but I've forgotten. I use Gromacs 4.5.1. Thank you very much for your help and time. -Justin -- == Justin A. Lemkul, Ph.D. Postdoctoral Fellow Department of Pharmaceutical Sciences School of Pharmacy Health Sciences Facility II, Room 601 University of Maryland, Baltimore 20 Penn St. Baltimore, MD 21201 jalem...@outerbanks.umaryland.edu | (410) 706-7441 == -- 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] Re: Limitations of simulations?
Alright, I think I have enough information now to play around with some settings! Thank you very much Justin, you have been very helpful :) -Jonathan -- View this message in context: http://gromacs.5086.x6.nabble.com/Limitations-of-simulations-tp5010100p5010147.html Sent from the GROMACS Users Forum mailing list archive at Nabble.com. -- gmx-users mailing listgmx-users@gromacs.org http://lists.gromacs.org/mailman/listinfo/gmx-users * Please search the archive at http://www.gromacs.org/Support/Mailing_Lists/Search before posting! * Please don't post (un)subscribe requests to the list. Use the www interface or send it to gmx-users-requ...@gromacs.org. * Can't post? Read http://www.gromacs.org/Support/Mailing_Lists
[gmx-users] Free Energy Calculations in Parallel
Dear Gromacs Users, Is it possible to run free energy calculations in parallel using mpirun? If not I am wondering what the fastest way to run a free energy calculation is. I am basing my free energy calculation on the tutorial by Justin Lemkul found at http://www.bevanlab.biochem.vt.edu/Pages/Personal/justin/gmx-tutorials/free_energy/01_theory.html I am trying to use the Bennet Acceptance Ration method (g_bar) to get the free energy of solvation for an ionic liquid. I have tried decoupling an ion pair from the system as well as individual ions. These simulations will run locally but when I try to run them in parallel I get the error: Fatal error: There is no domain decomposition for 8 nodes that is compatible with the given box and a minimum cell size of 2.26125 nm Thank you Quintin Sheridan -- 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] Free Energy Simulations in Parallel
Dear Gromacs Users, Is it possible to run free energy calculations in parallel using mpirun? If not, what is the fastest way to run free energy calculations. I am trying to us the Bennet's Accepetance Ratio (g_bar) to get the free energy of solvation for an ionic liquid based on the tutorial by Justin Lemkul. I hav tried to decouple an ion pair as well as individual ions. In either case the simulations run locally but when I try to run them in parrallel I get the error: Fatal error: There is no domain decomposition for 8 nodes that is compatible with the given box and a minimum cell size of 2.26125 nm Thank You Quintin Sheridan -- 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] Calculate interaction energy dynamically
On 7/26/13 12:37 PM, Davit Hakobyan wrote: Date: Fri, 26 Jul 2013 11:30:14 -0400 From: jalem...@vt.edu To: gmx-users@gromacs.org Subject: Re: [gmx-users] Calculate interaction energy dynamically On 7/26/13 5:21 AM, Davit Hakobyan wrote: From: davhak2...@hotmail.com To: gmx-users@gromacs.org Subject: RE: [gmx-users] Calculate interaction energy dynamically Date: Thu, 25 Jul 2013 09:13:14 + Date: Wed, 24 Jul 2013 09:35:26 -0400 From: jalem...@vt.edu To: gmx-users@gromacs.org Subject: Re: [gmx-users] Calculate interaction energy dynamically On 7/24/13 9:25 AM, Davit Hakobyan wrote: Thank you again for your time and help. Performing rerun on the original system passes without warnings: The subsequent g_energy also runs smoothly. One gets no warning also if passing the ordered *ordered.tpr and the new index (*ordered.ndx) files to the mdrun -rerun of the original system. This indeed leaves with the ordered trajectory file as the only problematic point which might cause the warnings. Good to know. Interestingly, the ordered *.trr file (product of trjorder command) is nearly twice as small as the original *.trr file. Since the coordinates in the ordered *.trr file look fine then probably it is the velocities which might be missing in it. Can the There's no might about it - your gmxcheck output from before shows that velocities are definitely absent, which should account for the difference in file size. warnings come from this issue ? In this case the g_energy should not be able to calculate the kinetic energy. If this is the case may one still rely on the numbers obtained by g_energy for the coulomb and vdw terms (despite the warning messages) ? This is a very weird issue. The .trr file clearly has step and time information (from your earlier post), but the resulting .edr file from mdrun -rerun does not have the correct step number (though it appears to have time). Your problem does not derive from anything to do with velocities, as far as I can tell. Proceed carefully and scrutinize the output to see if it indeed makes sense - you can verify energy values from snapshots that have not been ordered by creating index files that correspond to the same lipids, for instance. This is a very non-standard exercise, so it's probably never been tested and debugged. The energy calculation between the nearest neighbors at a given time frame based on the original trajectory file (with help of new index file) indeed gives the same values as for the ordered trajectory file. So one can probably rely on the energy output from the ordered trajectory despite the warnings. In case of more information about the warnings I will write in this thread. It seems that the problem appears when for the trjorder command one indicates *.trr instead of *.xtc as an output. When *.xtc is mentioned like trjorder ... -o output.xtc both the mdrun -rerun and g_energy subsequently pass without warnings. Interesting. That leads me to suspect that there is something wrong with the way trjorder writes .trr frames. It writes .xtc files correctly, but incorrect step information in .trr files. Thanks for following up. I'll look into it and probably file a bug report on it. What Gromacs version are you using? Maybe you said that before, but I've forgotten. I use Gromacs 4.5.1. If you happen to get a chance to run trjorder and mdrun -rerun with version 4.6.3, that would be very useful information. The 4.5.x series is no longer being actively developed, so debugging that version is not very useful. -Justin -- == Justin A. Lemkul, Ph.D. Postdoctoral Fellow Department of Pharmaceutical Sciences School of Pharmacy Health Sciences Facility II, Room 601 University of Maryland, Baltimore 20 Penn St. Baltimore, MD 21201 jalem...@outerbanks.umaryland.edu | (410) 706-7441 == -- 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] Free Energy Calculations in Parallel
On 7/26/13 4:06 PM, Quintin Sheridan wrote: Dear Gromacs Users, Is it possible to run free energy calculations in parallel using mpirun? If Yes, the free energy code is parallelized, either via MPI or threads. not I am wondering what the fastest way to run a free energy calculation is. I am basing my free energy calculation on the tutorial by Justin Lemkul found at http://www.bevanlab.biochem.vt.edu/Pages/Personal/justin/gmx-tutorials/free_energy/01_theory.html I am trying to use the Bennet Acceptance Ration method (g_bar) to get the free energy of solvation for an ionic liquid. I have tried decoupling an ion pair from the system as well as individual ions. These simulations will run locally but when I try to run them in parallel I get the error: Fatal error: There is no domain decomposition for 8 nodes that is compatible with the given box and a minimum cell size of 2.26125 nm Consult 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. This topic is discussed frequently, particularly in the context of free energy simulations, because there are unique factors that come into play in terms of the DD algorithm. Search the archive and you will surely find the detailed explanations. -Justin -- == Justin A. Lemkul, Ph.D. Postdoctoral Fellow Department of Pharmaceutical Sciences School of Pharmacy Health Sciences Facility II, Room 601 University of Maryland, Baltimore 20 Penn St. Baltimore, MD 21201 jalem...@outerbanks.umaryland.edu | (410) 706-7441 == -- 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: Aw: [gmx-users] unable to equilibrate pressure in npt
Thanks for the explanation and the suggestions. The values I mentioned in the last mail were obtained when using tau_p = 2. I have tried changing tau_p in the mdp file to 0.8 and 1.0. With tau_p = 1.0 the average pressure is around 0.9 but with tau_p = 0.8 the pressure is ~ -0.3 (time = 10ns). Is there something else I might be doing wrong? To add more information, I am using the same mdp file as given in the lysozyme tutorial except using amber03 ILDN and changing the barostat as I mentioned. Should I try increasing the coupling time to more than 2 ? Any suggestions would be highly appreciated. Amin I did such things in the past, and had similar issues.nbsp; First of, somone may have beter suggestions.nbsp; The coupling time will bring it down a bit in the .mdp file (0.2 Vs 2 picoseconds.nbsp; However, I have found that you will still see large fluctuations around a mean once equilibrated which can vary by 10 to 60, 70 Atomospheres depending on the pressure coupling algorythim used.nbsp; If it oscillates around this mean, without large up or down changes I assume it is eq#39;d.nbsp; Basically, if you look at it from the perspective of what pressure itself entails, it is the direct kenetic force of small molecules bumping into each other, so a simple protein movmeent in a small unit cell is turned by the computer into a molar ration, which looks crazy at the macroscopic level, but in reality you do not usually have the unit cell in such a deminsion.nbsp; This is of course unless your system is made of nothing other than small molecules, which take off a 0 to the oscillation range.nbsp; the pressure algorythms just keep it from going to insanly large limits, such as an ice cube because of 200 kcal/mol - enthalpy or boiling in the reverse, by pretending there is a universal pressure applied to the unit cell...so a generalized correction for the fact that the unit cell is not really representing infinaty as it is treated by the computer...based on real macroscopic pressures.nbsp; I may be wrong , but this is my assumption./div divnbsp;/div divStephan Watkins/div div name=quote style=margin:10px 5px 5px 10px; padding: 10px 0 10px 10px; border-left:2px solid #C3D9E5; word-wrap: break-word; -webkit-nbsp-mode: space; -webkit-line-break: after-white-space; div style=margin:0 0 10px 0;bGesendet:/bnbsp;Donnerstag, 25. Juli 2013 um 18:22 Uhrbr/ bVon:/bnbsp;a...@imtech.res.inbr/ bAn:/bnbsp;gmx-users@gromacs.orgbr/ bBetreff:/bnbsp;[gmx-users] unable to equilibrate pressure in npt/div div name=quoted-contentDear gromacs users,br/ I know similar issues have been raised many times on the list but I am unable tobr/ solve the problem so I am seeking your advice. I am trying to simulate a proteinbr/ in a dodecahedron box. The system size is ~30k atoms. I have followed thebr/ methodology given in the lysozyme tutorial i.e. minimization, nvt, npt andbr/ production. However the average pressure values after 5, 10, 20 ns of nptbr/ equilibration are not coming close to the reference pressure i.e 1 bar whenbr/ using parrinello rahman. I checked the mailing list and tried using Berendsonbr/ and after 10 ns I get average pressure close to 1 (0.85). However when I switchbr/ to parrinello rahman for production run the average pressure again goes far frombr/ 1. Can someone please help me with this? Here are the g_energy outputs frombr/ equilibration (parrinello rahman, 20ns) and production run (10ns afterbr/ Berendson).br/ br/ Pressure 0.0300644 0.38 134.243 1.57106 (bar)br/ br/ Pressure -0.0405509 0.79 134.204 0.151638 (bar)br/ br/ Can someone please help me with this?br/ br/ Regards.br/ Amin.br/ br/ __br/ #2360;#2370;#2325;#2381;#2359;#2381;#2350;#2332;#2368;#2357; #2346;#2381;#2352;#2380;#2342;#2381;#2351;#2379;#2327;#2367;#2325;#2368; #2360;#2306;#2360;#2381;#2341;#2366;#2344; (#2357;#2376;#2332;#2381;#2334;#2366;#2344;#2367;#2325; #2324;#2342;#2381;#2351;#2379;#2327;#2367;#2325; #2309;#2344;#2369;#2360;#2306;#2343;#2366;#2344; #2346;#2352;#2367;#2359;#2342;)br/ Institute of Microbial Technology (A CONSTITUENT ESTABLISHMENT OF CSIR)br/ #2360;#2376;#2325;#2381;#2335;#2352; 39 #2319;, #2330;#2339;#2381;#2337;#2368;#2327;#2338;#2364; / Sector 39-A, Chandigarhbr/ #2346;#2367;#2344; #2325;#2379;#2337;/PIN CODE :160036br/ #2342;#2370;#2352;#2349;#2366;#2359;/EPABX :0172 6665 201-202br/ --br/ gmx-users mailing list gmx-users@gromacs.orgbr/ a href=http://lists.gromacs.org/mailman/listinfo/gmx-users; target=_blankhttp://lists.gromacs.org/mailman/listinfo/gmx-users/abr/ * Please search the archive at a href=http://www.gromacs.org/Support/Mailing_Lists/Search; target=_blankhttp://www.gromacs.org/Support/Mailing_Lists/Search/a before posting!br/ * Please don#39;t post (un)subscribe requests to the list. Use thebr/ www interface or send it to
